source: palm/trunk/SOURCE/surface_data_output_mod.f90 @ 4502

Last change on this file since 4502 was 4502, checked in by schwenkel, 14 months ago

Implementation of ice microphysics

  • Property svn:keywords set to Id
File size: 266.7 KB
Line 
1!> @file surface_data_output_mod.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2020 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: surface_data_output_mod.f90 4502 2020-04-17 16:14:16Z schwenkel $
27! Implementation of ice microphysics
28!
29! 4500 2020-04-17 10:12:45Z suehring
30! - Correct output of ground/wall heat flux at USM surfaces
31! - Add conversion factor to heat and momentum-flux output
32!
33! 4495 2020-04-13 20:11:20Z raasch
34! restart data handling with MPI-IO added
35!
36! 4444 2020-03-05 15:59:50Z raasch
37! bugfix: cpp-directives for serial mode added
38!
39! 4360 2020-01-07 11:25:50Z suehring
40! Fix wrongly declared nc_stat variable in surface_data_output_mod
41!
42! 4205 2019-08-30 13:25:00Z suehring
43! - Correct x,y-coordinates of vertical surfaces in netcdf output
44! - Change definition of azimuth angle, reference is north 0 degree
45! - zenith angle is always defined, also for vertical surfaces where it is
46!   90 degree, while azimuth angle is only defined for vertical surfaces, not
47!   for horizontal ones
48!
49! 4182 2019-08-22 15:20:23Z scharf
50! Corrected "Former revisions" section
51!
52! 4129 2019-07-31 12:56:07Z gronemeier
53! - bugfix: corrected loop over horizontal default surfaces
54! - change default setting of to_vtk and to_netcdf
55!
56! 4029 2019-06-14 14:04:35Z raasch
57! netcdf variable NF90_NOFILL is used as argument instead of "1" in call to NF90_DEF_VAR_FILL
58!
59! 3881 2019-04-10 09:31:22Z suehring
60! Check for zero output timestep (not allowed in parallel NetCDF output mode)
61!
62! 3817 2019-03-26 13:53:57Z suehring
63! Correct output coordinates of vertical surface elements
64!
65! 3766 2019-02-26 16:23:41Z raasch
66! bugfix in surface_data_output_rrd_local (variable k removed)
67!
68! 3762 2019-02-25 16:54:16Z suehring
69! Remove unused variables and add preprocessor directives for variables that
70! are used only when netcdf4 is defined
71!
72! 3745 2019-02-15 18:57:56Z suehring
73! Output of waste_heat and innermost wall flux from indoor model
74!
75! 3744 2019-02-15 18:38:58Z suehring
76! Add azimuth and zenith to output file; set long-name attributes;
77! clean-up coding layout
78!
79! 3735 2019-02-12 09:52:40Z suehring
80! - Split initialization into initialization of arrays and further initialization
81!   in order to enable reading of restart data.
82! - Consider restarts in surface data averaging.
83! - Correct error message numbers
84!
85! 3731 2019-02-11 13:06:27Z suehring
86! Bugfix: add cpp options
87!
88! 3727 2019-02-08 14:52:10Z gronemeier
89! Enable NetCDF output for surface data (suehring, gronemeier)
90!
91! 3691 2019-01-23 09:57:04Z suehring
92! Add output of surface-parallel flow speed
93!
94! 3648 2019-01-02 16:35:46Z suehring
95! Rename module and subroutines
96! 3420 2018-10-24 17:30:08Z gronemeier
97! Initial implementation from Klaus Ketelsen and Matthias Suehring
98!
99!
100! Authors:
101! --------
102! @author Klaus Ketelsen, Matthias Suehring, Tobias Gronemeier
103!
104! Description:
105! ------------
106!> Generate output for surface data.
107!>
108!> @todo Create namelist file for post-processing tool.
109!------------------------------------------------------------------------------!
110
111MODULE surface_data_output_mod
112
113   USE kinds
114
115   USE arrays_3d,                                                                                  &
116       ONLY:  heatflux_output_conversion,                                                          &
117              momentumflux_output_conversion,                                                      &
118              waterflux_output_conversion,                                                         &
119              zu,                                                                                  &
120              zw
121
122   USE control_parameters,                                                                         &
123       ONLY:  coupling_char, data_output_during_spinup, end_time,                                  &
124              message_string,                                                                      &
125              restart_data_format_output, run_description_header, simulated_time_at_begin,         &
126              spinup_time, surface_output
127
128   USE grid_variables,                                                         &
129       ONLY: dx,dy
130
131   USE indices,                                                                &
132       ONLY: nxl, nxr, nys, nyn, nzb, nzt
133
134#if defined( __netcdf )
135   USE NETCDF
136#endif
137
138   USE netcdf_data_input_mod,                                                  &
139       ONLY:  init_model
140
141   USE netcdf_interface,                                                       &
142       ONLY:  nc_stat, netcdf_create_att, netcdf_create_dim,                   &
143              netcdf_create_file, netcdf_create_global_atts,                   &
144              netcdf_create_var, netcdf_data_format, netcdf_handle_error
145
146   USE pegrid
147
148   USE restart_data_mpi_io_mod,                                                                    &
149       ONLY:  rrd_mpi_io, wrd_mpi_io
150
151   USE surface_mod,                                                                                &
152       ONLY:  ind_pav_green,                                                                       &
153              ind_veg_wall,                                                                        &
154              ind_wat_win,                                                                         &
155              surf_def_h,                                                                          &
156              surf_def_v,                                                                          &
157              surf_lsm_h,                                                                          &
158              surf_lsm_v,                                                                          &
159              surf_usm_h,                                                                          &
160              surf_usm_v
161
162   IMPLICIT NONE
163
164   TYPE surf_out                      !< data structure which contains all surfaces elements of all types on subdomain
165
166      INTEGER(iwp) ::  ns             !< number of surface elements on subdomain
167      INTEGER(iwp) ::  ns_total       !< total number of surface elements
168      INTEGER(iwp) ::  npoints        !< number of points / vertices which define a surface element (on subdomain)
169      INTEGER(iwp) ::  npoints_total  !< total number of points / vertices which define a surface element
170
171      INTEGER(iwp), DIMENSION(:), ALLOCATABLE   ::  s        !< coordinate for NetCDF output, number of the surface element
172
173      REAL(wp) ::  fillvalue = -9999.0_wp !< fillvalue for surface elements which are not defined
174
175      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  azimuth  !< azimuth orientation coordinate for NetCDF output
176      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  es_utm   !< E-UTM coordinate for NetCDF output
177      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  ns_utm   !< E-UTM coordinate for NetCDF output
178      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  xs       !< x-coordinate for NetCDF output
179      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  ys       !< y-coordinate for NetCDF output
180      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  zs       !< z-coordinate for NetCDF output
181      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  zenith   !< zenith orientation coordinate for NetCDF output
182      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  var_out  !< output variables
183      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var_av   !< variables used for averaging
184      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  points   !< points  / vertices of a surface element
185      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  polygons !< polygon data of a surface element
186   END TYPE surf_out
187
188   CHARACTER(LEN=100), DIMENSION(300)       ::  data_output_surf = ' '  !< namelist variable which describes the output variables
189   CHARACTER(LEN=100), DIMENSION(0:1,300)   ::  dosurf = ' '            !< internal variable which describes the output variables
190                                                                        !<  and separates averaged from non-averaged output
191   CHARACTER(LEN=100), DIMENSION(0:1,300)   ::  dosurf_unit = ' '       !< internal variable which holds the unit of the given output variable
192
193   INTEGER(iwp) ::  average_count_surf = 0   !< number of ensemble members used for averaging
194   INTEGER(iwp) ::  dosurf_no(0:1) = 0       !< number of surface output quantities
195#if defined( __netcdf4_parallel )
196   INTEGER(iwp) ::  oldmode                  !< save old set-fill-mode of netcdf file (not needed, but required for routine call)
197
198   INTEGER(iwp), DIMENSION(0:1) ::  dosurf_time_count = 0 !< count of output time steps
199   INTEGER(iwp), DIMENSION(0:1) ::  id_dim_s_surf         !< netcdf ID for dimension s
200   INTEGER(iwp), DIMENSION(0:1) ::  id_dim_time_surf      !< netcdf ID for dimension time
201   INTEGER(iwp), DIMENSION(0:1) ::  id_set_surf           !< netcdf ID for file
202   INTEGER(iwp), DIMENSION(0:1) ::  id_var_azimuth_surf   !< netcdf ID for variable azimuth
203   INTEGER(iwp), DIMENSION(0:1) ::  id_var_etum_surf      !< netcdf ID for variable Es_UTM
204   INTEGER(iwp), DIMENSION(0:1) ::  id_var_nutm_surf      !< netcdf ID for variable Ns_UTM
205   INTEGER(iwp), DIMENSION(0:1) ::  id_var_time_surf      !< netcdf ID for variable time
206   INTEGER(iwp), DIMENSION(0:1) ::  id_var_s_surf         !< netcdf ID for variable s
207   INTEGER(iwp), DIMENSION(0:1) ::  id_var_xs_surf        !< netcdf ID for variable xs
208   INTEGER(iwp), DIMENSION(0:1) ::  id_var_ys_surf        !< netcdf ID for variable ys
209   INTEGER(iwp), DIMENSION(0:1) ::  id_var_zenith_surf    !< netcdf ID for variable zenith
210   INTEGER(iwp), DIMENSION(0:1) ::  id_var_zs_surf        !< netcdf ID for variable zs
211   INTEGER(iwp), DIMENSION(0:1) ::  ntdim_surf            !< number of output time steps
212
213   INTEGER(iwp), DIMENSION(0:1,300) ::  id_var_dosurf     !< netcdf ID for output variables
214#endif
215
216   LOGICAL :: first_output(0:1) = .FALSE.                 !< true if first output was already called
217   LOGICAL :: to_netcdf = .FALSE.                         !< flag indicating parallel NetCDF output
218   LOGICAL :: to_vtk = .FALSE.                            !< flag indicating binary surface-data output that can be further
219                                                          !< processed to VTK format
220
221   REAL(wp) ::  averaging_interval_surf  = 9999999.9_wp   !< averaging interval
222   REAL(wp) ::  dt_dosurf = 9999999.9_wp                  !< time interval for instantaneous data output
223   REAL(wp) ::  dt_dosurf_av = 9999999.9_wp               !< time interval for averaged data output
224   REAL(wp) ::  skip_time_dosurf = 0.0_wp                 !< skip time for instantaneous data output
225   REAL(wp) ::  skip_time_dosurf_av = 0.0_wp              !< skip time for averaged data output
226   REAL(wp) ::  time_dosurf = 0.0_wp                      !< internal counter variable to check for instantaneous data output
227   REAL(wp) ::  time_dosurf_av = 0.0_wp                   !< internal counter variable to check for averaged data output
228
229   TYPE(surf_out) ::  surfaces      !< variable which contains all required output information
230
231   SAVE
232
233   PRIVATE
234
235   INTERFACE  surface_data_output
236      MODULE PROCEDURE surface_data_output
237   END INTERFACE  surface_data_output
238
239   INTERFACE  surface_data_output_averaging
240      MODULE PROCEDURE surface_data_output_averaging
241   END INTERFACE  surface_data_output_averaging
242
243   INTERFACE  surface_data_output_check_parameters
244      MODULE PROCEDURE surface_data_output_check_parameters
245   END INTERFACE  surface_data_output_check_parameters
246
247   INTERFACE  surface_data_output_init
248      MODULE PROCEDURE surface_data_output_init
249   END INTERFACE  surface_data_output_init
250
251   INTERFACE  surface_data_output_init_arrays
252      MODULE PROCEDURE surface_data_output_init_arrays
253   END INTERFACE  surface_data_output_init_arrays
254
255   INTERFACE  surface_data_output_last_action
256      MODULE PROCEDURE surface_data_output_last_action
257   END INTERFACE  surface_data_output_last_action
258
259   INTERFACE  surface_data_output_parin
260      MODULE PROCEDURE surface_data_output_parin
261   END INTERFACE  surface_data_output_parin
262
263   INTERFACE  surface_data_output_rrd_global
264      MODULE PROCEDURE surface_data_output_rrd_global_ftn
265      MODULE PROCEDURE surface_data_output_rrd_global_mpi
266   END INTERFACE  surface_data_output_rrd_global
267
268   INTERFACE  surface_data_output_rrd_local
269      MODULE PROCEDURE surface_data_output_rrd_local
270   END INTERFACE  surface_data_output_rrd_local
271
272   INTERFACE  surface_data_output_wrd_global
273      MODULE PROCEDURE surface_data_output_wrd_global
274   END INTERFACE  surface_data_output_wrd_global
275
276   INTERFACE  surface_data_output_wrd_local
277      MODULE PROCEDURE surface_data_output_wrd_local
278   END INTERFACE  surface_data_output_wrd_local
279
280!
281!--Public subroutines
282   PUBLIC surface_data_output, surface_data_output_averaging,                  &
283          surface_data_output_check_parameters, surface_data_output_init,      &
284          surface_data_output_init_arrays, surface_data_output_last_action,    &
285          surface_data_output_parin, surface_data_output_rrd_global,           &
286          surface_data_output_rrd_local, surface_data_output_wrd_local,        &
287          surface_data_output_wrd_global
288!
289!--Public variables
290   PUBLIC average_count_surf, averaging_interval_surf, dt_dosurf, dt_dosurf_av,&
291          skip_time_dosurf, skip_time_dosurf_av, time_dosurf, time_dosurf_av
292
293 CONTAINS
294
295!------------------------------------------------------------------------------!
296! Description:
297! ------------
298!> This routine counts the number of surfaces on each core and allocates
299!> arrays.
300!------------------------------------------------------------------------------!
301   SUBROUTINE surface_data_output_init_arrays
302
303      IMPLICIT NONE
304
305!
306!--   Determine the number of surface elements on subdomain
307      surfaces%ns = surf_def_h(0)%ns + surf_lsm_h%ns + surf_usm_h%ns           & !horizontal upward-facing
308                  + surf_def_h(1)%ns                                           & !horizontal downard-facing
309                  + surf_def_v(0)%ns + surf_lsm_v(0)%ns + surf_usm_v(0)%ns     & !northward-facing
310                  + surf_def_v(1)%ns + surf_lsm_v(1)%ns + surf_usm_v(1)%ns     & !southward-facing
311                  + surf_def_v(2)%ns + surf_lsm_v(2)%ns + surf_usm_v(2)%ns     & !westward-facing
312                  + surf_def_v(3)%ns + surf_lsm_v(3)%ns + surf_usm_v(3)%ns       !eastward-facing
313!
314!--    Determine the total number of surfaces in the model domain
315#if defined( __parallel )
316       CALL MPI_ALLREDUCE( surfaces%ns, surfaces%ns_total, 1,                  &
317                           MPI_INTEGER, MPI_SUM, comm2d, ierr )
318#else
319       surfaces%ns_total = surfaces%ns
320#endif
321!
322!--   Allocate output variable and set to _FillValue attribute
323      ALLOCATE ( surfaces%var_out(1:surfaces%ns) )
324      surfaces%var_out = surfaces%fillvalue
325!
326!--   If there is an output of time average output variables, allocate the
327!--   required array.
328      IF ( dosurf_no(1) > 0 )  THEN
329         ALLOCATE ( surfaces%var_av(1:surfaces%ns,1:dosurf_no(1)) )
330         surfaces%var_av = 0.0_wp
331      ENDIF
332
333   END SUBROUTINE surface_data_output_init_arrays
334
335
336!------------------------------------------------------------------------------!
337! Description:
338! ------------
339!> Initialization surface-data output data structure: calculation of vertices
340!> and polygon data for the surface elements, allocation of required arrays.
341!------------------------------------------------------------------------------!
342   SUBROUTINE surface_data_output_init
343
344      IMPLICIT NONE
345
346#if defined( __netcdf4_parallel )
347      CHARACTER (LEN=100)  :: filename            !< name of output file
348      CHARACTER (LEN=80)   ::  time_average_text  !< string written to file attribute time_avg
349      CHARACTER (LEN=4000) ::  var_list           !< list of variables written to NetCDF file
350
351      INTEGER(iwp) ::  av                !< flag for averaged (=1) and non-averaged (=0) data
352#endif
353      INTEGER(iwp) ::  i                 !< grid index in x-direction, also running variable for counting non-average data output
354      INTEGER(iwp) ::  j                 !< grid index in y-direction, also running variable for counting average data output
355      INTEGER(iwp) ::  k                 !< grid index in z-direction
356      INTEGER(iwp) ::  l                 !< running index for surface-element orientation
357      INTEGER(iwp) ::  m                 !< running index for surface elements
358      INTEGER(iwp) ::  mm                !< local counting variable for surface elements
359      INTEGER(iwp) ::  npg               !< counter variable for all surface elements ( or polygons )
360      INTEGER(iwp) ::  point_index_count !< local counter variable for point index
361      INTEGER(iwp) ::  start_count       !< local start counter for the surface index
362
363      INTEGER(iwp), DIMENSION(0:numprocs-1) :: num_points_on_pe   !< array which contains the number of points on all mpi ranks
364      INTEGER(iwp), DIMENSION(0:numprocs-1) :: num_surfaces_on_pe !< array which contains the number of surfaces on all mpi ranks
365      INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:,:) ::  point_index !< dummy array used to check where the reference points for surface polygons are located
366
367      REAL(wp) ::  az    !< azimuth angle, indicated the vertical orientation of a surface element
368      REAL(wp) ::  off_x !< grid offset in x-direction between the stored grid index and the actual wall
369      REAL(wp) ::  off_y !< grid offset in y-direction between the stored grid index and the actual wall
370#if defined( __netcdf4_parallel )
371      REAL(wp), DIMENSION(:), ALLOCATABLE ::  netcdf_data_1d  !< dummy array to output 1D data into netcdf file
372#endif
373
374!
375!--   If output to VTK format is enabled, initialize point and polygon data.
376!--   In a first step, count the number of points which are defining
377!--   the surfaces and the polygons.
378      IF ( to_vtk )  THEN
379         ALLOCATE( point_index(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
380         point_index = -1
381!
382!--      Horizontal default surfaces
383         surfaces%npoints = 0
384         DO  l = 0, 1
385            DO  m = 1, surf_def_h(0)%ns
386!
387!--            Determine the indices of the respective grid cell inside the topography
388               i = surf_def_h(0)%i(m) + surf_def_h(0)%ioff
389               j = surf_def_h(0)%j(m) + surf_def_h(0)%joff
390               k = surf_def_h(0)%k(m) + surf_def_h(0)%koff
391!
392!--            Check if the vertices that define the surface element are already
393!--            defined, if not, increment the counter.
394               IF ( point_index(k,j,i) < 0 )  THEN
395                  surfaces%npoints   = surfaces%npoints + 1
396                  point_index(k,j,i) = surfaces%npoints - 1
397               ENDIF
398               IF ( point_index(k,j,i+1) < 0 )  THEN
399                  surfaces%npoints     = surfaces%npoints + 1
400                  point_index(k,j,i+1) = surfaces%npoints - 1
401               ENDIF
402               IF ( point_index(k,j+1,i+1) < 0 )  THEN
403                  surfaces%npoints       = surfaces%npoints + 1
404                  point_index(k,j+1,i+1) = surfaces%npoints - 1
405               ENDIF
406               IF ( point_index(k,j+1,i) < 0 )  THEN
407                  surfaces%npoints     = surfaces%npoints + 1
408                  point_index(k,j+1,i) = surfaces%npoints - 1
409               ENDIF
410            ENDDO
411         ENDDO
412         DO  m = 1, surf_lsm_h%ns
413            i = surf_lsm_h%i(m) + surf_lsm_h%ioff
414            j = surf_lsm_h%j(m) + surf_lsm_h%joff
415            k = surf_lsm_h%k(m) + surf_lsm_h%koff
416
417            IF ( point_index(k,j,i) < 0 )  THEN
418               surfaces%npoints   = surfaces%npoints + 1
419               point_index(k,j,i) = surfaces%npoints - 1
420            ENDIF
421            IF ( point_index(k,j,i+1) < 0 )  THEN
422               surfaces%npoints     = surfaces%npoints + 1
423               point_index(k,j,i+1) = surfaces%npoints - 1
424            ENDIF
425            IF ( point_index(k,j+1,i+1) < 0 )  THEN
426               surfaces%npoints       = surfaces%npoints + 1
427               point_index(k,j+1,i+1) = surfaces%npoints - 1
428            ENDIF
429            IF ( point_index(k,j+1,i) < 0 )  THEN
430               surfaces%npoints     = surfaces%npoints + 1
431               point_index(k,j+1,i) = surfaces%npoints - 1
432            ENDIF
433         ENDDO
434         DO  m = 1, surf_usm_h%ns
435            i = surf_usm_h%i(m) + surf_usm_h%ioff
436            j = surf_usm_h%j(m) + surf_usm_h%joff
437            k = surf_usm_h%k(m) + surf_usm_h%koff
438
439            IF ( point_index(k,j,i) < 0 )  THEN
440               surfaces%npoints   = surfaces%npoints + 1
441               point_index(k,j,i) = surfaces%npoints - 1
442            ENDIF
443            IF ( point_index(k,j,i+1) < 0 )  THEN
444               surfaces%npoints     = surfaces%npoints + 1
445               point_index(k,j,i+1) = surfaces%npoints - 1
446            ENDIF
447            IF ( point_index(k,j+1,i+1) < 0 )  THEN
448               surfaces%npoints       = surfaces%npoints + 1
449               point_index(k,j+1,i+1) = surfaces%npoints - 1
450            ENDIF
451            IF ( point_index(k,j+1,i) < 0 )  THEN
452               surfaces%npoints     = surfaces%npoints + 1
453               point_index(k,j+1,i) = surfaces%npoints - 1
454            ENDIF
455         ENDDO
456!
457!--      Vertical surfaces
458         DO  l = 0, 3
459            DO  m = 1, surf_def_v(l)%ns
460!
461!--            Determine the indices of the respective grid cell inside the
462!--            topography. Please note, j-index for north-facing surfaces
463!--            ( l==0 ) is identical to the reference j-index outside the grid
464!--            box. Equivalent for east-facing surfaces and i-index.
465               i = surf_def_v(l)%i(m) + MERGE( surf_def_v(l)%ioff, 0, l == 3 )
466               j = surf_def_v(l)%j(m) + MERGE( surf_def_v(l)%joff, 0, l == 1 )
467               k = surf_def_v(l)%k(m) + surf_def_v(l)%koff
468!
469!--            Lower left /front vertex
470               IF ( point_index(k,j,i) < 0 )  THEN
471                  surfaces%npoints   = surfaces%npoints + 1
472                  point_index(k,j,i) = surfaces%npoints - 1
473               ENDIF
474!
475!--            Upper / lower right index for north- and south-facing surfaces
476               IF ( l == 0  .OR.  l == 1 )  THEN
477                  IF ( point_index(k,j,i+1) < 0 )  THEN
478                     surfaces%npoints     = surfaces%npoints + 1
479                     point_index(k,j,i+1) = surfaces%npoints - 1
480                  ENDIF
481                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
482                     surfaces%npoints       = surfaces%npoints + 1
483                     point_index(k+1,j,i+1) = surfaces%npoints - 1
484                  ENDIF
485!
486!--            Upper / lower front index for east- and west-facing surfaces
487               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
488                  IF ( point_index(k,j+1,i) < 0 )  THEN
489                     surfaces%npoints     = surfaces%npoints + 1
490                     point_index(k,j+1,i) = surfaces%npoints - 1
491                  ENDIF
492                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
493                     surfaces%npoints       = surfaces%npoints + 1
494                     point_index(k+1,j+1,i) = surfaces%npoints - 1
495                  ENDIF
496               ENDIF
497!
498!--            Upper left / front vertex
499               IF ( point_index(k+1,j,i) < 0 )  THEN
500                  surfaces%npoints     = surfaces%npoints + 1
501                  point_index(k+1,j,i) = surfaces%npoints - 1
502               ENDIF
503            ENDDO
504            DO  m = 1, surf_lsm_v(l)%ns
505               i = surf_lsm_v(l)%i(m) + MERGE( surf_lsm_v(l)%ioff, 0, l == 3 )
506               j = surf_lsm_v(l)%j(m) + MERGE( surf_lsm_v(l)%joff, 0, l == 1 )
507               k = surf_lsm_v(l)%k(m) + surf_lsm_v(l)%koff
508!
509!--            Lower left /front vertex
510               IF ( point_index(k,j,i) < 0 )  THEN
511                  surfaces%npoints   = surfaces%npoints + 1
512                  point_index(k,j,i) = surfaces%npoints - 1
513               ENDIF
514!
515!--            Upper / lower right index for north- and south-facing surfaces
516               IF ( l == 0  .OR.  l == 1 )  THEN
517                  IF ( point_index(k,j,i+1) < 0 )  THEN
518                     surfaces%npoints     = surfaces%npoints + 1
519                     point_index(k,j,i+1) = surfaces%npoints - 1
520                  ENDIF
521                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
522                     surfaces%npoints       = surfaces%npoints + 1
523                     point_index(k+1,j,i+1) = surfaces%npoints - 1
524                  ENDIF
525!
526!--            Upper / lower front index for east- and west-facing surfaces
527               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
528                  IF ( point_index(k,j+1,i) < 0 )  THEN
529                     surfaces%npoints     = surfaces%npoints + 1
530                     point_index(k,j+1,i) = surfaces%npoints - 1
531                  ENDIF
532                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
533                     surfaces%npoints       = surfaces%npoints + 1
534                     point_index(k+1,j+1,i) = surfaces%npoints - 1
535                  ENDIF
536               ENDIF
537!
538!--            Upper left / front vertex
539               IF ( point_index(k+1,j,i) < 0 )  THEN
540                  surfaces%npoints     = surfaces%npoints + 1
541                  point_index(k+1,j,i) = surfaces%npoints - 1
542               ENDIF
543            ENDDO
544
545            DO  m = 1, surf_usm_v(l)%ns
546               i = surf_usm_v(l)%i(m) + MERGE( surf_usm_v(l)%ioff, 0, l == 3 )
547               j = surf_usm_v(l)%j(m) + MERGE( surf_usm_v(l)%joff, 0, l == 1 )
548               k = surf_usm_v(l)%k(m) + surf_usm_v(l)%koff
549!
550!--            Lower left /front vertex
551               IF ( point_index(k,j,i) < 0 )  THEN
552                  surfaces%npoints   = surfaces%npoints + 1
553                  point_index(k,j,i) = surfaces%npoints - 1
554               ENDIF
555!
556!--            Upper / lower right index for north- and south-facing surfaces
557               IF ( l == 0  .OR.  l == 1 )  THEN
558                  IF ( point_index(k,j,i+1) < 0 )  THEN
559                     surfaces%npoints     = surfaces%npoints + 1
560                     point_index(k,j,i+1) = surfaces%npoints - 1
561                  ENDIF
562                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
563                     surfaces%npoints       = surfaces%npoints + 1
564                     point_index(k+1,j,i+1) = surfaces%npoints - 1
565                  ENDIF
566!
567!--            Upper / lower front index for east- and west-facing surfaces
568               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
569                  IF ( point_index(k,j+1,i) < 0 )  THEN
570                     surfaces%npoints     = surfaces%npoints + 1
571                     point_index(k,j+1,i) = surfaces%npoints - 1
572                  ENDIF
573                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
574                     surfaces%npoints       = surfaces%npoints + 1
575                     point_index(k+1,j+1,i) = surfaces%npoints - 1
576                  ENDIF
577               ENDIF
578!
579!--            Upper left / front vertex
580               IF ( point_index(k+1,j,i) < 0 )  THEN
581                  surfaces%npoints     = surfaces%npoints + 1
582                  point_index(k+1,j,i) = surfaces%npoints - 1
583               ENDIF
584            ENDDO
585
586         ENDDO
587!
588!--      Allocate the number of points and polygons. Note, the number of
589!--      polygons is identical to the number of surfaces elements, whereas the
590!--      number of points (vertices), which define the polygons, can be larger.
591         ALLOCATE( surfaces%points(3,1:surfaces%npoints) )
592         ALLOCATE( surfaces%polygons(5,1:surfaces%ns)    )
593!
594!--      Note, PARAVIEW expects consecutively ordered points, in order to
595!--      unambiguously identify surfaces.
596!--      Hence, all PEs should know where they start counting, depending on the
597!--      number of points on the other PE's with lower MPI rank.
598#if defined( __parallel )
599         CALL MPI_ALLGATHER( surfaces%npoints, 1, MPI_INTEGER,                 &
600                             num_points_on_pe, 1, MPI_INTEGER, comm2d, ierr  )
601#else
602         num_points_on_pe = surfaces%npoints
603#endif
604
605!
606!--      After the number of vertices is counted, repeat the loops and define
607!--      the vertices. Start with the horizontal default surfaces.
608!--      First, however, determine the offset where couting of points should be
609!--      started, which is the sum of points of all PE's with lower MPI rank.
610         i                 = 0
611         point_index_count = 0
612         DO WHILE ( i < myid  .AND.  i <= SIZE( num_points_on_pe ) )
613            point_index_count = point_index_count + num_points_on_pe(i)
614            i                 = i + 1
615         ENDDO
616
617         surfaces%npoints = 0
618         point_index      = -1
619         npg              = 0
620
621         DO  l = 0, 1
622            DO  m = 1, surf_def_h(l)%ns
623!
624!--            Determine the indices of the respective grid cell inside the
625!--            topography.
626               i = surf_def_h(l)%i(m) + surf_def_h(l)%ioff
627               j = surf_def_h(l)%j(m) + surf_def_h(l)%joff
628               k = surf_def_h(l)%k(m) + surf_def_h(l)%koff
629!
630!--            Check if the vertices that define the surface element are
631!--            already defined, if not, increment the counter.
632               IF ( point_index(k,j,i) < 0 )  THEN
633                  surfaces%npoints   = surfaces%npoints + 1
634                  point_index(k,j,i) = point_index_count
635                  point_index_count  = point_index_count + 1
636                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
637                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
638                  surfaces%points(3,surfaces%npoints) = zw(k)
639               ENDIF
640               IF ( point_index(k,j,i+1) < 0 )  THEN
641                  surfaces%npoints     = surfaces%npoints + 1
642                  point_index(k,j,i+1) = point_index_count
643                  point_index_count    = point_index_count + 1
644                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
645                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
646                  surfaces%points(3,surfaces%npoints) = zw(k)
647               ENDIF
648               IF ( point_index(k,j+1,i+1) < 0 )  THEN
649                  surfaces%npoints       = surfaces%npoints + 1
650                  point_index(k,j+1,i+1) = point_index_count
651                  point_index_count      = point_index_count + 1
652                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
653                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
654                  surfaces%points(3,surfaces%npoints) = zw(k)
655               ENDIF
656               IF ( point_index(k,j+1,i) < 0 )  THEN
657                  surfaces%npoints     = surfaces%npoints + 1
658                  point_index(k,j+1,i) = point_index_count
659                  point_index_count    = point_index_count + 1
660                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
661                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
662                  surfaces%points(3,surfaces%npoints) = zw(k)
663               ENDIF
664
665               npg                        = npg + 1
666               surfaces%polygons(1,npg)   = 4
667               surfaces%polygons(2,npg)   = point_index(k,j,i)
668               surfaces%polygons(3,npg)   = point_index(k,j,i+1)
669               surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
670               surfaces%polygons(5,npg)   = point_index(k,j+1,i)
671            ENDDO
672         ENDDO
673         DO  m = 1, surf_lsm_h%ns
674            i = surf_lsm_h%i(m) + surf_lsm_h%ioff
675            j = surf_lsm_h%j(m) + surf_lsm_h%joff
676            k = surf_lsm_h%k(m) + surf_lsm_h%koff
677            IF ( point_index(k,j,i) < 0 )  THEN
678               surfaces%npoints   = surfaces%npoints + 1
679               point_index(k,j,i) = point_index_count
680               point_index_count  = point_index_count + 1
681               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
682               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
683               surfaces%points(3,surfaces%npoints) = zw(k)
684            ENDIF
685            IF ( point_index(k,j,i+1) < 0 )  THEN
686               surfaces%npoints     = surfaces%npoints + 1
687               point_index(k,j,i+1) = point_index_count
688               point_index_count    = point_index_count + 1
689               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
690               surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
691               surfaces%points(3,surfaces%npoints) = zw(k)
692            ENDIF
693            IF ( point_index(k,j+1,i+1) < 0 )  THEN
694               surfaces%npoints       = surfaces%npoints + 1
695               point_index(k,j+1,i+1) = point_index_count
696               point_index_count      = point_index_count + 1
697               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
698               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
699               surfaces%points(3,surfaces%npoints) = zw(k)
700            ENDIF
701            IF ( point_index(k,j+1,i) < 0 )  THEN
702               surfaces%npoints     = surfaces%npoints + 1
703               point_index(k,j+1,i) = point_index_count
704               point_index_count    = point_index_count + 1
705               surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
706               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
707               surfaces%points(3,surfaces%npoints) = zw(k)
708            ENDIF
709
710            npg                        = npg + 1
711            surfaces%polygons(1,npg)   = 4
712            surfaces%polygons(2,npg)   = point_index(k,j,i)
713            surfaces%polygons(3,npg)   = point_index(k,j,i+1)
714            surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
715            surfaces%polygons(5,npg)   = point_index(k,j+1,i)
716         ENDDO
717
718         DO  m = 1, surf_usm_h%ns
719            i = surf_usm_h%i(m) + surf_usm_h%ioff
720            j = surf_usm_h%j(m) + surf_usm_h%joff
721            k = surf_usm_h%k(m) + surf_usm_h%koff
722
723            IF ( point_index(k,j,i) < 0 )  THEN
724               surfaces%npoints   = surfaces%npoints + 1
725               point_index(k,j,i) = point_index_count
726               point_index_count  = point_index_count + 1
727               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
728               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
729               surfaces%points(3,surfaces%npoints) = zw(k)
730            ENDIF
731            IF ( point_index(k,j,i+1) < 0 )  THEN
732               surfaces%npoints     = surfaces%npoints + 1
733               point_index(k,j,i+1) = point_index_count
734               point_index_count    = point_index_count + 1
735               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
736               surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
737               surfaces%points(3,surfaces%npoints) = zw(k)
738            ENDIF
739            IF ( point_index(k,j+1,i+1) < 0 )  THEN
740               surfaces%npoints       = surfaces%npoints + 1
741               point_index(k,j+1,i+1) = point_index_count
742               point_index_count      = point_index_count + 1
743               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
744               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
745               surfaces%points(3,surfaces%npoints) = zw(k)
746            ENDIF
747            IF ( point_index(k,j+1,i) < 0 )  THEN
748               surfaces%npoints     = surfaces%npoints + 1
749               point_index(k,j+1,i) = point_index_count
750               point_index_count    = point_index_count + 1
751               surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
752               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
753               surfaces%points(3,surfaces%npoints) = zw(k)
754            ENDIF
755
756            npg                        = npg + 1
757            surfaces%polygons(1,npg)   = 4
758            surfaces%polygons(2,npg)   = point_index(k,j,i)
759            surfaces%polygons(3,npg)   = point_index(k,j,i+1)
760            surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
761            surfaces%polygons(5,npg)   = point_index(k,j+1,i)
762         ENDDO
763
764         DO  l = 0, 3
765            DO  m = 1, surf_def_v(l)%ns
766!
767!--            Determine the indices of the respective grid cell inside the
768!--            topography.
769!--            NOTE, j-index for north-facing surfaces ( l==0 ) is
770!--            identical to the reference j-index outside the grid box.
771!--            Equivalent for east-facing surfaces and i-index.
772               i = surf_def_v(l)%i(m) + MERGE( surf_def_v(l)%ioff, 0, l == 3 )
773               j = surf_def_v(l)%j(m) + MERGE( surf_def_v(l)%joff, 0, l == 1 )
774               k = surf_def_v(l)%k(m) + surf_def_v(l)%koff
775!
776!--            Lower left /front vertex
777               IF ( point_index(k,j,i) < 0 )  THEN
778                  surfaces%npoints   = surfaces%npoints + 1
779                  point_index(k,j,i) = point_index_count
780                  point_index_count  = point_index_count + 1
781                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
782                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
783                  surfaces%points(3,surfaces%npoints) = zw(k-1)
784               ENDIF
785!
786!--            Upper / lower right index for north- and south-facing surfaces
787               IF ( l == 0  .OR.  l == 1 )  THEN
788                  IF ( point_index(k,j,i+1) < 0 )  THEN
789                     surfaces%npoints     = surfaces%npoints + 1
790                     point_index(k,j,i+1) = point_index_count
791                     point_index_count    = point_index_count + 1
792                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
793                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
794                     surfaces%points(3,surfaces%npoints) = zw(k-1)
795                  ENDIF
796                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
797                     surfaces%npoints       = surfaces%npoints + 1
798                     point_index(k+1,j,i+1) = point_index_count
799                     point_index_count      = point_index_count + 1
800                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
801                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
802                     surfaces%points(3,surfaces%npoints) = zw(k)
803                  ENDIF
804!
805!--            Upper / lower front index for east- and west-facing surfaces
806               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
807                  IF ( point_index(k,j+1,i) < 0 )  THEN
808                     surfaces%npoints     = surfaces%npoints + 1
809                     point_index(k,j+1,i) = point_index_count
810                     point_index_count    = point_index_count + 1
811                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
812                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
813                     surfaces%points(3,surfaces%npoints) = zw(k-1)
814                  ENDIF
815                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
816                     surfaces%npoints       = surfaces%npoints + 1
817                     point_index(k+1,j+1,i) = point_index_count
818                     point_index_count      = point_index_count + 1
819                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
820                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
821                     surfaces%points(3,surfaces%npoints) = zw(k)
822                  ENDIF
823               ENDIF
824!
825!--            Upper left / front vertex
826               IF ( point_index(k+1,j,i) < 0 )  THEN
827                  surfaces%npoints     = surfaces%npoints + 1
828                  point_index(k+1,j,i) = point_index_count
829                  point_index_count    = point_index_count + 1
830                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
831                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
832                  surfaces%points(3,surfaces%npoints) = zw(k)
833               ENDIF
834
835               npg = npg + 1
836               IF ( l == 0  .OR.  l == 1 )  THEN
837                  surfaces%polygons(1,npg)   = 4
838                  surfaces%polygons(2,npg)   = point_index(k,j,i)
839                  surfaces%polygons(3,npg)   = point_index(k,j,i+1)
840                  surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
841                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
842               ELSE
843                  surfaces%polygons(1,npg)   = 4
844                  surfaces%polygons(2,npg)   = point_index(k,j,i)
845                  surfaces%polygons(3,npg)   = point_index(k,j+1,i)
846                  surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
847                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
848               ENDIF
849
850            ENDDO
851
852            DO  m = 1, surf_lsm_v(l)%ns
853               i = surf_lsm_v(l)%i(m) + MERGE( surf_lsm_v(l)%ioff, 0, l == 3 )
854               j = surf_lsm_v(l)%j(m) + MERGE( surf_lsm_v(l)%joff, 0, l == 1 )
855               k = surf_lsm_v(l)%k(m) + surf_lsm_v(l)%koff
856!
857!--            Lower left /front vertex
858               IF ( point_index(k,j,i) < 0 )  THEN
859                  surfaces%npoints   = surfaces%npoints + 1
860                  point_index(k,j,i) = point_index_count
861                  point_index_count  = point_index_count + 1
862                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
863                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
864                  surfaces%points(3,surfaces%npoints) = zw(k-1)
865               ENDIF
866!
867!--            Upper / lower right index for north- and south-facing surfaces
868               IF ( l == 0  .OR.  l == 1 )  THEN
869                  IF ( point_index(k,j,i+1) < 0 )  THEN
870                     surfaces%npoints     = surfaces%npoints + 1
871                     point_index(k,j,i+1) = point_index_count
872                     point_index_count    = point_index_count + 1
873                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
874                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
875                     surfaces%points(3,surfaces%npoints) = zw(k-1)
876                  ENDIF
877                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
878                     surfaces%npoints       = surfaces%npoints + 1
879                     point_index(k+1,j,i+1) = point_index_count
880                     point_index_count      = point_index_count + 1
881                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
882                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
883                     surfaces%points(3,surfaces%npoints) = zw(k)
884                  ENDIF
885!
886!--            Upper / lower front index for east- and west-facing surfaces
887               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
888                  IF ( point_index(k,j+1,i) < 0 )  THEN
889                     surfaces%npoints     = surfaces%npoints + 1
890                     point_index(k,j+1,i) = point_index_count
891                     point_index_count    = point_index_count + 1
892                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
893                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
894                     surfaces%points(3,surfaces%npoints) = zw(k-1)
895                  ENDIF
896                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
897                     surfaces%npoints       = surfaces%npoints + 1
898                     point_index(k+1,j+1,i) = point_index_count
899                     point_index_count      = point_index_count + 1
900                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
901                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
902                     surfaces%points(3,surfaces%npoints) = zw(k)
903                  ENDIF
904               ENDIF
905!
906!--            Upper left / front vertex
907               IF ( point_index(k+1,j,i) < 0 )  THEN
908                  surfaces%npoints     = surfaces%npoints + 1
909                  point_index(k+1,j,i) = point_index_count
910                  point_index_count    = point_index_count + 1
911                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
912                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
913                  surfaces%points(3,surfaces%npoints) = zw(k)
914               ENDIF
915
916               npg = npg + 1
917               IF ( l == 0  .OR.  l == 1 )  THEN
918                  surfaces%polygons(1,npg)   = 4
919                  surfaces%polygons(2,npg)   = point_index(k,j,i)
920                  surfaces%polygons(3,npg)   = point_index(k,j,i+1)
921                  surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
922                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
923               ELSE
924                  surfaces%polygons(1,npg)   = 4
925                  surfaces%polygons(2,npg)   = point_index(k,j,i)
926                  surfaces%polygons(3,npg)   = point_index(k,j+1,i)
927                  surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
928                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
929               ENDIF
930            ENDDO
931            DO  m = 1, surf_usm_v(l)%ns
932               i = surf_usm_v(l)%i(m) + MERGE( surf_usm_v(l)%ioff, 0, l == 3 )
933               j = surf_usm_v(l)%j(m) + MERGE( surf_usm_v(l)%joff, 0, l == 1 )
934               k = surf_usm_v(l)%k(m) + surf_usm_v(l)%koff
935!
936!--            Lower left /front vertex
937               IF ( point_index(k,j,i) < 0 )  THEN
938                  surfaces%npoints   = surfaces%npoints + 1
939                  point_index(k,j,i) = point_index_count
940                  point_index_count  = point_index_count + 1
941                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
942                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
943                  surfaces%points(3,surfaces%npoints) = zw(k-1)
944               ENDIF
945!
946!--            Upper / lower right index for north- and south-facing surfaces
947               IF ( l == 0  .OR.  l == 1 )  THEN
948                  IF ( point_index(k,j,i+1) < 0 )  THEN
949                     surfaces%npoints     = surfaces%npoints + 1
950                     point_index(k,j,i+1) = point_index_count
951                     point_index_count    = point_index_count + 1
952                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
953                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
954                     surfaces%points(3,surfaces%npoints) = zw(k-1)
955                  ENDIF
956                  IF ( point_index(k+1,j,i+1) < 0 )  THEN
957                     surfaces%npoints       = surfaces%npoints + 1
958                     point_index(k+1,j,i+1) = point_index_count
959                     point_index_count      = point_index_count + 1
960                     surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
961                     surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
962                     surfaces%points(3,surfaces%npoints) = zw(k)
963                  ENDIF
964!
965!--            Upper / lower front index for east- and west-facing surfaces
966               ELSEIF ( l == 2  .OR.  l == 3 )  THEN
967                  IF ( point_index(k,j+1,i) < 0 )  THEN
968                     surfaces%npoints     = surfaces%npoints + 1
969                     point_index(k,j+1,i) = point_index_count
970                     point_index_count    = point_index_count + 1
971                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
972                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
973                     surfaces%points(3,surfaces%npoints) = zw(k-1)
974                  ENDIF
975                  IF ( point_index(k+1,j+1,i) < 0 )  THEN
976                     surfaces%npoints       = surfaces%npoints + 1
977                     point_index(k+1,j+1,i) = point_index_count
978                     point_index_count      = point_index_count + 1
979                     surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
980                     surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
981                     surfaces%points(3,surfaces%npoints) = zw(k)
982                  ENDIF
983               ENDIF
984!
985!--            Upper left / front vertex
986               IF ( point_index(k+1,j,i) < 0 )  THEN
987                  surfaces%npoints     = surfaces%npoints + 1
988                  point_index(k+1,j,i) = point_index_count
989                  point_index_count    = point_index_count + 1
990                  surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
991                  surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
992                  surfaces%points(3,surfaces%npoints) = zw(k)
993               ENDIF
994
995               npg = npg + 1
996               IF ( l == 0  .OR.  l == 1 )  THEN
997                  surfaces%polygons(1,npg)   = 4
998                  surfaces%polygons(2,npg)   = point_index(k,j,i)
999                  surfaces%polygons(3,npg)   = point_index(k,j,i+1)
1000                  surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
1001                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
1002               ELSE
1003                  surfaces%polygons(1,npg)   = 4
1004                  surfaces%polygons(2,npg)   = point_index(k,j,i)
1005                  surfaces%polygons(3,npg)   = point_index(k,j+1,i)
1006                  surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
1007                  surfaces%polygons(5,npg)   = point_index(k+1,j,i)
1008               ENDIF
1009            ENDDO
1010
1011         ENDDO
1012!
1013!--      Deallocate temporary dummy variable
1014         DEALLOCATE ( point_index )
1015!
1016!--      Sum-up total number of vertices on domain. This
1017!--      will be needed for post-processing.
1018         surfaces%npoints_total = 0
1019#if defined( __parallel )
1020          CALL MPI_ALLREDUCE( surfaces%npoints, surfaces%npoints_total, 1,     &
1021                              MPI_INTEGER, MPI_SUM, comm2d, ierr )
1022#else
1023          surfaces%npoints_total = surfaces%npoints
1024#endif
1025       ENDIF
1026!
1027!--    If output to netcdf is enabled, set-up the coordinate arrays that
1028!--    unambiguously describe the position and orientation of each surface
1029!--    element.
1030       IF ( to_netcdf )  THEN
1031!
1032!--       Allocate local coordinate arrays
1033          ALLOCATE( surfaces%s(1:surfaces%ns)       )
1034          ALLOCATE( surfaces%xs(1:surfaces%ns)      )
1035          ALLOCATE( surfaces%ys(1:surfaces%ns)      )
1036          ALLOCATE( surfaces%zs(1:surfaces%ns)      )
1037          ALLOCATE( surfaces%azimuth(1:surfaces%ns) )
1038          ALLOCATE( surfaces%zenith(1:surfaces%ns)  )
1039          ALLOCATE( surfaces%es_utm(1:surfaces%ns)  )
1040          ALLOCATE( surfaces%ns_utm(1:surfaces%ns)  )
1041!
1042!--       Gather the number of surface on each processor, in order to number
1043!--       the surface elements in ascending order with respect to the total
1044!--       number of surfaces in the domain.
1045#if defined( __parallel )
1046          CALL MPI_ALLGATHER( surfaces%ns, 1, MPI_INTEGER,                     &
1047                              num_surfaces_on_pe, 1, MPI_INTEGER, comm2d, ierr  )
1048#else
1049          num_surfaces_on_pe = surfaces%ns
1050#endif
1051!
1052!--       First, however, determine the offset where couting of the surfaces
1053!--       should start (the sum of surfaces on all PE's with lower MPI rank).
1054          i           = 0
1055          start_count = 1
1056          DO WHILE ( i < myid  .AND.  i <= SIZE( num_surfaces_on_pe ) )
1057             start_count = start_count + num_surfaces_on_pe(i)
1058             i           = i + 1
1059          ENDDO
1060!
1061!--       Set coordinate arrays. For horizontal surfaces, azimuth
1062!--       angles are not defined (fill value). Zenith angle is 0 (180) for
1063!--       upward (downward)-facing surfaces.
1064          i  = start_count
1065          mm = 1
1066          DO  m = 1, surf_def_h(0)%ns
1067             surfaces%s(mm)       = i
1068             surfaces%xs(mm)      = ( surf_def_h(0)%i(m) + 0.5_wp ) * dx
1069             surfaces%ys(mm)      = ( surf_def_h(0)%j(m) + 0.5_wp ) * dy
1070             surfaces%zs(mm)      = zw(surf_def_h(0)%k(m)+surf_def_h(0)%koff)
1071             surfaces%azimuth(mm) = surfaces%fillvalue
1072             surfaces%zenith(mm)  = 0.0
1073             i                    = i + 1
1074             mm                   = mm + 1
1075          ENDDO
1076          DO  m = 1, surf_lsm_h%ns
1077             surfaces%s(mm)       = i
1078             surfaces%xs(mm)      = ( surf_lsm_h%i(m) + 0.5_wp ) * dx
1079             surfaces%ys(mm)      = ( surf_lsm_h%j(m) + 0.5_wp ) * dy
1080             surfaces%zs(mm)      = zw(surf_lsm_h%k(m)+surf_lsm_h%koff)
1081             surfaces%azimuth(mm) = surfaces%fillvalue
1082             surfaces%zenith(mm)  = 0.0
1083             i                    = i + 1
1084             mm                   = mm + 1
1085          ENDDO
1086          DO  m = 1, surf_usm_h%ns
1087             surfaces%s(mm)       = i
1088             surfaces%xs(mm)      = ( surf_usm_h%i(m) + 0.5_wp ) * dx
1089             surfaces%ys(mm)      = ( surf_usm_h%j(m) + 0.5_wp ) * dy
1090             surfaces%zs(mm)      = zw(surf_usm_h%k(m)+surf_usm_h%koff)
1091             surfaces%azimuth(mm) = surfaces%fillvalue
1092             surfaces%zenith(mm)  = 0.0
1093             i                    = i + 1
1094             mm                   = mm + 1
1095          ENDDO
1096          DO  m = 1, surf_def_h(1)%ns
1097             surfaces%s(mm)       = i
1098             surfaces%xs(mm)      = ( surf_def_h(1)%i(m) + 0.5_wp ) * dx
1099             surfaces%ys(mm)      = ( surf_def_h(1)%j(m) + 0.5_wp ) * dy
1100             surfaces%zs(mm)      = zw(surf_def_h(1)%k(m)+surf_def_h(1)%koff)
1101             surfaces%azimuth(mm) = surfaces%fillvalue
1102             surfaces%zenith(mm)  = 180.0
1103             i                    = i + 1
1104             mm                   = mm + 1
1105          ENDDO
1106!
1107!--       For vertical surfaces, zenith angles are not defined (fill value).
1108!--       Azimuth angle: northward (0), eastward (90), southward (180),
1109!--       westward (270).
1110!--       Note, for vertical surfaces, zenith angles are 90.0_wp.
1111          DO  l = 0, 3
1112             IF ( l == 0 )  THEN
1113                az    = 0.0_wp
1114                off_x = 0.5_wp
1115                off_y = 0.0_wp
1116             ELSEIF ( l == 1 )  THEN
1117                az    = 180.0_wp
1118                off_x = 0.5_wp
1119                off_y = 1.0_wp
1120             ELSEIF ( l == 2 )  THEN
1121                az    = 90.0_wp
1122                off_x = 0.0_wp
1123                off_y = 0.5_wp
1124             ELSEIF ( l == 3 )  THEN
1125                az    = 270.0_wp
1126                off_x = 1.0_wp
1127                off_y = 0.5_wp
1128             ENDIF
1129
1130             DO  m = 1, surf_def_v(l)%ns
1131                surfaces%s(mm)       = i
1132                surfaces%xs(mm)      = ( surf_def_v(l)%i(m) + off_x ) * dx
1133                surfaces%ys(mm)      = ( surf_def_v(l)%j(m) + off_y ) * dy
1134                surfaces%zs(mm)      = zu(surf_def_v(l)%k(m))
1135                surfaces%azimuth(mm) = az
1136                surfaces%zenith(mm)  = 90.0_wp
1137                i                    = i + 1
1138                mm                   = mm + 1
1139             ENDDO
1140             DO  m = 1, surf_lsm_v(l)%ns
1141                surfaces%s(mm)       = i
1142                surfaces%xs(mm)      = ( surf_lsm_v(l)%i(m) + off_x ) * dx
1143                surfaces%ys(mm)      = ( surf_lsm_v(l)%j(m) + off_y ) * dy
1144                surfaces%zs(mm)      = zu(surf_lsm_v(l)%k(m))
1145                surfaces%azimuth(mm) = az
1146                surfaces%zenith(mm)  = 90.0_wp
1147                i                    = i + 1
1148                mm                   = mm + 1
1149             ENDDO
1150             DO  m = 1, surf_usm_v(l)%ns
1151                surfaces%s(mm)       = i
1152                surfaces%xs(mm)      = ( surf_usm_v(l)%i(m) + off_x ) * dx
1153                surfaces%ys(mm)      = ( surf_usm_v(l)%j(m) + off_y ) * dy
1154                surfaces%zs(mm)      = zu(surf_usm_v(l)%k(m))
1155                surfaces%azimuth(mm) = az
1156                surfaces%zenith(mm)  = 90.0_wp
1157                i                    = i + 1
1158                mm                   = mm + 1
1159             ENDDO
1160          ENDDO
1161!
1162!--       Finally, define UTM coordinates, which are the x/y-coordinates
1163!--       plus the origin (lower-left coordinate of the model domain).
1164          surfaces%es_utm = surfaces%xs + init_model%origin_x
1165          surfaces%ns_utm = surfaces%ys + init_model%origin_y
1166!
1167!--       Initialize NetCDF data output. Please note, local start position for
1168!--       the surface elements in the NetCDF file is surfaces%s(1), while
1169!--       the number of surfaces on the subdomain is given by surfaces%ns.
1170#if defined( __netcdf4_parallel )
1171
1172!
1173!--       Calculate number of time steps to be output
1174          ntdim_surf(0) = dosurf_time_count(0) + CEILING(                      &
1175                        ( end_time - MAX(                                      &
1176                            MERGE( skip_time_dosurf,                           &
1177                                   skip_time_dosurf + spinup_time,             &
1178                                   data_output_during_spinup ),                &
1179                            simulated_time_at_begin )                          &
1180                        ) / dt_dosurf )
1181
1182          ntdim_surf(1) = dosurf_time_count(1) + CEILING(                      &
1183                        ( end_time - MAX(                                      &
1184                            MERGE( skip_time_dosurf_av,                        &
1185                                   skip_time_dosurf_av + spinup_time,          &
1186                                   data_output_during_spinup ),                &
1187                            simulated_time_at_begin )                          &
1188                        ) / dt_dosurf_av )
1189
1190!
1191!--       Create NetCDF4 files for parallel writing
1192          DO  av = 0, 1
1193!
1194!--          If there is no instantaneous data (av=0) or averaged data (av=1)
1195!--          requested, do not create the corresponding NetCDF file
1196             IF ( dosurf_no(av) == 0 ) CYCLE
1197
1198             IF ( av == 0 )  THEN
1199                filename = 'SURFACE_DATA_NETCDF' // TRIM( coupling_char )
1200             ELSE
1201                filename = 'SURFACE_DATA_AV_NETCDF' // TRIM( coupling_char )
1202             ENDIF
1203!
1204!--          Open file using netCDF4/HDF5 format, parallel
1205             nc_stat = NF90_CREATE( TRIM(filename),                            &
1206                                    IOR( NF90_NOCLOBBER,                       &
1207                                    IOR( NF90_NETCDF4, NF90_MPIIO ) ),         &
1208                                    id_set_surf(av),                           &
1209                                    COMM = comm2d, INFO = MPI_INFO_NULL )
1210             CALL netcdf_handle_error( 'surface_data_output_mod', 5550 )
1211
1212             !- Write some global attributes
1213             IF ( av == 0 )  THEN
1214                CALL netcdf_create_global_atts( id_set_surf(av),               &
1215                                                'surface-data',                &
1216                                                TRIM( run_description_header ),&
1217                                                5551 )
1218                time_average_text = ' '
1219             ELSE
1220                CALL netcdf_create_global_atts( id_set_surf(av),               &
1221                                                'surface-data_av',             &
1222                                                TRIM( run_description_header ),&
1223                                                5552 )
1224                WRITE ( time_average_text,'(F7.1,'' s avg'')' )  &
1225                   averaging_interval_surf
1226                nc_stat = NF90_PUT_ATT( id_set_surf(av), NF90_GLOBAL,          &
1227                                        'time_avg',                            &
1228                                        TRIM( time_average_text ) )
1229                CALL netcdf_handle_error( 'surface_data_output_mod', 5553 )
1230             ENDIF
1231
1232
1233!
1234!--          Define time coordinate for surface data.
1235!--          For parallel output the time dimension has to be limited
1236!--          (ntdim_surf), otherwise the performance drops significantly.
1237             CALL netcdf_create_dim( id_set_surf(av), 'time', ntdim_surf(av),  &
1238                                     id_dim_time_surf(av), 5554 )
1239
1240             CALL netcdf_create_var( id_set_surf(av),                          &
1241                                     (/ id_dim_time_surf(av) /),               &
1242                                     'time', NF90_DOUBLE,                      &
1243                                     id_var_time_surf(av),                     &
1244                                     'seconds since '//                        &
1245                                     TRIM(init_model%origin_time),             &
1246                                     'time', 5555, 5555, 5555 )
1247
1248             CALL netcdf_create_att( id_set_surf(av), id_var_time_surf(av),    &
1249                                     'standard_name', 'time', 5556)
1250
1251             CALL netcdf_create_att( id_set_surf(av), id_var_time_surf(av),    &
1252                                     'axis', 'T', 5557)
1253!
1254!--          Define spatial dimensions and coordinates:
1255!--          Define index of surface element
1256             CALL netcdf_create_dim( id_set_surf(av), 's', surfaces%ns_total,  &
1257                                     id_dim_s_surf(av), 5558 )
1258             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1259                                     's', NF90_DOUBLE, id_var_s_surf(av),      &
1260                                     '1', 'number of surface element',         &
1261                                     5559, 5559, 5559 )
1262!
1263!--          Define x coordinate
1264             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1265                                     'xs', NF90_DOUBLE, id_var_xs_surf(av),    &
1266                                     'meters',                                 &
1267                                     'distance to origin in x-direction',      &
1268                                     5561, 5561, 5561 )
1269!
1270!--           Define y coordinate
1271             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1272                                     'ys', NF90_DOUBLE, id_var_ys_surf(av),    &
1273                                     'meters',                                 &
1274                                     'distance to origin in y-direction',      &
1275                                     5562, 5562, 5562 )
1276!
1277!--          Define z coordinate
1278             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1279                                     'zs', NF90_DOUBLE, id_var_zs_surf(av),    &
1280                                     'meters', 'height', 5560, 5560, 5560 )
1281             CALL netcdf_create_att( id_set_surf(av), id_var_zs_surf(av),      &
1282                                     'standard_name', 'height', 5583 )
1283
1284!
1285!--          Define UTM coordinates
1286             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1287                                     'Es_UTM', NF90_DOUBLE,                    &
1288                                     id_var_etum_surf(av),                     &
1289                                     'meters', '', 5563, 5563, 5563 )
1290             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1291                                     'Ns_UTM', NF90_DOUBLE,                    &
1292                                     id_var_nutm_surf(av),                     &
1293                                     'meters', '', 5564, 5564, 5564 )
1294
1295!
1296!--          Define angles
1297             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1298                                     'azimuth', NF90_DOUBLE,                   &
1299                                     id_var_azimuth_surf(av),                  &
1300                                     'degree', 'azimuth angle',                &
1301                                     5577, 5578, 5579,                         &
1302                                     fill = .TRUE. )
1303             CALL netcdf_create_att( id_set_surf(av), id_var_azimuth_surf(av), &
1304                                     'standard_name', 'surface_azimuth_angle', &
1305                                     5584 )
1306
1307             CALL netcdf_create_var( id_set_surf(av), (/ id_dim_s_surf(av) /), &
1308                                     'zenith', NF90_DOUBLE,                    &
1309                                     id_var_zenith_surf(av),                   &
1310                                     'degree', '', 5580, 5581, 5582,           &
1311                                     fill = .TRUE. )
1312!
1313!--          Define the variables
1314             var_list = ';'
1315             i = 1
1316
1317             DO WHILE ( dosurf(av,i)(1:1) /= ' ' )
1318
1319                CALL netcdf_create_var( id_set_surf(av),(/ id_dim_s_surf(av),  &
1320                                        id_dim_time_surf(av) /), dosurf(av,i), &
1321                                        NF90_REAL4, id_var_dosurf(av,i),       &
1322                                        dosurf_unit(av,i), dosurf(av,i), 5565, &
1323                                        5565, 5565, .TRUE. )
1324!
1325!--                Set no fill for every variable to increase performance.
1326                nc_stat = NF90_DEF_VAR_FILL( id_set_surf(av),                  &
1327                                             id_var_dosurf(av,i),              &
1328                                             NF90_NOFILL, 0 )
1329                CALL netcdf_handle_error( 'surface_data_output_init', 5566 )
1330!
1331!--                Set collective io operations for parallel io
1332                nc_stat = NF90_VAR_PAR_ACCESS( id_set_surf(av),                &
1333                                               id_var_dosurf(av,i),            &
1334                                               NF90_COLLECTIVE )
1335                CALL netcdf_handle_error( 'surface_data_output_init', 5567 )
1336                var_list = TRIM( var_list ) // TRIM( dosurf(av,i) ) // ';'
1337
1338                i = i + 1
1339
1340             ENDDO
1341!
1342!--          Write the list of variables as global attribute (this is used by
1343!--          restart runs and by combine_plot_fields)
1344             nc_stat = NF90_PUT_ATT( id_set_surf(av), NF90_GLOBAL, 'VAR_LIST', &
1345                                     var_list )
1346             CALL netcdf_handle_error( 'surface_data_output_init', 5568 )
1347
1348!
1349!--          Set general no fill, otherwise the performance drops significantly
1350!--          for parallel output.
1351             nc_stat = NF90_SET_FILL( id_set_surf(av), NF90_NOFILL, oldmode )
1352             CALL netcdf_handle_error( 'surface_data_output_init', 5569 )
1353
1354!
1355!--          Leave netCDF define mode
1356             nc_stat = NF90_ENDDEF( id_set_surf(av) )
1357             CALL netcdf_handle_error( 'surface_data_output_init', 5570 )
1358
1359!
1360!--          These data are only written by PE0 for parallel output to increase
1361!--          the performance.
1362             IF ( myid == 0 )  THEN
1363!
1364!--             Write data for surface indices
1365                ALLOCATE( netcdf_data_1d(1:surfaces%ns_total) )
1366
1367                DO  i = 1, surfaces%ns_total
1368                   netcdf_data_1d(i) = i
1369                ENDDO
1370
1371                nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_s_surf(av),    &
1372                                        netcdf_data_1d, start = (/ 1 /),       &
1373                                        count = (/ surfaces%ns_total /) )
1374                CALL netcdf_handle_error( 'surface_data_output_init', 5571 )
1375
1376                DEALLOCATE( netcdf_data_1d )
1377
1378             ENDIF
1379
1380!
1381!--          Write surface positions to file
1382             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_xs_surf(av),      &
1383                                     surfaces%xs, start = (/ surfaces%s(1) /), &
1384                                     count = (/ surfaces%ns /) )
1385             CALL netcdf_handle_error( 'surface_data_output_init', 5572 )
1386
1387             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_ys_surf(av),      &
1388                                     surfaces%ys, start = (/ surfaces%s(1) /), &
1389                                     count = (/ surfaces%ns /) )
1390             CALL netcdf_handle_error( 'surface_data_output_init', 5573 )
1391
1392             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_zs_surf(av),      &
1393                                     surfaces%zs, start = (/ surfaces%s(1) /), &
1394                                     count = (/ surfaces%ns /) )
1395             CALL netcdf_handle_error( 'surface_data_output_init', 5574 )
1396
1397             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_etum_surf(av),    &
1398                                     surfaces%es_utm,                          &
1399                                     start = (/ surfaces%s(1) /),              &
1400                                     count = (/ surfaces%ns /) )
1401             CALL netcdf_handle_error( 'surface_data_output_init', 5575 )
1402
1403             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_nutm_surf(av),    &
1404                                     surfaces%ns_utm,                          &
1405                                     start = (/ surfaces%s(1) /),              &
1406                                     count = (/ surfaces%ns /) )
1407             CALL netcdf_handle_error( 'surface_data_output_init', 5576 )
1408
1409             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_azimuth_surf(av), &
1410                                     surfaces%azimuth,                         &
1411                                     start = (/ surfaces%s(1) /),              &
1412                                     count = (/ surfaces%ns /) )
1413             CALL netcdf_handle_error( 'surface_data_output_init', 5585 )
1414
1415             nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_zenith_surf(av),  &
1416                                     surfaces%zenith,                          &
1417                                     start = (/ surfaces%s(1) /),              &
1418                                     count = (/ surfaces%ns /) )
1419             CALL netcdf_handle_error( 'surface_data_output_init', 5586 )
1420
1421          ENDDO
1422#endif
1423
1424       ENDIF
1425
1426   END SUBROUTINE surface_data_output_init
1427
1428!------------------------------------------------------------------------------!
1429! Description:
1430! ------------
1431!> Routine for controlling the data output. Surface data is collected from
1432!> different types of surfaces (default, natural, urban) and different
1433!> orientation and written to one 1D-output array. Further, NetCDF routines
1434!> are called to write the surface data in the respective NetCDF files.
1435!------------------------------------------------------------------------------!
1436   SUBROUTINE surface_data_output( av )
1437
1438      USE control_parameters,                                                  &
1439          ONLY:  io_blocks, io_group, time_since_reference_point
1440
1441#if defined( __parallel )
1442      USE pegrid,                                                              &
1443          ONLY:  comm2d, ierr
1444#endif
1445
1446
1447      IMPLICIT NONE
1448
1449      CHARACTER(LEN=100) ::  trimvar = ' ' !< dummy for single output variable
1450
1451      INTEGER(iwp) ::  av     !< id indicating average or non-average data output
1452      INTEGER(iwp) ::  i      !< loop index
1453      INTEGER(iwp) ::  l      !< running index for surface orientation
1454      INTEGER(iwp) ::  m      !< running index for surface elements
1455      INTEGER(iwp) ::  n_out  !< counter variables for surface output
1456
1457!
1458!--   Return, if nothing to output
1459      IF ( dosurf_no(av) == 0 )  RETURN
1460!
1461!--   In case of VTK output, check if binary files are open and write coordinates.
1462      IF ( to_vtk )  THEN
1463
1464         CALL check_open( 25+av )
1465
1466         IF ( .NOT. first_output(av) )  THEN
1467            DO  i = 0, io_blocks-1
1468               IF ( i == io_group )  THEN
1469                  WRITE ( 25+av )  surfaces%npoints
1470                  WRITE ( 25+av )  surfaces%npoints_total
1471                  WRITE ( 25+av )  surfaces%ns
1472                  WRITE ( 25+av )  surfaces%ns_total
1473                  WRITE ( 25+av )  surfaces%points
1474                  WRITE ( 25+av )  surfaces%polygons
1475               ENDIF
1476#if defined( __parallel )
1477               CALL MPI_BARRIER( comm2d, ierr )
1478#endif
1479               first_output(av) = .TRUE.
1480            ENDDO
1481         ENDIF
1482      ENDIF
1483!
1484!--   In case of NetCDF output, check if enough time steps are available in file
1485!--   and update time axis.
1486      IF ( to_netcdf )  THEN
1487#if defined( __netcdf4_parallel )
1488         IF ( dosurf_time_count(av) + 1 > ntdim_surf(av) )  THEN
1489            WRITE ( message_string, * )                               &
1490               'Output of surface data is not given at t=',           &
1491               time_since_reference_point, 's because the maximum ',  &
1492               'number of output time levels is exceeded.'
1493            CALL message( 'surface_data_output', 'PA0539', 0, 1, 0, 6, 0 )
1494
1495            RETURN
1496
1497         ENDIF
1498!
1499!--      Update the netCDF time axis
1500!--      In case of parallel output, this is only done by PE0 to increase the
1501!--      performance.
1502         dosurf_time_count(av) = dosurf_time_count(av) + 1
1503         IF ( myid == 0 )  THEN
1504            nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_time_surf(av),  &
1505                                    (/ time_since_reference_point /),       &
1506                                    start = (/ dosurf_time_count(av) /),    &
1507                                    count = (/ 1 /) )
1508            CALL netcdf_handle_error( 'surface_data_output', 6666 )
1509         ENDIF
1510#endif
1511      ENDIF
1512
1513!
1514!--   Cycle through output quantities and write them to file.
1515      n_out = 0
1516      DO  WHILE ( dosurf(av,n_out+1)(1:1) /= ' ' )
1517
1518         n_out   = n_out + 1
1519         trimvar = TRIM( dosurf(av,n_out) )
1520!
1521!--      Set the output array to the _FillValue in case it is not
1522!--      defined for each type of surface.
1523         surfaces%var_out = surfaces%fillvalue
1524         SELECT CASE ( trimvar )
1525
1526            CASE ( 'us' )
1527!
1528!--            Output of instantaneous data
1529               IF ( av == 0 )  THEN
1530                  CALL surface_data_output_collect( surf_def_h(0)%us,          &
1531                                               surf_def_h(1)%us,               &
1532                                               surf_lsm_h%us,                  &
1533                                               surf_usm_h%us,                  &
1534                                               surf_def_v(0)%us,               &
1535                                               surf_lsm_v(0)%us,               &
1536                                               surf_usm_v(0)%us,               &
1537                                               surf_def_v(1)%us,               &
1538                                               surf_lsm_v(1)%us,               &
1539                                               surf_usm_v(1)%us,               &
1540                                               surf_def_v(2)%us,               &
1541                                               surf_lsm_v(2)%us,               &
1542                                               surf_usm_v(2)%us,               &
1543                                               surf_def_v(3)%us,               &
1544                                               surf_lsm_v(3)%us,               &
1545                                               surf_usm_v(3)%us )
1546               ELSE
1547!
1548!--               Output of averaged data
1549                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1550                                        REAL( average_count_surf, KIND=wp )
1551                  surfaces%var_av(:,n_out) = 0.0_wp
1552
1553               ENDIF
1554
1555            CASE ( 'ts' )
1556!
1557!--            Output of instantaneous data
1558               IF ( av == 0 )  THEN
1559                  CALL surface_data_output_collect( surf_def_h(0)%ts,          &
1560                                               surf_def_h(1)%ts,               &
1561                                               surf_lsm_h%ts,                  &
1562                                               surf_usm_h%ts,                  &
1563                                               surf_def_v(0)%ts,               &
1564                                               surf_lsm_v(0)%ts,               &
1565                                               surf_usm_v(0)%ts,               &
1566                                               surf_def_v(1)%ts,               &
1567                                               surf_lsm_v(1)%ts,               &
1568                                               surf_usm_v(1)%ts,               &
1569                                               surf_def_v(2)%ts,               &
1570                                               surf_lsm_v(2)%ts,               &
1571                                               surf_usm_v(2)%ts,               &
1572                                               surf_def_v(3)%ts,               &
1573                                               surf_lsm_v(3)%ts,               &
1574                                               surf_usm_v(3)%ts )
1575               ELSE
1576!
1577!--               Output of averaged data
1578                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1579                                        REAL( average_count_surf, KIND=wp )
1580                  surfaces%var_av(:,n_out) = 0.0_wp
1581
1582               ENDIF
1583
1584            CASE ( 'qs' )
1585!
1586!--            Output of instantaneous data
1587               IF ( av == 0 )  THEN
1588                  CALL surface_data_output_collect( surf_def_h(0)%qs,          &
1589                                               surf_def_h(1)%qs,               &
1590                                               surf_lsm_h%qs,                  &
1591                                               surf_usm_h%qs,                  &
1592                                               surf_def_v(0)%qs,               &
1593                                               surf_lsm_v(0)%qs,               &
1594                                               surf_usm_v(0)%qs,               &
1595                                               surf_def_v(1)%qs,               &
1596                                               surf_lsm_v(1)%qs,               &
1597                                               surf_usm_v(1)%qs,               &
1598                                               surf_def_v(2)%qs,               &
1599                                               surf_lsm_v(2)%qs,               &
1600                                               surf_usm_v(2)%qs,               &
1601                                               surf_def_v(3)%qs,               &
1602                                               surf_lsm_v(3)%qs,               &
1603                                               surf_usm_v(3)%qs )
1604               ELSE
1605!
1606!--               Output of averaged data
1607                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1608                                        REAL( average_count_surf, KIND=wp )
1609                  surfaces%var_av(:,n_out) = 0.0_wp
1610
1611               ENDIF
1612
1613            CASE ( 'ss' )
1614!
1615!--            Output of instantaneous data
1616               IF ( av == 0 )  THEN
1617                  CALL surface_data_output_collect( surf_def_h(0)%ss,          &
1618                                               surf_def_h(1)%ss,               &
1619                                               surf_lsm_h%ss,                  &
1620                                               surf_usm_h%ss,                  &
1621                                               surf_def_v(0)%ss,               &
1622                                               surf_lsm_v(0)%ss,               &
1623                                               surf_usm_v(0)%ss,               &
1624                                               surf_def_v(1)%ss,               &
1625                                               surf_lsm_v(1)%ss,               &
1626                                               surf_usm_v(1)%ss,               &
1627                                               surf_def_v(2)%ss,               &
1628                                               surf_lsm_v(2)%ss,               &
1629                                               surf_usm_v(2)%ss,               &
1630                                               surf_def_v(3)%ss,               &
1631                                               surf_lsm_v(3)%ss,               &
1632                                               surf_usm_v(3)%ss )
1633               ELSE
1634!
1635!--               Output of averaged data
1636                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1637                                        REAL( average_count_surf, KIND=wp )
1638                  surfaces%var_av(:,n_out) = 0.0_wp
1639
1640               ENDIF
1641
1642            CASE ( 'qcs' )
1643!
1644!--            Output of instantaneous data
1645               IF ( av == 0 )  THEN
1646                  CALL surface_data_output_collect( surf_def_h(0)%qcs,         &
1647                                               surf_def_h(1)%qcs,              &
1648                                               surf_lsm_h%qcs,                 &
1649                                               surf_usm_h%qcs,                 &
1650                                               surf_def_v(0)%qcs,              &
1651                                               surf_lsm_v(0)%qcs,              &
1652                                               surf_usm_v(0)%qcs,              &
1653                                               surf_def_v(1)%qcs,              &
1654                                               surf_lsm_v(1)%qcs,              &
1655                                               surf_usm_v(1)%qcs,              &
1656                                               surf_def_v(2)%qcs,              &
1657                                               surf_lsm_v(2)%qcs,              &
1658                                               surf_usm_v(2)%qcs,              &
1659                                               surf_def_v(3)%qcs,              &
1660                                               surf_lsm_v(3)%qcs,              &
1661                                               surf_usm_v(3)%qcs )
1662               ELSE
1663!
1664!--               Output of averaged data
1665                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1666                                        REAL( average_count_surf, KIND=wp )
1667                  surfaces%var_av(:,n_out) = 0.0_wp
1668
1669               ENDIF
1670
1671            CASE ( 'ncs' )
1672!
1673!--            Output of instantaneous data
1674               IF ( av == 0 )  THEN
1675                  CALL surface_data_output_collect( surf_def_h(0)%ncs,         &
1676                                               surf_def_h(1)%ncs,              &
1677                                               surf_lsm_h%ncs,                 &
1678                                               surf_usm_h%ncs,                 &
1679                                               surf_def_v(0)%ncs,              &
1680                                               surf_lsm_v(0)%ncs,              &
1681                                               surf_usm_v(0)%ncs,              &
1682                                               surf_def_v(1)%ncs,              &
1683                                               surf_lsm_v(1)%ncs,              &
1684                                               surf_usm_v(1)%ncs,              &
1685                                               surf_def_v(2)%ncs,              &
1686                                               surf_lsm_v(2)%ncs,              &
1687                                               surf_usm_v(2)%ncs,              &
1688                                               surf_def_v(3)%ncs,              &
1689                                               surf_lsm_v(3)%ncs,              &
1690                                               surf_usm_v(3)%ncs )
1691               ELSE
1692!
1693!--               Output of averaged data
1694                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1695                                        REAL( average_count_surf, KIND=wp )
1696                  surfaces%var_av(:,n_out) = 0.0_wp
1697
1698               ENDIF
1699
1700            CASE ( 'qis' )
1701!
1702!--            Output of instantaneous data
1703               IF ( av == 0 )  THEN
1704                  CALL surface_data_output_collect( surf_def_h(0)%qis,         &
1705                                               surf_def_h(1)%qis,              &
1706                                               surf_lsm_h%qis,                 &
1707                                               surf_usm_h%qis,                 &
1708                                               surf_def_v(0)%qis,              &
1709                                               surf_lsm_v(0)%qis,              &
1710                                               surf_usm_v(0)%qis,              &
1711                                               surf_def_v(1)%qis,              &
1712                                               surf_lsm_v(1)%qis,              &
1713                                               surf_usm_v(1)%qis,              &
1714                                               surf_def_v(2)%qis,              &
1715                                               surf_lsm_v(2)%qis,              &
1716                                               surf_usm_v(2)%qis,              &
1717                                               surf_def_v(3)%qis,              &
1718                                               surf_lsm_v(3)%qis,              &
1719                                               surf_usm_v(3)%qis )
1720               ELSE
1721!
1722!--               Output of averaged data
1723                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1724                                        REAL( average_count_surf, KIND=wp )
1725                  surfaces%var_av(:,n_out) = 0.0_wp
1726
1727               ENDIF
1728
1729            CASE ( 'nis' )
1730!
1731!--            Output of instantaneous data
1732               IF ( av == 0 )  THEN
1733                  CALL surface_data_output_collect( surf_def_h(0)%nis,         &
1734                                               surf_def_h(1)%nis,              &
1735                                               surf_lsm_h%nis,                 &
1736                                               surf_usm_h%nis,                 &
1737                                               surf_def_v(0)%nis,              &
1738                                               surf_lsm_v(0)%nis,              &
1739                                               surf_usm_v(0)%nis,              &
1740                                               surf_def_v(1)%nis,              &
1741                                               surf_lsm_v(1)%nis,              &
1742                                               surf_usm_v(1)%nis,              &
1743                                               surf_def_v(2)%nis,              &
1744                                               surf_lsm_v(2)%nis,              &
1745                                               surf_usm_v(2)%nis,              &
1746                                               surf_def_v(3)%nis,              &
1747                                               surf_lsm_v(3)%nis,              &
1748                                               surf_usm_v(3)%nis )
1749               ELSE
1750!
1751!--               Output of averaged data
1752                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1753                                        REAL( average_count_surf, KIND=wp )
1754                  surfaces%var_av(:,n_out) = 0.0_wp
1755
1756               ENDIF
1757
1758            CASE ( 'qrs' )
1759!
1760!--            Output of instantaneous data
1761               IF ( av == 0 )  THEN
1762                  CALL surface_data_output_collect( surf_def_h(0)%qrs,         &
1763                                               surf_def_h(1)%qrs,              &
1764                                               surf_lsm_h%qrs,                 &
1765                                               surf_usm_h%qrs,                 &
1766                                               surf_def_v(0)%qrs,              &
1767                                               surf_lsm_v(0)%qrs,              &
1768                                               surf_usm_v(0)%qrs,              &
1769                                               surf_def_v(1)%qrs,              &
1770                                               surf_lsm_v(1)%qrs,              &
1771                                               surf_usm_v(1)%qrs,              &
1772                                               surf_def_v(2)%qrs,              &
1773                                               surf_lsm_v(2)%qrs,              &
1774                                               surf_usm_v(2)%qrs,              &
1775                                               surf_def_v(3)%qrs,              &
1776                                               surf_lsm_v(3)%qrs,              &
1777                                               surf_usm_v(3)%qrs )
1778               ELSE
1779!
1780!--               Output of averaged data
1781                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1782                                        REAL( average_count_surf, KIND=wp )
1783                  surfaces%var_av(:,n_out) = 0.0_wp
1784
1785               ENDIF
1786
1787            CASE ( 'nrs' )
1788!
1789!--            Output of instantaneous data
1790               IF ( av == 0 )  THEN
1791                  CALL surface_data_output_collect( surf_def_h(0)%nrs,         &
1792                                               surf_def_h(1)%nrs,              &
1793                                               surf_lsm_h%nrs,                 &
1794                                               surf_usm_h%nrs,                 &
1795                                               surf_def_v(0)%nrs,              &
1796                                               surf_lsm_v(0)%nrs,              &
1797                                               surf_usm_v(0)%nrs,              &
1798                                               surf_def_v(1)%nrs,              &
1799                                               surf_lsm_v(1)%nrs,              &
1800                                               surf_usm_v(1)%nrs,              &
1801                                               surf_def_v(2)%nrs,              &
1802                                               surf_lsm_v(2)%nrs,              &
1803                                               surf_usm_v(2)%nrs,              &
1804                                               surf_def_v(3)%nrs,              &
1805                                               surf_lsm_v(3)%nrs,              &
1806                                               surf_usm_v(3)%nrs )
1807               ELSE
1808!
1809!--               Output of averaged data
1810                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1811                                        REAL( average_count_surf, KIND=wp )
1812                  surfaces%var_av(:,n_out) = 0.0_wp
1813
1814               ENDIF
1815
1816            CASE ( 'ol' )
1817!
1818!--            Output of instantaneous data
1819               IF ( av == 0 )  THEN
1820                  CALL surface_data_output_collect( surf_def_h(0)%ol,          &
1821                                               surf_def_h(1)%ol,               &
1822                                               surf_lsm_h%ol,                  &
1823                                               surf_usm_h%ol,                  &
1824                                               surf_def_v(0)%ol,               &
1825                                               surf_lsm_v(0)%ol,               &
1826                                               surf_usm_v(0)%ol,               &
1827                                               surf_def_v(1)%ol,               &
1828                                               surf_lsm_v(1)%ol,               &
1829                                               surf_usm_v(1)%ol,               &
1830                                               surf_def_v(2)%ol,               &
1831                                               surf_lsm_v(2)%ol,               &
1832                                               surf_usm_v(2)%ol,               &
1833                                               surf_def_v(3)%ol,               &
1834                                               surf_lsm_v(3)%ol,               &
1835                                               surf_usm_v(3)%ol )
1836               ELSE
1837!
1838!--               Output of averaged data
1839                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1840                                        REAL( average_count_surf, KIND=wp )
1841                  surfaces%var_av(:,n_out) = 0.0_wp
1842
1843               ENDIF
1844
1845            CASE ( 'z0' )
1846!
1847!--            Output of instantaneous data
1848               IF ( av == 0 )  THEN
1849                  CALL surface_data_output_collect( surf_def_h(0)%z0,          &
1850                                               surf_def_h(1)%z0,               &
1851                                               surf_lsm_h%z0,                  &
1852                                               surf_usm_h%z0,                  &
1853                                               surf_def_v(0)%z0,               &
1854                                               surf_lsm_v(0)%z0,               &
1855                                               surf_usm_v(0)%z0,               &
1856                                               surf_def_v(1)%z0,               &
1857                                               surf_lsm_v(1)%z0,               &
1858                                               surf_usm_v(1)%z0,               &
1859                                               surf_def_v(2)%z0,               &
1860                                               surf_lsm_v(2)%z0,               &
1861                                               surf_usm_v(2)%z0,               &
1862                                               surf_def_v(3)%z0,               &
1863                                               surf_lsm_v(3)%z0,               &
1864                                               surf_usm_v(3)%z0 )
1865               ELSE
1866!
1867!--               Output of averaged data
1868                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1869                                        REAL( average_count_surf, KIND=wp )
1870                  surfaces%var_av(:,n_out) = 0.0_wp
1871
1872               ENDIF
1873
1874            CASE ( 'z0h' )
1875!
1876!--            Output of instantaneous data
1877               IF ( av == 0 )  THEN
1878                  CALL surface_data_output_collect( surf_def_h(0)%z0h,         &
1879                                               surf_def_h(1)%z0h,              &
1880                                               surf_lsm_h%z0h,                 &
1881                                               surf_usm_h%z0h,                 &
1882                                               surf_def_v(0)%z0h,              &
1883                                               surf_lsm_v(0)%z0h,              &
1884                                               surf_usm_v(0)%z0h,              &
1885                                               surf_def_v(1)%z0h,              &
1886                                               surf_lsm_v(1)%z0h,              &
1887                                               surf_usm_v(1)%z0h,              &
1888                                               surf_def_v(2)%z0h,              &
1889                                               surf_lsm_v(2)%z0h,              &
1890                                               surf_usm_v(2)%z0h,              &
1891                                               surf_def_v(3)%z0h,              &
1892                                               surf_lsm_v(3)%z0h,              &
1893                                               surf_usm_v(3)%z0h )
1894               ELSE
1895!
1896!--               Output of averaged data
1897                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1898                                        REAL( average_count_surf, KIND=wp )
1899                  surfaces%var_av(:,n_out) = 0.0_wp
1900
1901               ENDIF
1902
1903            CASE ( 'z0q' )
1904!
1905!--            Output of instantaneous data
1906               IF ( av == 0 )  THEN
1907                  CALL surface_data_output_collect( surf_def_h(0)%z0q,         &
1908                                               surf_def_h(1)%z0q,              &
1909                                               surf_lsm_h%z0q,                 &
1910                                               surf_usm_h%z0q,                 &
1911                                               surf_def_v(0)%z0q,              &
1912                                               surf_lsm_v(0)%z0q,              &
1913                                               surf_usm_v(0)%z0q,              &
1914                                               surf_def_v(1)%z0q,              &
1915                                               surf_lsm_v(1)%z0q,              &
1916                                               surf_usm_v(1)%z0q,              &
1917                                               surf_def_v(2)%z0q,              &
1918                                               surf_lsm_v(2)%z0q,              &
1919                                               surf_usm_v(2)%z0q,              &
1920                                               surf_def_v(3)%z0q,              &
1921                                               surf_lsm_v(3)%z0q,              &
1922                                               surf_usm_v(3)%z0q )
1923               ELSE
1924!
1925!--               Output of averaged data
1926                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1927                                        REAL( average_count_surf, KIND=wp )
1928                  surfaces%var_av(:,n_out) = 0.0_wp
1929
1930               ENDIF
1931
1932            CASE ( 'theta1' )
1933!
1934!--            Output of instantaneous data
1935               IF ( av == 0 )  THEN
1936                  CALL surface_data_output_collect( surf_def_h(0)%pt1,         &
1937                                               surf_def_h(1)%pt1,              &
1938                                               surf_lsm_h%pt1,                 &
1939                                               surf_usm_h%pt1,                 &
1940                                               surf_def_v(0)%pt1,              &
1941                                               surf_lsm_v(0)%pt1,              &
1942                                               surf_usm_v(0)%pt1,              &
1943                                               surf_def_v(1)%pt1,              &
1944                                               surf_lsm_v(1)%pt1,              &
1945                                               surf_usm_v(1)%pt1,              &
1946                                               surf_def_v(2)%pt1,              &
1947                                               surf_lsm_v(2)%pt1,              &
1948                                               surf_usm_v(2)%pt1,              &
1949                                               surf_def_v(3)%pt1,              &
1950                                               surf_lsm_v(3)%pt1,              &
1951                                               surf_usm_v(3)%pt1 )
1952               ELSE
1953!
1954!--               Output of averaged data
1955                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1956                                        REAL( average_count_surf, KIND=wp )
1957                  surfaces%var_av(:,n_out) = 0.0_wp
1958
1959               ENDIF
1960
1961            CASE ( 'qv1' )
1962!
1963!--            Output of instantaneous data
1964               IF ( av == 0 )  THEN
1965                  CALL surface_data_output_collect( surf_def_h(0)%qv1,         &
1966                                               surf_def_h(1)%qv1,              &
1967                                               surf_lsm_h%qv1,                 &
1968                                               surf_usm_h%qv1,                 &
1969                                               surf_def_v(0)%qv1,              &
1970                                               surf_lsm_v(0)%qv1,              &
1971                                               surf_usm_v(0)%qv1,              &
1972                                               surf_def_v(1)%qv1,              &
1973                                               surf_lsm_v(1)%qv1,              &
1974                                               surf_usm_v(1)%qv1,              &
1975                                               surf_def_v(2)%qv1,              &
1976                                               surf_lsm_v(2)%qv1,              &
1977                                               surf_usm_v(2)%qv1,              &
1978                                               surf_def_v(3)%qv1,              &
1979                                               surf_lsm_v(3)%qv1,              &
1980                                               surf_usm_v(3)%qv1 )
1981               ELSE
1982!
1983!--               Output of averaged data
1984                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1985                                        REAL( average_count_surf, KIND=wp )
1986                  surfaces%var_av(:,n_out) = 0.0_wp
1987
1988               ENDIF
1989
1990            CASE ( 'thetav1' )
1991!
1992!--            Output of instantaneous data
1993               IF ( av == 0 )  THEN
1994                  CALL surface_data_output_collect( surf_def_h(0)%vpt1,        &
1995                                               surf_def_h(1)%vpt1,             &
1996                                               surf_lsm_h%vpt1,                &
1997                                               surf_usm_h%vpt1,                &
1998                                               surf_def_v(0)%vpt1,             &
1999                                               surf_lsm_v(0)%vpt1,             &
2000                                               surf_usm_v(0)%vpt1,             &
2001                                               surf_def_v(1)%vpt1,             &
2002                                               surf_lsm_v(1)%vpt1,             &
2003                                               surf_usm_v(1)%vpt1,             &
2004                                               surf_def_v(2)%vpt1,             &
2005                                               surf_lsm_v(2)%vpt1,             &
2006                                               surf_usm_v(2)%vpt1,             &
2007                                               surf_def_v(3)%vpt1,             &
2008                                               surf_lsm_v(3)%vpt1,             &
2009                                               surf_usm_v(3)%vpt1 )
2010               ELSE
2011!
2012!--               Output of averaged data
2013                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2014                                        REAL( average_count_surf, KIND=wp )
2015                  surfaces%var_av(:,n_out) = 0.0_wp
2016
2017               ENDIF
2018
2019            CASE ( 'usws' )
2020!
2021!--            Output of instantaneous data
2022               IF ( av == 0 )  THEN
2023                  CALL surface_data_output_collect( surf_def_h(0)%usws,        &
2024                                               surf_def_h(1)%usws,             &
2025                                               surf_lsm_h%usws,                &
2026                                               surf_usm_h%usws,                &
2027                                               surf_def_v(0)%usws,             &
2028                                               surf_lsm_v(0)%usws,             &
2029                                               surf_usm_v(0)%usws,             &
2030                                               surf_def_v(1)%usws,             &
2031                                               surf_lsm_v(1)%usws,             &
2032                                               surf_usm_v(1)%usws,             &
2033                                               surf_def_v(2)%usws,             &
2034                                               surf_lsm_v(2)%usws,             &
2035                                               surf_usm_v(2)%usws,             &
2036                                               surf_def_v(3)%usws,             &
2037                                               surf_lsm_v(3)%usws,             &
2038                                               surf_usm_v(3)%usws,             &
2039                                               momentumflux_output_conversion )
2040               ELSE
2041!
2042!--               Output of averaged data
2043                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2044                                        REAL( average_count_surf, KIND=wp )
2045                  surfaces%var_av(:,n_out) = 0.0_wp
2046
2047               ENDIF
2048
2049            CASE ( 'vsws' )
2050!
2051!--            Output of instantaneous data
2052               IF ( av == 0 )  THEN
2053                  CALL surface_data_output_collect( surf_def_h(0)%vsws,        &
2054                                               surf_def_h(1)%vsws,             &
2055                                               surf_lsm_h%vsws,                &
2056                                               surf_usm_h%vsws,                &
2057                                               surf_def_v(0)%vsws,             &
2058                                               surf_lsm_v(0)%vsws,             &
2059                                               surf_usm_v(0)%vsws,             &
2060                                               surf_def_v(1)%vsws,             &
2061                                               surf_lsm_v(1)%vsws,             &
2062                                               surf_usm_v(1)%vsws,             &
2063                                               surf_def_v(2)%vsws,             &
2064                                               surf_lsm_v(2)%vsws,             &
2065                                               surf_usm_v(2)%vsws,             &
2066                                               surf_def_v(3)%vsws,             &
2067                                               surf_lsm_v(3)%vsws,             &
2068                                               surf_usm_v(3)%vsws,             &
2069                                               momentumflux_output_conversion )
2070               ELSE
2071!
2072!--               Output of averaged data
2073                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2074                                        REAL( average_count_surf, KIND=wp )
2075                  surfaces%var_av(:,n_out) = 0.0_wp
2076
2077               ENDIF
2078
2079            CASE ( 'shf' )
2080!
2081!--            Output of instantaneous data
2082               IF ( av == 0 )  THEN
2083                  CALL surface_data_output_collect( surf_def_h(0)%shf,         &
2084                                               surf_def_h(1)%shf,              &
2085                                               surf_lsm_h%shf,                 &
2086                                               surf_usm_h%shf,                 &
2087                                               surf_def_v(0)%shf,              &
2088                                               surf_lsm_v(0)%shf,              &
2089                                               surf_usm_v(0)%shf,              &
2090                                               surf_def_v(1)%shf,              &
2091                                               surf_lsm_v(1)%shf,              &
2092                                               surf_usm_v(1)%shf,              &
2093                                               surf_def_v(2)%shf,              &
2094                                               surf_lsm_v(2)%shf,              &
2095                                               surf_usm_v(2)%shf,              &
2096                                               surf_def_v(3)%shf,              &
2097                                               surf_lsm_v(3)%shf,              &
2098                                               surf_usm_v(3)%shf,              &
2099                                               heatflux_output_conversion )
2100               ELSE
2101!
2102!--               Output of averaged data
2103                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2104                                        REAL( average_count_surf, KIND=wp )
2105                  surfaces%var_av(:,n_out) = 0.0_wp
2106               ENDIF
2107
2108            CASE ( 'qsws' )
2109!
2110!--            Output of instantaneous data
2111               IF ( av == 0 )  THEN
2112                  CALL surface_data_output_collect( surf_def_h(0)%qsws,        &
2113                                               surf_def_h(1)%qsws,             &
2114                                               surf_lsm_h%qsws,                &
2115                                               surf_usm_h%qsws,                &
2116                                               surf_def_v(0)%qsws,             &
2117                                               surf_lsm_v(0)%qsws,             &
2118                                               surf_usm_v(0)%qsws,             &
2119                                               surf_def_v(1)%qsws,             &
2120                                               surf_lsm_v(1)%qsws,             &
2121                                               surf_usm_v(1)%qsws,             &
2122                                               surf_def_v(2)%qsws,             &
2123                                               surf_lsm_v(2)%qsws,             &
2124                                               surf_usm_v(2)%qsws,             &
2125                                               surf_def_v(3)%qsws,             &
2126                                               surf_lsm_v(3)%qsws,             &
2127                                               surf_usm_v(3)%qsws,             &
2128                                               waterflux_output_conversion )
2129               ELSE
2130!
2131!--               Output of averaged data
2132                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2133                                        REAL( average_count_surf, KIND=wp )
2134                  surfaces%var_av(:,n_out) = 0.0_wp
2135
2136               ENDIF
2137
2138            CASE ( 'ssws' )
2139!
2140!--            Output of instantaneous data
2141               IF ( av == 0 )  THEN
2142                  CALL surface_data_output_collect( surf_def_h(0)%ssws,        &
2143                                               surf_def_h(1)%ssws,             &
2144                                               surf_lsm_h%ssws,                &
2145                                               surf_usm_h%ssws,                &
2146                                               surf_def_v(0)%ssws,             &
2147                                               surf_lsm_v(0)%ssws,             &
2148                                               surf_usm_v(0)%ssws,             &
2149                                               surf_def_v(1)%ssws,             &
2150                                               surf_lsm_v(1)%ssws,             &
2151                                               surf_usm_v(1)%ssws,             &
2152                                               surf_def_v(2)%ssws,             &
2153                                               surf_lsm_v(2)%ssws,             &
2154                                               surf_usm_v(2)%ssws,             &
2155                                               surf_def_v(3)%ssws,             &
2156                                               surf_lsm_v(3)%ssws,             &
2157                                               surf_usm_v(3)%ssws )
2158               ELSE
2159!
2160!--               Output of averaged data
2161                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2162                                        REAL( average_count_surf, KIND=wp )
2163                  surfaces%var_av(:,n_out) = 0.0_wp
2164
2165               ENDIF
2166
2167            CASE ( 'qcsws' )
2168!
2169!--            Output of instantaneous data
2170               IF ( av == 0 )  THEN
2171                  CALL surface_data_output_collect( surf_def_h(0)%qcsws,       &
2172                                               surf_def_h(1)%qcsws,            &
2173                                               surf_lsm_h%qcsws,               &
2174                                               surf_usm_h%qcsws,               &
2175                                               surf_def_v(0)%qcsws,            &
2176                                               surf_lsm_v(0)%qcsws,            &
2177                                               surf_usm_v(0)%qcsws,            &
2178                                               surf_def_v(1)%qcsws,            &
2179                                               surf_lsm_v(1)%qcsws,            &
2180                                               surf_usm_v(1)%qcsws,            &
2181                                               surf_def_v(2)%qcsws,            &
2182                                               surf_lsm_v(2)%qcsws,            &
2183                                               surf_usm_v(2)%qcsws,            &
2184                                               surf_def_v(3)%qcsws,            &
2185                                               surf_lsm_v(3)%qcsws,            &
2186                                               surf_usm_v(3)%qcsws )
2187               ELSE
2188!
2189!--               Output of averaged data
2190                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2191                                        REAL( average_count_surf, KIND=wp )
2192                  surfaces%var_av(:,n_out) = 0.0_wp
2193
2194               ENDIF
2195
2196            CASE ( 'ncsws' )
2197!
2198!--            Output of instantaneous data
2199               IF ( av == 0 )  THEN
2200                  CALL surface_data_output_collect( surf_def_h(0)%ncsws,       &
2201                                               surf_def_h(1)%ncsws,            &
2202                                               surf_lsm_h%ncsws,               &
2203                                               surf_usm_h%ncsws,               &
2204                                               surf_def_v(0)%ncsws,            &
2205                                               surf_lsm_v(0)%ncsws,            &
2206                                               surf_usm_v(0)%ncsws,            &
2207                                               surf_def_v(1)%ncsws,            &
2208                                               surf_lsm_v(1)%ncsws,            &
2209                                               surf_usm_v(1)%ncsws,            &
2210                                               surf_def_v(2)%ncsws,            &
2211                                               surf_lsm_v(2)%ncsws,            &
2212                                               surf_usm_v(2)%ncsws,            &
2213                                               surf_def_v(3)%ncsws,            &
2214                                               surf_lsm_v(3)%ncsws,            &
2215                                               surf_usm_v(3)%ncsws )
2216               ELSE
2217!
2218!--               Output of averaged data
2219                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2220                                        REAL( average_count_surf, KIND=wp )
2221                  surfaces%var_av(:,n_out) = 0.0_wp
2222
2223               ENDIF
2224
2225
2226            CASE ( 'qisws' )
2227!
2228!--            Output of instantaneous data
2229               IF ( av == 0 )  THEN
2230                  CALL surface_data_output_collect( surf_def_h(0)%qisws,       &
2231                                               surf_def_h(1)%qisws,            &
2232                                               surf_lsm_h%qisws,               &
2233                                               surf_usm_h%qisws,               &
2234                                               surf_def_v(0)%qisws,            &
2235                                               surf_lsm_v(0)%qisws,            &
2236                                               surf_usm_v(0)%qisws,            &
2237                                               surf_def_v(1)%qisws,            &
2238                                               surf_lsm_v(1)%qisws,            &
2239                                               surf_usm_v(1)%qisws,            &
2240                                               surf_def_v(2)%qisws,            &
2241                                               surf_lsm_v(2)%qisws,            &
2242                                               surf_usm_v(2)%qisws,            &
2243                                               surf_def_v(3)%qisws,            &
2244                                               surf_lsm_v(3)%qisws,            &
2245                                               surf_usm_v(3)%qisws )
2246               ELSE
2247!
2248!--               Output of averaged data
2249                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2250                                        REAL( average_count_surf, KIND=wp )
2251                  surfaces%var_av(:,n_out) = 0.0_wp
2252
2253               ENDIF
2254
2255            CASE ( 'nisws' )
2256!
2257!--            Output of instantaneous data
2258               IF ( av == 0 )  THEN
2259                  CALL surface_data_output_collect( surf_def_h(0)%nisws,       &
2260                                               surf_def_h(1)%nisws,            &
2261                                               surf_lsm_h%nisws,               &
2262                                               surf_usm_h%nisws,               &
2263                                               surf_def_v(0)%nisws,            &
2264                                               surf_lsm_v(0)%nisws,            &
2265                                               surf_usm_v(0)%nisws,            &
2266                                               surf_def_v(1)%nisws,            &
2267                                               surf_lsm_v(1)%nisws,            &
2268                                               surf_usm_v(1)%nisws,            &
2269                                               surf_def_v(2)%nisws,            &
2270                                               surf_lsm_v(2)%nisws,            &
2271                                               surf_usm_v(2)%nisws,            &
2272                                               surf_def_v(3)%nisws,            &
2273                                               surf_lsm_v(3)%nisws,            &
2274                                               surf_usm_v(3)%nisws )
2275               ELSE
2276!
2277!--               Output of averaged data
2278                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2279                                        REAL( average_count_surf, KIND=wp )
2280                  surfaces%var_av(:,n_out) = 0.0_wp
2281
2282               ENDIF
2283
2284            CASE ( 'qrsws' )
2285!
2286!--            Output of instantaneous data
2287               IF ( av == 0 )  THEN
2288                  CALL surface_data_output_collect( surf_def_h(0)%qrsws,       &
2289                                               surf_def_h(1)%qrsws,            &
2290                                               surf_lsm_h%qrsws,               &
2291                                               surf_usm_h%qrsws,               &
2292                                               surf_def_v(0)%qrsws,            &
2293                                               surf_lsm_v(0)%qrsws,            &
2294                                               surf_usm_v(0)%qrsws,            &
2295                                               surf_def_v(1)%qrsws,            &
2296                                               surf_lsm_v(1)%qrsws,            &
2297                                               surf_usm_v(1)%qrsws,            &
2298                                               surf_def_v(2)%qrsws,            &
2299                                               surf_lsm_v(2)%qrsws,            &
2300                                               surf_usm_v(2)%qrsws,            &
2301                                               surf_def_v(3)%qrsws,            &
2302                                               surf_lsm_v(3)%qrsws,            &
2303                                               surf_usm_v(3)%qrsws )
2304               ELSE
2305!
2306!--               Output of averaged data
2307                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2308                                        REAL( average_count_surf, KIND=wp )
2309                  surfaces%var_av(:,n_out) = 0.0_wp
2310
2311               ENDIF
2312
2313            CASE ( 'nrsws' )
2314!
2315!--            Output of instantaneous data
2316               IF ( av == 0 )  THEN
2317                  CALL surface_data_output_collect( surf_def_h(0)%nrsws,       &
2318                                               surf_def_h(1)%nrsws,            &
2319                                               surf_lsm_h%nrsws,               &
2320                                               surf_usm_h%nrsws,               &
2321                                               surf_def_v(0)%nrsws,            &
2322                                               surf_lsm_v(0)%nrsws,            &
2323                                               surf_usm_v(0)%nrsws,            &
2324                                               surf_def_v(1)%nrsws,            &
2325                                               surf_lsm_v(1)%nrsws,            &
2326                                               surf_usm_v(1)%nrsws,            &
2327                                               surf_def_v(2)%nrsws,            &
2328                                               surf_lsm_v(2)%nrsws,            &
2329                                               surf_usm_v(2)%nrsws,            &
2330                                               surf_def_v(3)%nrsws,            &
2331                                               surf_lsm_v(3)%nrsws,            &
2332                                               surf_usm_v(3)%nrsws )
2333               ELSE
2334!
2335!--               Output of averaged data
2336                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2337                                        REAL( average_count_surf, KIND=wp )
2338                  surfaces%var_av(:,n_out) = 0.0_wp
2339
2340               ENDIF
2341
2342            CASE ( 'sasws' )
2343!
2344!--            Output of instantaneous data
2345               IF ( av == 0 )  THEN
2346                  CALL surface_data_output_collect( surf_def_h(0)%sasws,       &
2347                                               surf_def_h(1)%sasws,            &
2348                                               surf_lsm_h%sasws,               &
2349                                               surf_usm_h%sasws,               &
2350                                               surf_def_v(0)%sasws,            &
2351                                               surf_lsm_v(0)%sasws,            &
2352                                               surf_usm_v(0)%sasws,            &
2353                                               surf_def_v(1)%sasws,            &
2354                                               surf_lsm_v(1)%sasws,            &
2355                                               surf_usm_v(1)%sasws,            &
2356                                               surf_def_v(2)%sasws,            &
2357                                               surf_lsm_v(2)%sasws,            &
2358                                               surf_usm_v(2)%sasws,            &
2359                                               surf_def_v(3)%sasws,            &
2360                                               surf_lsm_v(3)%sasws,            &
2361                                               surf_usm_v(3)%sasws )
2362               ELSE
2363!
2364!--               Output of averaged data
2365                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2366                                        REAL( average_count_surf, KIND=wp )
2367                  surfaces%var_av(:,n_out) = 0.0_wp
2368
2369               ENDIF
2370
2371            CASE ( 'q_surface' )
2372!
2373!--            Output of instantaneous data
2374               IF ( av == 0 )  THEN
2375                  CALL surface_data_output_collect( surf_def_h(0)%q_surface,   &
2376                                               surf_def_h(1)%q_surface,        &
2377                                               surf_lsm_h%q_surface,           &
2378                                               surf_usm_h%q_surface,           &
2379                                               surf_def_v(0)%q_surface,        &
2380                                               surf_lsm_v(0)%q_surface,        &
2381                                               surf_usm_v(0)%q_surface,        &
2382                                               surf_def_v(1)%q_surface,        &
2383                                               surf_lsm_v(1)%q_surface,        &
2384                                               surf_usm_v(1)%q_surface,        &
2385                                               surf_def_v(2)%q_surface,        &
2386                                               surf_lsm_v(2)%q_surface,        &
2387                                               surf_usm_v(2)%q_surface,        &
2388                                               surf_def_v(3)%q_surface,        &
2389                                               surf_lsm_v(3)%q_surface,        &
2390                                               surf_usm_v(3)%q_surface )
2391               ELSE
2392!
2393!--               Output of averaged data
2394                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2395                                        REAL( average_count_surf, KIND=wp )
2396                  surfaces%var_av(:,n_out) = 0.0_wp
2397
2398               ENDIF
2399
2400            CASE ( 'theta_surface' )
2401!
2402!--            Output of instantaneous data
2403               IF ( av == 0 )  THEN
2404                  CALL surface_data_output_collect( surf_def_h(0)%pt_surface,  &
2405                                               surf_def_h(1)%pt_surface,       &
2406                                               surf_lsm_h%pt_surface,          &
2407                                               surf_usm_h%pt_surface,          &
2408                                               surf_def_v(0)%pt_surface,       &
2409                                               surf_lsm_v(0)%pt_surface,       &
2410                                               surf_usm_v(0)%pt_surface,       &
2411                                               surf_def_v(1)%pt_surface,       &
2412                                               surf_lsm_v(1)%pt_surface,       &
2413                                               surf_usm_v(1)%pt_surface,       &
2414                                               surf_def_v(2)%pt_surface,       &
2415                                               surf_lsm_v(2)%pt_surface,       &
2416                                               surf_usm_v(2)%pt_surface,       &
2417                                               surf_def_v(3)%pt_surface,       &
2418                                               surf_lsm_v(3)%pt_surface,       &
2419                                               surf_usm_v(3)%pt_surface )
2420               ELSE
2421!
2422!--               Output of averaged data
2423                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2424                                        REAL( average_count_surf, KIND=wp )
2425                  surfaces%var_av(:,n_out) = 0.0_wp
2426
2427               ENDIF
2428
2429            CASE ( 'thetav_surface' )
2430!
2431!--            Output of instantaneous data
2432               IF ( av == 0 )  THEN
2433                  CALL surface_data_output_collect( surf_def_h(0)%vpt_surface, &
2434                                               surf_def_h(1)%vpt_surface,      &
2435                                               surf_lsm_h%vpt_surface,         &
2436                                               surf_usm_h%vpt_surface,         &
2437                                               surf_def_v(0)%vpt_surface,      &
2438                                               surf_lsm_v(0)%vpt_surface,      &
2439                                               surf_usm_v(0)%vpt_surface,      &
2440                                               surf_def_v(1)%vpt_surface,      &
2441                                               surf_lsm_v(1)%vpt_surface,      &
2442                                               surf_usm_v(1)%vpt_surface,      &
2443                                               surf_def_v(2)%vpt_surface,      &
2444                                               surf_lsm_v(2)%vpt_surface,      &
2445                                               surf_usm_v(2)%vpt_surface,      &
2446                                               surf_def_v(3)%vpt_surface,      &
2447                                               surf_lsm_v(3)%vpt_surface,      &
2448                                               surf_usm_v(3)%vpt_surface)
2449               ELSE
2450!
2451!--               Output of averaged data
2452                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2453                                        REAL( average_count_surf, KIND=wp )
2454                  surfaces%var_av(:,n_out) = 0.0_wp
2455
2456               ENDIF
2457
2458            CASE ( 'rad_net' )
2459!
2460!--            Output of instantaneous data
2461               IF ( av == 0 )  THEN
2462                  CALL surface_data_output_collect( surf_def_h(0)%rad_net,     &
2463                                               surf_def_h(1)%rad_net,          &
2464                                               surf_lsm_h%rad_net,             &
2465                                               surf_usm_h%rad_net,             &
2466                                               surf_def_v(0)%rad_net,          &
2467                                               surf_lsm_v(0)%rad_net,          &
2468                                               surf_usm_v(0)%rad_net,          &
2469                                               surf_def_v(1)%rad_net,          &
2470                                               surf_lsm_v(1)%rad_net,          &
2471                                               surf_usm_v(1)%rad_net,          &
2472                                               surf_def_v(2)%rad_net,          &
2473                                               surf_lsm_v(2)%rad_net,          &
2474                                               surf_usm_v(2)%rad_net,          &
2475                                               surf_def_v(3)%rad_net,          &
2476                                               surf_lsm_v(3)%rad_net,          &
2477                                               surf_usm_v(3)%rad_net )
2478               ELSE
2479!
2480!--               Output of averaged data
2481                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2482                                        REAL( average_count_surf, KIND=wp )
2483                  surfaces%var_av(:,n_out) = 0.0_wp
2484
2485               ENDIF
2486
2487            CASE ( 'rad_lw_in' )
2488!
2489!--            Output of instantaneous data
2490               IF ( av == 0 )  THEN
2491                  CALL surface_data_output_collect( surf_def_h(0)%rad_lw_in,   &
2492                                               surf_def_h(1)%rad_lw_in,        &
2493                                               surf_lsm_h%rad_lw_in,           &
2494                                               surf_usm_h%rad_lw_in,           &
2495                                               surf_def_v(0)%rad_lw_in,        &
2496                                               surf_lsm_v(0)%rad_lw_in,        &
2497                                               surf_usm_v(0)%rad_lw_in,        &
2498                                               surf_def_v(1)%rad_lw_in,        &
2499                                               surf_lsm_v(1)%rad_lw_in,        &
2500                                               surf_usm_v(1)%rad_lw_in,        &
2501                                               surf_def_v(2)%rad_lw_in,        &
2502                                               surf_lsm_v(2)%rad_lw_in,        &
2503                                               surf_usm_v(2)%rad_lw_in,        &
2504                                               surf_def_v(3)%rad_lw_in,        &
2505                                               surf_lsm_v(3)%rad_lw_in,        &
2506                                               surf_usm_v(3)%rad_lw_in )
2507               ELSE
2508!
2509!--               Output of averaged data
2510                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2511                                        REAL( average_count_surf, KIND=wp )
2512                  surfaces%var_av(:,n_out) = 0.0_wp
2513
2514               ENDIF
2515
2516            CASE ( 'rad_lw_out' )
2517!
2518!--            Output of instantaneous data
2519               IF ( av == 0 )  THEN
2520                  CALL surface_data_output_collect( surf_def_h(0)%rad_lw_out,  &
2521                                               surf_def_h(1)%rad_lw_out,       &
2522                                               surf_lsm_h%rad_lw_out,          &
2523                                               surf_usm_h%rad_lw_out,          &
2524                                               surf_def_v(0)%rad_lw_out,       &
2525                                               surf_lsm_v(0)%rad_lw_out,       &
2526                                               surf_usm_v(0)%rad_lw_out,       &
2527                                               surf_def_v(1)%rad_lw_out,       &
2528                                               surf_lsm_v(1)%rad_lw_out,       &
2529                                               surf_usm_v(1)%rad_lw_out,       &
2530                                               surf_def_v(2)%rad_lw_out,       &
2531                                               surf_lsm_v(2)%rad_lw_out,       &
2532                                               surf_usm_v(2)%rad_lw_out,       &
2533                                               surf_def_v(3)%rad_lw_out,       &
2534                                               surf_lsm_v(3)%rad_lw_out,       &
2535                                               surf_usm_v(3)%rad_lw_out )
2536               ELSE
2537!
2538!--               Output of averaged data
2539                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2540                                        REAL( average_count_surf, KIND=wp )
2541                  surfaces%var_av(:,n_out) = 0.0_wp
2542
2543               ENDIF
2544
2545            CASE ( 'rad_sw_in' )
2546!
2547!--            Output of instantaneous data
2548               IF ( av == 0 )  THEN
2549                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_in,   &
2550                                               surf_def_h(1)%rad_sw_in,        &
2551                                               surf_lsm_h%rad_sw_in,           &
2552                                               surf_usm_h%rad_sw_in,           &
2553                                               surf_def_v(0)%rad_sw_in,        &
2554                                               surf_lsm_v(0)%rad_sw_in,        &
2555                                               surf_usm_v(0)%rad_sw_in,        &
2556                                               surf_def_v(1)%rad_sw_in,        &
2557                                               surf_lsm_v(1)%rad_sw_in,        &
2558                                               surf_usm_v(1)%rad_sw_in,        &
2559                                               surf_def_v(2)%rad_sw_in,        &
2560                                               surf_lsm_v(2)%rad_sw_in,        &
2561                                               surf_usm_v(2)%rad_sw_in,        &
2562                                               surf_def_v(3)%rad_sw_in,        &
2563                                               surf_lsm_v(3)%rad_sw_in,        &
2564                                               surf_usm_v(3)%rad_sw_in )
2565               ELSE
2566!
2567!--               Output of averaged data
2568                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2569                                        REAL( average_count_surf, KIND=wp )
2570                  surfaces%var_av(:,n_out) = 0.0_wp
2571
2572               ENDIF
2573
2574            CASE ( 'rad_sw_out' )
2575!
2576!--            Output of instantaneous data
2577               IF ( av == 0 )  THEN
2578                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_out,  &
2579                                               surf_def_h(1)%rad_sw_out,       &
2580                                               surf_lsm_h%rad_sw_out,          &
2581                                               surf_usm_h%rad_sw_out,          &
2582                                               surf_def_v(0)%rad_sw_out,       &
2583                                               surf_lsm_v(0)%rad_sw_out,       &
2584                                               surf_usm_v(0)%rad_sw_out,       &
2585                                               surf_def_v(1)%rad_sw_out,       &
2586                                               surf_lsm_v(1)%rad_sw_out,       &
2587                                               surf_usm_v(1)%rad_sw_out,       &
2588                                               surf_def_v(2)%rad_sw_out,       &
2589                                               surf_lsm_v(2)%rad_sw_out,       &
2590                                               surf_usm_v(2)%rad_sw_out,       &
2591                                               surf_def_v(3)%rad_sw_out,       &
2592                                               surf_lsm_v(3)%rad_sw_out,       &
2593                                               surf_usm_v(3)%rad_sw_out )
2594               ELSE
2595!
2596!--               Output of averaged data
2597                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2598                                        REAL( average_count_surf, KIND=wp )
2599                  surfaces%var_av(:,n_out) = 0.0_wp
2600
2601               ENDIF
2602
2603            CASE ( 'ghf' )
2604!
2605!--            Output of instantaneous data
2606               IF ( av == 0 )  THEN
2607!
2608!--               Sum up ground / wall heat flux. Note, for urban surfaces the
2609!--               wall heat flux is aggregated from the different green, window
2610!--               and wall tiles.
2611                  DO  m = 1, surf_usm_h%ns
2612                     surf_usm_h%ghf(m) = surf_usm_h%frac(m,ind_veg_wall)       &
2613                                       * surf_usm_h%wghf_eb(m) +               &
2614                                         surf_usm_h%frac(m,ind_pav_green)      &
2615                                       * surf_usm_h%wghf_eb_green(m) +         &
2616                                         surf_usm_h%frac(m,ind_wat_win)        &
2617                                       * surf_usm_h%wghf_eb_window(m)
2618                  ENDDO
2619                  DO  l = 0, 3
2620                     DO  m = 1, surf_usm_v(l)%ns
2621                        surf_usm_v(l)%ghf(m) = surf_usm_v(l)%frac(m,ind_veg_wall) &
2622                                             * surf_usm_v(l)%wghf_eb(m) +          &
2623                                               surf_usm_v(l)%frac(m,ind_pav_green) &
2624                                             * surf_usm_v(l)%wghf_eb_green(m) +    &
2625                                               surf_usm_v(l)%frac(m,ind_wat_win)   &
2626                                             * surf_usm_v(l)%wghf_eb_window(m)
2627                     ENDDO
2628                  ENDDO
2629
2630                  CALL surface_data_output_collect( surf_def_h(0)%ghf,         &
2631                                               surf_def_h(1)%ghf,              &
2632                                               surf_lsm_h%ghf,                 &
2633                                               surf_usm_h%ghf,                 &
2634                                               surf_def_v(0)%ghf,              &
2635                                               surf_lsm_v(0)%ghf,              &
2636                                               surf_usm_v(0)%ghf,              &
2637                                               surf_def_v(1)%ghf,              &
2638                                               surf_lsm_v(1)%ghf,              &
2639                                               surf_usm_v(1)%ghf,              &
2640                                               surf_def_v(2)%ghf,              &
2641                                               surf_lsm_v(2)%ghf,              &
2642                                               surf_usm_v(2)%ghf,              &
2643                                               surf_def_v(3)%ghf,              &
2644                                               surf_lsm_v(3)%ghf,              &
2645                                               surf_usm_v(3)%ghf )
2646               ELSE
2647!
2648!--               Output of averaged data
2649                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2650                                        REAL( average_count_surf, KIND=wp )
2651                  surfaces%var_av(:,n_out) = 0.0_wp
2652
2653               ENDIF
2654
2655            CASE ( 'r_a' )
2656!
2657!--            Output of instantaneous data
2658               IF ( av == 0 )  THEN
2659                  CALL surface_data_output_collect( surf_def_h(0)%r_a,         &
2660                                               surf_def_h(1)%r_a,              &
2661                                               surf_lsm_h%r_a,                 &
2662                                               surf_usm_h%r_a,                 &
2663                                               surf_def_v(0)%r_a,              &
2664                                               surf_lsm_v(0)%r_a,              &
2665                                               surf_usm_v(0)%r_a,              &
2666                                               surf_def_v(1)%r_a,              &
2667                                               surf_lsm_v(1)%r_a,              &
2668                                               surf_usm_v(1)%r_a,              &
2669                                               surf_def_v(2)%r_a,              &
2670                                               surf_lsm_v(2)%r_a,              &
2671                                               surf_usm_v(2)%r_a,              &
2672                                               surf_def_v(3)%r_a,              &
2673                                               surf_lsm_v(3)%r_a,              &
2674                                               surf_usm_v(3)%r_a )
2675                                                                      ELSE
2676!
2677!--               Output of averaged data
2678                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2679                                        REAL( average_count_surf, KIND=wp )
2680                  surfaces%var_av(:,n_out) = 0.0_wp
2681
2682               ENDIF
2683
2684            CASE ( 'r_soil' )
2685!
2686!--            Output of instantaneous data
2687               IF ( av == 0 )  THEN
2688                  CALL surface_data_output_collect( surf_def_h(0)%r_soil,      &
2689                                               surf_def_h(1)%r_soil,           &
2690                                               surf_lsm_h%r_soil,              &
2691                                               surf_usm_h%r_soil,              &
2692                                               surf_def_v(0)%r_soil,           &
2693                                               surf_lsm_v(0)%r_soil,           &
2694                                               surf_usm_v(0)%r_soil,           &
2695                                               surf_def_v(1)%r_soil,           &
2696                                               surf_lsm_v(1)%r_soil,           &
2697                                               surf_usm_v(1)%r_soil,           &
2698                                               surf_def_v(2)%r_soil,           &
2699                                               surf_lsm_v(2)%r_soil,           &
2700                                               surf_usm_v(2)%r_soil,           &
2701                                               surf_def_v(3)%r_soil,           &
2702                                               surf_lsm_v(3)%r_soil,           &
2703                                               surf_usm_v(3)%r_soil )
2704               ELSE
2705!
2706!--               Output of averaged data
2707                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2708                                        REAL( average_count_surf, KIND=wp )
2709                  surfaces%var_av(:,n_out) = 0.0_wp
2710
2711               ENDIF
2712
2713            CASE ( 'r_canopy' )
2714!
2715!--            Output of instantaneous data
2716               IF ( av == 0 )  THEN
2717                  CALL surface_data_output_collect( surf_def_h(0)%r_canopy,    &
2718                                               surf_def_h(1)%r_canopy,         &
2719                                               surf_lsm_h%r_canopy,            &
2720                                               surf_usm_h%r_canopy,            &
2721                                               surf_def_v(0)%r_canopy,         &
2722                                               surf_lsm_v(0)%r_canopy,         &
2723                                               surf_usm_v(0)%r_canopy,         &
2724                                               surf_def_v(1)%r_canopy,         &
2725                                               surf_lsm_v(1)%r_canopy,         &
2726                                               surf_usm_v(1)%r_canopy,         &
2727                                               surf_def_v(2)%r_canopy,         &
2728                                               surf_lsm_v(2)%r_canopy,         &
2729                                               surf_usm_v(2)%r_canopy,         &
2730                                               surf_def_v(3)%r_canopy,         &
2731                                               surf_lsm_v(3)%r_canopy,         &
2732                                               surf_usm_v(3)%r_canopy )
2733               ELSE
2734!
2735!--               Output of averaged data
2736                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2737                                        REAL( average_count_surf, KIND=wp )
2738                  surfaces%var_av(:,n_out) = 0.0_wp
2739
2740               ENDIF
2741
2742            CASE ( 'r_s' )
2743!
2744!--            Output of instantaneous data
2745               IF ( av == 0 )  THEN
2746                  CALL surface_data_output_collect( surf_def_h(0)%r_s,         &
2747                                               surf_def_h(1)%r_s,              &
2748                                               surf_lsm_h%r_s,                 &
2749                                               surf_usm_h%r_s,                 &
2750                                               surf_def_v(0)%r_s,              &
2751                                               surf_lsm_v(0)%r_s,              &
2752                                               surf_usm_v(0)%r_s,              &
2753                                               surf_def_v(1)%r_s,              &
2754                                               surf_lsm_v(1)%r_s,              &
2755                                               surf_usm_v(1)%r_s,              &
2756                                               surf_def_v(2)%r_s,              &
2757                                               surf_lsm_v(2)%r_s,              &
2758                                               surf_usm_v(2)%r_s,              &
2759                                               surf_def_v(3)%r_s,              &
2760                                               surf_lsm_v(3)%r_s,              &
2761                                               surf_usm_v(3)%r_s )
2762               ELSE
2763!
2764!--               Output of averaged data
2765                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2766                                        REAL( average_count_surf, KIND=wp )
2767                  surfaces%var_av(:,n_out) = 0.0_wp
2768
2769               ENDIF
2770
2771            CASE ( 'rad_sw_dir' )
2772!
2773!--            Output of instantaneous data
2774               IF ( av == 0 )  THEN
2775                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_dir,  &
2776                                               surf_def_h(1)%rad_sw_dir,       &
2777                                               surf_lsm_h%rad_sw_dir,          &
2778                                               surf_usm_h%rad_sw_dir,          &
2779                                               surf_def_v(0)%rad_sw_dir,       &
2780                                               surf_lsm_v(0)%rad_sw_dir,       &
2781                                               surf_usm_v(0)%rad_sw_dir,       &
2782                                               surf_def_v(1)%rad_sw_dir,       &
2783                                               surf_lsm_v(1)%rad_sw_dir,       &
2784                                               surf_usm_v(1)%rad_sw_dir,       &
2785                                               surf_def_v(2)%rad_sw_dir,       &
2786                                               surf_lsm_v(2)%rad_sw_dir,       &
2787                                               surf_usm_v(2)%rad_sw_dir,       &
2788                                               surf_def_v(3)%rad_sw_dir,       &
2789                                               surf_lsm_v(3)%rad_sw_dir,       &
2790                                               surf_usm_v(3)%rad_sw_dir )
2791               ELSE
2792!
2793!--               Output of averaged data
2794                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2795                                        REAL( average_count_surf, KIND=wp )
2796                  surfaces%var_av(:,n_out) = 0.0_wp
2797
2798               ENDIF
2799
2800            CASE ( 'rad_sw_dif' )
2801!
2802!--            Output of instantaneous data
2803               IF ( av == 0 )  THEN
2804                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_dif,  &
2805                                               surf_def_h(1)%rad_sw_dif,       &
2806                                               surf_lsm_h%rad_sw_dif,          &
2807                                               surf_usm_h%rad_sw_dif,          &
2808                                               surf_def_v(0)%rad_sw_dif,       &
2809                                               surf_lsm_v(0)%rad_sw_dif,       &
2810                                               surf_usm_v(0)%rad_sw_dif,       &
2811                                               surf_def_v(1)%rad_sw_dif,       &
2812                                               surf_lsm_v(1)%rad_sw_dif,       &
2813                                               surf_usm_v(1)%rad_sw_dif,       &
2814                                               surf_def_v(2)%rad_sw_dif,       &
2815                                               surf_lsm_v(2)%rad_sw_dif,       &
2816                                               surf_usm_v(2)%rad_sw_dif,       &
2817                                               surf_def_v(3)%rad_sw_dif,       &
2818                                               surf_lsm_v(3)%rad_sw_dif,       &
2819                                               surf_usm_v(3)%rad_sw_dif )
2820               ELSE
2821!
2822!--               Output of averaged data
2823                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2824                                        REAL( average_count_surf, KIND=wp )
2825                  surfaces%var_av(:,n_out) = 0.0_wp
2826
2827               ENDIF
2828
2829            CASE ( 'rad_sw_ref' )
2830!
2831!--            Output of instantaneous data
2832               IF ( av == 0 )  THEN
2833                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_ref,  &
2834                                               surf_def_h(1)%rad_sw_ref,       &
2835                                               surf_lsm_h%rad_sw_ref,          &
2836                                               surf_usm_h%rad_sw_ref,          &
2837                                               surf_def_v(0)%rad_sw_ref,       &
2838                                               surf_lsm_v(0)%rad_sw_ref,       &
2839                                               surf_usm_v(0)%rad_sw_ref,       &
2840                                               surf_def_v(1)%rad_sw_ref,       &
2841                                               surf_lsm_v(1)%rad_sw_ref,       &
2842                                               surf_usm_v(1)%rad_sw_ref,       &
2843                                               surf_def_v(2)%rad_sw_ref,       &
2844                                               surf_lsm_v(2)%rad_sw_ref,       &
2845                                               surf_usm_v(2)%rad_sw_ref,       &
2846                                               surf_def_v(3)%rad_sw_ref,       &
2847                                               surf_lsm_v(3)%rad_sw_ref,       &
2848                                               surf_usm_v(3)%rad_sw_ref )
2849               ELSE
2850!
2851!--               Output of averaged data
2852                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2853                                        REAL( average_count_surf, KIND=wp )
2854                  surfaces%var_av(:,n_out) = 0.0_wp
2855
2856               ENDIF
2857
2858            CASE ( 'rad_sw_res' )
2859!
2860!--            Output of instantaneous data
2861               IF ( av == 0 )  THEN
2862                  CALL surface_data_output_collect( surf_def_h(0)%rad_sw_res,  &
2863                                               surf_def_h(1)%rad_sw_res,       &
2864                                               surf_lsm_h%rad_sw_res,          &
2865                                               surf_usm_h%rad_sw_res,          &
2866                                               surf_def_v(0)%rad_sw_res,       &
2867                                               surf_lsm_v(0)%rad_sw_res,       &
2868                                               surf_usm_v(0)%rad_sw_res,       &
2869                                               surf_def_v(1)%rad_sw_res,       &
2870                                               surf_lsm_v(1)%rad_sw_res,       &
2871                                               surf_usm_v(1)%rad_sw_res,       &
2872                                               surf_def_v(2)%rad_sw_res,       &
2873                                               surf_lsm_v(2)%rad_sw_res,       &
2874                                               surf_usm_v(2)%rad_sw_res,       &
2875                                               surf_def_v(3)%rad_sw_res,       &
2876                                               surf_lsm_v(3)%rad_sw_res,       &
2877                                               surf_usm_v(3)%rad_sw_res )
2878               ELSE
2879!
2880!--               Output of averaged data
2881                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2882                                        REAL( average_count_surf, KIND=wp )
2883                  surfaces%var_av(:,n_out) = 0.0_wp
2884
2885               ENDIF
2886
2887            CASE ( 'rad_lw_dif' )
2888!
2889!--            Output of instantaneous data
2890               IF ( av == 0 )  THEN
2891                  CALL surface_data_output_collect( surf_def_h(0)%rad_lw_dif,  &
2892                                               surf_def_h(1)%rad_lw_dif,       &
2893                                               surf_lsm_h%rad_lw_dif,          &
2894                                               surf_usm_h%rad_lw_dif,          &
2895                                               surf_def_v(0)%rad_lw_dif,       &
2896                                               surf_lsm_v(0)%rad_lw_dif,       &
2897                                               surf_usm_v(0)%rad_lw_dif,       &
2898                                               surf_def_v(1)%rad_lw_dif,       &
2899                                               surf_lsm_v(1)%rad_lw_dif,       &
2900                                               surf_usm_v(1)%rad_lw_dif,       &
2901                                               surf_def_v(2)%rad_lw_dif,       &
2902                                               surf_lsm_v(2)%rad_lw_dif,       &
2903                                               surf_usm_v(2)%rad_lw_dif,       &
2904                                               surf_def_v(3)%rad_lw_dif,       &
2905                                               surf_lsm_v(3)%rad_lw_dif,       &
2906                                               surf_usm_v(3)%rad_lw_dif )
2907               ELSE
2908!
2909!--               Output of averaged data
2910                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2911                                        REAL( average_count_surf, KIND=wp )
2912                  surfaces%var_av(:,n_out) = 0.0_wp
2913
2914               ENDIF
2915
2916            CASE ( 'rad_lw_ref' )
2917!
2918!--            Output of instantaneous data
2919               IF ( av == 0 )  THEN
2920                  CALL surface_data_output_collect( surf_def_h(0)%rad_lw_ref,  &
2921                                               surf_def_h(1)%rad_lw_ref,       &
2922                                               surf_lsm_h%rad_lw_ref,          &
2923                                               surf_usm_h%rad_lw_ref,          &
2924                                               surf_def_v(0)%rad_lw_ref,       &
2925                                               surf_lsm_v(0)%rad_lw_ref,       &
2926                                               surf_usm_v(0)%rad_lw_ref,       &
2927                                               surf_def_v(1)%rad_lw_ref,       &
2928                                               surf_lsm_v(1)%rad_lw_ref,       &
2929                                               surf_usm_v(1)%rad_lw_ref,       &
2930                                               surf_def_v(2)%rad_lw_ref,       &
2931                                               surf_lsm_v(2)%rad_lw_ref,       &
2932                                               surf_usm_v(2)%rad_lw_ref,       &
2933                                               surf_def_v(3)%rad_lw_ref,       &
2934                                               surf_lsm_v(3)%rad_lw_ref,       &
2935                                               surf_usm_v(3)%rad_lw_ref )
2936               ELSE
2937!
2938!--               Output of averaged data
2939                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2940                                        REAL( average_count_surf, KIND=wp )
2941                  surfaces%var_av(:,n_out) = 0.0_wp
2942
2943               ENDIF
2944
2945            CASE ( 'rad_lw_res' )
2946!
2947!--            Output of instantaneous data
2948               IF ( av == 0 )  THEN
2949                  CALL surface_data_output_collect( surf_def_h(0)%rad_lw_res,  &
2950                                               surf_def_h(1)%rad_lw_res,       &
2951                                               surf_lsm_h%rad_lw_res,          &
2952                                               surf_usm_h%rad_lw_res,          &
2953                                               surf_def_v(0)%rad_lw_res,       &
2954                                               surf_lsm_v(0)%rad_lw_res,       &
2955                                               surf_usm_v(0)%rad_lw_res,       &
2956                                               surf_def_v(1)%rad_lw_res,       &
2957                                               surf_lsm_v(1)%rad_lw_res,       &
2958                                               surf_usm_v(1)%rad_lw_res,       &
2959                                               surf_def_v(2)%rad_lw_res,       &
2960                                               surf_lsm_v(2)%rad_lw_res,       &
2961                                               surf_usm_v(2)%rad_lw_res,       &
2962                                               surf_def_v(3)%rad_lw_res,       &
2963                                               surf_lsm_v(3)%rad_lw_res,       &
2964                                               surf_usm_v(3)%rad_lw_res )
2965               ELSE
2966!
2967!--               Output of averaged data
2968                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2969                                        REAL( average_count_surf, KIND=wp )
2970                  surfaces%var_av(:,n_out) = 0.0_wp
2971
2972               ENDIF
2973
2974            CASE ( 'uvw1' )
2975!
2976!--            Output of instantaneous data
2977               IF ( av == 0 )  THEN
2978                  CALL surface_data_output_collect( surf_def_h(0)%uvw_abs,     &
2979                                               surf_def_h(1)%uvw_abs,          &
2980                                               surf_lsm_h%uvw_abs,             &
2981                                               surf_usm_h%uvw_abs,             &
2982                                               surf_def_v(0)%uvw_abs,          &
2983                                               surf_lsm_v(0)%uvw_abs,          &
2984                                               surf_usm_v(0)%uvw_abs,          &
2985                                               surf_def_v(1)%uvw_abs,          &
2986                                               surf_lsm_v(1)%uvw_abs,          &
2987                                               surf_usm_v(1)%uvw_abs,          &
2988                                               surf_def_v(2)%uvw_abs,          &
2989                                               surf_lsm_v(2)%uvw_abs,          &
2990                                               surf_usm_v(2)%uvw_abs,          &
2991                                               surf_def_v(3)%uvw_abs,          &
2992                                               surf_lsm_v(3)%uvw_abs,          &
2993                                               surf_usm_v(3)%uvw_abs )
2994               ELSE
2995!
2996!--               Output of averaged data
2997                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2998                                        REAL( average_count_surf, KIND=wp )
2999                  surfaces%var_av(:,n_out) = 0.0_wp
3000
3001               ENDIF
3002!
3003!--         Waste heat from indoor model
3004            CASE ( 'waste_heat' )
3005!
3006!--            Output of instantaneous data
3007               IF ( av == 0 )  THEN
3008                  CALL surface_data_output_collect( surf_def_h(0)%waste_heat,  &
3009                                               surf_def_h(1)%waste_heat,       &
3010                                               surf_lsm_h%waste_heat,          &
3011                                               surf_usm_h%waste_heat,          &
3012                                               surf_def_v(0)%waste_heat,       &
3013                                               surf_lsm_v(0)%waste_heat,       &
3014                                               surf_usm_v(0)%waste_heat,       &
3015                                               surf_def_v(1)%waste_heat,       &
3016                                               surf_lsm_v(1)%waste_heat,       &
3017                                               surf_usm_v(1)%waste_heat,       &
3018                                               surf_def_v(2)%waste_heat,       &
3019                                               surf_lsm_v(2)%waste_heat,       &
3020                                               surf_usm_v(2)%waste_heat,       &
3021                                               surf_def_v(3)%waste_heat,       &
3022                                               surf_lsm_v(3)%waste_heat,       &
3023                                               surf_usm_v(3)%waste_heat )
3024               ELSE
3025!
3026!--               Output of averaged data
3027                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
3028                                        REAL( average_count_surf, KIND=wp )
3029                  surfaces%var_av(:,n_out) = 0.0_wp
3030
3031               ENDIF
3032!
3033!--         Innermost building wall flux from indoor model
3034            CASE ( 'im_hf' )
3035!
3036!--            Output of instantaneous data
3037               IF ( av == 0 )  THEN
3038                  CALL surface_data_output_collect( surf_def_h(0)%iwghf_eb,    &
3039                                               surf_def_h(1)%iwghf_eb,         &
3040                                               surf_lsm_h%iwghf_eb,            &
3041                                               surf_usm_h%iwghf_eb,            &
3042                                               surf_def_v(0)%iwghf_eb,         &
3043                                               surf_lsm_v(0)%iwghf_eb,         &
3044                                               surf_usm_v(0)%iwghf_eb,         &
3045                                               surf_def_v(1)%iwghf_eb,         &
3046                                               surf_lsm_v(1)%iwghf_eb,         &
3047                                               surf_usm_v(1)%iwghf_eb,         &
3048                                               surf_def_v(2)%iwghf_eb,         &
3049                                               surf_lsm_v(2)%iwghf_eb,         &
3050                                               surf_usm_v(2)%iwghf_eb,         &
3051                                               surf_def_v(3)%iwghf_eb,         &
3052                                               surf_lsm_v(3)%iwghf_eb,         &
3053                                               surf_usm_v(3)%iwghf_eb )
3054               ELSE
3055!
3056!--               Output of averaged data
3057                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
3058                                        REAL( average_count_surf, KIND=wp )
3059                  surfaces%var_av(:,n_out) = 0.0_wp
3060
3061               ENDIF
3062!
3063!--            Add further variables:
3064!--            'css', 'cssws', 'qsws_liq', 'qsws_soil', 'qsws_veg'
3065
3066         END SELECT
3067!
3068!--      Write to binary file:
3069!--      - surfaces%points ( 3, 1-npoints )
3070!--      - surfaces%polygons ( 5, 1-ns )
3071!--      - surfaces%var_out ( 1-ns, time )
3072!--      - Dimension: 1-nsurfaces, 1-npoints - can be ordered consecutively
3073!--      - Distinguish between average and non-average data
3074         IF ( to_vtk )  THEN
3075            DO  i = 0, io_blocks-1
3076               IF ( i == io_group )  THEN
3077                  WRITE ( 25+av )  LEN_TRIM( 'time' )
3078                  WRITE ( 25+av )  'time'
3079                  WRITE ( 25+av )  time_since_reference_point
3080                  WRITE ( 25+av )  LEN_TRIM( trimvar )
3081                  WRITE ( 25+av )  TRIM( trimvar )
3082                  WRITE ( 25+av )  surfaces%var_out
3083               ENDIF
3084#if defined( __parallel )
3085               CALL MPI_BARRIER( comm2d, ierr )
3086#endif
3087            ENDDO
3088         ENDIF
3089
3090         IF ( to_netcdf )  THEN
3091#if defined( __netcdf4_parallel )
3092!
3093!--         Write output array to file
3094            nc_stat = NF90_PUT_VAR( id_set_surf(av), id_var_dosurf(av,n_out),  &
3095                                    surfaces%var_out,                          &
3096                                    start = (/ surfaces%s(1),                  &
3097                                               dosurf_time_count(av) /),       &
3098                                    count = (/ surfaces%ns, 1 /) )
3099            CALL netcdf_handle_error( 'surface_data_output', 6667 )
3100#endif
3101         ENDIF
3102
3103      ENDDO
3104
3105!
3106!--   If averaged output was written to NetCDF file, set the counter to zero
3107      IF ( av == 1 )  average_count_surf = 0
3108
3109   END SUBROUTINE surface_data_output
3110
3111!------------------------------------------------------------------------------!
3112! Description:
3113! ------------
3114!> Routine for controlling the data averaging.
3115!------------------------------------------------------------------------------!
3116   SUBROUTINE surface_data_output_averaging
3117
3118      IMPLICIT NONE
3119
3120      CHARACTER(LEN=100) ::  trimvar !< dummy variable for current output variable
3121
3122      INTEGER(iwp) ::  l      !< running index for surface orientation
3123      INTEGER(iwp) ::  m      !< running index for surface elements
3124      INTEGER(iwp) ::  n_out  !< counter variables for surface output
3125
3126      n_out = 0
3127      DO  WHILE ( dosurf(1,n_out+1)(1:1) /= ' ' )
3128
3129         n_out   = n_out + 1
3130         trimvar = TRIM( dosurf(1,n_out) )
3131
3132         SELECT CASE ( trimvar )
3133
3134            CASE ( 'us' )
3135               CALL surface_data_output_sum_up( surf_def_h(0)%us,              &
3136                                           surf_def_h(1)%us,                   &
3137                                           surf_lsm_h%us,                      &
3138                                           surf_usm_h%us,                      &
3139                                           surf_def_v(0)%us,                   &
3140                                           surf_lsm_v(0)%us,                   &
3141                                           surf_usm_v(0)%us,                   &
3142                                           surf_def_v(1)%us,                   &
3143                                           surf_lsm_v(1)%us,                   &
3144                                           surf_usm_v(1)%us,                   &
3145                                           surf_def_v(2)%us,                   &
3146                                           surf_lsm_v(2)%us,                   &
3147                                           surf_usm_v(2)%us,                   &
3148                                           surf_def_v(3)%us,                   &
3149                                           surf_lsm_v(3)%us,                   &
3150                                           surf_usm_v(3)%us, n_out )
3151
3152            CASE ( 'ts' )
3153               CALL surface_data_output_sum_up( surf_def_h(0)%ts,              &
3154                                           surf_def_h(1)%ts,                   &
3155                                           surf_lsm_h%ts,                      &
3156                                           surf_usm_h%ts,                      &
3157                                           surf_def_v(0)%ts,                   &
3158                                           surf_lsm_v(0)%ts,                   &
3159                                           surf_usm_v(0)%ts,                   &
3160                                           surf_def_v(1)%ts,                   &
3161                                           surf_lsm_v(1)%ts,                   &
3162                                           surf_usm_v(1)%ts,                   &
3163                                           surf_def_v(2)%ts,                   &
3164                                           surf_lsm_v(2)%ts,                   &
3165                                           surf_usm_v(2)%ts,                   &
3166                                           surf_def_v(3)%ts,                   &
3167                                           surf_lsm_v(3)%ts,                   &
3168                                           surf_usm_v(3)%ts, n_out )
3169
3170            CASE ( 'qs' )
3171               CALL surface_data_output_sum_up( surf_def_h(0)%qs,              &
3172                                           surf_def_h(1)%qs,                   &
3173                                           surf_lsm_h%qs,                      &
3174                                           surf_usm_h%qs,                      &
3175                                           surf_def_v(0)%qs,                   &
3176                                           surf_lsm_v(0)%qs,                   &
3177                                           surf_usm_v(0)%qs,                   &
3178                                           surf_def_v(1)%qs,                   &
3179                                           surf_lsm_v(1)%qs,                   &
3180                                           surf_usm_v(1)%qs,                   &
3181                                           surf_def_v(2)%qs,                   &
3182                                           surf_lsm_v(2)%qs,                   &
3183                                           surf_usm_v(2)%qs,                   &
3184                                           surf_def_v(3)%qs,                   &
3185                                           surf_lsm_v(3)%qs,                   &
3186                                           surf_usm_v(3)%qs, n_out )
3187
3188            CASE ( 'ss' )
3189               CALL surface_data_output_sum_up( surf_def_h(0)%ss,              &
3190                                           surf_def_h(1)%ss,                   &
3191                                           surf_lsm_h%ss,                      &
3192                                           surf_usm_h%ss,                      &
3193                                           surf_def_v(0)%ss,                   &
3194                                           surf_lsm_v(0)%ss,                   &
3195                                           surf_usm_v(0)%ss,                   &
3196                                           surf_def_v(1)%ss,                   &
3197                                           surf_lsm_v(1)%ss,                   &
3198                                           surf_usm_v(1)%ss,                   &
3199                                           surf_def_v(2)%ss,                   &
3200                                           surf_lsm_v(2)%ss,                   &
3201                                           surf_usm_v(2)%ss,                   &
3202                                           surf_def_v(3)%ss,                   &
3203                                           surf_lsm_v(3)%ss,                   &
3204                                           surf_usm_v(3)%ss, n_out )
3205
3206            CASE ( 'qcs' )
3207               CALL surface_data_output_sum_up( surf_def_h(0)%qcs,             &
3208                                           surf_def_h(1)%qcs,                  &
3209                                           surf_lsm_h%qcs,                     &
3210                                           surf_usm_h%qcs,                     &
3211                                           surf_def_v(0)%qcs,                  &
3212                                           surf_lsm_v(0)%qcs,                  &
3213                                           surf_usm_v(0)%qcs,                  &
3214                                           surf_def_v(1)%qcs,                  &
3215                                           surf_lsm_v(1)%qcs,                  &
3216                                           surf_usm_v(1)%qcs,                  &
3217                                           surf_def_v(2)%qcs,                  &
3218                                           surf_lsm_v(2)%qcs,                  &
3219                                           surf_usm_v(2)%qcs,                  &
3220                                           surf_def_v(3)%qcs,                  &
3221                                           surf_lsm_v(3)%qcs,                  &
3222                                           surf_usm_v(3)%qcs, n_out )
3223
3224            CASE ( 'ncs' )
3225               CALL surface_data_output_sum_up( surf_def_h(0)%ncs,             &
3226                                           surf_def_h(1)%ncs,                  &
3227                                           surf_lsm_h%ncs,                     &
3228                                           surf_usm_h%ncs,                     &
3229                                           surf_def_v(0)%ncs,                  &
3230                                           surf_lsm_v(0)%ncs,                  &
3231                                           surf_usm_v(0)%ncs,                  &
3232                                           surf_def_v(1)%ncs,                  &
3233                                           surf_lsm_v(1)%ncs,                  &
3234                                           surf_usm_v(1)%ncs,                  &
3235                                           surf_def_v(2)%ncs,                  &
3236                                           surf_lsm_v(2)%ncs,                  &
3237                                           surf_usm_v(2)%ncs,                  &
3238                                           surf_def_v(3)%ncs,                  &
3239                                           surf_lsm_v(3)%ncs,                  &
3240                                           surf_usm_v(3)%ncs, n_out )
3241
3242            CASE ( 'qis' )
3243               CALL surface_data_output_sum_up( surf_def_h(0)%qis,             &
3244                                           surf_def_h(1)%qis,                  &
3245                                           surf_lsm_h%qis,                     &
3246                                           surf_usm_h%qis,                     &
3247                                           surf_def_v(0)%qis,                  &
3248                                           surf_lsm_v(0)%qis,                  &
3249                                           surf_usm_v(0)%qis,                  &
3250                                           surf_def_v(1)%qis,                  &
3251                                           surf_lsm_v(1)%qis,                  &
3252                                           surf_usm_v(1)%qis,                  &
3253                                           surf_def_v(2)%qis,                  &
3254                                           surf_lsm_v(2)%qis,                  &
3255                                           surf_usm_v(2)%qis,                  &
3256                                           surf_def_v(3)%qis,                  &
3257                                           surf_lsm_v(3)%qis,                  &
3258                                           surf_usm_v(3)%qrs, n_out )
3259
3260            CASE ( 'nis' )
3261               CALL surface_data_output_sum_up( surf_def_h(0)%nis,             &
3262                                           surf_def_h(1)%nis,                  &
3263                                           surf_lsm_h%nis,                     &
3264                                           surf_usm_h%nis,                     &
3265                                           surf_def_v(0)%nis,                  &
3266                                           surf_lsm_v(0)%nis,                  &
3267                                           surf_usm_v(0)%nis,                  &
3268                                           surf_def_v(1)%nis,                  &
3269                                           surf_lsm_v(1)%nis,                  &
3270                                           surf_usm_v(1)%nis,                  &
3271                                           surf_def_v(2)%nis,                  &
3272                                           surf_lsm_v(2)%nis,                  &
3273                                           surf_usm_v(2)%nis,                  &
3274                                           surf_def_v(3)%nis,                  &
3275                                           surf_lsm_v(3)%nis,                  &
3276                                           surf_usm_v(3)%nis, n_out )
3277
3278            CASE ( 'qrs' )
3279               CALL surface_data_output_sum_up( surf_def_h(0)%qrs,             &
3280                                           surf_def_h(1)%qrs,                  &
3281                                           surf_lsm_h%qrs,                     &
3282                                           surf_usm_h%qrs,                     &
3283                                           surf_def_v(0)%qrs,                  &
3284                                           surf_lsm_v(0)%qrs,                  &
3285                                           surf_usm_v(0)%qrs,                  &
3286                                           surf_def_v(1)%qrs,                  &
3287                                           surf_lsm_v(1)%qrs,                  &
3288                                           surf_usm_v(1)%qrs,                  &
3289                                           surf_def_v(2)%qrs,                  &
3290                                           surf_lsm_v(2)%qrs,                  &
3291                                           surf_usm_v(2)%qrs,                  &
3292                                           surf_def_v(3)%qrs,                  &
3293                                           surf_lsm_v(3)%qrs,                  &
3294                                           surf_usm_v(3)%qrs, n_out )
3295
3296            CASE ( 'nrs' )
3297               CALL surface_data_output_sum_up( surf_def_h(0)%nrs,             &
3298                                           surf_def_h(1)%nrs,                  &
3299                                           surf_lsm_h%nrs,                     &
3300                                           surf_usm_h%nrs,                     &
3301                                           surf_def_v(0)%nrs,                  &
3302                                           surf_lsm_v(0)%nrs,                  &
3303                                           surf_usm_v(0)%nrs,                  &
3304                                           surf_def_v(1)%nrs,                  &
3305                                           surf_lsm_v(1)%nrs,                  &
3306                                           surf_usm_v(1)%nrs,                  &
3307                                           surf_def_v(2)%nrs,                  &
3308                                           surf_lsm_v(2)%nrs,                  &
3309                                           surf_usm_v(2)%nrs,                  &
3310                                           surf_def_v(3)%nrs,                  &
3311                                           surf_lsm_v(3)%nrs,                  &
3312                                           surf_usm_v(3)%nrs, n_out )
3313
3314            CASE ( 'ol' )
3315               CALL surface_data_output_sum_up( surf_def_h(0)%ol,              &
3316                                           surf_def_h(1)%ol,                   &
3317                                           surf_lsm_h%ol,                      &
3318                                           surf_usm_h%ol,                      &
3319                                           surf_def_v(0)%ol,                   &
3320                                           surf_lsm_v(0)%ol,                   &
3321                                           surf_usm_v(0)%ol,                   &
3322                                           surf_def_v(1)%ol,                   &
3323                                           surf_lsm_v(1)%ol,                   &
3324                                           surf_usm_v(1)%ol,                   &
3325                                           surf_def_v(2)%ol,                   &
3326                                           surf_lsm_v(2)%ol,                   &
3327                                           surf_usm_v(2)%ol,                   &
3328                                           surf_def_v(3)%ol,                   &
3329                                           surf_lsm_v(3)%ol,                   &
3330                                           surf_usm_v(3)%ol, n_out )
3331
3332            CASE ( 'z0' )
3333               CALL surface_data_output_sum_up( surf_def_h(0)%z0,              &
3334                                           surf_def_h(1)%z0,                   &
3335                                           surf_lsm_h%z0,                      &
3336                                           surf_usm_h%z0,                      &
3337                                           surf_def_v(0)%z0,                   &
3338                                           surf_lsm_v(0)%z0,                   &
3339                                           surf_usm_v(0)%z0,                   &
3340                                           surf_def_v(1)%z0,                   &
3341                                           surf_lsm_v(1)%z0,                   &
3342                                           surf_usm_v(1)%z0,                   &
3343                                           surf_def_v(2)%z0,                   &
3344                                           surf_lsm_v(2)%z0,                   &
3345                                           surf_usm_v(2)%z0,                   &
3346                                           surf_def_v(3)%z0,                   &
3347                                           surf_lsm_v(3)%z0,                   &
3348                                           surf_usm_v(3)%z0, n_out )
3349
3350            CASE ( 'z0h' )
3351               CALL surface_data_output_sum_up( surf_def_h(0)%z0h,             &
3352                                           surf_def_h(1)%z0h,                  &
3353                                           surf_lsm_h%z0h,                     &
3354                                           surf_usm_h%z0h,                     &
3355                                           surf_def_v(0)%z0h,                  &
3356                                           surf_lsm_v(0)%z0h,                  &
3357                                           surf_usm_v(0)%z0h,                  &
3358                                           surf_def_v(1)%z0h,                  &
3359                                           surf_lsm_v(1)%z0h,                  &
3360                                           surf_usm_v(1)%z0h,                  &
3361                                           surf_def_v(2)%z0h,                  &
3362                                           surf_lsm_v(2)%z0h,                  &
3363                                           surf_usm_v(2)%z0h,                  &
3364                                           surf_def_v(3)%z0h,                  &
3365                                           surf_lsm_v(3)%z0h,                  &
3366                                           surf_usm_v(3)%z0h, n_out )
3367
3368            CASE ( 'z0q' )
3369               CALL surface_data_output_sum_up( surf_def_h(0)%z0q,             &
3370                                           surf_def_h(1)%z0q,                  &
3371                                           surf_lsm_h%z0q,                     &
3372                                           surf_usm_h%z0q,                     &
3373                                           surf_def_v(0)%z0q,                  &
3374                                           surf_lsm_v(0)%z0q,                  &
3375                                           surf_usm_v(0)%z0q,                  &
3376                                           surf_def_v(1)%z0q,                  &
3377                                           surf_lsm_v(1)%z0q,                  &
3378                                           surf_usm_v(1)%z0q,                  &
3379                                           surf_def_v(2)%z0q,                  &
3380                                           surf_lsm_v(2)%z0q,                  &
3381                                           surf_usm_v(2)%z0q,                  &
3382                                           surf_def_v(3)%z0q,                  &
3383                                           surf_lsm_v(3)%z0q,                  &
3384                                           surf_usm_v(3)%z0q, n_out )
3385
3386            CASE ( 'theta1' )
3387               CALL surface_data_output_sum_up( surf_def_h(0)%pt1,             &
3388                                           surf_def_h(1)%pt1,                  &
3389                                           surf_lsm_h%pt1,                     &
3390                                           surf_usm_h%pt1,                     &
3391                                           surf_def_v(0)%pt1,                  &
3392                                           surf_lsm_v(0)%pt1,                  &
3393                                           surf_usm_v(0)%pt1,                  &
3394                                           surf_def_v(1)%pt1,                  &
3395                                           surf_lsm_v(1)%pt1,                  &
3396                                           surf_usm_v(1)%pt1,                  &
3397                                           surf_def_v(2)%pt1,                  &
3398                                           surf_lsm_v(2)%pt1,                  &
3399                                           surf_usm_v(2)%pt1,                  &
3400                                           surf_def_v(3)%pt1,                  &
3401                                           surf_lsm_v(3)%pt1,                  &
3402                                           surf_usm_v(3)%pt1, n_out )
3403
3404            CASE ( 'qv1' )
3405               CALL surface_data_output_sum_up( surf_def_h(0)%qv1,             &
3406                                           surf_def_h(1)%qv1,                  &
3407                                           surf_lsm_h%qv1,                     &
3408                                           surf_usm_h%qv1,                     &
3409                                           surf_def_v(0)%qv1,                  &
3410                                           surf_lsm_v(0)%qv1,                  &
3411                                           surf_usm_v(0)%qv1,                  &
3412                                           surf_def_v(1)%qv1,                  &
3413                                           surf_lsm_v(1)%qv1,                  &
3414                                           surf_usm_v(1)%qv1,                  &
3415                                           surf_def_v(2)%qv1,                  &
3416                                           surf_lsm_v(2)%qv1,                  &
3417                                           surf_usm_v(2)%qv1,                  &
3418                                           surf_def_v(3)%qv1,                  &
3419                                           surf_lsm_v(3)%qv1,                  &
3420                                           surf_usm_v(3)%qv1, n_out )
3421
3422            CASE ( 'thetav1' )
3423               CALL surface_data_output_sum_up( surf_def_h(0)%vpt1,            &
3424                                           surf_def_h(1)%vpt1,                 &
3425                                           surf_lsm_h%vpt1,                    &
3426                                           surf_usm_h%vpt1,                    &
3427                                           surf_def_v(0)%vpt1,                 &
3428                                           surf_lsm_v(0)%vpt1,                 &
3429                                           surf_usm_v(0)%vpt1,                 &
3430                                           surf_def_v(1)%vpt1,                 &
3431                                           surf_lsm_v(1)%vpt1,                 &
3432                                           surf_usm_v(1)%vpt1,                 &
3433                                           surf_def_v(2)%vpt1,                 &
3434                                           surf_lsm_v(2)%vpt1,                 &
3435                                           surf_usm_v(2)%vpt1,                 &
3436                                           surf_def_v(3)%vpt1,                 &
3437                                           surf_lsm_v(3)%vpt1,                 &
3438                                           surf_usm_v(3)%vpt1, n_out )
3439
3440            CASE ( 'usws' )
3441               CALL surface_data_output_sum_up( surf_def_h(0)%usws,            &
3442                                           surf_def_h(1)%usws,                 &
3443                                           surf_lsm_h%usws,                    &
3444                                           surf_usm_h%usws,                    &
3445                                           surf_def_v(0)%usws,                 &
3446                                           surf_lsm_v(0)%usws,                 &
3447                                           surf_usm_v(0)%usws,                 &
3448                                           surf_def_v(1)%usws,                 &
3449                                           surf_lsm_v(1)%usws,                 &
3450                                           surf_usm_v(1)%usws,                 &
3451                                           surf_def_v(2)%usws,                 &
3452                                           surf_lsm_v(2)%usws,                 &
3453                                           surf_usm_v(2)%usws,                 &
3454                                           surf_def_v(3)%usws,                 &
3455                                           surf_lsm_v(3)%usws,                 &
3456                                           surf_usm_v(3)%usws, n_out,          &
3457                                           momentumflux_output_conversion )
3458
3459            CASE ( 'vsws' )
3460               CALL surface_data_output_sum_up( surf_def_h(0)%vsws,            &
3461                                           surf_def_h(1)%vsws,                 &
3462                                           surf_lsm_h%vsws,                    &
3463                                           surf_usm_h%vsws,                    &
3464                                           surf_def_v(0)%vsws,                 &
3465                                           surf_lsm_v(0)%vsws,                 &
3466                                           surf_usm_v(0)%vsws,                 &
3467                                           surf_def_v(1)%vsws,                 &
3468                                           surf_lsm_v(1)%vsws,                 &
3469                                           surf_usm_v(1)%vsws,                 &
3470                                           surf_def_v(2)%vsws,                 &
3471                                           surf_lsm_v(2)%vsws,                 &
3472                                           surf_usm_v(2)%vsws,                 &
3473                                           surf_def_v(3)%vsws,                 &
3474                                           surf_lsm_v(3)%vsws,                 &
3475                                           surf_usm_v(3)%vsws, n_out,          &
3476                                           momentumflux_output_conversion )
3477
3478            CASE ( 'shf' )
3479               CALL surface_data_output_sum_up( surf_def_h(0)%shf,             &
3480                                           surf_def_h(1)%shf,                  &
3481                                           surf_lsm_h%shf,                     &
3482                                           surf_usm_h%shf,                     &
3483                                           surf_def_v(0)%shf,                  &
3484                                           surf_lsm_v(0)%shf,                  &
3485                                           surf_usm_v(0)%shf,                  &
3486                                           surf_def_v(1)%shf,                  &
3487                                           surf_lsm_v(1)%shf,                  &
3488                                           surf_usm_v(1)%shf,                  &
3489                                           surf_def_v(2)%shf,                  &
3490                                           surf_lsm_v(2)%shf,                  &
3491                                           surf_usm_v(2)%shf,                  &
3492                                           surf_def_v(3)%shf,                  &
3493                                           surf_lsm_v(3)%shf,                  &
3494                                           surf_usm_v(3)%shf, n_out,          &
3495                                           heatflux_output_conversion )
3496
3497            CASE ( 'qsws' )
3498               CALL surface_data_output_sum_up( surf_def_h(0)%qsws,            &
3499                                           surf_def_h(1)%qsws,                 &
3500                                           surf_lsm_h%qsws,                    &
3501                                           surf_usm_h%qsws,                    &
3502                                           surf_def_v(0)%qsws,                 &
3503                                           surf_lsm_v(0)%qsws,                 &
3504                                           surf_usm_v(0)%qsws,                 &
3505                                           surf_def_v(1)%qsws,                 &
3506                                           surf_lsm_v(1)%qsws,                 &
3507                                           surf_usm_v(1)%qsws,                 &
3508                                           surf_def_v(2)%qsws,                 &
3509                                           surf_lsm_v(2)%qsws,                 &
3510                                           surf_usm_v(2)%qsws,                 &
3511                                           surf_def_v(3)%qsws,                 &
3512                                           surf_lsm_v(3)%qsws,                 &
3513                                           surf_usm_v(3)%qsws, n_out,          &
3514                                           waterflux_output_conversion )
3515
3516            CASE ( 'ssws' )
3517               CALL surface_data_output_sum_up( surf_def_h(0)%ssws,            &
3518                                           surf_def_h(1)%ssws,                 &
3519                                           surf_lsm_h%ssws,                    &
3520                                           surf_usm_h%ssws,                    &
3521                                           surf_def_v(0)%ssws,                 &
3522                                           surf_lsm_v(0)%ssws,                 &
3523                                           surf_usm_v(0)%ssws,                 &
3524                                           surf_def_v(1)%ssws,                 &
3525                                           surf_lsm_v(1)%ssws,                 &
3526                                           surf_usm_v(1)%ssws,                 &
3527                                           surf_def_v(2)%ssws,                 &
3528                                           surf_lsm_v(2)%ssws,                 &
3529                                           surf_usm_v(2)%ssws,                 &
3530                                           surf_def_v(3)%ssws,                 &
3531                                           surf_lsm_v(3)%ssws,                 &
3532                                           surf_usm_v(3)%ssws, n_out )
3533
3534            CASE ( 'qcsws' )
3535               CALL surface_data_output_sum_up( surf_def_h(0)%qcsws,           &
3536                                           surf_def_h(1)%qcsws,                &
3537                                           surf_lsm_h%qcsws,                   &
3538                                           surf_usm_h%qcsws,                   &
3539                                           surf_def_v(0)%qcsws,                &
3540                                           surf_lsm_v(0)%qcsws,                &
3541                                           surf_usm_v(0)%qcsws,                &
3542                                           surf_def_v(1)%qcsws,                &
3543                                           surf_lsm_v(1)%qcsws,                &
3544                                           surf_usm_v(1)%qcsws,                &
3545                                           surf_def_v(2)%qcsws,                &
3546                                           surf_lsm_v(2)%qcsws,                &
3547                                           surf_usm_v(2)%qcsws,                &
3548                                           surf_def_v(3)%qcsws,                &
3549                                           surf_lsm_v(3)%qcsws,                &
3550                                           surf_usm_v(3)%qcsws, n_out )
3551
3552            CASE ( 'ncsws' )
3553               CALL surface_data_output_sum_up( surf_def_h(0)%ncsws,           &
3554                                           surf_def_h(1)%ncsws,                &
3555                                           surf_lsm_h%ncsws,                   &
3556                                           surf_usm_h%ncsws,                   &
3557                                           surf_def_v(0)%ncsws,                &
3558                                           surf_lsm_v(0)%ncsws,                &
3559                                           surf_usm_v(0)%ncsws,                &
3560                                           surf_def_v(1)%ncsws,                &
3561                                           surf_lsm_v(1)%ncsws,                &
3562                                           surf_usm_v(1)%ncsws,                &
3563                                           surf_def_v(2)%ncsws,                &
3564                                           surf_lsm_v(2)%ncsws,                &
3565                                           surf_usm_v(2)%ncsws,                &
3566                                           surf_def_v(3)%ncsws,                &
3567                                           surf_lsm_v(3)%ncsws,                &
3568                                           surf_usm_v(3)%ncsws, n_out )
3569
3570            CASE ( 'qisws' )
3571               CALL surface_data_output_sum_up( surf_def_h(0)%qisws,           &
3572                                           surf_def_h(1)%qisws,                &
3573                                           surf_lsm_h%qisws,                   &
3574                                           surf_usm_h%qisws,                   &
3575                                           surf_def_v(0)%qisws,                &
3576                                           surf_lsm_v(0)%qisws,                &
3577                                           surf_usm_v(0)%qisws,                &
3578                                           surf_def_v(1)%qisws,                &
3579                                           surf_lsm_v(1)%qisws,                &
3580                                           surf_usm_v(1)%qisws,                &
3581                                           surf_def_v(2)%qisws,                &
3582                                           surf_lsm_v(2)%qisws,                &
3583                                           surf_usm_v(2)%qisws,                &
3584                                           surf_def_v(3)%qisws,                &
3585                                           surf_lsm_v(3)%qisws,                &
3586                                           surf_usm_v(3)%qisws, n_out )
3587
3588            CASE ( 'nisws' )
3589               CALL surface_data_output_sum_up( surf_def_h(0)%nisws,           &
3590                                           surf_def_h(1)%nisws,                &
3591                                           surf_lsm_h%nisws,                   &
3592                                           surf_usm_h%nisws,                   &
3593                                           surf_def_v(0)%nisws,                &
3594                                           surf_lsm_v(0)%nisws,                &
3595                                           surf_usm_v(0)%nisws,                &
3596                                           surf_def_v(1)%nisws,                &
3597                                           surf_lsm_v(1)%nisws,                &
3598                                           surf_usm_v(1)%nisws,                &
3599                                           surf_def_v(2)%nisws,                &
3600                                           surf_lsm_v(2)%nisws,                &
3601                                           surf_usm_v(2)%nisws,                &
3602                                           surf_def_v(3)%nisws,                &
3603                                           surf_lsm_v(3)%nisws,                &
3604                                           surf_usm_v(3)%nisws, n_out )
3605
3606            CASE ( 'qrsws' )
3607               CALL surface_data_output_sum_up( surf_def_h(0)%qrsws,           &
3608                                           surf_def_h(1)%qrsws,                &
3609                                           surf_lsm_h%qrsws,                   &
3610                                           surf_usm_h%qrsws,                   &
3611                                           surf_def_v(0)%qrsws,                &
3612                                           surf_lsm_v(0)%qrsws,                &
3613                                           surf_usm_v(0)%qrsws,                &
3614                                           surf_def_v(1)%qrsws,                &
3615                                           surf_lsm_v(1)%qrsws,                &
3616                                           surf_usm_v(1)%qrsws,                &
3617                                           surf_def_v(2)%qrsws,                &
3618                                           surf_lsm_v(2)%qrsws,                &
3619                                           surf_usm_v(2)%qrsws,                &
3620                                           surf_def_v(3)%qrsws,                &
3621                                           surf_lsm_v(3)%qrsws,                &
3622                                           surf_usm_v(3)%qrsws, n_out )
3623
3624            CASE ( 'nrsws' )
3625               CALL surface_data_output_sum_up( surf_def_h(0)%nrsws,           &
3626                                           surf_def_h(1)%nrsws,                &
3627                                           surf_lsm_h%nrsws,                   &
3628                                           surf_usm_h%nrsws,                   &
3629                                           surf_def_v(0)%nrsws,                &
3630                                           surf_lsm_v(0)%nrsws,                &
3631                                           surf_usm_v(0)%nrsws,                &
3632                                           surf_def_v(1)%nrsws,                &
3633                                           surf_lsm_v(1)%nrsws,                &
3634                                           surf_usm_v(1)%nrsws,                &
3635                                           surf_def_v(2)%nrsws,                &
3636                                           surf_lsm_v(2)%nrsws,                &
3637                                           surf_usm_v(2)%nrsws,                &
3638                                           surf_def_v(3)%nrsws,                &
3639                                           surf_lsm_v(3)%nrsws,                &
3640                                           surf_usm_v(3)%nrsws, n_out )
3641
3642            CASE ( 'sasws' )
3643               CALL surface_data_output_sum_up( surf_def_h(0)%sasws,           &
3644                                           surf_def_h(1)%sasws,                &
3645                                           surf_lsm_h%sasws,                   &
3646                                           surf_usm_h%sasws,                   &
3647                                           surf_def_v(0)%sasws,                &
3648                                           surf_lsm_v(0)%sasws,                &
3649                                           surf_usm_v(0)%sasws,                &
3650                                           surf_def_v(1)%sasws,                &
3651                                           surf_lsm_v(1)%sasws,                &
3652                                           surf_usm_v(1)%sasws,                &
3653                                           surf_def_v(2)%sasws,                &
3654                                           surf_lsm_v(2)%sasws,                &
3655                                           surf_usm_v(2)%sasws,                &
3656                                           surf_def_v(3)%sasws,                &
3657                                           surf_lsm_v(3)%sasws,                &
3658                                           surf_usm_v(3)%sasws, n_out )
3659
3660            CASE ( 'q_surface' )
3661               CALL surface_data_output_sum_up( surf_def_h(0)%q_surface,       &
3662                                           surf_def_h(1)%q_surface,            &
3663                                           surf_lsm_h%q_surface,               &
3664                                           surf_usm_h%q_surface,               &
3665                                           surf_def_v(0)%q_surface,            &
3666                                           surf_lsm_v(0)%q_surface,            &
3667                                           surf_usm_v(0)%q_surface,            &
3668                                           surf_def_v(1)%q_surface,            &
3669                                           surf_lsm_v(1)%q_surface,            &
3670                                           surf_usm_v(1)%q_surface,            &
3671                                           surf_def_v(2)%q_surface,            &
3672                                           surf_lsm_v(2)%q_surface,            &
3673                                           surf_usm_v(2)%q_surface,            &
3674                                           surf_def_v(3)%q_surface,            &
3675                                           surf_lsm_v(3)%q_surface,            &
3676                                           surf_usm_v(3)%q_surface, n_out )
3677
3678
3679            CASE ( 'theta_surface' )
3680               CALL surface_data_output_sum_up( surf_def_h(0)%pt_surface,      &
3681                                           surf_def_h(1)%pt_surface,           &
3682                                           surf_lsm_h%pt_surface,              &
3683                                           surf_usm_h%pt_surface,              &
3684                                           surf_def_v(0)%pt_surface,           &
3685                                           surf_lsm_v(0)%pt_surface,           &
3686                                           surf_usm_v(0)%pt_surface,           &
3687                                           surf_def_v(1)%pt_surface,           &
3688                                           surf_lsm_v(1)%pt_surface,           &
3689                                           surf_usm_v(1)%pt_surface,           &
3690                                           surf_def_v(2)%pt_surface,           &
3691                                           surf_lsm_v(2)%pt_surface,           &
3692                                           surf_usm_v(2)%pt_surface,           &
3693                                           surf_def_v(3)%pt_surface,           &
3694                                           surf_lsm_v(3)%pt_surface,           &
3695                                           surf_usm_v(3)%pt_surface, n_out )
3696
3697            CASE ( 'thetav_surface' )
3698               CALL surface_data_output_sum_up( surf_def_h(0)%vpt_surface,     &
3699                                           surf_def_h(1)%vpt_surface,          &
3700                                           surf_lsm_h%vpt_surface,             &
3701                                           surf_usm_h%vpt_surface,             &
3702                                           surf_def_v(0)%vpt_surface,          &
3703                                           surf_lsm_v(0)%vpt_surface,          &
3704                                           surf_usm_v(0)%vpt_surface,          &
3705                                           surf_def_v(1)%vpt_surface,          &
3706                                           surf_lsm_v(1)%vpt_surface,          &
3707                                           surf_usm_v(1)%vpt_surface,          &
3708                                           surf_def_v(2)%vpt_surface,          &
3709                                           surf_lsm_v(2)%vpt_surface,          &
3710                                           surf_usm_v(2)%vpt_surface,          &
3711                                           surf_def_v(3)%vpt_surface,          &
3712                                           surf_lsm_v(3)%vpt_surface,          &
3713                                           surf_usm_v(3)%vpt_surface, n_out )
3714
3715            CASE ( 'rad_net' )
3716               CALL surface_data_output_sum_up( surf_def_h(0)%rad_net,         &
3717                                           surf_def_h(1)%rad_net,              &
3718                                           surf_lsm_h%rad_net,                 &
3719                                           surf_usm_h%rad_net,                 &
3720                                           surf_def_v(0)%rad_net,              &
3721                                           surf_lsm_v(0)%rad_net,              &
3722                                           surf_usm_v(0)%rad_net,              &
3723                                           surf_def_v(1)%rad_net,              &
3724                                           surf_lsm_v(1)%rad_net,              &
3725                                           surf_usm_v(1)%rad_net,              &
3726                                           surf_def_v(2)%rad_net,              &
3727                                           surf_lsm_v(2)%rad_net,              &
3728                                           surf_usm_v(2)%rad_net,              &
3729                                           surf_def_v(3)%rad_net,              &
3730                                           surf_lsm_v(3)%rad_net,              &
3731                                           surf_usm_v(3)%rad_net, n_out )
3732
3733            CASE ( 'rad_lw_in' )
3734               CALL surface_data_output_sum_up( surf_def_h(0)%rad_lw_in,       &
3735                                           surf_def_h(1)%rad_lw_in,            &
3736                                           surf_lsm_h%rad_lw_in,               &
3737                                           surf_usm_h%rad_lw_in,               &
3738                                           surf_def_v(0)%rad_lw_in,            &
3739                                           surf_lsm_v(0)%rad_lw_in,            &
3740                                           surf_usm_v(0)%rad_lw_in,            &
3741                                           surf_def_v(1)%rad_lw_in,            &
3742                                           surf_lsm_v(1)%rad_lw_in,            &
3743                                           surf_usm_v(1)%rad_lw_in,            &
3744                                           surf_def_v(2)%rad_lw_in,            &
3745                                           surf_lsm_v(2)%rad_lw_in,            &
3746                                           surf_usm_v(2)%rad_lw_in,            &
3747                                           surf_def_v(3)%rad_lw_in,            &
3748                                           surf_lsm_v(3)%rad_lw_in,            &
3749                                           surf_usm_v(3)%rad_lw_in, n_out )
3750
3751            CASE ( 'rad_lw_out' )
3752               CALL surface_data_output_sum_up( surf_def_h(0)%rad_lw_out,      &
3753                                           surf_def_h(1)%rad_lw_out,           &
3754                                           surf_lsm_h%rad_lw_out,              &
3755                                           surf_usm_h%rad_lw_out,              &
3756                                           surf_def_v(0)%rad_lw_out,           &
3757                                           surf_lsm_v(0)%rad_lw_out,           &
3758                                           surf_usm_v(0)%rad_lw_out,           &
3759                                           surf_def_v(1)%rad_lw_out,           &
3760                                           surf_lsm_v(1)%rad_lw_out,           &
3761                                           surf_usm_v(1)%rad_lw_out,           &
3762                                           surf_def_v(2)%rad_lw_out,           &
3763                                           surf_lsm_v(2)%rad_lw_out,           &
3764                                           surf_usm_v(2)%rad_lw_out,           &
3765                                           surf_def_v(3)%rad_lw_out,           &
3766                                           surf_lsm_v(3)%rad_lw_out,           &
3767                                           surf_usm_v(3)%rad_lw_out, n_out )
3768
3769            CASE ( 'rad_sw_in' )
3770               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_in,       &
3771                                           surf_def_h(1)%rad_sw_in,            &
3772                                           surf_lsm_h%rad_sw_in,               &
3773                                           surf_usm_h%rad_sw_in,               &
3774                                           surf_def_v(0)%rad_sw_in,            &
3775                                           surf_lsm_v(0)%rad_sw_in,            &
3776                                           surf_usm_v(0)%rad_sw_in,            &
3777                                           surf_def_v(1)%rad_sw_in,            &
3778                                           surf_lsm_v(1)%rad_sw_in,            &
3779                                           surf_usm_v(1)%rad_sw_in,            &
3780                                           surf_def_v(2)%rad_sw_in,            &
3781                                           surf_lsm_v(2)%rad_sw_in,            &
3782                                           surf_usm_v(2)%rad_sw_in,            &
3783                                           surf_def_v(3)%rad_sw_in,            &
3784                                           surf_lsm_v(3)%rad_sw_in,            &
3785                                           surf_usm_v(3)%rad_sw_in, n_out )
3786
3787            CASE ( 'rad_sw_out' )
3788               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_out,      &
3789                                           surf_def_h(1)%rad_sw_out,           &
3790                                           surf_lsm_h%rad_sw_out,              &
3791                                           surf_usm_h%rad_sw_out,              &
3792                                           surf_def_v(0)%rad_sw_out,           &
3793                                           surf_lsm_v(0)%rad_sw_out,           &
3794                                           surf_usm_v(0)%rad_sw_out,           &
3795                                           surf_def_v(1)%rad_sw_out,           &
3796                                           surf_lsm_v(1)%rad_sw_out,           &
3797                                           surf_usm_v(1)%rad_sw_out,           &
3798                                           surf_def_v(2)%rad_sw_out,           &
3799                                           surf_lsm_v(2)%rad_sw_out,           &
3800                                           surf_usm_v(2)%rad_sw_out,           &
3801                                           surf_def_v(3)%rad_sw_out,           &
3802                                           surf_lsm_v(3)%rad_sw_out,           &
3803                                           surf_usm_v(3)%rad_sw_out, n_out )
3804
3805            CASE ( 'ghf' )
3806!
3807!--            Sum up ground / wall heat flux. Note, for urban surfaces the
3808!--            wall heat flux is aggregated from the different green, window
3809!--            and wall tiles.
3810               DO  m = 1, surf_usm_h%ns
3811                  surf_usm_h%ghf(m) = surf_usm_h%frac(m,ind_veg_wall)          &
3812                                    * surf_usm_h%wghf_eb(m) +                  &
3813                                      surf_usm_h%frac(m,ind_pav_green)         &
3814                                    * surf_usm_h%wghf_eb_green(m) +            &
3815                                      surf_usm_h%frac(m,ind_wat_win)           &
3816                                    * surf_usm_h%wghf_eb_window(m)
3817               ENDDO
3818               DO  l = 0, 3
3819                  DO  m = 1, surf_usm_v(l)%ns
3820                     surf_usm_v(l)%ghf(m) = surf_usm_v(l)%frac(m,ind_veg_wall) &
3821                                          * surf_usm_v(l)%wghf_eb(m) +         &
3822                                            surf_usm_v(l)%frac(m,ind_pav_green)&
3823                                          * surf_usm_v(l)%wghf_eb_green(m) +   &
3824                                            surf_usm_v(l)%frac(m,ind_wat_win)  &
3825                                          * surf_usm_v(l)%wghf_eb_window(m)
3826                  ENDDO
3827               ENDDO
3828
3829               CALL surface_data_output_sum_up( surf_def_h(0)%ghf,             &
3830                                           surf_def_h(1)%ghf,                  &
3831                                           surf_lsm_h%ghf,                     &
3832                                           surf_usm_h%ghf,                     &
3833                                           surf_def_v(0)%ghf,                  &
3834                                           surf_lsm_v(0)%ghf,                  &
3835                                           surf_usm_v(0)%ghf,                  &
3836                                           surf_def_v(1)%ghf,                  &
3837                                           surf_lsm_v(1)%ghf,                  &
3838                                           surf_usm_v(1)%ghf,                  &
3839                                           surf_def_v(2)%ghf,                  &
3840                                           surf_lsm_v(2)%ghf,                  &
3841                                           surf_usm_v(2)%ghf,                  &
3842                                           surf_def_v(3)%ghf,                  &
3843                                           surf_lsm_v(3)%ghf,                  &
3844                                           surf_usm_v(3)%ghf, n_out )
3845
3846            CASE ( 'r_a' )
3847               CALL surface_data_output_sum_up( surf_def_h(0)%r_a,             &
3848                                           surf_def_h(1)%r_a,                  &
3849                                           surf_lsm_h%r_a,                     &
3850                                           surf_usm_h%r_a,                     &
3851                                           surf_def_v(0)%r_a,                  &
3852                                           surf_lsm_v(0)%r_a,                  &
3853                                           surf_usm_v(0)%r_a,                  &
3854                                           surf_def_v(1)%r_a,                  &
3855                                           surf_lsm_v(1)%r_a,                  &
3856                                           surf_usm_v(1)%r_a,                  &
3857                                           surf_def_v(2)%r_a,                  &
3858                                           surf_lsm_v(2)%r_a,                  &
3859                                           surf_usm_v(2)%r_a,                  &
3860                                           surf_def_v(3)%r_a,                  &
3861                                           surf_lsm_v(3)%r_a,                  &
3862                                           surf_usm_v(3)%r_a, n_out )
3863
3864            CASE ( 'r_soil' )
3865               CALL surface_data_output_sum_up( surf_def_h(0)%r_soil,          &
3866                                           surf_def_h(1)%r_soil,               &
3867                                           surf_lsm_h%r_soil,                  &
3868                                           surf_usm_h%r_soil,                  &
3869                                           surf_def_v(0)%r_soil,               &
3870                                           surf_lsm_v(0)%r_soil,               &
3871                                           surf_usm_v(0)%r_soil,               &
3872                                           surf_def_v(1)%r_soil,               &
3873                                           surf_lsm_v(1)%r_soil,               &
3874                                           surf_usm_v(1)%r_soil,               &
3875                                           surf_def_v(2)%r_soil,               &
3876                                           surf_lsm_v(2)%r_soil,               &
3877                                           surf_usm_v(2)%r_soil,               &
3878                                           surf_def_v(3)%r_soil,               &
3879                                           surf_lsm_v(3)%r_soil,               &
3880                                           surf_usm_v(3)%r_soil, n_out )
3881
3882            CASE ( 'r_canopy' )
3883               CALL surface_data_output_sum_up( surf_def_h(0)%r_canopy,        &
3884                                           surf_def_h(1)%r_canopy,             &
3885                                           surf_lsm_h%r_canopy,                &
3886                                           surf_usm_h%r_canopy,                &
3887                                           surf_def_v(0)%r_canopy,             &
3888                                           surf_lsm_v(0)%r_canopy,             &
3889                                           surf_usm_v(0)%r_canopy,             &
3890                                           surf_def_v(1)%r_canopy,             &
3891                                           surf_lsm_v(1)%r_canopy,             &
3892                                           surf_usm_v(1)%r_canopy,             &
3893                                           surf_def_v(2)%r_canopy,             &
3894                                           surf_lsm_v(2)%r_canopy,             &
3895                                           surf_usm_v(2)%r_canopy,             &
3896                                           surf_def_v(3)%r_canopy,             &
3897                                           surf_lsm_v(3)%r_canopy,             &
3898                                           surf_usm_v(3)%r_canopy, n_out )
3899
3900            CASE ( 'r_s' )
3901               CALL surface_data_output_sum_up( surf_def_h(0)%r_s,             &
3902                                           surf_def_h(1)%r_s,                  &
3903                                           surf_lsm_h%r_s,                     &
3904                                           surf_usm_h%r_s,                     &
3905                                           surf_def_v(0)%r_s,                  &
3906                                           surf_lsm_v(0)%r_s,                  &
3907                                           surf_usm_v(0)%r_s,                  &
3908                                           surf_def_v(1)%r_s,                  &
3909                                           surf_lsm_v(1)%r_s,                  &
3910                                           surf_usm_v(1)%r_s,                  &
3911                                           surf_def_v(2)%r_s,                  &
3912                                           surf_lsm_v(2)%r_s,                  &
3913                                           surf_usm_v(2)%r_s,                  &
3914                                           surf_def_v(3)%r_s,                  &
3915                                           surf_lsm_v(3)%r_s,                  &
3916                                           surf_usm_v(3)%r_s, n_out )
3917
3918
3919            CASE ( 'rad_sw_dir' )
3920               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_dir,      &
3921                                           surf_def_h(1)%rad_sw_dir,           &
3922                                           surf_lsm_h%rad_sw_dir,              &
3923                                           surf_usm_h%rad_sw_dir,              &
3924                                           surf_def_v(0)%rad_sw_dir,           &
3925                                           surf_lsm_v(0)%rad_sw_dir,           &
3926                                           surf_usm_v(0)%rad_sw_dir,           &
3927                                           surf_def_v(1)%rad_sw_dir,           &
3928                                           surf_lsm_v(1)%rad_sw_dir,           &
3929                                           surf_usm_v(1)%rad_sw_dir,           &
3930                                           surf_def_v(2)%rad_sw_dir,           &
3931                                           surf_lsm_v(2)%rad_sw_dir,           &
3932                                           surf_usm_v(2)%rad_sw_dir,           &
3933                                           surf_def_v(3)%rad_sw_dir,           &
3934                                           surf_lsm_v(3)%rad_sw_dir,           &
3935                                           surf_usm_v(3)%rad_sw_dir, n_out )
3936            CASE ( 'rad_sw_dif' )
3937               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_dif,      &
3938                                           surf_def_h(1)%rad_sw_dif,           &
3939                                           surf_lsm_h%rad_sw_dif,              &
3940                                           surf_usm_h%rad_sw_dif,              &
3941                                           surf_def_v(0)%rad_sw_dif,           &
3942                                           surf_lsm_v(0)%rad_sw_dif,           &
3943                                           surf_usm_v(0)%rad_sw_dif,           &
3944                                           surf_def_v(1)%rad_sw_dif,           &
3945                                           surf_lsm_v(1)%rad_sw_dif,           &
3946                                           surf_usm_v(1)%rad_sw_dif,           &
3947                                           surf_def_v(2)%rad_sw_dif,           &
3948                                           surf_lsm_v(2)%rad_sw_dif,           &
3949                                           surf_usm_v(2)%rad_sw_dif,           &
3950                                           surf_def_v(3)%rad_sw_dif,           &
3951                                           surf_lsm_v(3)%rad_sw_dif,           &
3952                                           surf_usm_v(3)%rad_sw_dif, n_out )
3953
3954            CASE ( 'rad_sw_ref' )
3955               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_ref,      &
3956                                           surf_def_h(1)%rad_sw_ref,           &
3957                                           surf_lsm_h%rad_sw_ref,              &
3958                                           surf_usm_h%rad_sw_ref,              &
3959                                           surf_def_v(0)%rad_sw_ref,           &
3960                                           surf_lsm_v(0)%rad_sw_ref,           &
3961                                           surf_usm_v(0)%rad_sw_ref,           &
3962                                           surf_def_v(1)%rad_sw_ref,           &
3963                                           surf_lsm_v(1)%rad_sw_ref,           &
3964                                           surf_usm_v(1)%rad_sw_ref,           &
3965                                           surf_def_v(2)%rad_sw_ref,           &
3966                                           surf_lsm_v(2)%rad_sw_ref,           &
3967                                           surf_usm_v(2)%rad_sw_ref,           &
3968                                           surf_def_v(3)%rad_sw_ref,           &
3969                                           surf_lsm_v(3)%rad_sw_ref,           &
3970                                           surf_usm_v(3)%rad_sw_ref, n_out )
3971
3972            CASE ( 'rad_sw_res' )
3973               CALL surface_data_output_sum_up( surf_def_h(0)%rad_sw_res,      &
3974                                           surf_def_h(1)%rad_sw_res,           &
3975                                           surf_lsm_h%rad_sw_res,              &
3976                                           surf_usm_h%rad_sw_res,              &
3977                                           surf_def_v(0)%rad_sw_res,           &
3978                                           surf_lsm_v(0)%rad_sw_res,           &
3979                                           surf_usm_v(0)%rad_sw_res,           &
3980                                           surf_def_v(1)%rad_sw_res,           &
3981                                           surf_lsm_v(1)%rad_sw_res,           &
3982                                           surf_usm_v(1)%rad_sw_res,           &
3983                                           surf_def_v(2)%rad_sw_res,           &
3984                                           surf_lsm_v(2)%rad_sw_res,           &
3985                                           surf_usm_v(2)%rad_sw_res,           &
3986                                           surf_def_v(3)%rad_sw_res,           &
3987                                           surf_lsm_v(3)%rad_sw_res,           &
3988                                           surf_usm_v(3)%rad_sw_res, n_out )
3989
3990            CASE ( 'rad_lw_dif' )
3991               CALL surface_data_output_sum_up( surf_def_h(0)%rad_lw_dif,      &
3992                                           surf_def_h(1)%rad_lw_dif,           &
3993                                           surf_lsm_h%rad_lw_dif,              &
3994                                           surf_usm_h%rad_lw_dif,              &
3995                                           surf_def_v(0)%rad_lw_dif,           &
3996                                           surf_lsm_v(0)%rad_lw_dif,           &
3997                                           surf_usm_v(0)%rad_lw_dif,           &
3998                                           surf_def_v(1)%rad_lw_dif,           &
3999                                           surf_lsm_v(1)%rad_lw_dif,           &
4000                                           surf_usm_v(1)%rad_lw_dif,           &
4001                                           surf_def_v(2)%rad_lw_dif,           &
4002                                           surf_lsm_v(2)%rad_lw_dif,           &
4003                                           surf_usm_v(2)%rad_lw_dif,           &
4004                                           surf_def_v(3)%rad_lw_dif,           &
4005                                           surf_lsm_v(3)%rad_lw_dif,           &
4006                                           surf_usm_v(3)%rad_lw_dif, n_out )
4007
4008            CASE ( 'rad_lw_ref' )
4009               CALL surface_data_output_sum_up( surf_def_h(0)%rad_lw_ref,      &
4010                                           surf_def_h(1)%rad_lw_ref,           &
4011                                           surf_lsm_h%rad_lw_ref,              &
4012                                           surf_usm_h%rad_lw_ref,              &
4013                                           surf_def_v(0)%rad_lw_ref,           &
4014                                           surf_lsm_v(0)%rad_lw_ref,           &
4015                                           surf_usm_v(0)%rad_lw_ref,           &
4016                                           surf_def_v(1)%rad_lw_ref,           &
4017                                           surf_lsm_v(1)%rad_lw_ref,           &
4018                                           surf_usm_v(1)%rad_lw_ref,           &
4019                                           surf_def_v(2)%rad_lw_ref,           &
4020                                           surf_lsm_v(2)%rad_lw_ref,           &
4021                                           surf_usm_v(2)%rad_lw_ref,           &
4022                                           surf_def_v(3)%rad_lw_ref,           &
4023                                           surf_lsm_v(3)%rad_lw_ref,           &
4024                                           surf_usm_v(3)%rad_lw_ref, n_out )
4025
4026            CASE ( 'rad_lw_res' )
4027               CALL surface_data_output_sum_up( surf_def_h(0)%rad_lw_res,      &
4028                                           surf_def_h(1)%rad_lw_res,           &
4029                                           surf_lsm_h%rad_lw_res,              &
4030                                           surf_usm_h%rad_lw_res,              &
4031                                           surf_def_v(0)%rad_lw_res,           &
4032                                           surf_lsm_v(0)%rad_lw_res,           &
4033                                           surf_usm_v(0)%rad_lw_res,           &
4034                                           surf_def_v(1)%rad_lw_res,           &
4035                                           surf_lsm_v(1)%rad_lw_res,           &
4036                                           surf_usm_v(1)%rad_lw_res,           &
4037                                           surf_def_v(2)%rad_lw_res,           &
4038                                           surf_lsm_v(2)%rad_lw_res,           &
4039                                           surf_usm_v(2)%rad_lw_res,           &
4040                                           surf_def_v(3)%rad_lw_res,           &
4041                                           surf_lsm_v(3)%rad_lw_res,           &
4042                                           surf_usm_v(3)%rad_lw_res, n_out )
4043
4044            CASE ( 'uvw1' )
4045               CALL surface_data_output_sum_up( surf_def_h(0)%uvw_abs,         &
4046                                           surf_def_h(1)%uvw_abs,              &
4047                                           surf_lsm_h%uvw_abs,                 &
4048                                           surf_usm_h%uvw_abs,                 &
4049                                           surf_def_v(0)%uvw_abs,              &
4050                                           surf_lsm_v(0)%uvw_abs,              &
4051                                           surf_usm_v(0)%uvw_abs,              &
4052                                           surf_def_v(1)%uvw_abs,              &
4053                                           surf_lsm_v(1)%uvw_abs,              &
4054                                           surf_usm_v(1)%uvw_abs,              &
4055                                           surf_def_v(2)%uvw_abs,              &
4056                                           surf_lsm_v(2)%uvw_abs,              &
4057                                           surf_usm_v(2)%uvw_abs,              &
4058                                           surf_def_v(3)%uvw_abs,              &
4059                                           surf_lsm_v(3)%uvw_abs,              &
4060                                           surf_usm_v(3)%uvw_abs, n_out )
4061
4062            CASE ( 'waste_heat' )
4063               CALL surface_data_output_sum_up( surf_def_h(0)%waste_heat,      &
4064                                           surf_def_h(1)%waste_heat,           &
4065                                           surf_lsm_h%waste_heat,              &
4066                                           surf_usm_h%waste_heat,              &
4067                                           surf_def_v(0)%waste_heat,           &
4068                                           surf_lsm_v(0)%waste_heat,           &
4069                                           surf_usm_v(0)%waste_heat,           &
4070                                           surf_def_v(1)%waste_heat,           &
4071                                           surf_lsm_v(1)%waste_heat,           &
4072                                           surf_usm_v(1)%waste_heat,           &
4073                                           surf_def_v(2)%waste_heat,           &
4074                                           surf_lsm_v(2)%waste_heat,           &
4075                                           surf_usm_v(2)%waste_heat,           &
4076                                           surf_def_v(3)%waste_heat,           &
4077                                           surf_lsm_v(3)%waste_heat,           &
4078                                           surf_usm_v(3)%waste_heat, n_out )
4079
4080            CASE ( 'im_hf' )
4081               CALL surface_data_output_sum_up( surf_def_h(0)%iwghf_eb,        &
4082                                           surf_def_h(1)%iwghf_eb,             &
4083                                           surf_lsm_h%iwghf_eb,                &
4084                                           surf_usm_h%iwghf_eb,                &
4085                                           surf_def_v(0)%iwghf_eb,             &
4086                                           surf_lsm_v(0)%iwghf_eb,             &
4087                                           surf_usm_v(0)%iwghf_eb,             &
4088                                           surf_def_v(1)%iwghf_eb,             &
4089                                           surf_lsm_v(1)%iwghf_eb,             &
4090                                           surf_usm_v(1)%iwghf_eb,             &
4091                                           surf_def_v(2)%iwghf_eb,             &
4092                                           surf_lsm_v(2)%iwghf_eb,             &
4093                                           surf_usm_v(2)%iwghf_eb,             &
4094                                           surf_def_v(3)%iwghf_eb,             &
4095                                           surf_lsm_v(3)%iwghf_eb,             &
4096                                           surf_usm_v(3)%iwghf_eb, n_out )
4097
4098         END SELECT
4099      ENDDO
4100
4101
4102   END SUBROUTINE surface_data_output_averaging
4103
4104!------------------------------------------------------------------------------!
4105! Description:
4106! ------------
4107!> Sum-up the surface data for average output variables.
4108!------------------------------------------------------------------------------!
4109   SUBROUTINE surface_data_output_sum_up( var_def_h0, var_def_h1,              &
4110                                     var_lsm_h,  var_usm_h,                    &
4111                                     var_def_v0, var_lsm_v0, var_usm_v0,       &
4112                                     var_def_v1, var_lsm_v1, var_usm_v1,       &
4113                                     var_def_v2, var_lsm_v2, var_usm_v2,       &
4114                                     var_def_v3, var_lsm_v3, var_usm_v3, n_out,&
4115                                     fac )
4116
4117      IMPLICIT NONE
4118
4119      INTEGER(iwp) ::  k          !< height index of surface element
4120      INTEGER(iwp) ::  m          !< running index for surface elements
4121      INTEGER(iwp) ::  n_out      !< index for output variable
4122      INTEGER(iwp) ::  n_surf     !< running index for surface elements
4123
4124      REAL(wp), DIMENSION(:), OPTIONAL                ::  fac               !< passed output conversion factor for heatflux output
4125      REAL(wp), DIMENSION(nzb:nzt+1)                  ::  conversion_factor !< effective array for output conversion factor
4126
4127      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h0 !< output variable at upward-facing default-type surfaces
4128      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h1 !< output variable at downward-facing default-type surfaces
4129      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_h  !< output variable at upward-facing natural-type surfaces
4130      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_h  !< output variable at upward-facing urban-type surfaces
4131      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v0 !< output variable at northward-facing default-type surfaces
4132      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v1 !< output variable at southward-facing default-type surfaces
4133      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v2 !< output variable at eastward-facing default-type surfaces
4134      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v3 !< output variable at westward-facing default-type surfaces
4135      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v0 !< output variable at northward-facing natural-type surfaces
4136      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v1 !< output variable at southward-facing natural-type surfaces
4137      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v2 !< output variable at eastward-facing natural-type surfaces
4138      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v3 !< output variable at westward-facing natural-type surfaces
4139      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v0 !< output variable at northward-facing urban-type surfaces
4140      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v1 !< output variable at southward-facing urban-type surfaces
4141      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v2 !< output variable at eastward-facing urban-type surfaces
4142      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v3 !< output variable at westward-facing urban-type surfaces
4143
4144!
4145!--   Set conversion factor to one if not present
4146      IF ( .NOT. PRESENT( fac ) )  THEN
4147         conversion_factor = 1.0_wp
4148      ELSE
4149         conversion_factor = fac
4150      ENDIF
4151!
4152!--   Set counter variable to zero before the variable is written to
4153!--   the output array.
4154      n_surf = 0
4155
4156!
4157!--   Write the horizontal surfaces.
4158!--   Before each the variable is written to the output data structure, first
4159!--   check if the variable for the respective surface type is defined.
4160!--   If a variable is not defined, skip the block and increment the counter
4161!--   variable by the number of surface elements of this type. Usually this
4162!--   is zere, however, there might be the situation that e.g. urban surfaces
4163!--   are defined but the respective variable is not allocated for this surface
4164!--   type. To write the data on the exact position, increment the counter.
4165      IF ( ALLOCATED( var_def_h0 ) )  THEN
4166         DO  m = 1, surf_def_h(0)%ns
4167            n_surf                        = n_surf + 1
4168            k                             = surf_def_h(0)%k(m)
4169            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4170                                          + var_def_h0(m) * conversion_factor(k)
4171         ENDDO
4172      ELSE
4173         n_surf = n_surf + surf_def_h(0)%ns
4174      ENDIF
4175      IF ( ALLOCATED( var_def_h1 ) )  THEN
4176         DO  m = 1, surf_def_h(1)%ns
4177            n_surf                   = n_surf + 1
4178            k                             = surf_def_h(1)%k(m)
4179            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4180                                          + var_def_h1(m) * conversion_factor(k)
4181         ENDDO
4182      ELSE
4183         n_surf = n_surf + surf_def_h(1)%ns
4184      ENDIF
4185      IF ( ALLOCATED( var_lsm_h ) )  THEN
4186         DO  m = 1, surf_lsm_h%ns
4187            n_surf                        = n_surf + 1
4188            k                             = surf_lsm_h%k(m)
4189            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4190                                          + var_lsm_h(m) * conversion_factor(k)
4191         ENDDO
4192      ELSE
4193         n_surf = n_surf + surf_lsm_h%ns
4194      ENDIF
4195      IF ( ALLOCATED( var_usm_h ) )  THEN
4196         DO  m = 1, surf_usm_h%ns
4197            n_surf                        = n_surf + 1
4198            k                             = surf_usm_h%k(m)
4199            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4200                                          + var_usm_h(m) * conversion_factor(k)
4201         ENDDO
4202      ELSE
4203         n_surf = n_surf + surf_usm_h%ns
4204      ENDIF
4205!
4206!--   Write northward-facing
4207      IF ( ALLOCATED( var_def_v0 ) )  THEN
4208         DO  m = 1, surf_def_v(0)%ns
4209            n_surf                        = n_surf + 1
4210            k                             = surf_def_v(0)%k(m)
4211            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4212                                          + var_def_v0(m) * conversion_factor(k)
4213         ENDDO
4214      ELSE
4215         n_surf = n_surf + surf_def_v(0)%ns
4216      ENDIF
4217      IF ( ALLOCATED( var_lsm_v0 ) )  THEN
4218         DO  m = 1, surf_lsm_v(0)%ns
4219            n_surf                        = n_surf + 1
4220            k                             = surf_lsm_v(0)%k(m)
4221            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4222                                          + var_lsm_v0(m) * conversion_factor(k)
4223         ENDDO
4224      ELSE
4225         n_surf = n_surf + surf_lsm_v(0)%ns
4226      ENDIF
4227      IF ( ALLOCATED( var_usm_v0 ) )  THEN
4228         DO  m = 1, surf_usm_v(0)%ns
4229            n_surf                        = n_surf + 1
4230            k                             = surf_usm_v(0)%k(m)
4231            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4232                                          + var_usm_v0(m) * conversion_factor(k)
4233         ENDDO
4234      ELSE
4235         n_surf = n_surf + surf_usm_v(0)%ns
4236      ENDIF
4237!
4238!--   Write southward-facing
4239      IF ( ALLOCATED( var_def_v1 ) )  THEN
4240         DO  m = 1, surf_def_v(1)%ns
4241            n_surf                        = n_surf + 1
4242            k                             = surf_def_v(1)%k(m)
4243            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4244                                          + var_def_v1(m) * conversion_factor(k)
4245         ENDDO
4246      ELSE
4247         n_surf = n_surf + surf_def_v(1)%ns
4248      ENDIF
4249      IF ( ALLOCATED( var_lsm_v1 ) )  THEN
4250         DO  m = 1, surf_lsm_v(1)%ns
4251            n_surf                        = n_surf + 1
4252            k                             = surf_lsm_v(1)%k(m)
4253            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4254                                          + var_lsm_v1(m) * conversion_factor(k)
4255         ENDDO
4256      ELSE
4257         n_surf = n_surf + surf_lsm_v(1)%ns
4258      ENDIF
4259      IF ( ALLOCATED( var_usm_v1 ) )  THEN
4260         DO  m = 1, surf_usm_v(1)%ns
4261            n_surf                        = n_surf + 1
4262            k                             = surf_usm_v(1)%k(m)
4263            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4264                                          + var_usm_v1(m) * conversion_factor(k)
4265         ENDDO
4266      ELSE
4267         n_surf = n_surf + surf_usm_v(1)%ns
4268      ENDIF
4269!
4270!--   Write eastward-facing
4271      IF ( ALLOCATED( var_def_v2 ) )  THEN
4272         DO  m = 1, surf_def_v(2)%ns
4273            n_surf                        = n_surf + 1
4274            k                             = surf_def_v(2)%k(m)
4275            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4276                                          + var_def_v2(m) * conversion_factor(k)
4277         ENDDO
4278      ELSE
4279         n_surf = n_surf + surf_def_v(2)%ns
4280      ENDIF
4281      IF ( ALLOCATED( var_lsm_v2 ) )  THEN
4282         DO  m = 1, surf_lsm_v(2)%ns
4283            n_surf                        = n_surf + 1
4284            k                             = surf_lsm_v(2)%k(m)
4285            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4286                                          + var_lsm_v2(m) * conversion_factor(k)
4287         ENDDO
4288      ELSE
4289         n_surf = n_surf + surf_lsm_v(2)%ns
4290      ENDIF
4291      IF ( ALLOCATED( var_usm_v2 ) )  THEN
4292         DO  m = 1, surf_usm_v(2)%ns
4293            n_surf                        = n_surf + 1
4294            k                             = surf_usm_v(2)%k(m)
4295            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4296                                          + var_usm_v2(m) * conversion_factor(k)
4297         ENDDO
4298      ELSE
4299         n_surf = n_surf + surf_usm_v(2)%ns
4300      ENDIF
4301!
4302!--   Write westward-facing
4303      IF ( ALLOCATED( var_def_v3 ) )  THEN
4304         DO  m = 1, surf_def_v(3)%ns
4305            n_surf                        = n_surf + 1
4306            k                             = surf_def_v(3)%k(m)
4307            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4308                                          + var_def_v3(m) * conversion_factor(k)
4309         ENDDO
4310      ELSE
4311         n_surf = n_surf + surf_def_v(3)%ns
4312      ENDIF
4313      IF ( ALLOCATED( var_lsm_v3 ) )  THEN
4314         DO  m = 1, surf_lsm_v(3)%ns
4315            n_surf                        = n_surf + 1
4316            k                             = surf_lsm_v(3)%k(m)
4317            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4318                                          + var_lsm_v3(m) * conversion_factor(k)
4319         ENDDO
4320      ELSE
4321         n_surf = n_surf + surf_lsm_v(3)%ns
4322      ENDIF
4323      IF ( ALLOCATED( var_usm_v3 ) )  THEN
4324         DO  m = 1, surf_usm_v(3)%ns
4325            n_surf                        = n_surf + 1
4326            k                             = surf_usm_v(3)%k(m)
4327            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
4328                                          + var_usm_v3(m) * conversion_factor(k)
4329         ENDDO
4330      ELSE
4331         n_surf = n_surf + surf_usm_v(3)%ns
4332      ENDIF
4333
4334   END SUBROUTINE surface_data_output_sum_up
4335
4336!------------------------------------------------------------------------------!
4337! Description:
4338! ------------
4339!> Collect the surface data from different types and different orientation.
4340!------------------------------------------------------------------------------!
4341   SUBROUTINE surface_data_output_collect( var_def_h0, var_def_h1,             &
4342                                      var_lsm_h,  var_usm_h,                   &</