source: palm/trunk/SOURCE/surface_output_mod.f90 @ 3552

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

Surface output revised and some bugs are fixed + new post-processing tool to convert binary surface output to Paraview readable VTK files

  • Property svn:keywords set to Id
File size: 183.9 KB
Line 
1!> @file surface_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-2018 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: surface_output_mod.f90 3494 2018-11-06 14:51:27Z suehring $
27! Bugfix in gathering surface data from different types and orientation.
28! Output of total number of surfaces and vertices added.
29! String length is output, for more convinient post-processing.
30! Last actions added.
31!
32! 3483 2018-11-02 14:19:26Z raasch
33! bugfix: missed directives for MPI added
34!
35! 3421 2018-10-24 18:39:32Z gronemeier
36! Add output variables
37! Write data into binary file
38!
39! 3420 2018-10-24 17:30:08Z gronemeier
40! Initial implementation from Klaus Ketelsen and Matthias Suehring
41!
42!
43! Authors:
44! --------
45! @author Klaus Ketelsen, Matthias Suehring
46!
47! Description:
48! ------------
49!> Generate output for surface data.
50!>
51!> @todo Create namelist file for post-processing tool.
52!------------------------------------------------------------------------------!
53
54MODULE surface_output_mod
55
56   USE kinds
57
58   USE arrays_3d,                                                              &
59       ONLY:  zu, zw
60       
61   USE control_parameters,                                                     &
62       ONLY:  surface_data_output
63       
64   USE grid_variables,                                                         &
65       ONLY: dx,dy
66       
67   USE indices,                                                                &
68       ONLY: nxl, nxr, nys, nyn, nzb, nzt
69       
70   USE pegrid
71       
72   USE surface_mod,                                                            &
73       ONLY:  surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v,                  &
74              surf_usm_h, surf_usm_v
75
76   IMPLICIT NONE
77
78   TYPE surf_out                      !< data structure which contains all surfaces elements of all types on subdomain
79   
80      INTEGER(iwp) ::  ns             !< number of surface elements on subdomain
81      INTEGER(iwp) ::  ns_total       !< total number of surface elements
82      INTEGER(iwp) ::  npoints        !< number of points / vertices which define a surface element (on subdomain)
83      INTEGER(iwp) ::  npoints_total  !< total number of points / vertices which define a surface element
84     
85      REAL(wp) ::  fillvalue = -9999.9_wp !< fillvalue for surface elements which are not defined
86     
87      REAL(wp), DIMENSION(:), ALLOCATABLE   ::  var_out  !< output variables
88      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  var_av   !< variables used for averaging
89      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  points   !< points  / vertices of a surface element
90      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  polygons !< polygon data of a surface element
91   END TYPE surf_out
92   
93   CHARACTER(LEN=100), DIMENSION(300)       ::  data_output_surf = ' '  !< namelist variable which describes the output variables
94   CHARACTER(LEN=100), DIMENSION(0:1,300)   ::  dosurf = ' '            !< internal variable which describes the output variables and separates
95                                                                        !< averaged from non-averaged output
96   
97   INTEGER(iwp) ::  average_count_surf = 0                !< number of ensemble members used for averaging
98   INTEGER(iwp) ::  dosurf_no(0:1) = 0                    !< number of surface output quantities
99
100   LOGICAL :: first_output(0:1) = .FALSE.                 !< true if first output was already called
101
102   REAL(wp) ::  averaging_interval_surf  = 9999999.9_wp   !< averaging interval
103   REAL(wp) ::  dt_dosurf = 9999999.9_wp                  !< time interval for instantaneous data output
104   REAL(wp) ::  dt_dosurf_av = 9999999.9_wp               !< time interval for averaged data output
105   REAL(wp) ::  skip_time_dosurf = 0.0_wp                 !< skip time for instantaneous data output
106   REAL(wp) ::  skip_time_dosurf_av = 0.0_wp              !< skip time for averaged data output
107   REAL(wp) ::  time_dosurf = 0.0_wp                      !< internal counter variable to check for instantaneous data output
108   REAL(wp) ::  time_dosurf_av = 0.0_wp                   !< internal counter variable to check for averaged data output
109
110   
111   TYPE(surf_out) ::  surfaces      !< variable which contains all required output information
112   
113   SAVE
114
115   PRIVATE
116
117   INTERFACE  surface_output
118      MODULE PROCEDURE surface_output
119   END INTERFACE  surface_output
120   
121   INTERFACE  surface_output_averaging
122      MODULE PROCEDURE surface_output_averaging
123   END INTERFACE  surface_output_averaging
124   
125   INTERFACE  surface_output_check_parameters
126      MODULE PROCEDURE surface_output_check_parameters
127   END INTERFACE  surface_output_check_parameters
128   
129   INTERFACE  surface_output_init
130      MODULE PROCEDURE surface_output_init
131   END INTERFACE  surface_output_init
132   
133   INTERFACE  surface_output_last_action
134      MODULE PROCEDURE surface_output_last_action
135   END INTERFACE  surface_output_last_action
136     
137   INTERFACE  surface_output_parin
138      MODULE PROCEDURE surface_output_parin
139   END INTERFACE  surface_output_parin
140
141!
142!--Public subroutines
143   PUBLIC surface_output, surface_output_averaging,                            &
144          surface_output_check_parameters, surface_output_init,                &
145          surface_output_last_action, surface_output_parin
146!
147!--Public variables
148   PUBLIC average_count_surf, averaging_interval_surf, dt_dosurf, dt_dosurf_av,&
149          skip_time_dosurf, skip_time_dosurf_av, time_dosurf, time_dosurf_av
150
151 CONTAINS
152
153!------------------------------------------------------------------------------!
154! Description:
155! ------------
156!> Initialization surface-data output data structure: calculation of vertices
157!> and polygon data for the surface elements, allocation of required arrays.
158!------------------------------------------------------------------------------!
159   SUBROUTINE surface_output_init
160   
161      IMPLICIT NONE
162
163      INTEGER(iwp) ::  i                 !< grid index in x-direction, also running variable for counting non-average data output
164      INTEGER(iwp) ::  ilen              !< string length
165      INTEGER(iwp) ::  j                 !< grid index in y-direction, also running variable for counting average data output
166      INTEGER(iwp) ::  k                 !< grid index in z-direction
167      INTEGER(iwp) ::  l                 !< running index for surface-element orientation
168      INTEGER(iwp) ::  m                 !< running index for surface elements
169      INTEGER(iwp) ::  n_out             !< running index for number of output variables   
170      INTEGER(iwp) ::  npg               !< counter variable for all surface elements ( or polygons )
171      INTEGER(iwp) ::  point_index_count !< local counter variable for point index
172     
173      INTEGER(iwp), DIMENSION(0:numprocs-1) :: num_points_on_pe   !< array which contains the number of points on all mpi ranks
174      INTEGER(iwp), ALLOCATABLE, DIMENSION(:,:,:) ::  point_index !< dummy array used to check where the reference points for surface polygons are located
175           
176!
177!--   Determine the number of surface elements on subdomain
178      surfaces%ns = surf_def_h(0)%ns + surf_lsm_h%ns + surf_usm_h%ns           & !horizontal upward-facing
179                  + surf_def_h(1)%ns                                           & !horizontal downard-facing   
180                  + surf_def_v(0)%ns + surf_lsm_v(0)%ns + surf_usm_v(0)%ns     & !northward-facing
181                  + surf_def_v(1)%ns + surf_lsm_v(1)%ns + surf_usm_v(1)%ns     & !southward-facing   
182                  + surf_def_v(2)%ns + surf_lsm_v(2)%ns + surf_usm_v(2)%ns     & !westward-facing   
183                  + surf_def_v(3)%ns + surf_lsm_v(3)%ns + surf_usm_v(3)%ns       !eastward-facing   
184!
185!--   Allocate output variable and set to _FillValue attribute
186      ALLOCATE ( surfaces%var_out(1:surfaces%ns) )
187      surfaces%var_out = surfaces%fillvalue
188!
189!--   If there is an output of time average output variables, allocate the
190!--   required array.
191      IF ( dosurf_no(1) > 0 )  THEN
192         ALLOCATE ( surfaces%var_av(1:surfaces%ns,1:dosurf_no(1)) )
193         surfaces%var_av = 0.0_wp
194      ENDIF
195!
196!--   Initialize points and polygon data.
197!--   In a first step, count the number of points which are defining
198!--   the surfaces and the polygons.
199      ALLOCATE( point_index(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) )
200      point_index = -1
201!
202!--   Horizontal default surfaces
203      surfaces%npoints = 0
204      DO  l = 0, 1
205         DO  m = 1, surf_def_h(0)%ns
206!
207!--         Determine the indices of the respective grid cell inside the topography
208            i = surf_def_h(0)%i(m) + surf_def_h(0)%ioff
209            j = surf_def_h(0)%j(m) + surf_def_h(0)%joff
210            k = surf_def_h(0)%k(m) + surf_def_h(0)%koff
211!
212!--         Check if the vertices that define the surface element are already
213!--         defined, if not, increment the counter.
214            IF ( point_index(k,j,i) < 0 )  THEN
215               surfaces%npoints   = surfaces%npoints + 1   
216               point_index(k,j,i) = surfaces%npoints - 1
217            ENDIF
218            IF ( point_index(k,j,i+1) < 0 )  THEN
219               surfaces%npoints     = surfaces%npoints + 1   
220               point_index(k,j,i+1) = surfaces%npoints - 1
221            ENDIF
222            IF ( point_index(k,j+1,i+1) < 0 )  THEN
223               surfaces%npoints       = surfaces%npoints + 1   
224               point_index(k,j+1,i+1) = surfaces%npoints - 1
225            ENDIF
226            IF ( point_index(k,j+1,i) < 0 )  THEN
227               surfaces%npoints     = surfaces%npoints + 1   
228               point_index(k,j+1,i) = surfaces%npoints - 1
229            ENDIF           
230         ENDDO
231      ENDDO
232      DO  m = 1, surf_lsm_h%ns
233         i = surf_lsm_h%i(m) + surf_lsm_h%ioff
234         j = surf_lsm_h%j(m) + surf_lsm_h%joff
235         k = surf_lsm_h%k(m) + surf_lsm_h%koff
236
237         IF ( point_index(k,j,i) < 0 )  THEN
238            surfaces%npoints   = surfaces%npoints + 1   
239            point_index(k,j,i) = surfaces%npoints - 1
240         ENDIF
241         IF ( point_index(k,j,i+1) < 0 )  THEN
242            surfaces%npoints     = surfaces%npoints + 1   
243            point_index(k,j,i+1) = surfaces%npoints - 1
244         ENDIF
245         IF ( point_index(k,j+1,i+1) < 0 )  THEN
246            surfaces%npoints       = surfaces%npoints + 1   
247            point_index(k,j+1,i+1) = surfaces%npoints - 1
248         ENDIF
249         IF ( point_index(k,j+1,i) < 0 )  THEN
250            surfaces%npoints     = surfaces%npoints + 1   
251            point_index(k,j+1,i) = surfaces%npoints - 1
252         ENDIF   
253      ENDDO
254      DO  m = 1, surf_usm_h%ns
255         i = surf_usm_h%i(m) + surf_usm_h%ioff
256         j = surf_usm_h%j(m) + surf_usm_h%joff
257         k = surf_usm_h%k(m) + surf_usm_h%koff
258
259         IF ( point_index(k,j,i) < 0 )  THEN
260            surfaces%npoints   = surfaces%npoints + 1   
261            point_index(k,j,i) = surfaces%npoints - 1
262         ENDIF
263         IF ( point_index(k,j,i+1) < 0 )  THEN
264            surfaces%npoints     = surfaces%npoints + 1   
265            point_index(k,j,i+1) = surfaces%npoints - 1
266         ENDIF
267         IF ( point_index(k,j+1,i+1) < 0 )  THEN
268            surfaces%npoints       = surfaces%npoints + 1   
269            point_index(k,j+1,i+1) = surfaces%npoints - 1
270         ENDIF
271         IF ( point_index(k,j+1,i) < 0 )  THEN
272            surfaces%npoints     = surfaces%npoints + 1   
273            point_index(k,j+1,i) = surfaces%npoints - 1
274         ENDIF     
275      ENDDO
276!
277!--   Vertical surfaces
278      DO  l = 0, 3
279         DO  m = 1, surf_def_v(l)%ns
280!
281!--         Determine the indices of the respective grid cell inside the topography.
282!--         Please note, j-index for north-facing surfaces ( l==0 ) is
283!--         identical to the reference j-index outside the grid box.
284!--         Equivalent for east-facing surfaces and i-index.
285            i = surf_def_v(l)%i(m) + MERGE( surf_def_v(l)%ioff, 0, l == 3 )
286            j = surf_def_v(l)%j(m) + MERGE( surf_def_v(l)%joff, 0, l == 1 )
287            k = surf_def_v(l)%k(m) + surf_def_v(l)%koff
288!
289!--         Lower left /front vertex
290            IF ( point_index(k,j,i) < 0 )  THEN
291               surfaces%npoints   = surfaces%npoints + 1   
292               point_index(k,j,i) = surfaces%npoints - 1
293            ENDIF 
294!
295!--         Upper / lower right index for north- and south-facing surfaces
296            IF ( l == 0  .OR.  l == 1 )  THEN
297               IF ( point_index(k,j,i+1) < 0 )  THEN
298                  surfaces%npoints     = surfaces%npoints + 1   
299                  point_index(k,j,i+1) = surfaces%npoints - 1
300               ENDIF
301               IF ( point_index(k+1,j,i+1) < 0 )  THEN
302                  surfaces%npoints       = surfaces%npoints + 1   
303                  point_index(k+1,j,i+1) = surfaces%npoints - 1
304               ENDIF
305!
306!--         Upper / lower front index for east- and west-facing surfaces
307            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
308               IF ( point_index(k,j+1,i) < 0 )  THEN
309                  surfaces%npoints     = surfaces%npoints + 1   
310                  point_index(k,j+1,i) = surfaces%npoints - 1
311               ENDIF
312               IF ( point_index(k+1,j+1,i) < 0 )  THEN
313                  surfaces%npoints       = surfaces%npoints + 1   
314                  point_index(k+1,j+1,i) = surfaces%npoints - 1
315               ENDIF
316            ENDIF
317!
318!--         Upper left / front vertex
319            IF ( point_index(k+1,j,i) < 0 )  THEN
320               surfaces%npoints     = surfaces%npoints + 1   
321               point_index(k+1,j,i) = surfaces%npoints - 1
322            ENDIF
323         ENDDO
324         DO  m = 1, surf_lsm_v(l)%ns
325            i = surf_lsm_v(l)%i(m) + MERGE( surf_lsm_v(l)%ioff, 0, l == 3 )
326            j = surf_lsm_v(l)%j(m) + MERGE( surf_lsm_v(l)%joff, 0, l == 1 )
327            k = surf_lsm_v(l)%k(m) + surf_lsm_v(l)%koff
328!
329!--         Lower left /front vertex
330            IF ( point_index(k,j,i) < 0 )  THEN
331               surfaces%npoints   = surfaces%npoints + 1   
332               point_index(k,j,i) = surfaces%npoints - 1
333            ENDIF 
334!
335!--         Upper / lower right index for north- and south-facing surfaces
336            IF ( l == 0  .OR.  l == 1 )  THEN
337               IF ( point_index(k,j,i+1) < 0 )  THEN
338                  surfaces%npoints     = surfaces%npoints + 1   
339                  point_index(k,j,i+1) = surfaces%npoints - 1
340               ENDIF
341               IF ( point_index(k+1,j,i+1) < 0 )  THEN
342                  surfaces%npoints       = surfaces%npoints + 1   
343                  point_index(k+1,j,i+1) = surfaces%npoints - 1
344               ENDIF
345!
346!--         Upper / lower front index for east- and west-facing surfaces
347            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
348               IF ( point_index(k,j+1,i) < 0 )  THEN
349                  surfaces%npoints     = surfaces%npoints + 1   
350                  point_index(k,j+1,i) = surfaces%npoints - 1
351               ENDIF
352               IF ( point_index(k+1,j+1,i) < 0 )  THEN
353                  surfaces%npoints       = surfaces%npoints + 1   
354                  point_index(k+1,j+1,i) = surfaces%npoints - 1
355               ENDIF
356            ENDIF
357!
358!--         Upper left / front vertex
359            IF ( point_index(k+1,j,i) < 0 )  THEN
360               surfaces%npoints     = surfaces%npoints + 1   
361               point_index(k+1,j,i) = surfaces%npoints - 1
362            ENDIF
363         ENDDO
364
365         DO  m = 1, surf_usm_v(l)%ns
366            i = surf_usm_v(l)%i(m) + MERGE( surf_usm_v(l)%ioff, 0, l == 3 )
367            j = surf_usm_v(l)%j(m) + MERGE( surf_usm_v(l)%joff, 0, l == 1 )
368            k = surf_usm_v(l)%k(m) + surf_usm_v(l)%koff
369!
370!--         Lower left /front vertex
371            IF ( point_index(k,j,i) < 0 )  THEN
372               surfaces%npoints   = surfaces%npoints + 1   
373               point_index(k,j,i) = surfaces%npoints - 1
374            ENDIF 
375!
376!--         Upper / lower right index for north- and south-facing surfaces
377            IF ( l == 0  .OR.  l == 1 )  THEN
378               IF ( point_index(k,j,i+1) < 0 )  THEN
379                  surfaces%npoints     = surfaces%npoints + 1   
380                  point_index(k,j,i+1) = surfaces%npoints - 1
381               ENDIF
382               IF ( point_index(k+1,j,i+1) < 0 )  THEN
383                  surfaces%npoints       = surfaces%npoints + 1   
384                  point_index(k+1,j,i+1) = surfaces%npoints - 1
385               ENDIF
386!
387!--         Upper / lower front index for east- and west-facing surfaces
388            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
389               IF ( point_index(k,j+1,i) < 0 )  THEN
390                  surfaces%npoints     = surfaces%npoints + 1   
391                  point_index(k,j+1,i) = surfaces%npoints - 1
392               ENDIF
393               IF ( point_index(k+1,j+1,i) < 0 )  THEN
394                  surfaces%npoints       = surfaces%npoints + 1   
395                  point_index(k+1,j+1,i) = surfaces%npoints - 1
396               ENDIF
397            ENDIF
398!
399!--         Upper left / front vertex
400            IF ( point_index(k+1,j,i) < 0 )  THEN
401               surfaces%npoints     = surfaces%npoints + 1   
402               point_index(k+1,j,i) = surfaces%npoints - 1
403            ENDIF
404         ENDDO
405
406      ENDDO
407!
408!--   Allocate the number of points and polygons. Note, the number of polygons
409!--   is identical to the number of surfaces elements, whereas the number
410!--   of points (vertices) which define the polygons can be larger.
411      ALLOCATE( surfaces%points(3,1:surfaces%npoints) )
412      ALLOCATE( surfaces%polygons(5,1:surfaces%ns)    )
413!
414!--   Note, PARAVIEW expects consecutively ordered points which can be
415!--   unambiguously identified.
416!--   Hence, all PEs know where they should start counting, depending on the
417!--   number of points on the other PE's with lower MPI rank.
418#if defined( __parallel )
419      CALL MPI_ALLGATHER( surfaces%npoints, 1, MPI_INTEGER,                    &
420                          num_points_on_pe, 1, MPI_INTEGER, comm2d, ierr  )
421#else
422      num_points_on_pe = surfaces%npoints
423#endif
424
425!
426!--   After the number of vertices is counted, repeat the loops and define the
427!--   vertices. Start with the horizontal default surfaces
428!--   First, however, determine the offset where couting of points should be
429!--   started, which is the sum of points of all PE's with lower MPI rank.
430      i                 = 0
431      point_index_count = 0
432      DO WHILE ( i < myid  .AND.  i <= SIZE( num_points_on_pe ) )
433         point_index_count = point_index_count + num_points_on_pe(i)
434         i                 = i + 1
435      ENDDO
436         
437     
438      surfaces%npoints = 0
439      point_index      = -1
440      npg              = 0
441     
442      DO  l = 0, 1
443         DO  m = 1, surf_def_h(0)%ns
444!
445!--         Determine the indices of the respective grid cell inside the topography
446            i = surf_def_h(0)%i(m) + surf_def_h(0)%ioff
447            j = surf_def_h(0)%j(m) + surf_def_h(0)%joff
448            k = surf_def_h(0)%k(m) + surf_def_h(0)%koff
449!
450!--         Check if the vertices that define the surface element are already
451!--         defined, if not, increment the counter.
452            IF ( point_index(k,j,i) < 0 )  THEN
453               surfaces%npoints   = surfaces%npoints + 1   
454               point_index(k,j,i) = point_index_count
455               point_index_count  = point_index_count + 1               
456               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
457               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
458               surfaces%points(3,surfaces%npoints) = zw(k)
459            ENDIF
460            IF ( point_index(k,j,i+1) < 0 )  THEN
461               surfaces%npoints     = surfaces%npoints + 1   
462               point_index(k,j,i+1) = point_index_count
463               point_index_count    = point_index_count + 1 
464               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
465               surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
466               surfaces%points(3,surfaces%npoints) = zw(k)
467            ENDIF
468            IF ( point_index(k,j+1,i+1) < 0 )  THEN
469               surfaces%npoints       = surfaces%npoints + 1   
470               point_index(k,j+1,i+1) = point_index_count
471               point_index_count      = point_index_count + 1 
472               surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
473               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
474               surfaces%points(3,surfaces%npoints) = zw(k)
475            ENDIF
476            IF ( point_index(k,j+1,i) < 0 )  THEN
477               surfaces%npoints     = surfaces%npoints + 1   
478               point_index(k,j+1,i) = point_index_count
479               point_index_count    = point_index_count + 1 
480               surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
481               surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
482               surfaces%points(3,surfaces%npoints) = zw(k)
483            ENDIF
484           
485            npg                        = npg + 1
486            surfaces%polygons(1,npg)   = 4
487            surfaces%polygons(2,npg)   = point_index(k,j,i)
488            surfaces%polygons(3,npg)   = point_index(k,j,i+1)
489            surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
490            surfaces%polygons(5,npg)   = point_index(k,j+1,i)
491         ENDDO
492      ENDDO
493      DO  m = 1, surf_lsm_h%ns
494         i = surf_lsm_h%i(m) + surf_lsm_h%ioff
495         j = surf_lsm_h%j(m) + surf_lsm_h%joff
496         k = surf_lsm_h%k(m) + surf_lsm_h%koff
497         IF ( point_index(k,j,i) < 0 )  THEN
498            surfaces%npoints   = surfaces%npoints + 1   
499            point_index(k,j,i) = point_index_count
500            point_index_count  = point_index_count + 1 
501            surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
502            surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
503            surfaces%points(3,surfaces%npoints) = zw(k)
504         ENDIF
505         IF ( point_index(k,j,i+1) < 0 )  THEN
506            surfaces%npoints     = surfaces%npoints + 1   
507            point_index(k,j,i+1) = point_index_count
508            point_index_count    = point_index_count + 1 
509            surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
510            surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
511            surfaces%points(3,surfaces%npoints) = zw(k)
512         ENDIF
513         IF ( point_index(k,j+1,i+1) < 0 )  THEN
514            surfaces%npoints       = surfaces%npoints + 1   
515            point_index(k,j+1,i+1) = point_index_count
516            point_index_count      = point_index_count + 1 
517            surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
518            surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
519            surfaces%points(3,surfaces%npoints) = zw(k)
520         ENDIF
521         IF ( point_index(k,j+1,i) < 0 )  THEN
522            surfaces%npoints     = surfaces%npoints + 1   
523            point_index(k,j+1,i) = point_index_count
524            point_index_count    = point_index_count + 1 
525            surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
526            surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
527            surfaces%points(3,surfaces%npoints) = zw(k)
528         ENDIF
529         
530         npg                        = npg + 1
531         surfaces%polygons(1,npg)   = 4
532         surfaces%polygons(2,npg)   = point_index(k,j,i)
533         surfaces%polygons(3,npg)   = point_index(k,j,i+1)
534         surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
535         surfaces%polygons(5,npg)   = point_index(k,j+1,i) 
536      ENDDO
537
538      DO  m = 1, surf_usm_h%ns
539         i = surf_usm_h%i(m) + surf_usm_h%ioff
540         j = surf_usm_h%j(m) + surf_usm_h%joff
541         k = surf_usm_h%k(m) + surf_usm_h%koff
542
543         IF ( point_index(k,j,i) < 0 )  THEN
544            surfaces%npoints   = surfaces%npoints + 1   
545            point_index(k,j,i) = point_index_count
546            point_index_count  = point_index_count + 1 
547            surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
548            surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
549            surfaces%points(3,surfaces%npoints) = zw(k)
550         ENDIF
551         IF ( point_index(k,j,i+1) < 0 )  THEN
552            surfaces%npoints     = surfaces%npoints + 1   
553            point_index(k,j,i+1) = point_index_count
554            point_index_count    = point_index_count + 1 
555            surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
556            surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
557            surfaces%points(3,surfaces%npoints) = zw(k)
558         ENDIF
559         IF ( point_index(k,j+1,i+1) < 0 )  THEN
560            surfaces%npoints       = surfaces%npoints + 1   
561            point_index(k,j+1,i+1) = point_index_count
562            point_index_count      = point_index_count + 1 
563            surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
564            surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
565            surfaces%points(3,surfaces%npoints) = zw(k)
566         ENDIF
567         IF ( point_index(k,j+1,i) < 0 )  THEN
568            surfaces%npoints     = surfaces%npoints + 1   
569            point_index(k,j+1,i) = point_index_count
570            point_index_count    = point_index_count + 1 
571            surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
572            surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
573            surfaces%points(3,surfaces%npoints) = zw(k)
574         ENDIF
575         
576         npg                        = npg + 1
577         surfaces%polygons(1,npg)   = 4
578         surfaces%polygons(2,npg)   = point_index(k,j,i)
579         surfaces%polygons(3,npg)   = point_index(k,j,i+1)
580         surfaces%polygons(4,npg)   = point_index(k,j+1,i+1)
581         surfaces%polygons(5,npg)   = point_index(k,j+1,i)   
582      ENDDO
583
584      DO  l = 0, 3
585         DO  m = 1, surf_def_v(l)%ns
586!
587!--         Determine the indices of the respective grid cell inside the topography.
588!--         Please note, j-index for north-facing surfaces ( l==0 ) is
589!--         identical to the reference j-index outside the grid box.
590!--         Equivalent for east-facing surfaces and i-index.
591            i = surf_def_v(l)%i(m) + MERGE( surf_def_v(l)%ioff, 0, l == 3 )
592            j = surf_def_v(l)%j(m) + MERGE( surf_def_v(l)%joff, 0, l == 1 )
593            k = surf_def_v(l)%k(m) + surf_def_v(l)%koff
594!
595!--         Lower left /front vertex
596            IF ( point_index(k,j,i) < 0 )  THEN
597               surfaces%npoints   = surfaces%npoints + 1   
598               point_index(k,j,i) = point_index_count
599               point_index_count  = point_index_count + 1 
600               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
601               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
602               surfaces%points(3,surfaces%npoints) = zw(k-1)
603            ENDIF 
604!
605!--         Upper / lower right index for north- and south-facing surfaces
606            IF ( l == 0  .OR.  l == 1 )  THEN
607               IF ( point_index(k,j,i+1) < 0 )  THEN
608                  surfaces%npoints     = surfaces%npoints + 1   
609                  point_index(k,j,i+1) = point_index_count
610                  point_index_count    = point_index_count + 1
611                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
612                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
613                  surfaces%points(3,surfaces%npoints) = zw(k-1)
614               ENDIF
615               IF ( point_index(k+1,j,i+1) < 0 )  THEN
616                  surfaces%npoints       = surfaces%npoints + 1   
617                  point_index(k+1,j,i+1) = point_index_count
618                  point_index_count      = point_index_count + 1
619                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
620                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
621                  surfaces%points(3,surfaces%npoints) = zw(k)
622               ENDIF
623!
624!--         Upper / lower front index for east- and west-facing surfaces
625            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
626               IF ( point_index(k,j+1,i) < 0 )  THEN
627                  surfaces%npoints     = surfaces%npoints + 1   
628                  point_index(k,j+1,i) = point_index_count
629                  point_index_count    = point_index_count + 1
630                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
631                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
632                  surfaces%points(3,surfaces%npoints) = zw(k-1)
633               ENDIF
634               IF ( point_index(k+1,j+1,i) < 0 )  THEN
635                  surfaces%npoints       = surfaces%npoints + 1   
636                  point_index(k+1,j+1,i) = point_index_count
637                  point_index_count      = point_index_count + 1
638                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
639                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
640                  surfaces%points(3,surfaces%npoints) = zw(k)
641               ENDIF
642            ENDIF
643!
644!--         Upper left / front vertex
645            IF ( point_index(k+1,j,i) < 0 )  THEN
646               surfaces%npoints     = surfaces%npoints + 1   
647               point_index(k+1,j,i) = point_index_count
648               point_index_count    = point_index_count + 1
649               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
650               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
651               surfaces%points(3,surfaces%npoints) = zw(k)
652            ENDIF
653           
654            npg = npg + 1           
655            IF ( l == 0  .OR.  l == 1 )  THEN
656               surfaces%polygons(1,npg)   = 4
657               surfaces%polygons(2,npg)   = point_index(k,j,i)
658               surfaces%polygons(3,npg)   = point_index(k,j,i+1)
659               surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
660               surfaces%polygons(5,npg)   = point_index(k+1,j,i) 
661            ELSE
662               surfaces%polygons(1,npg)   = 4
663               surfaces%polygons(2,npg)   = point_index(k,j,i)
664               surfaces%polygons(3,npg)   = point_index(k,j+1,i)
665               surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
666               surfaces%polygons(5,npg)   = point_index(k+1,j,i)
667            ENDIF
668           
669         ENDDO
670
671         DO  m = 1, surf_lsm_v(l)%ns
672            i = surf_lsm_v(l)%i(m) + MERGE( surf_lsm_v(l)%ioff, 0, l == 3 )
673            j = surf_lsm_v(l)%j(m) + MERGE( surf_lsm_v(l)%joff, 0, l == 1 )
674            k = surf_lsm_v(l)%k(m) + surf_lsm_v(l)%koff
675!
676!--         Lower left /front vertex
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-1)
684            ENDIF 
685!
686!--         Upper / lower right index for north- and south-facing surfaces
687            IF ( l == 0  .OR.  l == 1 )  THEN
688               IF ( point_index(k,j,i+1) < 0 )  THEN
689                  surfaces%npoints     = surfaces%npoints + 1   
690                  point_index(k,j,i+1) = point_index_count
691                  point_index_count    = point_index_count + 1
692                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
693                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
694                  surfaces%points(3,surfaces%npoints) = zw(k-1)
695               ENDIF
696               IF ( point_index(k+1,j,i+1) < 0 )  THEN
697                  surfaces%npoints       = surfaces%npoints + 1   
698                  point_index(k+1,j,i+1) = point_index_count
699                  point_index_count      = point_index_count + 1
700                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
701                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
702                  surfaces%points(3,surfaces%npoints) = zw(k)
703               ENDIF
704!
705!--         Upper / lower front index for east- and west-facing surfaces
706            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
707               IF ( point_index(k,j+1,i) < 0 )  THEN
708                  surfaces%npoints     = surfaces%npoints + 1   
709                  point_index(k,j+1,i) = point_index_count
710                  point_index_count    = point_index_count + 1
711                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
712                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
713                  surfaces%points(3,surfaces%npoints) = zw(k-1)
714               ENDIF
715               IF ( point_index(k+1,j+1,i) < 0 )  THEN
716                  surfaces%npoints       = surfaces%npoints + 1   
717                  point_index(k+1,j+1,i) = point_index_count
718                  point_index_count      = point_index_count + 1
719                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
720                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
721                  surfaces%points(3,surfaces%npoints) = zw(k)
722               ENDIF
723            ENDIF
724!
725!--         Upper left / front vertex
726            IF ( point_index(k+1,j,i) < 0 )  THEN
727               surfaces%npoints     = surfaces%npoints + 1   
728               point_index(k+1,j,i) = point_index_count
729               point_index_count    = point_index_count + 1
730               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
731               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
732               surfaces%points(3,surfaces%npoints) = zw(k)
733            ENDIF
734           
735            npg = npg + 1           
736            IF ( l == 0  .OR.  l == 1 )  THEN
737               surfaces%polygons(1,npg)   = 4
738               surfaces%polygons(2,npg)   = point_index(k,j,i)
739               surfaces%polygons(3,npg)   = point_index(k,j,i+1)
740               surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
741               surfaces%polygons(5,npg)   = point_index(k+1,j,i) 
742            ELSE
743               surfaces%polygons(1,npg)   = 4
744               surfaces%polygons(2,npg)   = point_index(k,j,i)
745               surfaces%polygons(3,npg)   = point_index(k,j+1,i)
746               surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
747               surfaces%polygons(5,npg)   = point_index(k+1,j,i)
748            ENDIF
749         ENDDO
750         DO  m = 1, surf_usm_v(l)%ns
751            i = surf_usm_v(l)%i(m) + MERGE( surf_usm_v(l)%ioff, 0, l == 3 )
752            j = surf_usm_v(l)%j(m) + MERGE( surf_usm_v(l)%joff, 0, l == 1 )
753            k = surf_usm_v(l)%k(m) + surf_usm_v(l)%koff
754!
755!--         Lower left /front vertex
756            IF ( point_index(k,j,i) < 0 )  THEN
757               surfaces%npoints   = surfaces%npoints + 1   
758               point_index(k,j,i) = point_index_count
759               point_index_count  = point_index_count + 1
760               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
761               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
762               surfaces%points(3,surfaces%npoints) = zw(k-1)
763            ENDIF 
764!
765!--         Upper / lower right index for north- and south-facing surfaces
766            IF ( l == 0  .OR.  l == 1 )  THEN
767               IF ( point_index(k,j,i+1) < 0 )  THEN
768                  surfaces%npoints     = surfaces%npoints + 1   
769                  point_index(k,j,i+1) = point_index_count
770                  point_index_count    = point_index_count + 1
771                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
772                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
773                  surfaces%points(3,surfaces%npoints) = zw(k-1)
774               ENDIF
775               IF ( point_index(k+1,j,i+1) < 0 )  THEN
776                  surfaces%npoints       = surfaces%npoints + 1   
777                  point_index(k+1,j,i+1) = point_index_count
778                  point_index_count      = point_index_count + 1
779                  surfaces%points(1,surfaces%npoints) = ( i + 1 - 0.5_wp ) * dx
780                  surfaces%points(2,surfaces%npoints) = ( j     - 0.5_wp ) * dy
781                  surfaces%points(3,surfaces%npoints) = zw(k)
782               ENDIF
783!
784!--         Upper / lower front index for east- and west-facing surfaces
785            ELSEIF ( l == 2  .OR.  l == 3 )  THEN
786               IF ( point_index(k,j+1,i) < 0 )  THEN
787                  surfaces%npoints     = surfaces%npoints + 1   
788                  point_index(k,j+1,i) = point_index_count
789                  point_index_count    = point_index_count + 1
790                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
791                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
792                  surfaces%points(3,surfaces%npoints) = zw(k-1)
793               ENDIF
794               IF ( point_index(k+1,j+1,i) < 0 )  THEN
795                  surfaces%npoints       = surfaces%npoints + 1   
796                  point_index(k+1,j+1,i) = point_index_count
797                  point_index_count      = point_index_count + 1
798                  surfaces%points(1,surfaces%npoints) = ( i     - 0.5_wp ) * dx
799                  surfaces%points(2,surfaces%npoints) = ( j + 1 - 0.5_wp ) * dy
800                  surfaces%points(3,surfaces%npoints) = zw(k)
801               ENDIF
802            ENDIF
803!
804!--         Upper left / front vertex
805            IF ( point_index(k+1,j,i) < 0 )  THEN
806               surfaces%npoints     = surfaces%npoints + 1   
807               point_index(k+1,j,i) = point_index_count
808               point_index_count    = point_index_count + 1
809               surfaces%points(1,surfaces%npoints) = ( i - 0.5_wp ) * dx
810               surfaces%points(2,surfaces%npoints) = ( j - 0.5_wp ) * dy
811               surfaces%points(3,surfaces%npoints) = zw(k)
812            ENDIF
813           
814            npg = npg + 1           
815            IF ( l == 0  .OR.  l == 1 )  THEN
816               surfaces%polygons(1,npg)   = 4
817               surfaces%polygons(2,npg)   = point_index(k,j,i)
818               surfaces%polygons(3,npg)   = point_index(k,j,i+1)
819               surfaces%polygons(4,npg)   = point_index(k+1,j,i+1)
820               surfaces%polygons(5,npg)   = point_index(k+1,j,i) 
821            ELSE
822               surfaces%polygons(1,npg)   = 4
823               surfaces%polygons(2,npg)   = point_index(k,j,i)
824               surfaces%polygons(3,npg)   = point_index(k,j+1,i)
825               surfaces%polygons(4,npg)   = point_index(k+1,j+1,i)
826               surfaces%polygons(5,npg)   = point_index(k+1,j,i)
827            ENDIF
828         ENDDO
829
830      ENDDO
831!
832!--   Deallocate temporary dummy variable
833      DEALLOCATE ( point_index )
834!
835!--   Sum-up total number of surface elements and vertices on domain. This
836!--   will be needed for post-processing.
837      surfaces%npoints_total = 0
838#if defined( __parallel )
839       CALL MPI_ALLREDUCE( surfaces%npoints, surfaces%npoints_total, 1,        &
840                           MPI_INTEGER, MPI_SUM, comm2d, ierr )
841       CALL MPI_ALLREDUCE( surfaces%ns, surfaces%ns_total, 1,                  &
842                           MPI_INTEGER, MPI_SUM, comm2d, ierr )
843#else
844       surfaces%npoints_total = surfaces%npoints
845       surfaces%ns_total      = surfaces%ns
846#endif
847
848   END SUBROUTINE surface_output_init
849   
850!------------------------------------------------------------------------------!
851! Description:
852! ------------
853!> Routine for controlling the data output. Surface data is collected from
854!> different types of surfaces (default, natural, urban) and different
855!> orientation and written to one 1D-output array. Further, NetCDF routines
856!> are called to write the surface data in the respective NetCDF files.
857!------------------------------------------------------------------------------!   
858   SUBROUTINE surface_output( av )
859   
860      USE control_parameters,                                                  &
861          ONLY:  io_blocks, io_group, simulated_time
862
863      USE pegrid,                                                              &
864          ONLY:  comm2d, ierr
865
866   
867      IMPLICIT NONE
868
869      CHARACTER(LEN=100) ::  trimvar = ' ' !< dummy for single output variable
870     
871      INTEGER(iwp) ::  av     !< id indicating average or non-average data output
872      INTEGER(iwp) ::  i      !< loop index
873      INTEGER(iwp) ::  l      !< running index for surface-element orientation
874      INTEGER(iwp) ::  m      !< running index for surface elements
875      INTEGER(iwp) ::  n_out  !< counter variables for surface output
876
877!
878!--   Return, if nothing to output
879      IF ( dosurf_no(av) == 0 )  RETURN
880!
881!--   Open files
882      CALL check_open( 25+av )
883!
884!--   Write coordinates
885      IF ( .NOT. first_output(av) )  THEN
886         DO  i = 0, io_blocks-1
887            IF ( i == io_group )  THEN
888               WRITE ( 25+av )  surfaces%npoints
889               WRITE ( 25+av )  surfaces%npoints_total
890               WRITE ( 25+av )  surfaces%ns
891               WRITE ( 25+av )  surfaces%ns_total
892               WRITE ( 25+av )  surfaces%points
893               WRITE ( 25+av )  surfaces%polygons
894            ENDIF
895#if defined( __parallel )
896            CALL MPI_BARRIER( comm2d, ierr )
897#endif
898            first_output(av) = .TRUE.
899         ENDDO
900      ENDIF
901
902      n_out = 0
903      DO  WHILE ( dosurf(av,n_out+1)(1:1) /= ' ' )
904
905         n_out   = n_out + 1
906         trimvar = TRIM( dosurf(av,n_out) )
907!
908!--      Set the output array to the _FillValue in case it is not
909!--      defined for each type of surface.
910         surfaces%var_out = surfaces%fillvalue
911         SELECT CASE ( trimvar )
912
913            CASE ( 'us' )
914!
915!--            Output of instantaneous data
916               IF ( av == 0 )  THEN
917                  CALL surface_output_collect( surf_def_h(0)%us,               &
918                                               surf_def_h(1)%us,               &
919                                               surf_lsm_h%us,                  &
920                                               surf_usm_h%us,                  &
921                                               surf_def_v(0)%us,               &
922                                               surf_lsm_v(0)%us,               &
923                                               surf_usm_v(0)%us,               &
924                                               surf_def_v(1)%us,               &
925                                               surf_lsm_v(1)%us,               &
926                                               surf_usm_v(1)%us,               &
927                                               surf_def_v(2)%us,               &
928                                               surf_lsm_v(2)%us,               &
929                                               surf_usm_v(2)%us,               &
930                                               surf_def_v(3)%us,               &
931                                               surf_lsm_v(3)%us,               &
932                                               surf_usm_v(3)%us ) 
933               ELSE
934!
935!--               Output of averaged data
936                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
937                                        REAL( average_count_surf, KIND=wp )
938                  surfaces%var_av(:,n_out) = 0.0_wp
939                                                       
940               ENDIF
941
942            CASE ( 'ts' )
943!
944!--            Output of instantaneous data
945               IF ( av == 0 )  THEN
946                  CALL surface_output_collect( surf_def_h(0)%ts,               &
947                                               surf_def_h(1)%ts,               &
948                                               surf_lsm_h%ts,                  &
949                                               surf_usm_h%ts,                  &
950                                               surf_def_v(0)%ts,               &
951                                               surf_lsm_v(0)%ts,               &
952                                               surf_usm_v(0)%ts,               &
953                                               surf_def_v(1)%ts,               &
954                                               surf_lsm_v(1)%ts,               &
955                                               surf_usm_v(1)%ts,               &
956                                               surf_def_v(2)%ts,               &
957                                               surf_lsm_v(2)%ts,               &
958                                               surf_usm_v(2)%ts,               &
959                                               surf_def_v(3)%ts,               &
960                                               surf_lsm_v(3)%ts,               &
961                                               surf_usm_v(3)%ts ) 
962               ELSE
963!
964!--               Output of averaged data
965                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
966                                        REAL( average_count_surf, KIND=wp )
967                  surfaces%var_av(:,n_out) = 0.0_wp
968                                                       
969               ENDIF
970
971            CASE ( 'qs' )
972!
973!--            Output of instantaneous data
974               IF ( av == 0 )  THEN
975                  CALL surface_output_collect( surf_def_h(0)%qs,               &
976                                               surf_def_h(1)%qs,               &
977                                               surf_lsm_h%qs,                  &
978                                               surf_usm_h%qs,                  &
979                                               surf_def_v(0)%qs,               &
980                                               surf_lsm_v(0)%qs,               &
981                                               surf_usm_v(0)%qs,               &
982                                               surf_def_v(1)%qs,               &
983                                               surf_lsm_v(1)%qs,               &
984                                               surf_usm_v(1)%qs,               &
985                                               surf_def_v(2)%qs,               &
986                                               surf_lsm_v(2)%qs,               &
987                                               surf_usm_v(2)%qs,               &
988                                               surf_def_v(3)%qs,               &
989                                               surf_lsm_v(3)%qs,               &
990                                               surf_usm_v(3)%qs ) 
991               ELSE
992!
993!--               Output of averaged data
994                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
995                                        REAL( average_count_surf, KIND=wp )
996                  surfaces%var_av(:,n_out) = 0.0_wp
997                                                       
998               ENDIF
999
1000            CASE ( 'ss' )
1001!
1002!--            Output of instantaneous data
1003               IF ( av == 0 )  THEN
1004                  CALL surface_output_collect( surf_def_h(0)%ss,               &
1005                                               surf_def_h(1)%ss,               &
1006                                               surf_lsm_h%ss,                  &
1007                                               surf_usm_h%ss,                  &
1008                                               surf_def_v(0)%ss,               &
1009                                               surf_lsm_v(0)%ss,               &
1010                                               surf_usm_v(0)%ss,               &
1011                                               surf_def_v(1)%ss,               &
1012                                               surf_lsm_v(1)%ss,               &
1013                                               surf_usm_v(1)%ss,               &
1014                                               surf_def_v(2)%ss,               &
1015                                               surf_lsm_v(2)%ss,               &
1016                                               surf_usm_v(2)%ss,               &
1017                                               surf_def_v(3)%ss,               &
1018                                               surf_lsm_v(3)%ss,               &
1019                                               surf_usm_v(3)%ss ) 
1020               ELSE
1021!
1022!--               Output of averaged data
1023                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1024                                        REAL( average_count_surf, KIND=wp )
1025                  surfaces%var_av(:,n_out) = 0.0_wp
1026                                                       
1027               ENDIF
1028
1029            CASE ( 'qcs' )
1030!
1031!--            Output of instantaneous data
1032               IF ( av == 0 )  THEN
1033                  CALL surface_output_collect( surf_def_h(0)%qcs,              &
1034                                               surf_def_h(1)%qcs,              &
1035                                               surf_lsm_h%qcs,                 &
1036                                               surf_usm_h%qcs,                 &
1037                                               surf_def_v(0)%qcs,              &
1038                                               surf_lsm_v(0)%qcs,              &
1039                                               surf_usm_v(0)%qcs,              &
1040                                               surf_def_v(1)%qcs,              &
1041                                               surf_lsm_v(1)%qcs,              &
1042                                               surf_usm_v(1)%qcs,              &
1043                                               surf_def_v(2)%qcs,              &
1044                                               surf_lsm_v(2)%qcs,              &
1045                                               surf_usm_v(2)%qcs,              &
1046                                               surf_def_v(3)%qcs,              &
1047                                               surf_lsm_v(3)%qcs,              &
1048                                               surf_usm_v(3)%qcs ) 
1049               ELSE
1050!
1051!--               Output of averaged data
1052                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1053                                        REAL( average_count_surf, KIND=wp )
1054                  surfaces%var_av(:,n_out) = 0.0_wp
1055                                                       
1056               ENDIF
1057
1058            CASE ( 'ncs' )
1059!
1060!--            Output of instantaneous data
1061               IF ( av == 0 )  THEN
1062                  CALL surface_output_collect( surf_def_h(0)%ncs,              &
1063                                               surf_def_h(1)%ncs,              &
1064                                               surf_lsm_h%ncs,                 &
1065                                               surf_usm_h%ncs,                 &
1066                                               surf_def_v(0)%ncs,              &
1067                                               surf_lsm_v(0)%ncs,              &
1068                                               surf_usm_v(0)%ncs,              &
1069                                               surf_def_v(1)%ncs,              &
1070                                               surf_lsm_v(1)%ncs,              &
1071                                               surf_usm_v(1)%ncs,              &
1072                                               surf_def_v(2)%ncs,              &
1073                                               surf_lsm_v(2)%ncs,              &
1074                                               surf_usm_v(2)%ncs,              &
1075                                               surf_def_v(3)%ncs,              &
1076                                               surf_lsm_v(3)%ncs,              &
1077                                               surf_usm_v(3)%ncs ) 
1078               ELSE
1079!
1080!--               Output of averaged data
1081                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1082                                        REAL( average_count_surf, KIND=wp )
1083                  surfaces%var_av(:,n_out) = 0.0_wp
1084                                                       
1085               ENDIF
1086
1087            CASE ( 'qrs' )
1088!
1089!--            Output of instantaneous data
1090               IF ( av == 0 )  THEN
1091                  CALL surface_output_collect( surf_def_h(0)%qrs,              &
1092                                               surf_def_h(1)%qrs,              &
1093                                               surf_lsm_h%qrs,                 &
1094                                               surf_usm_h%qrs,                 &
1095                                               surf_def_v(0)%qrs,              &
1096                                               surf_lsm_v(0)%qrs,              &
1097                                               surf_usm_v(0)%qrs,              &
1098                                               surf_def_v(1)%qrs,              &
1099                                               surf_lsm_v(1)%qrs,              &
1100                                               surf_usm_v(1)%qrs,              &
1101                                               surf_def_v(2)%qrs,              &
1102                                               surf_lsm_v(2)%qrs,              &
1103                                               surf_usm_v(2)%qrs,              &
1104                                               surf_def_v(3)%qrs,              &
1105                                               surf_lsm_v(3)%qrs,              &
1106                                               surf_usm_v(3)%qrs ) 
1107               ELSE
1108!
1109!--               Output of averaged data
1110                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1111                                        REAL( average_count_surf, KIND=wp )
1112                  surfaces%var_av(:,n_out) = 0.0_wp
1113                                                       
1114               ENDIF
1115
1116            CASE ( 'nrs' )
1117!
1118!--            Output of instantaneous data
1119               IF ( av == 0 )  THEN
1120                  CALL surface_output_collect( surf_def_h(0)%nrs,              &
1121                                               surf_def_h(1)%nrs,              &
1122                                               surf_lsm_h%nrs,                 &
1123                                               surf_usm_h%nrs,                 &
1124                                               surf_def_v(0)%nrs,              &
1125                                               surf_lsm_v(0)%nrs,              &
1126                                               surf_usm_v(0)%nrs,              &
1127                                               surf_def_v(1)%nrs,              &
1128                                               surf_lsm_v(1)%nrs,              &
1129                                               surf_usm_v(1)%nrs,              &
1130                                               surf_def_v(2)%nrs,              &
1131                                               surf_lsm_v(2)%nrs,              &
1132                                               surf_usm_v(2)%nrs,              &
1133                                               surf_def_v(3)%nrs,              &
1134                                               surf_lsm_v(3)%nrs,              &
1135                                               surf_usm_v(3)%nrs ) 
1136               ELSE
1137!
1138!--               Output of averaged data
1139                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1140                                        REAL( average_count_surf, KIND=wp )
1141                  surfaces%var_av(:,n_out) = 0.0_wp
1142                                                       
1143               ENDIF
1144
1145            CASE ( 'ol' )
1146!
1147!--            Output of instantaneous data
1148               IF ( av == 0 )  THEN
1149                  CALL surface_output_collect( surf_def_h(0)%ol,               &
1150                                               surf_def_h(1)%ol,               &
1151                                               surf_lsm_h%ol,                  &
1152                                               surf_usm_h%ol,                  &
1153                                               surf_def_v(0)%ol,               &
1154                                               surf_lsm_v(0)%ol,               &
1155                                               surf_usm_v(0)%ol,               &
1156                                               surf_def_v(1)%ol,               &
1157                                               surf_lsm_v(1)%ol,               &
1158                                               surf_usm_v(1)%ol,               &
1159                                               surf_def_v(2)%ol,               &
1160                                               surf_lsm_v(2)%ol,               &
1161                                               surf_usm_v(2)%ol,               &
1162                                               surf_def_v(3)%ol,               &
1163                                               surf_lsm_v(3)%ol,               &
1164                                               surf_usm_v(3)%ol ) 
1165               ELSE
1166!
1167!--               Output of averaged data
1168                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1169                                        REAL( average_count_surf, KIND=wp )
1170                  surfaces%var_av(:,n_out) = 0.0_wp
1171                                                       
1172               ENDIF
1173
1174            CASE ( 'z0' )
1175!
1176!--            Output of instantaneous data
1177               IF ( av == 0 )  THEN
1178                  CALL surface_output_collect( surf_def_h(0)%z0,               &
1179                                               surf_def_h(1)%z0,               &
1180                                               surf_lsm_h%z0,                  &
1181                                               surf_usm_h%z0,                  &
1182                                               surf_def_v(0)%z0,               &
1183                                               surf_lsm_v(0)%z0,               &
1184                                               surf_usm_v(0)%z0,               &
1185                                               surf_def_v(1)%z0,               &
1186                                               surf_lsm_v(1)%z0,               &
1187                                               surf_usm_v(1)%z0,               &
1188                                               surf_def_v(2)%z0,               &
1189                                               surf_lsm_v(2)%z0,               &
1190                                               surf_usm_v(2)%z0,               &
1191                                               surf_def_v(3)%z0,               &
1192                                               surf_lsm_v(3)%z0,               &
1193                                               surf_usm_v(3)%z0 )
1194               ELSE
1195!
1196!--               Output of averaged data
1197                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1198                                        REAL( average_count_surf, KIND=wp )
1199                  surfaces%var_av(:,n_out) = 0.0_wp
1200                                                       
1201               ENDIF
1202
1203            CASE ( 'z0h' )
1204!
1205!--            Output of instantaneous data
1206               IF ( av == 0 )  THEN
1207                  CALL surface_output_collect( surf_def_h(0)%z0h,              &
1208                                               surf_def_h(1)%z0h,              &
1209                                               surf_lsm_h%z0h,                 &
1210                                               surf_usm_h%z0h,                 &
1211                                               surf_def_v(0)%z0h,              &
1212                                               surf_lsm_v(0)%z0h,              &
1213                                               surf_usm_v(0)%z0h,              &
1214                                               surf_def_v(1)%z0h,              &
1215                                               surf_lsm_v(1)%z0h,              &
1216                                               surf_usm_v(1)%z0h,              &
1217                                               surf_def_v(2)%z0h,              &
1218                                               surf_lsm_v(2)%z0h,              &
1219                                               surf_usm_v(2)%z0h,              &
1220                                               surf_def_v(3)%z0h,              &
1221                                               surf_lsm_v(3)%z0h,              &
1222                                               surf_usm_v(3)%z0h )
1223               ELSE
1224!
1225!--               Output of averaged data
1226                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1227                                        REAL( average_count_surf, KIND=wp )
1228                  surfaces%var_av(:,n_out) = 0.0_wp
1229                                                       
1230               ENDIF
1231
1232            CASE ( 'z0q' )
1233!
1234!--            Output of instantaneous data
1235               IF ( av == 0 )  THEN
1236                  CALL surface_output_collect( surf_def_h(0)%z0q,              &
1237                                               surf_def_h(1)%z0q,              &
1238                                               surf_lsm_h%z0q,                 &
1239                                               surf_usm_h%z0q,                 &
1240                                               surf_def_v(0)%z0q,              &
1241                                               surf_lsm_v(0)%z0q,              &
1242                                               surf_usm_v(0)%z0q,              &
1243                                               surf_def_v(1)%z0q,              &
1244                                               surf_lsm_v(1)%z0q,              &
1245                                               surf_usm_v(1)%z0q,              &
1246                                               surf_def_v(2)%z0q,              &
1247                                               surf_lsm_v(2)%z0q,              &
1248                                               surf_usm_v(2)%z0q,              &
1249                                               surf_def_v(3)%z0q,              &
1250                                               surf_lsm_v(3)%z0q,              &
1251                                               surf_usm_v(3)%z0q ) 
1252               ELSE
1253!
1254!--               Output of averaged data
1255                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1256                                        REAL( average_count_surf, KIND=wp )
1257                  surfaces%var_av(:,n_out) = 0.0_wp
1258
1259               ENDIF
1260
1261            CASE ( 'theta1' )
1262!
1263!--            Output of instantaneous data
1264               IF ( av == 0 )  THEN
1265                  CALL surface_output_collect( surf_def_h(0)%pt1,              &
1266                                               surf_def_h(1)%pt1,              &
1267                                               surf_lsm_h%pt1,                 &
1268                                               surf_usm_h%pt1,                 &
1269                                               surf_def_v(0)%pt1,              &
1270                                               surf_lsm_v(0)%pt1,              &
1271                                               surf_usm_v(0)%pt1,              &
1272                                               surf_def_v(1)%pt1,              &
1273                                               surf_lsm_v(1)%pt1,              &
1274                                               surf_usm_v(1)%pt1,              &
1275                                               surf_def_v(2)%pt1,              &
1276                                               surf_lsm_v(2)%pt1,              &
1277                                               surf_usm_v(2)%pt1,              &
1278                                               surf_def_v(3)%pt1,              &
1279                                               surf_lsm_v(3)%pt1,              &
1280                                               surf_usm_v(3)%pt1 )
1281               ELSE
1282!
1283!--               Output of averaged data
1284                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1285                                        REAL( average_count_surf, KIND=wp )
1286                  surfaces%var_av(:,n_out) = 0.0_wp
1287                                                       
1288               ENDIF
1289
1290            CASE ( 'qv1' )
1291!
1292!--            Output of instantaneous data
1293               IF ( av == 0 )  THEN
1294                  CALL surface_output_collect( surf_def_h(0)%qv1,              &
1295                                               surf_def_h(1)%qv1,              &
1296                                               surf_lsm_h%qv1,                 &
1297                                               surf_usm_h%qv1,                 &
1298                                               surf_def_v(0)%qv1,              &
1299                                               surf_lsm_v(0)%qv1,              &
1300                                               surf_usm_v(0)%qv1,              &
1301                                               surf_def_v(1)%qv1,              &
1302                                               surf_lsm_v(1)%qv1,              &
1303                                               surf_usm_v(1)%qv1,              &
1304                                               surf_def_v(2)%qv1,              &
1305                                               surf_lsm_v(2)%qv1,              &
1306                                               surf_usm_v(2)%qv1,              &
1307                                               surf_def_v(3)%qv1,              &
1308                                               surf_lsm_v(3)%qv1,              &
1309                                               surf_usm_v(3)%qv1 )
1310               ELSE
1311!
1312!--               Output of averaged data
1313                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1314                                        REAL( average_count_surf, KIND=wp )
1315                  surfaces%var_av(:,n_out) = 0.0_wp
1316                                                       
1317               ENDIF
1318
1319            CASE ( 'thetav1' )
1320!
1321!--            Output of instantaneous data
1322               IF ( av == 0 )  THEN
1323                  CALL surface_output_collect( surf_def_h(0)%vpt1,             &
1324                                               surf_def_h(1)%vpt1,             &
1325                                               surf_lsm_h%vpt1,                &
1326                                               surf_usm_h%vpt1,                &
1327                                               surf_def_v(0)%vpt1,             &
1328                                               surf_lsm_v(0)%vpt1,             &
1329                                               surf_usm_v(0)%vpt1,             &
1330                                               surf_def_v(1)%vpt1,             &
1331                                               surf_lsm_v(1)%vpt1,             &
1332                                               surf_usm_v(1)%vpt1,             &
1333                                               surf_def_v(2)%vpt1,             &
1334                                               surf_lsm_v(2)%vpt1,             &
1335                                               surf_usm_v(2)%vpt1,             &
1336                                               surf_def_v(3)%vpt1,             &
1337                                               surf_lsm_v(3)%vpt1,             &
1338                                               surf_usm_v(3)%vpt1 )
1339               ELSE
1340!
1341!--               Output of averaged data
1342                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1343                                        REAL( average_count_surf, KIND=wp )
1344                  surfaces%var_av(:,n_out) = 0.0_wp
1345                                                       
1346               ENDIF
1347
1348            CASE ( 'usws' )
1349!
1350!--            Output of instantaneous data
1351               IF ( av == 0 )  THEN
1352                  CALL surface_output_collect( surf_def_h(0)%usws,             &
1353                                               surf_def_h(1)%usws,             &
1354                                               surf_lsm_h%usws,                &
1355                                               surf_usm_h%usws,                &
1356                                               surf_def_v(0)%usws,             &
1357                                               surf_lsm_v(0)%usws,             &
1358                                               surf_usm_v(0)%usws,             &
1359                                               surf_def_v(1)%usws,             &
1360                                               surf_lsm_v(1)%usws,             &
1361                                               surf_usm_v(1)%usws,             &
1362                                               surf_def_v(2)%usws,             &
1363                                               surf_lsm_v(2)%usws,             &
1364                                               surf_usm_v(2)%usws,             &
1365                                               surf_def_v(3)%usws,             &
1366                                               surf_lsm_v(3)%usws,             &
1367                                               surf_usm_v(3)%usws )
1368               ELSE
1369!
1370!--               Output of averaged data
1371                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1372                                        REAL( average_count_surf, KIND=wp )
1373                  surfaces%var_av(:,n_out) = 0.0_wp
1374                                                       
1375               ENDIF
1376
1377            CASE ( 'vsws' )
1378!
1379!--            Output of instantaneous data
1380               IF ( av == 0 )  THEN
1381                  CALL surface_output_collect( surf_def_h(0)%vsws,             &
1382                                               surf_def_h(1)%vsws,             &
1383                                               surf_lsm_h%vsws,                &
1384                                               surf_usm_h%vsws,                &
1385                                               surf_def_v(0)%vsws,             &
1386                                               surf_lsm_v(0)%vsws,             &
1387                                               surf_usm_v(0)%vsws,             &
1388                                               surf_def_v(1)%vsws,             &
1389                                               surf_lsm_v(1)%vsws,             &
1390                                               surf_usm_v(1)%vsws,             &
1391                                               surf_def_v(2)%vsws,             &
1392                                               surf_lsm_v(2)%vsws,             &
1393                                               surf_usm_v(2)%vsws,             &
1394                                               surf_def_v(3)%vsws,             &
1395                                               surf_lsm_v(3)%vsws,             &
1396                                               surf_usm_v(3)%vsws )
1397               ELSE
1398!
1399!--               Output of averaged data
1400                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1401                                        REAL( average_count_surf, KIND=wp )
1402                  surfaces%var_av(:,n_out) = 0.0_wp
1403                                                       
1404               ENDIF
1405
1406            CASE ( 'shf' )
1407!
1408!--            Output of instantaneous data
1409               IF ( av == 0 )  THEN
1410                  CALL surface_output_collect( surf_def_h(0)%shf,              &
1411                                               surf_def_h(1)%shf,              &
1412                                               surf_lsm_h%shf,                 &
1413                                               surf_usm_h%shf,                 &
1414                                               surf_def_v(0)%shf,              &
1415                                               surf_lsm_v(0)%shf,              &
1416                                               surf_usm_v(0)%shf,              &
1417                                               surf_def_v(1)%shf,              &
1418                                               surf_lsm_v(1)%shf,              &
1419                                               surf_usm_v(1)%shf,              &
1420                                               surf_def_v(2)%shf,              &
1421                                               surf_lsm_v(2)%shf,              &
1422                                               surf_usm_v(2)%shf,              &
1423                                               surf_def_v(3)%shf,              &
1424                                               surf_lsm_v(3)%shf,              &
1425                                               surf_usm_v(3)%shf )
1426               ELSE
1427!
1428!--               Output of averaged data
1429                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1430                                        REAL( average_count_surf, KIND=wp )
1431                  surfaces%var_av(:,n_out) = 0.0_wp
1432               ENDIF
1433
1434            CASE ( 'qsws' )
1435!
1436!--            Output of instantaneous data
1437               IF ( av == 0 )  THEN
1438                  CALL surface_output_collect( surf_def_h(0)%qsws,             &
1439                                               surf_def_h(1)%qsws,             &
1440                                               surf_lsm_h%qsws,                &
1441                                               surf_usm_h%qsws,                &
1442                                               surf_def_v(0)%qsws,             &
1443                                               surf_lsm_v(0)%qsws,             &
1444                                               surf_usm_v(0)%qsws,             &
1445                                               surf_def_v(1)%qsws,             &
1446                                               surf_lsm_v(1)%qsws,             &
1447                                               surf_usm_v(1)%qsws,             &
1448                                               surf_def_v(2)%qsws,             &
1449                                               surf_lsm_v(2)%qsws,             &
1450                                               surf_usm_v(2)%qsws,             &
1451                                               surf_def_v(3)%qsws,             &
1452                                               surf_lsm_v(3)%qsws,             &
1453                                               surf_usm_v(3)%qsws )
1454               ELSE
1455!
1456!--               Output of averaged data
1457                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1458                                        REAL( average_count_surf, KIND=wp )
1459                  surfaces%var_av(:,n_out) = 0.0_wp
1460                                                       
1461               ENDIF
1462
1463            CASE ( 'ssws' )
1464!
1465!--            Output of instantaneous data
1466               IF ( av == 0 )  THEN
1467                  CALL surface_output_collect( surf_def_h(0)%ssws,             &
1468                                               surf_def_h(1)%ssws,             &
1469                                               surf_lsm_h%ssws,                &
1470                                               surf_usm_h%ssws,                &
1471                                               surf_def_v(0)%ssws,             &
1472                                               surf_lsm_v(0)%ssws,             &
1473                                               surf_usm_v(0)%ssws,             &
1474                                               surf_def_v(1)%ssws,             &
1475                                               surf_lsm_v(1)%ssws,             &
1476                                               surf_usm_v(1)%ssws,             &
1477                                               surf_def_v(2)%ssws,             &
1478                                               surf_lsm_v(2)%ssws,             &
1479                                               surf_usm_v(2)%ssws,             &
1480                                               surf_def_v(3)%ssws,             &
1481                                               surf_lsm_v(3)%ssws,             &
1482                                               surf_usm_v(3)%ssws )
1483               ELSE
1484!
1485!--               Output of averaged data
1486                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1487                                        REAL( average_count_surf, KIND=wp )
1488                  surfaces%var_av(:,n_out) = 0.0_wp
1489                                                       
1490               ENDIF
1491
1492            CASE ( 'qcsws' )
1493!
1494!--            Output of instantaneous data
1495               IF ( av == 0 )  THEN
1496                  CALL surface_output_collect( surf_def_h(0)%qcsws,            &
1497                                               surf_def_h(1)%qcsws,            &
1498                                               surf_lsm_h%qcsws,               &
1499                                               surf_usm_h%qcsws,               &
1500                                               surf_def_v(0)%qcsws,            &
1501                                               surf_lsm_v(0)%qcsws,            &
1502                                               surf_usm_v(0)%qcsws,            &
1503                                               surf_def_v(1)%qcsws,            &
1504                                               surf_lsm_v(1)%qcsws,            &
1505                                               surf_usm_v(1)%qcsws,            &
1506                                               surf_def_v(2)%qcsws,            &
1507                                               surf_lsm_v(2)%qcsws,            &
1508                                               surf_usm_v(2)%qcsws,            &
1509                                               surf_def_v(3)%qcsws,            &
1510                                               surf_lsm_v(3)%qcsws,            &
1511                                               surf_usm_v(3)%qcsws )
1512               ELSE
1513!
1514!--               Output of averaged data
1515                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1516                                        REAL( average_count_surf, KIND=wp )
1517                  surfaces%var_av(:,n_out) = 0.0_wp
1518                                                       
1519               ENDIF
1520
1521            CASE ( 'ncsws' )
1522!
1523!--            Output of instantaneous data
1524               IF ( av == 0 )  THEN
1525                  CALL surface_output_collect( surf_def_h(0)%ncsws,            &
1526                                               surf_def_h(1)%ncsws,            &
1527                                               surf_lsm_h%ncsws,               &
1528                                               surf_usm_h%ncsws,               &
1529                                               surf_def_v(0)%ncsws,            &
1530                                               surf_lsm_v(0)%ncsws,            &
1531                                               surf_usm_v(0)%ncsws,            &
1532                                               surf_def_v(1)%ncsws,            &
1533                                               surf_lsm_v(1)%ncsws,            &
1534                                               surf_usm_v(1)%ncsws,            &
1535                                               surf_def_v(2)%ncsws,            &
1536                                               surf_lsm_v(2)%ncsws,            &
1537                                               surf_usm_v(2)%ncsws,            &
1538                                               surf_def_v(3)%ncsws,            &
1539                                               surf_lsm_v(3)%ncsws,            &
1540                                               surf_usm_v(3)%ncsws )
1541               ELSE
1542!
1543!--               Output of averaged data
1544                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1545                                        REAL( average_count_surf, KIND=wp )
1546                  surfaces%var_av(:,n_out) = 0.0_wp
1547                                                       
1548               ENDIF
1549
1550            CASE ( 'qrsws' )
1551!
1552!--            Output of instantaneous data
1553               IF ( av == 0 )  THEN
1554                  CALL surface_output_collect( surf_def_h(0)%qrsws,            &
1555                                               surf_def_h(1)%qrsws,            &
1556                                               surf_lsm_h%qrsws,               &
1557                                               surf_usm_h%qrsws,               &
1558                                               surf_def_v(0)%qrsws,            &
1559                                               surf_lsm_v(0)%qrsws,            &
1560                                               surf_usm_v(0)%qrsws,            &
1561                                               surf_def_v(1)%qrsws,            &
1562                                               surf_lsm_v(1)%qrsws,            &
1563                                               surf_usm_v(1)%qrsws,            &
1564                                               surf_def_v(2)%qrsws,            &
1565                                               surf_lsm_v(2)%qrsws,            &
1566                                               surf_usm_v(2)%qrsws,            &
1567                                               surf_def_v(3)%qrsws,            &
1568                                               surf_lsm_v(3)%qrsws,            &
1569                                               surf_usm_v(3)%qrsws )
1570               ELSE
1571!
1572!--               Output of averaged data
1573                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1574                                        REAL( average_count_surf, KIND=wp )
1575                  surfaces%var_av(:,n_out) = 0.0_wp
1576                                                       
1577               ENDIF
1578
1579            CASE ( 'nrsws' )
1580!
1581!--            Output of instantaneous data
1582               IF ( av == 0 )  THEN
1583                  CALL surface_output_collect( surf_def_h(0)%nrsws,            &
1584                                               surf_def_h(1)%nrsws,            &
1585                                               surf_lsm_h%nrsws,               &
1586                                               surf_usm_h%nrsws,               &
1587                                               surf_def_v(0)%nrsws,            &
1588                                               surf_lsm_v(0)%nrsws,            &
1589                                               surf_usm_v(0)%nrsws,            &
1590                                               surf_def_v(1)%nrsws,            &
1591                                               surf_lsm_v(1)%nrsws,            &
1592                                               surf_usm_v(1)%nrsws,            &
1593                                               surf_def_v(2)%nrsws,            &
1594                                               surf_lsm_v(2)%nrsws,            &
1595                                               surf_usm_v(2)%nrsws,            &
1596                                               surf_def_v(3)%nrsws,            &
1597                                               surf_lsm_v(3)%nrsws,            &
1598                                               surf_usm_v(3)%nrsws )
1599               ELSE
1600!
1601!--               Output of averaged data
1602                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1603                                        REAL( average_count_surf, KIND=wp )
1604                  surfaces%var_av(:,n_out) = 0.0_wp
1605                                                       
1606               ENDIF
1607
1608            CASE ( 'sasws' )
1609!
1610!--            Output of instantaneous data
1611               IF ( av == 0 )  THEN
1612                  CALL surface_output_collect( surf_def_h(0)%sasws,            &
1613                                               surf_def_h(1)%sasws,            &
1614                                               surf_lsm_h%sasws,               &
1615                                               surf_usm_h%sasws,               &
1616                                               surf_def_v(0)%sasws,            &
1617                                               surf_lsm_v(0)%sasws,            &
1618                                               surf_usm_v(0)%sasws,            &
1619                                               surf_def_v(1)%sasws,            &
1620                                               surf_lsm_v(1)%sasws,            &
1621                                               surf_usm_v(1)%sasws,            &
1622                                               surf_def_v(2)%sasws,            &
1623                                               surf_lsm_v(2)%sasws,            &
1624                                               surf_usm_v(2)%sasws,            &
1625                                               surf_def_v(3)%sasws,            &
1626                                               surf_lsm_v(3)%sasws,            &
1627                                               surf_usm_v(3)%sasws )
1628               ELSE
1629!
1630!--               Output of averaged data
1631                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1632                                        REAL( average_count_surf, KIND=wp )
1633                  surfaces%var_av(:,n_out) = 0.0_wp
1634                                                       
1635               ENDIF
1636
1637            CASE ( 'q_surface' )
1638!
1639!--            Output of instantaneous data
1640               IF ( av == 0 )  THEN
1641                  CALL surface_output_collect( surf_def_h(0)%q_surface,        &
1642                                               surf_def_h(1)%q_surface,        &
1643                                               surf_lsm_h%q_surface,           &
1644                                               surf_usm_h%q_surface,           &
1645                                               surf_def_v(0)%q_surface,        &
1646                                               surf_lsm_v(0)%q_surface,        &
1647                                               surf_usm_v(0)%q_surface,        &
1648                                               surf_def_v(1)%q_surface,        &
1649                                               surf_lsm_v(1)%q_surface,        &
1650                                               surf_usm_v(1)%q_surface,        &
1651                                               surf_def_v(2)%q_surface,        &
1652                                               surf_lsm_v(2)%q_surface,        &
1653                                               surf_usm_v(2)%q_surface,        &
1654                                               surf_def_v(3)%q_surface,        &
1655                                               surf_lsm_v(3)%q_surface,        &
1656                                               surf_usm_v(3)%q_surface )
1657               ELSE
1658!
1659!--               Output of averaged data
1660                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1661                                        REAL( average_count_surf, KIND=wp )
1662                  surfaces%var_av(:,n_out) = 0.0_wp
1663                                                       
1664               ENDIF
1665
1666            CASE ( 'theta_surface' )
1667!
1668!--            Output of instantaneous data
1669               IF ( av == 0 )  THEN
1670                  CALL surface_output_collect( surf_def_h(0)%pt_surface,       &
1671                                               surf_def_h(1)%pt_surface,       &
1672                                               surf_lsm_h%pt_surface,          &
1673                                               surf_usm_h%pt_surface,          &
1674                                               surf_def_v(0)%pt_surface,       &
1675                                               surf_lsm_v(0)%pt_surface,       &
1676                                               surf_usm_v(0)%pt_surface,       &
1677                                               surf_def_v(1)%pt_surface,       &
1678                                               surf_lsm_v(1)%pt_surface,       &
1679                                               surf_usm_v(1)%pt_surface,       &
1680                                               surf_def_v(2)%pt_surface,       &
1681                                               surf_lsm_v(2)%pt_surface,       &
1682                                               surf_usm_v(2)%pt_surface,       &
1683                                               surf_def_v(3)%pt_surface,       &
1684                                               surf_lsm_v(3)%pt_surface,       &
1685                                               surf_usm_v(3)%pt_surface )
1686               ELSE
1687!
1688!--               Output of averaged data
1689                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1690                                        REAL( average_count_surf, KIND=wp )
1691                  surfaces%var_av(:,n_out) = 0.0_wp
1692                                                       
1693               ENDIF
1694
1695            CASE ( 'thetav_surface' )
1696!
1697!--            Output of instantaneous data
1698               IF ( av == 0 )  THEN
1699                  CALL surface_output_collect( surf_def_h(0)%vpt_surface,      &
1700                                               surf_def_h(1)%vpt_surface,      &
1701                                               surf_lsm_h%vpt_surface,         &
1702                                               surf_usm_h%vpt_surface,         &
1703                                               surf_def_v(0)%vpt_surface,      &
1704                                               surf_lsm_v(0)%vpt_surface,      &
1705                                               surf_usm_v(0)%vpt_surface,      &
1706                                               surf_def_v(1)%vpt_surface,      &
1707                                               surf_lsm_v(1)%vpt_surface,      &
1708                                               surf_usm_v(1)%vpt_surface,      &
1709                                               surf_def_v(2)%vpt_surface,      &
1710                                               surf_lsm_v(2)%vpt_surface,      &
1711                                               surf_usm_v(2)%vpt_surface,      &
1712                                               surf_def_v(3)%vpt_surface,      &
1713                                               surf_lsm_v(3)%vpt_surface,      &
1714                                               surf_usm_v(3)%vpt_surface)
1715               ELSE
1716!
1717!--               Output of averaged data
1718                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1719                                        REAL( average_count_surf, KIND=wp )
1720                  surfaces%var_av(:,n_out) = 0.0_wp
1721                                                       
1722               ENDIF
1723
1724            CASE ( 'rad_net' )
1725!
1726!--            Output of instantaneous data
1727               IF ( av == 0 )  THEN
1728                  CALL surface_output_collect( surf_def_h(0)%rad_net,          &
1729                                               surf_def_h(1)%rad_net,          &
1730                                               surf_lsm_h%rad_net,             &
1731                                               surf_usm_h%rad_net,             &
1732                                               surf_def_v(0)%rad_net,          &
1733                                               surf_lsm_v(0)%rad_net,          &
1734                                               surf_usm_v(0)%rad_net,          &
1735                                               surf_def_v(1)%rad_net,          &
1736                                               surf_lsm_v(1)%rad_net,          &
1737                                               surf_usm_v(1)%rad_net,          &
1738                                               surf_def_v(2)%rad_net,          &
1739                                               surf_lsm_v(2)%rad_net,          &
1740                                               surf_usm_v(2)%rad_net,          &
1741                                               surf_def_v(3)%rad_net,          &
1742                                               surf_lsm_v(3)%rad_net,          &
1743                                               surf_usm_v(3)%rad_net ) 
1744               ELSE
1745!
1746!--               Output of averaged data
1747                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1748                                        REAL( average_count_surf, KIND=wp )
1749                  surfaces%var_av(:,n_out) = 0.0_wp
1750                                                       
1751               ENDIF
1752
1753            CASE ( 'rad_lw_in' )
1754!
1755!--            Output of instantaneous data
1756               IF ( av == 0 )  THEN
1757                  CALL surface_output_collect( surf_def_h(0)%rad_lw_in,        &
1758                                               surf_def_h(1)%rad_lw_in,        &
1759                                               surf_lsm_h%rad_lw_in,           &
1760                                               surf_usm_h%rad_lw_in,           &
1761                                               surf_def_v(0)%rad_lw_in,        &
1762                                               surf_lsm_v(0)%rad_lw_in,        &
1763                                               surf_usm_v(0)%rad_lw_in,        &
1764                                               surf_def_v(1)%rad_lw_in,        &
1765                                               surf_lsm_v(1)%rad_lw_in,        &
1766                                               surf_usm_v(1)%rad_lw_in,        &
1767                                               surf_def_v(2)%rad_lw_in,        &
1768                                               surf_lsm_v(2)%rad_lw_in,        &
1769                                               surf_usm_v(2)%rad_lw_in,        &
1770                                               surf_def_v(3)%rad_lw_in,        &
1771                                               surf_lsm_v(3)%rad_lw_in,        &
1772                                               surf_usm_v(3)%rad_lw_in )
1773               ELSE
1774!
1775!--               Output of averaged data
1776                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1777                                        REAL( average_count_surf, KIND=wp )
1778                  surfaces%var_av(:,n_out) = 0.0_wp
1779                                                       
1780               ENDIF
1781
1782            CASE ( 'rad_lw_out' )
1783!
1784!--            Output of instantaneous data
1785               IF ( av == 0 )  THEN
1786                  CALL surface_output_collect( surf_def_h(0)%rad_lw_out,       &
1787                                               surf_def_h(1)%rad_lw_out,       &
1788                                               surf_lsm_h%rad_lw_out,          &
1789                                               surf_usm_h%rad_lw_out,          &
1790                                               surf_def_v(0)%rad_lw_out,       &
1791                                               surf_lsm_v(0)%rad_lw_out,       &
1792                                               surf_usm_v(0)%rad_lw_out,       &
1793                                               surf_def_v(1)%rad_lw_out,       &
1794                                               surf_lsm_v(1)%rad_lw_out,       &
1795                                               surf_usm_v(1)%rad_lw_out,       &
1796                                               surf_def_v(2)%rad_lw_out,       &
1797                                               surf_lsm_v(2)%rad_lw_out,       &
1798                                               surf_usm_v(2)%rad_lw_out,       &
1799                                               surf_def_v(3)%rad_lw_out,       &
1800                                               surf_lsm_v(3)%rad_lw_out,       &
1801                                               surf_usm_v(3)%rad_lw_out )
1802               ELSE
1803!
1804!--               Output of averaged data
1805                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1806                                        REAL( average_count_surf, KIND=wp )
1807                  surfaces%var_av(:,n_out) = 0.0_wp
1808                                                       
1809               ENDIF
1810
1811            CASE ( 'rad_sw_in' )
1812!
1813!--            Output of instantaneous data
1814               IF ( av == 0 )  THEN
1815                  CALL surface_output_collect( surf_def_h(0)%rad_sw_in,        &
1816                                               surf_def_h(1)%rad_sw_in,        &
1817                                               surf_lsm_h%rad_sw_in,           &
1818                                               surf_usm_h%rad_sw_in,           &
1819                                               surf_def_v(0)%rad_sw_in,        &
1820                                               surf_lsm_v(0)%rad_sw_in,        &
1821                                               surf_usm_v(0)%rad_sw_in,        &
1822                                               surf_def_v(1)%rad_sw_in,        &
1823                                               surf_lsm_v(1)%rad_sw_in,        &
1824                                               surf_usm_v(1)%rad_sw_in,        &
1825                                               surf_def_v(2)%rad_sw_in,        &
1826                                               surf_lsm_v(2)%rad_sw_in,        &
1827                                               surf_usm_v(2)%rad_sw_in,        &
1828                                               surf_def_v(3)%rad_sw_in,        &
1829                                               surf_lsm_v(3)%rad_sw_in,        &
1830                                               surf_usm_v(3)%rad_sw_in )
1831               ELSE
1832!
1833!--               Output of averaged data
1834                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1835                                        REAL( average_count_surf, KIND=wp )
1836                  surfaces%var_av(:,n_out) = 0.0_wp
1837                                                       
1838               ENDIF
1839
1840            CASE ( 'rad_sw_out' )
1841!
1842!--            Output of instantaneous data
1843               IF ( av == 0 )  THEN
1844                  CALL surface_output_collect( surf_def_h(0)%rad_sw_out,       &
1845                                               surf_def_h(1)%rad_sw_out,       &
1846                                               surf_lsm_h%rad_sw_out,          &
1847                                               surf_usm_h%rad_sw_out,          &
1848                                               surf_def_v(0)%rad_sw_out,       &
1849                                               surf_lsm_v(0)%rad_sw_out,       &
1850                                               surf_usm_v(0)%rad_sw_out,       &
1851                                               surf_def_v(1)%rad_sw_out,       &
1852                                               surf_lsm_v(1)%rad_sw_out,       &
1853                                               surf_usm_v(1)%rad_sw_out,       &
1854                                               surf_def_v(2)%rad_sw_out,       &
1855                                               surf_lsm_v(2)%rad_sw_out,       &
1856                                               surf_usm_v(2)%rad_sw_out,       &
1857                                               surf_def_v(3)%rad_sw_out,       &
1858                                               surf_lsm_v(3)%rad_sw_out,       &
1859                                               surf_usm_v(3)%rad_sw_out )
1860               ELSE
1861!
1862!--               Output of averaged data
1863                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1864                                        REAL( average_count_surf, KIND=wp )
1865                  surfaces%var_av(:,n_out) = 0.0_wp
1866                                                       
1867               ENDIF
1868
1869            CASE ( 'ghf' )
1870!
1871!--            Output of instantaneous data
1872               IF ( av == 0 )  THEN
1873                  CALL surface_output_collect( surf_def_h(0)%ghf,              &
1874                                               surf_def_h(1)%ghf,              &
1875                                               surf_lsm_h%ghf,                 &
1876                                               surf_usm_h%ghf,                 &
1877                                               surf_def_v(0)%ghf,              &
1878                                               surf_lsm_v(0)%ghf,              &
1879                                               surf_usm_v(0)%ghf,              &
1880                                               surf_def_v(1)%ghf,              &
1881                                               surf_lsm_v(1)%ghf,              &
1882                                               surf_usm_v(1)%ghf,              &
1883                                               surf_def_v(2)%ghf,              &
1884                                               surf_lsm_v(2)%ghf,              &
1885                                               surf_usm_v(2)%ghf,              &
1886                                               surf_def_v(3)%ghf,              &
1887                                               surf_lsm_v(3)%ghf,              &
1888                                               surf_usm_v(3)%ghf )     
1889                                                                        ELSE
1890!
1891!--               Output of averaged data
1892                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1893                                        REAL( average_count_surf, KIND=wp )
1894                  surfaces%var_av(:,n_out) = 0.0_wp
1895                                                       
1896               ENDIF
1897                     
1898            CASE ( 'r_a' )                                                     
1899!
1900!--            Output of instantaneous data
1901               IF ( av == 0 )  THEN
1902                  CALL surface_output_collect( surf_def_h(0)%r_a,              &
1903                                               surf_def_h(1)%r_a,              &
1904                                               surf_lsm_h%r_a,                 &
1905                                               surf_usm_h%r_a,                 &
1906                                               surf_def_v(0)%r_a,              &
1907                                               surf_lsm_v(0)%r_a,              &
1908                                               surf_usm_v(0)%r_a,              &
1909                                               surf_def_v(1)%r_a,              &
1910                                               surf_lsm_v(1)%r_a,              &
1911                                               surf_usm_v(1)%r_a,              &
1912                                               surf_def_v(2)%r_a,              &
1913                                               surf_lsm_v(2)%r_a,              &
1914                                               surf_usm_v(2)%r_a,              &
1915                                               surf_def_v(3)%r_a,              &
1916                                               surf_lsm_v(3)%r_a,              &
1917                                               surf_usm_v(3)%r_a )     
1918                                                                      ELSE
1919!
1920!--               Output of averaged data
1921                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1922                                        REAL( average_count_surf, KIND=wp )
1923                  surfaces%var_av(:,n_out) = 0.0_wp
1924                                                       
1925               ENDIF
1926                       
1927            CASE ( 'r_soil' )                                                 
1928!
1929!--            Output of instantaneous data
1930               IF ( av == 0 )  THEN
1931                  CALL surface_output_collect( surf_def_h(0)%r_soil,           &
1932                                               surf_def_h(1)%r_soil,           &
1933                                               surf_lsm_h%r_soil,              &
1934                                               surf_usm_h%r_soil,              &
1935                                               surf_def_v(0)%r_soil,           &
1936                                               surf_lsm_v(0)%r_soil,           &
1937                                               surf_usm_v(0)%r_soil,           &
1938                                               surf_def_v(1)%r_soil,           &
1939                                               surf_lsm_v(1)%r_soil,           &
1940                                               surf_usm_v(1)%r_soil,           &
1941                                               surf_def_v(2)%r_soil,           &
1942                                               surf_lsm_v(2)%r_soil,           &
1943                                               surf_usm_v(2)%r_soil,           &
1944                                               surf_def_v(3)%r_soil,           &
1945                                               surf_lsm_v(3)%r_soil,           &
1946                                               surf_usm_v(3)%r_soil ) 
1947               ELSE
1948!
1949!--               Output of averaged data
1950                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1951                                        REAL( average_count_surf, KIND=wp )
1952                  surfaces%var_av(:,n_out) = 0.0_wp
1953                                                       
1954               ENDIF
1955
1956            CASE ( 'r_canopy' )
1957!
1958!--            Output of instantaneous data
1959               IF ( av == 0 )  THEN
1960                  CALL surface_output_collect( surf_def_h(0)%r_canopy,         &
1961                                               surf_def_h(1)%r_canopy,         &
1962                                               surf_lsm_h%r_canopy,            &
1963                                               surf_usm_h%r_canopy,            &
1964                                               surf_def_v(0)%r_canopy,         &
1965                                               surf_lsm_v(0)%r_canopy,         &
1966                                               surf_usm_v(0)%r_canopy,         &
1967                                               surf_def_v(1)%r_canopy,         &
1968                                               surf_lsm_v(1)%r_canopy,         &
1969                                               surf_usm_v(1)%r_canopy,         &
1970                                               surf_def_v(2)%r_canopy,         &
1971                                               surf_lsm_v(2)%r_canopy,         &
1972                                               surf_usm_v(2)%r_canopy,         &
1973                                               surf_def_v(3)%r_canopy,         &
1974                                               surf_lsm_v(3)%r_canopy,         &
1975                                               surf_usm_v(3)%r_canopy ) 
1976               ELSE
1977!
1978!--               Output of averaged data
1979                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
1980                                        REAL( average_count_surf, KIND=wp )
1981                  surfaces%var_av(:,n_out) = 0.0_wp
1982                                                       
1983               ENDIF
1984
1985            CASE ( 'r_s' )
1986!
1987!--            Output of instantaneous data
1988               IF ( av == 0 )  THEN
1989                  CALL surface_output_collect( surf_def_h(0)%r_s,              &
1990                                               surf_def_h(1)%r_s,              &
1991                                               surf_lsm_h%r_s,                 &
1992                                               surf_usm_h%r_s,                 &
1993                                               surf_def_v(0)%r_s,              &
1994                                               surf_lsm_v(0)%r_s,              &
1995                                               surf_usm_v(0)%r_s,              &
1996                                               surf_def_v(1)%r_s,              &
1997                                               surf_lsm_v(1)%r_s,              &
1998                                               surf_usm_v(1)%r_s,              &
1999                                               surf_def_v(2)%r_s,              &
2000                                               surf_lsm_v(2)%r_s,              &
2001                                               surf_usm_v(2)%r_s,              &
2002                                               surf_def_v(3)%r_s,              &
2003                                               surf_lsm_v(3)%r_s,              &
2004                                               surf_usm_v(3)%r_s ) 
2005               ELSE
2006!
2007!--               Output of averaged data
2008                  surfaces%var_out(:) = surfaces%var_av(:,n_out) /             &
2009                                        REAL( average_count_surf, KIND=wp )
2010                  surfaces%var_av(:,n_out) = 0.0_wp
2011                                                       
2012               ENDIF
2013
2014!
2015!--            Add further variables:
2016!--            'css', 'cssws', 'qsws_liq', 'qsws_soil', 'qsws_veg'
2017
2018         END SELECT
2019!
2020!--      Write to binary file:
2021!--      - surfaces%points ( 3, 1-npoints )
2022!--      - surfaces%polygons ( 5, 1-ns )
2023!--      - surfaces%var_out ( 1-ns, time )
2024!--      - Dimension: 1-nsurfaces, 1-npoints - can be ordered consecutively
2025!--      - Distinguish between average and non-average data
2026         DO  i = 0, io_blocks-1
2027            IF ( i == io_group )  THEN
2028               WRITE ( 25+av )  LEN_TRIM( 'time' )
2029               WRITE ( 25+av )  'time'
2030               WRITE ( 25+av )  simulated_time
2031               WRITE ( 25+av )  LEN_TRIM( trimvar )
2032               WRITE ( 25+av )  TRIM( trimvar )
2033               WRITE ( 25+av )  surfaces%var_out
2034            ENDIF
2035#if defined( __parallel )
2036            CALL MPI_BARRIER( comm2d, ierr )
2037#endif
2038         ENDDO
2039
2040      ENDDO
2041     
2042!
2043!--   If averaged output was written to NetCDF file, set the counter to zero
2044      IF ( av == 1 )  average_count_surf = 0
2045             
2046   END SUBROUTINE surface_output
2047   
2048!------------------------------------------------------------------------------!
2049! Description:
2050! ------------
2051!> Routine for controlling the data averaging.
2052!------------------------------------------------------------------------------!     
2053   SUBROUTINE surface_output_averaging
2054   
2055      IMPLICIT NONE
2056     
2057      CHARACTER(LEN=100) ::  trimvar !< dummy variable for current output variable
2058       
2059      INTEGER(iwp) ::  l      !< running index for surface-element orientation
2060      INTEGER(iwp) ::  m      !< running index for surface elements
2061      INTEGER(iwp) ::  n_out  !< counter variables for surface output
2062     
2063      n_out = 0
2064      DO  WHILE ( dosurf(1,n_out+1)(1:1) /= ' ' )
2065         
2066         n_out   = n_out + 1
2067         trimvar = TRIM( dosurf(1,n_out) )
2068
2069         SELECT CASE ( trimvar )
2070
2071            CASE ( 'us' )
2072               CALL surface_output_sum_up( surf_def_h(0)%us,                   &
2073                                           surf_def_h(1)%us,                   &
2074                                           surf_lsm_h%us,                      &
2075                                           surf_usm_h%us,                      &
2076                                           surf_def_v(0)%us,                   &
2077                                           surf_lsm_v(0)%us,                   &
2078                                           surf_usm_v(0)%us,                   &
2079                                           surf_def_v(1)%us,                   &
2080                                           surf_lsm_v(1)%us,                   &
2081                                           surf_usm_v(1)%us,                   &
2082                                           surf_def_v(2)%us,                   &
2083                                           surf_lsm_v(2)%us,                   &
2084                                           surf_usm_v(2)%us,                   &
2085                                           surf_def_v(3)%us,                   &
2086                                           surf_lsm_v(3)%us,                   &
2087                                           surf_usm_v(3)%us, n_out ) 
2088
2089            CASE ( 'ts' )
2090               CALL surface_output_sum_up( surf_def_h(0)%ts,                   &
2091                                           surf_def_h(1)%ts,                   &
2092                                           surf_lsm_h%ts,                      &
2093                                           surf_usm_h%ts,                      &
2094                                           surf_def_v(0)%ts,                   &
2095                                           surf_lsm_v(0)%ts,                   &
2096                                           surf_usm_v(0)%ts,                   &
2097                                           surf_def_v(1)%ts,                   &
2098                                           surf_lsm_v(1)%ts,                   &
2099                                           surf_usm_v(1)%ts,                   &
2100                                           surf_def_v(2)%ts,                   &
2101                                           surf_lsm_v(2)%ts,                   &
2102                                           surf_usm_v(2)%ts,                   &
2103                                           surf_def_v(3)%ts,                   &
2104                                           surf_lsm_v(3)%ts,                   &
2105                                           surf_usm_v(3)%ts, n_out ) 
2106
2107            CASE ( 'qs' )
2108               CALL surface_output_sum_up( surf_def_h(0)%qs,                   &
2109                                           surf_def_h(1)%qs,                   &
2110                                           surf_lsm_h%qs,                      &
2111                                           surf_usm_h%qs,                      &
2112                                           surf_def_v(0)%qs,                   &
2113                                           surf_lsm_v(0)%qs,                   &
2114                                           surf_usm_v(0)%qs,                   &
2115                                           surf_def_v(1)%qs,                   &
2116                                           surf_lsm_v(1)%qs,                   &
2117                                           surf_usm_v(1)%qs,                   &
2118                                           surf_def_v(2)%qs,                   &
2119                                           surf_lsm_v(2)%qs,                   &
2120                                           surf_usm_v(2)%qs,                   &
2121                                           surf_def_v(3)%qs,                   &
2122                                           surf_lsm_v(3)%qs,                   &
2123                                           surf_usm_v(3)%qs, n_out ) 
2124
2125            CASE ( 'ss' )
2126               CALL surface_output_sum_up( surf_def_h(0)%ss,                   &
2127                                           surf_def_h(1)%ss,                   &
2128                                           surf_lsm_h%ss,                      &
2129                                           surf_usm_h%ss,                      &
2130                                           surf_def_v(0)%ss,                   &
2131                                           surf_lsm_v(0)%ss,                   &
2132                                           surf_usm_v(0)%ss,                   &
2133                                           surf_def_v(1)%ss,                   &
2134                                           surf_lsm_v(1)%ss,                   &
2135                                           surf_usm_v(1)%ss,                   &
2136                                           surf_def_v(2)%ss,                   &
2137                                           surf_lsm_v(2)%ss,                   &
2138                                           surf_usm_v(2)%ss,                   &
2139                                           surf_def_v(3)%ss,                   &
2140                                           surf_lsm_v(3)%ss,                   &
2141                                           surf_usm_v(3)%ss, n_out ) 
2142
2143            CASE ( 'qcs' )
2144               CALL surface_output_sum_up( surf_def_h(0)%qcs,                  &
2145                                           surf_def_h(1)%qcs,                  &
2146                                           surf_lsm_h%qcs,                     &
2147                                           surf_usm_h%qcs,                     &
2148                                           surf_def_v(0)%qcs,                  &
2149                                           surf_lsm_v(0)%qcs,                  &
2150                                           surf_usm_v(0)%qcs,                  &
2151                                           surf_def_v(1)%qcs,                  &
2152                                           surf_lsm_v(1)%qcs,                  &
2153                                           surf_usm_v(1)%qcs,                  &
2154                                           surf_def_v(2)%qcs,                  &
2155                                           surf_lsm_v(2)%qcs,                  &
2156                                           surf_usm_v(2)%qcs,                  &
2157                                           surf_def_v(3)%qcs,                  &
2158                                           surf_lsm_v(3)%qcs,                  &
2159                                           surf_usm_v(3)%qcs, n_out ) 
2160
2161            CASE ( 'ncs' )
2162               CALL surface_output_sum_up( surf_def_h(0)%ncs,                  &
2163                                           surf_def_h(1)%ncs,                  &
2164                                           surf_lsm_h%ncs,                     &
2165                                           surf_usm_h%ncs,                     &
2166                                           surf_def_v(0)%ncs,                  &
2167                                           surf_lsm_v(0)%ncs,                  &
2168                                           surf_usm_v(0)%ncs,                  &
2169                                           surf_def_v(1)%ncs,                  &
2170                                           surf_lsm_v(1)%ncs,                  &
2171                                           surf_usm_v(1)%ncs,                  &
2172                                           surf_def_v(2)%ncs,                  &
2173                                           surf_lsm_v(2)%ncs,                  &
2174                                           surf_usm_v(2)%ncs,                  &
2175                                           surf_def_v(3)%ncs,                  &
2176                                           surf_lsm_v(3)%ncs,                  &
2177                                           surf_usm_v(3)%ncs, n_out ) 
2178
2179            CASE ( 'qrs' )
2180               CALL surface_output_sum_up( surf_def_h(0)%qrs,                  &
2181                                           surf_def_h(1)%qrs,                  &
2182                                           surf_lsm_h%qrs,                     &
2183                                           surf_usm_h%qrs,                     &
2184                                           surf_def_v(0)%qrs,                  &
2185                                           surf_lsm_v(0)%qrs,                  &
2186                                           surf_usm_v(0)%qrs,                  &
2187                                           surf_def_v(1)%qrs,                  &
2188                                           surf_lsm_v(1)%qrs,                  &
2189                                           surf_usm_v(1)%qrs,                  &
2190                                           surf_def_v(2)%qrs,                  &
2191                                           surf_lsm_v(2)%qrs,                  &
2192                                           surf_usm_v(2)%qrs,                  &
2193                                           surf_def_v(3)%qrs,                  &
2194                                           surf_lsm_v(3)%qrs,                  &
2195                                           surf_usm_v(3)%qrs, n_out ) 
2196
2197            CASE ( 'nrs' )
2198               CALL surface_output_sum_up( surf_def_h(0)%nrs,                  &
2199                                           surf_def_h(1)%nrs,                  &
2200                                           surf_lsm_h%nrs,                     &
2201                                           surf_usm_h%nrs,                     &
2202                                           surf_def_v(0)%nrs,                  &
2203                                           surf_lsm_v(0)%nrs,                  &
2204                                           surf_usm_v(0)%nrs,                  &
2205                                           surf_def_v(1)%nrs,                  &
2206                                           surf_lsm_v(1)%nrs,                  &
2207                                           surf_usm_v(1)%nrs,                  &
2208                                           surf_def_v(2)%nrs,                  &
2209                                           surf_lsm_v(2)%nrs,                  &
2210                                           surf_usm_v(2)%nrs,                  &
2211                                           surf_def_v(3)%nrs,                  &
2212                                           surf_lsm_v(3)%nrs,                  &
2213                                           surf_usm_v(3)%nrs, n_out ) 
2214
2215            CASE ( 'ol' )
2216               CALL surface_output_sum_up( surf_def_h(0)%ol,                   &
2217                                           surf_def_h(1)%ol,                   &
2218                                           surf_lsm_h%ol,                      &
2219                                           surf_usm_h%ol,                      &
2220                                           surf_def_v(0)%ol,                   &
2221                                           surf_lsm_v(0)%ol,                   &
2222                                           surf_usm_v(0)%ol,                   &
2223                                           surf_def_v(1)%ol,                   &
2224                                           surf_lsm_v(1)%ol,                   &
2225                                           surf_usm_v(1)%ol,                   &
2226                                           surf_def_v(2)%ol,                   &
2227                                           surf_lsm_v(2)%ol,                   &
2228                                           surf_usm_v(2)%ol,                   &
2229                                           surf_def_v(3)%ol,                   &
2230                                           surf_lsm_v(3)%ol,                   &
2231                                           surf_usm_v(3)%ol, n_out ) 
2232
2233            CASE ( 'z0' )
2234               CALL surface_output_sum_up( surf_def_h(0)%z0,                   &
2235                                           surf_def_h(1)%z0,                   &
2236                                           surf_lsm_h%z0,                      &
2237                                           surf_usm_h%z0,                      &
2238                                           surf_def_v(0)%z0,                   &
2239                                           surf_lsm_v(0)%z0,                   &
2240                                           surf_usm_v(0)%z0,                   &
2241                                           surf_def_v(1)%z0,                   &
2242                                           surf_lsm_v(1)%z0,                   &
2243                                           surf_usm_v(1)%z0,                   &
2244                                           surf_def_v(2)%z0,                   &
2245                                           surf_lsm_v(2)%z0,                   &
2246                                           surf_usm_v(2)%z0,                   &
2247                                           surf_def_v(3)%z0,                   &
2248                                           surf_lsm_v(3)%z0,                   &
2249                                           surf_usm_v(3)%z0, n_out ) 
2250
2251            CASE ( 'z0h' )
2252               CALL surface_output_sum_up( surf_def_h(0)%z0h,                  &
2253                                           surf_def_h(1)%z0h,                  &
2254                                           surf_lsm_h%z0h,                     &
2255                                           surf_usm_h%z0h,                     &
2256                                           surf_def_v(0)%z0h,                  &
2257                                           surf_lsm_v(0)%z0h,                  &
2258                                           surf_usm_v(0)%z0h,                  &
2259                                           surf_def_v(1)%z0h,                  &
2260                                           surf_lsm_v(1)%z0h,                  &
2261                                           surf_usm_v(1)%z0h,                  &
2262                                           surf_def_v(2)%z0h,                  &
2263                                           surf_lsm_v(2)%z0h,                  &
2264                                           surf_usm_v(2)%z0h,                  &
2265                                           surf_def_v(3)%z0h,                  &
2266                                           surf_lsm_v(3)%z0h,                  &
2267                                           surf_usm_v(3)%z0h, n_out )         
2268
2269            CASE ( 'z0q' )
2270               CALL surface_output_sum_up( surf_def_h(0)%z0q,                  &
2271                                           surf_def_h(1)%z0q,                  &
2272                                           surf_lsm_h%z0q,                     &
2273                                           surf_usm_h%z0q,                     &
2274                                           surf_def_v(0)%z0q,                  &
2275                                           surf_lsm_v(0)%z0q,                  &
2276                                           surf_usm_v(0)%z0q,                  &
2277                                           surf_def_v(1)%z0q,                  &
2278                                           surf_lsm_v(1)%z0q,                  &
2279                                           surf_usm_v(1)%z0q,                  &
2280                                           surf_def_v(2)%z0q,                  &
2281                                           surf_lsm_v(2)%z0q,                  &
2282                                           surf_usm_v(2)%z0q,                  &
2283                                           surf_def_v(3)%z0q,                  &
2284                                           surf_lsm_v(3)%z0q,                  &
2285                                           surf_usm_v(3)%z0q, n_out ) 
2286
2287            CASE ( 'theta1' )                                                 
2288               CALL surface_output_sum_up( surf_def_h(0)%pt1,                  &
2289                                           surf_def_h(1)%pt1,                  &
2290                                           surf_lsm_h%pt1,                     &
2291                                           surf_usm_h%pt1,                     &
2292                                           surf_def_v(0)%pt1,                  &
2293                                           surf_lsm_v(0)%pt1,                  &
2294                                           surf_usm_v(0)%pt1,                  &
2295                                           surf_def_v(1)%pt1,                  &
2296                                           surf_lsm_v(1)%pt1,                  &
2297                                           surf_usm_v(1)%pt1,                  &
2298                                           surf_def_v(2)%pt1,                  &
2299                                           surf_lsm_v(2)%pt1,                  &
2300                                           surf_usm_v(2)%pt1,                  &
2301                                           surf_def_v(3)%pt1,                  &
2302                                           surf_lsm_v(3)%pt1,                  &
2303                                           surf_usm_v(3)%pt1, n_out )         
2304                                                                               
2305            CASE ( 'qv1' )                                                     
2306               CALL surface_output_sum_up( surf_def_h(0)%qv1,                  &
2307                                           surf_def_h(1)%qv1,                  &
2308                                           surf_lsm_h%qv1,                     &
2309                                           surf_usm_h%qv1,                     &
2310                                           surf_def_v(0)%qv1,                  &
2311                                           surf_lsm_v(0)%qv1,                  &
2312                                           surf_usm_v(0)%qv1,                  &
2313                                           surf_def_v(1)%qv1,                  &
2314                                           surf_lsm_v(1)%qv1,                  &
2315                                           surf_usm_v(1)%qv1,                  &
2316                                           surf_def_v(2)%qv1,                  &
2317                                           surf_lsm_v(2)%qv1,                  &
2318                                           surf_usm_v(2)%qv1,                  &
2319                                           surf_def_v(3)%qv1,                  &
2320                                           surf_lsm_v(3)%qv1,                  &
2321                                           surf_usm_v(3)%qv1, n_out ) 
2322
2323            CASE ( 'thetav1' )
2324               CALL surface_output_sum_up( surf_def_h(0)%vpt1,                 &
2325                                           surf_def_h(1)%vpt1,                 &
2326                                           surf_lsm_h%vpt1,                    &
2327                                           surf_usm_h%vpt1,                    &
2328                                           surf_def_v(0)%vpt1,                 &
2329                                           surf_lsm_v(0)%vpt1,                 &
2330                                           surf_usm_v(0)%vpt1,                 &
2331                                           surf_def_v(1)%vpt1,                 &
2332                                           surf_lsm_v(1)%vpt1,                 &
2333                                           surf_usm_v(1)%vpt1,                 &
2334                                           surf_def_v(2)%vpt1,                 &
2335                                           surf_lsm_v(2)%vpt1,                 &
2336                                           surf_usm_v(2)%vpt1,                 &
2337                                           surf_def_v(3)%vpt1,                 &
2338                                           surf_lsm_v(3)%vpt1,                 &
2339                                           surf_usm_v(3)%vpt1, n_out ) 
2340
2341            CASE ( 'usws' )
2342               CALL surface_output_sum_up( surf_def_h(0)%usws,                 &
2343                                           surf_def_h(1)%usws,                 &
2344                                           surf_lsm_h%usws,                    &
2345                                           surf_usm_h%usws,                    &
2346                                           surf_def_v(0)%usws,                 &
2347                                           surf_lsm_v(0)%usws,                 &
2348                                           surf_usm_v(0)%usws,                 &
2349                                           surf_def_v(1)%usws,                 &
2350                                           surf_lsm_v(1)%usws,                 &
2351                                           surf_usm_v(1)%usws,                 &
2352                                           surf_def_v(2)%usws,                 &
2353                                           surf_lsm_v(2)%usws,                 &
2354                                           surf_usm_v(2)%usws,                 &
2355                                           surf_def_v(3)%usws,                 &
2356                                           surf_lsm_v(3)%usws,                 &
2357                                           surf_usm_v(3)%usws, n_out ) 
2358
2359            CASE ( 'vsws' )
2360               CALL surface_output_sum_up( surf_def_h(0)%vsws,                 &
2361                                           surf_def_h(1)%vsws,                 &
2362                                           surf_lsm_h%vsws,                    &
2363                                           surf_usm_h%vsws,                    &
2364                                           surf_def_v(0)%vsws,                 &
2365                                           surf_lsm_v(0)%vsws,                 &
2366                                           surf_usm_v(0)%vsws,                 &
2367                                           surf_def_v(1)%vsws,                 &
2368                                           surf_lsm_v(1)%vsws,                 &
2369                                           surf_usm_v(1)%vsws,                 &
2370                                           surf_def_v(2)%vsws,                 &
2371                                           surf_lsm_v(2)%vsws,                 &
2372                                           surf_usm_v(2)%vsws,                 &
2373                                           surf_def_v(3)%vsws,                 &
2374                                           surf_lsm_v(3)%vsws,                 &
2375                                           surf_usm_v(3)%vsws, n_out ) 
2376
2377            CASE ( 'shf' )
2378               CALL surface_output_sum_up( surf_def_h(0)%shf,                  &
2379                                           surf_def_h(1)%shf,                  &
2380                                           surf_lsm_h%shf,                     &
2381                                           surf_usm_h%shf,                     & 
2382                                           surf_def_v(0)%shf,                  &
2383                                           surf_lsm_v(0)%shf,                  &
2384                                           surf_usm_v(0)%shf,                  &
2385                                           surf_def_v(1)%shf,                  &
2386                                           surf_lsm_v(1)%shf,                  &
2387                                           surf_usm_v(1)%shf,                  &
2388                                           surf_def_v(2)%shf,                  &
2389                                           surf_lsm_v(2)%shf,                  &
2390                                           surf_usm_v(2)%shf,                  &
2391                                           surf_def_v(3)%shf,                  &
2392                                           surf_lsm_v(3)%shf,                  &
2393                                           surf_usm_v(3)%shf, n_out )
2394
2395            CASE ( 'qsws' )
2396               CALL surface_output_sum_up( surf_def_h(0)%qsws,                 &
2397                                           surf_def_h(1)%qsws,                 &
2398                                           surf_lsm_h%qsws,                    &
2399                                           surf_usm_h%qsws,                    &
2400                                           surf_def_v(0)%qsws,                 &
2401                                           surf_lsm_v(0)%qsws,                 &
2402                                           surf_usm_v(0)%qsws,                 &
2403                                           surf_def_v(1)%qsws,                 &
2404                                           surf_lsm_v(1)%qsws,                 &
2405                                           surf_usm_v(1)%qsws,                 &
2406                                           surf_def_v(2)%qsws,                 &
2407                                           surf_lsm_v(2)%qsws,                 &
2408                                           surf_usm_v(2)%qsws,                 &
2409                                           surf_def_v(3)%qsws,                 &
2410                                           surf_lsm_v(3)%qsws,                 &
2411                                           surf_usm_v(3)%qsws, n_out ) 
2412
2413            CASE ( 'ssws' )
2414               CALL surface_output_sum_up( surf_def_h(0)%ssws,                 &
2415                                           surf_def_h(1)%ssws,                 &
2416                                           surf_lsm_h%ssws,                    &
2417                                           surf_usm_h%ssws,                    &
2418                                           surf_def_v(0)%ssws,                 &
2419                                           surf_lsm_v(0)%ssws,                 &
2420                                           surf_usm_v(0)%ssws,                 &
2421                                           surf_def_v(1)%ssws,                 &
2422                                           surf_lsm_v(1)%ssws,                 &
2423                                           surf_usm_v(1)%ssws,                 &
2424                                           surf_def_v(2)%ssws,                 &
2425                                           surf_lsm_v(2)%ssws,                 &
2426                                           surf_usm_v(2)%ssws,                 &
2427                                           surf_def_v(3)%ssws,                 &
2428                                           surf_lsm_v(3)%ssws,                 &
2429                                           surf_usm_v(3)%ssws, n_out ) 
2430
2431            CASE ( 'qcsws' )
2432               CALL surface_output_sum_up( surf_def_h(0)%qcsws,                &
2433                                           surf_def_h(1)%qcsws,                &
2434                                           surf_lsm_h%qcsws,                   &
2435                                           surf_usm_h%qcsws,                   &
2436                                           surf_def_v(0)%qcsws,                &
2437                                           surf_lsm_v(0)%qcsws,                &
2438                                           surf_usm_v(0)%qcsws,                &
2439                                           surf_def_v(1)%qcsws,                &
2440                                           surf_lsm_v(1)%qcsws,                &
2441                                           surf_usm_v(1)%qcsws,                &
2442                                           surf_def_v(2)%qcsws,                &
2443                                           surf_lsm_v(2)%qcsws,                &
2444                                           surf_usm_v(2)%qcsws,                &
2445                                           surf_def_v(3)%qcsws,                &
2446                                           surf_lsm_v(3)%qcsws,                &
2447                                           surf_usm_v(3)%qcsws, n_out ) 
2448
2449            CASE ( 'ncsws' )
2450               CALL surface_output_sum_up( surf_def_h(0)%ncsws,                &
2451                                           surf_def_h(1)%ncsws,                &
2452                                           surf_lsm_h%ncsws,                   &
2453                                           surf_usm_h%ncsws,                   &
2454                                           surf_def_v(0)%ncsws,                &
2455                                           surf_lsm_v(0)%ncsws,                &
2456                                           surf_usm_v(0)%ncsws,                &
2457                                           surf_def_v(1)%ncsws,                &
2458                                           surf_lsm_v(1)%ncsws,                &
2459                                           surf_usm_v(1)%ncsws,                &
2460                                           surf_def_v(2)%ncsws,                &
2461                                           surf_lsm_v(2)%ncsws,                &
2462                                           surf_usm_v(2)%ncsws,                &
2463                                           surf_def_v(3)%ncsws,                &
2464                                           surf_lsm_v(3)%ncsws,                &
2465                                           surf_usm_v(3)%ncsws, n_out ) 
2466
2467            CASE ( 'qrsws' )
2468               CALL surface_output_sum_up( surf_def_h(0)%qrsws,                &
2469                                           surf_def_h(1)%qrsws,                &
2470                                           surf_lsm_h%qrsws,                   &
2471                                           surf_usm_h%qrsws,                   &
2472                                           surf_def_v(0)%qrsws,                &
2473                                           surf_lsm_v(0)%qrsws,                &
2474                                           surf_usm_v(0)%qrsws,                &
2475                                           surf_def_v(1)%qrsws,                &
2476                                           surf_lsm_v(1)%qrsws,                &
2477                                           surf_usm_v(1)%qrsws,                &
2478                                           surf_def_v(2)%qrsws,                &
2479                                           surf_lsm_v(2)%qrsws,                &
2480                                           surf_usm_v(2)%qrsws,                &
2481                                           surf_def_v(3)%qrsws,                &
2482                                           surf_lsm_v(3)%qrsws,                &
2483                                           surf_usm_v(3)%qrsws, n_out ) 
2484
2485            CASE ( 'nrsws' )
2486               CALL surface_output_sum_up( surf_def_h(0)%nrsws,                &
2487                                           surf_def_h(1)%nrsws,                &
2488                                           surf_lsm_h%nrsws,                   &
2489                                           surf_usm_h%nrsws,                   &
2490                                           surf_def_v(0)%nrsws,                &
2491                                           surf_lsm_v(0)%nrsws,                &
2492                                           surf_usm_v(0)%nrsws,                &
2493                                           surf_def_v(1)%nrsws,                &
2494                                           surf_lsm_v(1)%nrsws,                &
2495                                           surf_usm_v(1)%nrsws,                &
2496                                           surf_def_v(2)%nrsws,                &
2497                                           surf_lsm_v(2)%nrsws,                &
2498                                           surf_usm_v(2)%nrsws,                &
2499                                           surf_def_v(3)%nrsws,                &
2500                                           surf_lsm_v(3)%nrsws,                &
2501                                           surf_usm_v(3)%nrsws, n_out ) 
2502
2503            CASE ( 'sasws' )
2504               CALL surface_output_sum_up( surf_def_h(0)%sasws,                &
2505                                           surf_def_h(1)%sasws,                &
2506                                           surf_lsm_h%sasws,                   &
2507                                           surf_usm_h%sasws,                   &
2508                                           surf_def_v(0)%sasws,                &
2509                                           surf_lsm_v(0)%sasws,                &
2510                                           surf_usm_v(0)%sasws,                &
2511                                           surf_def_v(1)%sasws,                &
2512                                           surf_lsm_v(1)%sasws,                &
2513                                           surf_usm_v(1)%sasws,                &
2514                                           surf_def_v(2)%sasws,                &
2515                                           surf_lsm_v(2)%sasws,                &
2516                                           surf_usm_v(2)%sasws,                &
2517                                           surf_def_v(3)%sasws,                &
2518                                           surf_lsm_v(3)%sasws,                &
2519                                           surf_usm_v(3)%sasws, n_out ) 
2520
2521            CASE ( 'q_surface' )
2522               CALL surface_output_sum_up( surf_def_h(0)%q_surface,            &
2523                                           surf_def_h(1)%q_surface,            &
2524                                           surf_lsm_h%q_surface,               &
2525                                           surf_usm_h%q_surface,               &
2526                                           surf_def_v(0)%q_surface,            &
2527                                           surf_lsm_v(0)%q_surface,            &
2528                                           surf_usm_v(0)%q_surface,            &
2529                                           surf_def_v(1)%q_surface,            &
2530                                           surf_lsm_v(1)%q_surface,            &
2531                                           surf_usm_v(1)%q_surface,            &
2532                                           surf_def_v(2)%q_surface,            &
2533                                           surf_lsm_v(2)%q_surface,            &
2534                                           surf_usm_v(2)%q_surface,            &
2535                                           surf_def_v(3)%q_surface,            &
2536                                           surf_lsm_v(3)%q_surface,            &
2537                                           surf_usm_v(3)%q_surface, n_out ) 
2538
2539            CASE ( 'theta_surface' )
2540               CALL surface_output_sum_up( surf_def_h(0)%pt_surface,           &
2541                                           surf_def_h(1)%pt_surface,           &
2542                                           surf_lsm_h%pt_surface,              &
2543                                           surf_usm_h%pt_surface,              &
2544                                           surf_def_v(0)%pt_surface,           &
2545                                           surf_lsm_v(0)%pt_surface,           &
2546                                           surf_usm_v(0)%pt_surface,           &
2547                                           surf_def_v(1)%pt_surface,           &
2548                                           surf_lsm_v(1)%pt_surface,           &
2549                                           surf_usm_v(1)%pt_surface,           &
2550                                           surf_def_v(2)%pt_surface,           &
2551                                           surf_lsm_v(2)%pt_surface,           &
2552                                           surf_usm_v(2)%pt_surface,           &
2553                                           surf_def_v(3)%pt_surface,           &
2554                                           surf_lsm_v(3)%pt_surface,           &
2555                                           surf_usm_v(3)%pt_surface, n_out ) 
2556
2557            CASE ( 'thetav_surface' )
2558               CALL surface_output_sum_up( surf_def_h(0)%vpt_surface,          &
2559                                           surf_def_h(1)%vpt_surface,          &
2560                                           surf_lsm_h%vpt_surface,             &
2561                                           surf_usm_h%vpt_surface,             &
2562                                           surf_def_v(0)%vpt_surface,          &
2563                                           surf_lsm_v(0)%vpt_surface,          &
2564                                           surf_usm_v(0)%vpt_surface,          &
2565                                           surf_def_v(1)%vpt_surface,          &
2566                                           surf_lsm_v(1)%vpt_surface,          &
2567                                           surf_usm_v(1)%vpt_surface,          &
2568                                           surf_def_v(2)%vpt_surface,          &
2569                                           surf_lsm_v(2)%vpt_surface,          &
2570                                           surf_usm_v(2)%vpt_surface,          &
2571                                           surf_def_v(3)%vpt_surface,          &
2572                                           surf_lsm_v(3)%vpt_surface,          &
2573                                           surf_usm_v(3)%vpt_surface, n_out ) 
2574
2575            CASE ( 'rad_net' )
2576               CALL surface_output_sum_up( surf_def_h(0)%rad_net,              &
2577                                           surf_def_h(1)%rad_net,              &
2578                                           surf_lsm_h%rad_net,                 &
2579                                           surf_usm_h%rad_net,                 &
2580                                           surf_def_v(0)%rad_net,              &
2581                                           surf_lsm_v(0)%rad_net,              &
2582                                           surf_usm_v(0)%rad_net,              &
2583                                           surf_def_v(1)%rad_net,              &
2584                                           surf_lsm_v(1)%rad_net,              &
2585                                           surf_usm_v(1)%rad_net,              &
2586                                           surf_def_v(2)%rad_net,              &
2587                                           surf_lsm_v(2)%rad_net,              &
2588                                           surf_usm_v(2)%rad_net,              &
2589                                           surf_def_v(3)%rad_net,              &
2590                                           surf_lsm_v(3)%rad_net,              &
2591                                           surf_usm_v(3)%rad_net, n_out ) 
2592
2593            CASE ( 'rad_lw_in' )
2594               CALL surface_output_sum_up( surf_def_h(0)%rad_lw_in,            &
2595                                           surf_def_h(1)%rad_lw_in,            &
2596                                           surf_lsm_h%rad_lw_in,               &
2597                                           surf_usm_h%rad_lw_in,               &
2598                                           surf_def_v(0)%rad_lw_in,            &
2599                                           surf_lsm_v(0)%rad_lw_in,            &
2600                                           surf_usm_v(0)%rad_lw_in,            &
2601                                           surf_def_v(1)%rad_lw_in,            &
2602                                           surf_lsm_v(1)%rad_lw_in,            &
2603                                           surf_usm_v(1)%rad_lw_in,            &
2604                                           surf_def_v(2)%rad_lw_in,            &
2605                                           surf_lsm_v(2)%rad_lw_in,            &
2606                                           surf_usm_v(2)%rad_lw_in,            &
2607                                           surf_def_v(3)%rad_lw_in,            &
2608                                           surf_lsm_v(3)%rad_lw_in,            &
2609                                           surf_usm_v(3)%rad_lw_in, n_out ) 
2610
2611            CASE ( 'rad_lw_out' )
2612               CALL surface_output_sum_up( surf_def_h(0)%rad_lw_out,           &
2613                                           surf_def_h(1)%rad_lw_out,           &
2614                                           surf_lsm_h%rad_lw_out,              &
2615                                           surf_usm_h%rad_lw_out,              &
2616                                           surf_def_v(0)%rad_lw_out,           &
2617                                           surf_lsm_v(0)%rad_lw_out,           &
2618                                           surf_usm_v(0)%rad_lw_out,           &
2619                                           surf_def_v(1)%rad_lw_out,           &
2620                                           surf_lsm_v(1)%rad_lw_out,           &
2621                                           surf_usm_v(1)%rad_lw_out,           &
2622                                           surf_def_v(2)%rad_lw_out,           &
2623                                           surf_lsm_v(2)%rad_lw_out,           &
2624                                           surf_usm_v(2)%rad_lw_out,           &
2625                                           surf_def_v(3)%rad_lw_out,           &
2626                                           surf_lsm_v(3)%rad_lw_out,           &
2627                                           surf_usm_v(3)%rad_lw_out, n_out ) 
2628
2629            CASE ( 'rad_sw_in' )
2630               CALL surface_output_sum_up( surf_def_h(0)%rad_sw_in,            &
2631                                           surf_def_h(1)%rad_sw_in,            &
2632                                           surf_lsm_h%rad_sw_in,               &
2633                                           surf_usm_h%rad_sw_in,               &
2634                                           surf_def_v(0)%rad_sw_in,            &
2635                                           surf_lsm_v(0)%rad_sw_in,            &
2636                                           surf_usm_v(0)%rad_sw_in,            &
2637                                           surf_def_v(1)%rad_sw_in,            &
2638                                           surf_lsm_v(1)%rad_sw_in,            &
2639                                           surf_usm_v(1)%rad_sw_in,            &
2640                                           surf_def_v(2)%rad_sw_in,            &
2641                                           surf_lsm_v(2)%rad_sw_in,            &
2642                                           surf_usm_v(2)%rad_sw_in,            &
2643                                           surf_def_v(3)%rad_sw_in,            &
2644                                           surf_lsm_v(3)%rad_sw_in,            &
2645                                           surf_usm_v(3)%rad_sw_in, n_out ) 
2646
2647            CASE ( 'rad_sw_out' )
2648               CALL surface_output_sum_up( surf_def_h(0)%rad_sw_out,           &
2649                                           surf_def_h(1)%rad_sw_out,           &
2650                                           surf_lsm_h%rad_sw_out,              &
2651                                           surf_usm_h%rad_sw_out,              &
2652                                           surf_def_v(0)%rad_sw_out,           &
2653                                           surf_lsm_v(0)%rad_sw_out,           &
2654                                           surf_usm_v(0)%rad_sw_out,           &
2655                                           surf_def_v(1)%rad_sw_out,           &
2656                                           surf_lsm_v(1)%rad_sw_out,           &
2657                                           surf_usm_v(1)%rad_sw_out,           &
2658                                           surf_def_v(2)%rad_sw_out,           &
2659                                           surf_lsm_v(2)%rad_sw_out,           &
2660                                           surf_usm_v(2)%rad_sw_out,           &
2661                                           surf_def_v(3)%rad_sw_out,           &
2662                                           surf_lsm_v(3)%rad_sw_out,           &
2663                                           surf_usm_v(3)%rad_sw_out, n_out ) 
2664
2665            CASE ( 'ghf' )
2666               CALL surface_output_sum_up( surf_def_h(0)%ghf,                  &
2667                                           surf_def_h(1)%ghf,                  &
2668                                           surf_lsm_h%ghf,                     &
2669                                           surf_usm_h%ghf,                     &
2670                                           surf_def_v(0)%ghf,                  &
2671                                           surf_lsm_v(0)%ghf,                  &
2672                                           surf_usm_v(0)%ghf,                  &
2673                                           surf_def_v(1)%ghf,                  &
2674                                           surf_lsm_v(1)%ghf,                  &
2675                                           surf_usm_v(1)%ghf,                  &
2676                                           surf_def_v(2)%ghf,                  &
2677                                           surf_lsm_v(2)%ghf,                  &
2678                                           surf_usm_v(2)%ghf,                  &
2679                                           surf_def_v(3)%ghf,                  &
2680                                           surf_lsm_v(3)%ghf,                  &
2681                                           surf_usm_v(3)%ghf, n_out )         
2682                                                                               
2683            CASE ( 'r_a' )                                                     
2684               CALL surface_output_sum_up( surf_def_h(0)%r_a,                  &
2685                                           surf_def_h(1)%r_a,                  &
2686                                           surf_lsm_h%r_a,                     &
2687                                           surf_usm_h%r_a,                     &
2688                                           surf_def_v(0)%r_a,                  &
2689                                           surf_lsm_v(0)%r_a,                  &
2690                                           surf_usm_v(0)%r_a,                  &
2691                                           surf_def_v(1)%r_a,                  &
2692                                           surf_lsm_v(1)%r_a,                  &
2693                                           surf_usm_v(1)%r_a,                  &
2694                                           surf_def_v(2)%r_a,                  &
2695                                           surf_lsm_v(2)%r_a,                  &
2696                                           surf_usm_v(2)%r_a,                  &
2697                                           surf_def_v(3)%r_a,                  &
2698                                           surf_lsm_v(3)%r_a,                  &
2699                                           surf_usm_v(3)%r_a, n_out )         
2700                                                                               
2701            CASE ( 'r_soil' )                                                 
2702               CALL surface_output_sum_up( surf_def_h(0)%r_soil,               &
2703                                           surf_def_h(1)%r_soil,               &
2704                                           surf_lsm_h%r_soil,                  &
2705                                           surf_usm_h%r_soil,                  &
2706                                           surf_def_v(0)%r_soil,               &
2707                                           surf_lsm_v(0)%r_soil,               &
2708                                           surf_usm_v(0)%r_soil,               &
2709                                           surf_def_v(1)%r_soil,               &
2710                                           surf_lsm_v(1)%r_soil,               &
2711                                           surf_usm_v(1)%r_soil,               &
2712                                           surf_def_v(2)%r_soil,               &
2713                                           surf_lsm_v(2)%r_soil,               &
2714                                           surf_usm_v(2)%r_soil,               &
2715                                           surf_def_v(3)%r_soil,               &
2716                                           surf_lsm_v(3)%r_soil,               &
2717                                           surf_usm_v(3)%r_soil, n_out ) 
2718
2719            CASE ( 'r_canopy' )
2720               CALL surface_output_sum_up( surf_def_h(0)%r_canopy,             &
2721                                           surf_def_h(1)%r_canopy,             &
2722                                           surf_lsm_h%r_canopy,                &
2723                                           surf_usm_h%r_canopy,                &
2724                                           surf_def_v(0)%r_canopy,             &
2725                                           surf_lsm_v(0)%r_canopy,             &
2726                                           surf_usm_v(0)%r_canopy,             &
2727                                           surf_def_v(1)%r_canopy,             &
2728                                           surf_lsm_v(1)%r_canopy,             &
2729                                           surf_usm_v(1)%r_canopy,             &
2730                                           surf_def_v(2)%r_canopy,             &
2731                                           surf_lsm_v(2)%r_canopy,             &
2732                                           surf_usm_v(2)%r_canopy,             &
2733                                           surf_def_v(3)%r_canopy,             &
2734                                           surf_lsm_v(3)%r_canopy,             &
2735                                           surf_usm_v(3)%r_canopy, n_out ) 
2736
2737            CASE ( 'r_s' )
2738               CALL surface_output_sum_up( surf_def_h(0)%r_s,                  &
2739                                           surf_def_h(1)%r_s,                  &
2740                                           surf_lsm_h%r_s,                     &
2741                                           surf_usm_h%r_s,                     &
2742                                           surf_def_v(0)%r_s,                  &
2743                                           surf_lsm_v(0)%r_s,                  &
2744                                           surf_usm_v(0)%r_s,                  &
2745                                           surf_def_v(1)%r_s,                  &
2746                                           surf_lsm_v(1)%r_s,                  &
2747                                           surf_usm_v(1)%r_s,                  &
2748                                           surf_def_v(2)%r_s,                  &
2749                                           surf_lsm_v(2)%r_s,                  &
2750                                           surf_usm_v(2)%r_s,                  &
2751                                           surf_def_v(3)%r_s,                  &
2752                                           surf_lsm_v(3)%r_s,                  &
2753                                           surf_usm_v(3)%r_s, n_out ) 
2754
2755         END SELECT
2756      ENDDO
2757     
2758   
2759   END SUBROUTINE surface_output_averaging
2760   
2761!------------------------------------------------------------------------------!
2762! Description:
2763! ------------
2764!> Sum-up the surface data for average output variables. 
2765!------------------------------------------------------------------------------!   
2766   SUBROUTINE surface_output_sum_up( var_def_h0, var_def_h1,                   &
2767                                     var_lsm_h,  var_usm_h,                    &
2768                                     var_def_v0, var_lsm_v0, var_usm_v0,       &
2769                                     var_def_v1, var_lsm_v1, var_usm_v1,       &
2770                                     var_def_v2, var_lsm_v2, var_usm_v2,       &
2771                                     var_def_v3, var_lsm_v3, var_usm_v3, n_out )
2772   
2773      IMPLICIT NONE
2774
2775      INTEGER(iwp) ::  m          !< running index for surface elements
2776      INTEGER(iwp) ::  n_out      !< index for output variable
2777      INTEGER(iwp) ::  n_surf     !< running index for surface elements
2778     
2779      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h0 !< output variable at upward-facing default-type surfaces
2780      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h1 !< output variable at downward-facing default-type surfaces
2781      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_h  !< output variable at upward-facing natural-type surfaces
2782      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_h  !< output variable at upward-facing urban-type surfaces
2783      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v0 !< output variable at northward-facing default-type surfaces
2784      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v1 !< output variable at southward-facing default-type surfaces
2785      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v2 !< output variable at eastward-facing default-type surfaces
2786      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v3 !< output variable at westward-facing default-type surfaces
2787      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v0 !< output variable at northward-facing natural-type surfaces
2788      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v1 !< output variable at southward-facing natural-type surfaces
2789      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v2 !< output variable at eastward-facing natural-type surfaces
2790      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v3 !< output variable at westward-facing natural-type surfaces
2791      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v0 !< output variable at northward-facing urban-type surfaces
2792      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v1 !< output variable at southward-facing urban-type surfaces
2793      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v2 !< output variable at eastward-facing urban-type surfaces
2794      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v3 !< output variable at westward-facing urban-type surfaces
2795     
2796!     
2797!--   Set counter variable to zero before the variable is written to
2798!--   the output array.
2799      n_surf = 0
2800     
2801!
2802!--   Write the horizontal surfaces.
2803!--   Before each the variable is written to the output data structure, first
2804!--   check if the variable for the respective surface type is defined.
2805!--   If a variable is not defined, skip the block and increment the counter
2806!--   variable by the number of surface elements of this type. Usually this
2807!--   is zere, however, there might be the situation that e.g. urban surfaces
2808!--   are defined but the respective variable is not allocated for this surface
2809!--   type. To write the data on the exact position, increment the counter.
2810      IF ( ALLOCATED( var_def_h0 ) )  THEN
2811         DO  m = 1, surf_def_h(0)%ns 
2812            n_surf                        = n_surf + 1
2813            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2814                                          + var_def_h0(m)
2815         ENDDO
2816      ELSE
2817         n_surf = n_surf + surf_def_h(0)%ns 
2818      ENDIF
2819      IF ( ALLOCATED( var_def_h1 ) )  THEN
2820         DO  m = 1, surf_def_h(1)%ns 
2821            n_surf                   = n_surf + 1
2822            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2823                                          + var_def_h1(m)
2824         ENDDO
2825      ELSE
2826         n_surf = n_surf + surf_def_h(1)%ns
2827      ENDIF
2828      IF ( ALLOCATED( var_lsm_h ) )  THEN
2829         DO  m = 1, surf_lsm_h%ns 
2830            n_surf                        = n_surf + 1
2831            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2832                                          + var_lsm_h(m)
2833         ENDDO
2834      ELSE
2835         n_surf = n_surf + surf_lsm_h%ns
2836      ENDIF
2837      IF ( ALLOCATED( var_usm_h ) )  THEN
2838         DO  m = 1, surf_usm_h%ns 
2839            n_surf                        = n_surf + 1
2840            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2841                                          + var_usm_h(m)
2842         ENDDO
2843      ELSE
2844         n_surf = n_surf + surf_usm_h%ns
2845      ENDIF
2846!
2847!--   Write northward-facing
2848      IF ( ALLOCATED( var_def_v0 ) )  THEN
2849         DO  m = 1, surf_def_v(0)%ns 
2850            n_surf                        = n_surf + 1
2851            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2852                                          + var_def_v0(m)
2853         ENDDO
2854      ELSE
2855         n_surf = n_surf + surf_def_v(0)%ns
2856      ENDIF
2857      IF ( ALLOCATED( var_lsm_v0 ) )  THEN
2858         DO  m = 1, surf_lsm_v(0)%ns 
2859            n_surf                        = n_surf + 1
2860            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2861                                          + var_lsm_v0(m)
2862         ENDDO
2863      ELSE
2864         n_surf = n_surf + surf_lsm_v(0)%ns
2865      ENDIF
2866      IF ( ALLOCATED( var_usm_v0 ) )  THEN
2867         DO  m = 1, surf_usm_v(0)%ns 
2868            n_surf                        = n_surf + 1
2869            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2870                                          + var_usm_v0(m)
2871         ENDDO
2872      ELSE
2873         n_surf = n_surf + surf_usm_v(0)%ns
2874      ENDIF
2875!
2876!--   Write southward-facing
2877      IF ( ALLOCATED( var_def_v1 ) )  THEN
2878         DO  m = 1, surf_def_v(1)%ns 
2879            n_surf                        = n_surf + 1
2880            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2881                                          + var_def_v1(m)
2882         ENDDO
2883      ELSE
2884         n_surf = n_surf + surf_def_v(1)%ns
2885      ENDIF
2886      IF ( ALLOCATED( var_lsm_v1 ) )  THEN
2887         DO  m = 1, surf_lsm_v(1)%ns 
2888            n_surf                        = n_surf + 1
2889            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2890                                          + var_lsm_v1(m)
2891         ENDDO
2892      ELSE
2893         n_surf = n_surf + surf_lsm_v(1)%ns
2894      ENDIF
2895      IF ( ALLOCATED( var_usm_v1 ) )  THEN
2896         DO  m = 1, surf_usm_v(1)%ns 
2897            n_surf                        = n_surf + 1
2898            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2899                                          + var_usm_v1(m)
2900         ENDDO
2901      ELSE
2902         n_surf = n_surf + surf_usm_v(1)%ns
2903      ENDIF
2904!
2905!--   Write eastward-facing
2906      IF ( ALLOCATED( var_def_v2 ) )  THEN
2907         DO  m = 1, surf_def_v(2)%ns 
2908            n_surf                        = n_surf + 1
2909            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2910                                          + var_def_v2(m)
2911         ENDDO
2912      ELSE
2913         n_surf = n_surf + surf_def_v(2)%ns
2914      ENDIF
2915      IF ( ALLOCATED( var_lsm_v2 ) )  THEN
2916         DO  m = 1, surf_lsm_v(2)%ns 
2917            n_surf                        = n_surf + 1
2918            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2919                                          + var_lsm_v2(m)
2920         ENDDO
2921      ELSE
2922         n_surf = n_surf + surf_lsm_v(2)%ns
2923      ENDIF
2924      IF ( ALLOCATED( var_usm_v2 ) )  THEN
2925         DO  m = 1, surf_usm_v(2)%ns 
2926            n_surf                        = n_surf + 1
2927            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2928                                          + var_usm_v2(m)
2929         ENDDO
2930      ELSE
2931         n_surf = n_surf + surf_usm_v(2)%ns
2932      ENDIF
2933!
2934!--   Write westward-facing
2935      IF ( ALLOCATED( var_def_v3 ) )  THEN
2936         DO  m = 1, surf_def_v(3)%ns 
2937            n_surf                        = n_surf + 1
2938            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2939                                          + var_def_v3(m)
2940         ENDDO
2941      ELSE
2942         n_surf = n_surf + surf_def_v(3)%ns
2943      ENDIF
2944      IF ( ALLOCATED( var_lsm_v3 ) )  THEN
2945         DO  m = 1, surf_lsm_v(3)%ns 
2946            n_surf                        = n_surf + 1
2947            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2948                                          + var_lsm_v3(m)
2949         ENDDO
2950      ELSE
2951         n_surf = n_surf + surf_lsm_v(3)%ns
2952      ENDIF
2953      IF ( ALLOCATED( var_usm_v3 ) )  THEN
2954         DO  m = 1, surf_usm_v(3)%ns 
2955            n_surf                        = n_surf + 1
2956            surfaces%var_av(n_surf,n_out) = surfaces%var_av(n_surf,n_out)      &
2957                                          + var_usm_v3(m)
2958         ENDDO
2959      ELSE
2960         n_surf = n_surf + surf_usm_v(3)%ns
2961      ENDIF
2962   
2963   END SUBROUTINE surface_output_sum_up
2964   
2965!------------------------------------------------------------------------------!
2966! Description:
2967! ------------
2968!> Collect the surface data from different types and different orientation.
2969!------------------------------------------------------------------------------!
2970   SUBROUTINE surface_output_collect( var_def_h0, var_def_h1,                  &
2971                                      var_lsm_h,  var_usm_h,                   &
2972                                      var_def_v0, var_lsm_v0, var_usm_v0,      &
2973                                      var_def_v1, var_lsm_v1, var_usm_v1,      &
2974                                      var_def_v2, var_lsm_v2, var_usm_v2,      &
2975                                      var_def_v3, var_lsm_v3, var_usm_v3 )                                               
2976   
2977      IMPLICIT NONE
2978
2979      INTEGER(iwp) ::  m      !< running index for surface elements
2980      INTEGER(iwp) ::  n_surf !< running index for surface elements
2981     
2982      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h0 !< output variable at upward-facing default-type surfaces
2983      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_h1 !< output variable at downward-facing default-type surfaces
2984      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_h  !< output variable at upward-facing natural-type surfaces
2985      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_h  !< output variable at upward-facing urban-type surfaces
2986      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v0 !< output variable at northward-facing default-type surfaces
2987      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v1 !< output variable at southward-facing default-type surfaces
2988      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v2 !< output variable at eastward-facing default-type surfaces
2989      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_def_v3 !< output variable at westward-facing default-type surfaces
2990      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v0 !< output variable at northward-facing natural-type surfaces
2991      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v1 !< output variable at southward-facing natural-type surfaces
2992      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v2 !< output variable at eastward-facing natural-type surfaces
2993      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_lsm_v3 !< output variable at westward-facing natural-type surfaces
2994      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v0 !< output variable at northward-facing urban-type surfaces
2995      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v1 !< output variable at southward-facing urban-type surfaces
2996      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v2 !< output variable at eastward-facing urban-type surfaces
2997      REAL(wp), DIMENSION(:), ALLOCATABLE, INTENT(IN) ::  var_usm_v3 !< output variable at westward-facing urban-type surfaces
2998     
2999!     
3000!--   Set counter variable to zero before the variable is written to
3001!--   the output array.
3002      n_surf = 0
3003     
3004!
3005!--   Write the horizontal surfaces.
3006!--   Before each the variable is written to the output data structure, first
3007!--   check if the variable for the respective surface type is defined.
3008!--   If a variable is not defined, skip the block and increment the counter
3009!--   variable by the number of surface elements of this type. Usually this
3010!--   is zere, however, there might be the situation that e.g. urban surfaces
3011!--   are defined but the respective variable is not allocated for this surface
3012!--   type. To write the data on the exact position, increment the counter.
3013      IF ( ALLOCATED( var_def_h0 ) )  THEN
3014         DO  m = 1, surf_def_h(0)%ns 
3015            n_surf                   = n_surf + 1
3016            surfaces%var_out(n_surf) = var_def_h0(m)
3017         ENDDO
3018      ELSE
3019         n_surf = n_surf + surf_def_h(0)%ns 
3020      ENDIF
3021      IF ( ALLOCATED( var_def_h1 ) )  THEN
3022         DO  m = 1, surf_def_h(1)%ns 
3023            n_surf                   = n_surf + 1
3024            surfaces%var_out(n_surf) = var_def_h1(m)
3025         ENDDO
3026      ELSE
3027         n_surf = n_surf + surf_def_h(1)%ns
3028      ENDIF
3029      IF ( ALLOCATED( var_lsm_h ) )  THEN
3030         DO  m = 1, surf_lsm_h%ns 
3031            n_surf                   = n_surf + 1
3032            surfaces%var_out(n_surf) = var_lsm_h(m)
3033         ENDDO
3034      ELSE
3035         n_surf = n_surf + surf_lsm_h%ns
3036      ENDIF
3037      IF ( ALLOCATED( var_usm_h ) )  THEN
3038         DO  m = 1, surf_usm_h%ns 
3039            n_surf                   = n_surf + 1
3040            surfaces%var_out(n_surf) = var_usm_h(m)
3041         ENDDO
3042      ELSE
3043         n_surf = n_surf + surf_usm_h%ns
3044      ENDIF
3045!
3046!--   Write northward-facing
3047      IF ( ALLOCATED( var_def_v0 ) )  THEN
3048         DO  m = 1, surf_def_v(0)%ns 
3049            n_surf                   = n_surf + 1
3050            surfaces%var_out(n_surf) = var_def_v0(m)
3051         ENDDO
3052      ELSE
3053         n_surf = n_surf + surf_def_v(0)%ns
3054      ENDIF
3055      IF ( ALLOCATED( var_lsm_v0 ) )  THEN
3056         DO  m = 1, surf_lsm_v(0)%ns 
3057            n_surf                   = n_surf + 1
3058            surfaces%var_out(n_surf) = var_lsm_v0(m)
3059         ENDDO
3060      ELSE
3061         n_surf = n_surf + surf_lsm_v(0)%ns
3062      ENDIF
3063      IF ( ALLOCATED( var_usm_v0 ) )  THEN
3064         DO  m = 1, surf_usm_v(0)%ns 
3065            n_surf                   = n_surf + 1
3066            surfaces%var_out(n_surf) = var_usm_v0(m)
3067         ENDDO
3068      ELSE
3069         n_surf = n_surf + surf_usm_v(0)%ns
3070      ENDIF
3071!
3072!--   Write southward-facing
3073      IF ( ALLOCATED( var_def_v1 ) )  THEN
3074         DO  m = 1, surf_def_v(1)%ns 
3075            n_surf                   = n_surf + 1
3076            surfaces%var_out(n_surf) = var_def_v1(m)
3077         ENDDO
3078      ELSE
3079         n_surf = n_surf + surf_def_v(1)%ns
3080      ENDIF
3081      IF ( ALLOCATED( var_lsm_v1 ) )  THEN
3082         DO  m = 1, surf_lsm_v(1)%ns 
3083            n_surf                   = n_surf + 1
3084            surfaces%var_out(n_surf) = var_lsm_v1(m)
3085         ENDDO
3086      ELSE
3087         n_surf = n_surf + surf_lsm_v(1)%ns
3088      ENDIF
3089      IF ( ALLOCATED( var_usm_v1 ) )  THEN
3090         DO  m = 1, surf_usm_v(1)%ns 
3091            n_surf                   = n_surf + 1
3092            surfaces%var_out(n_surf) = var_usm_v1(m)
3093         ENDDO
3094      ELSE
3095         n_surf = n_surf + surf_usm_v(1)%ns
3096      ENDIF
3097!
3098!--   Write eastward-facing
3099      IF ( ALLOCATED( var_def_v2 ) )  THEN
3100         DO  m = 1, surf_def_v(2)%ns 
3101            n_surf                   = n_surf + 1
3102            surfaces%var_out(n_surf) = var_def_v2(m)
3103         ENDDO
3104      ELSE
3105         n_surf = n_surf + surf_def_v(2)%ns
3106      ENDIF
3107      IF ( ALLOCATED( var_lsm_v2 ) )  THEN
3108         DO  m = 1, surf_lsm_v(2)%ns 
3109            n_surf                   = n_surf + 1
3110            surfaces%var_out(n_surf) = var_lsm_v2(m)
3111         ENDDO
3112      ELSE
3113         n_surf = n_surf + surf_lsm_v(2)%ns
3114      ENDIF
3115      IF ( ALLOCATED( var_usm_v2 ) )  THEN
3116         DO  m = 1, surf_usm_v(2)%ns 
3117            n_surf                   = n_surf + 1
3118            surfaces%var_out(n_surf) = var_usm_v2(m)
3119         ENDDO
3120      ELSE
3121         n_surf = n_surf + surf_usm_v(2)%ns
3122      ENDIF
3123!
3124!--   Write westward-facing
3125      IF ( ALLOCATED( var_def_v3 ) )  THEN
3126         DO  m = 1, surf_def_v(3)%ns 
3127            n_surf                   = n_surf + 1
3128            surfaces%var_out(n_surf) = var_def_v3(m)
3129         ENDDO
3130      ELSE
3131         n_surf = n_surf + surf_def_v(3)%ns
3132      ENDIF
3133      IF ( ALLOCATED( var_lsm_v3 ) )  THEN
3134         DO  m = 1, surf_lsm_v(3)%ns 
3135            n_surf                   = n_surf + 1
3136            surfaces%var_out(n_surf) = var_lsm_v3(m)
3137         ENDDO
3138      ELSE
3139         n_surf = n_surf + surf_lsm_v(3)%ns
3140      ENDIF
3141      IF ( ALLOCATED( var_usm_v3 ) )  THEN
3142         DO  m = 1, surf_usm_v(3)%ns 
3143            n_surf                   = n_surf + 1
3144            surfaces%var_out(n_surf) = var_usm_v3(m)
3145         ENDDO
3146      ELSE
3147         n_surf = n_surf + surf_usm_v(3)%ns
3148      ENDIF
3149   
3150   END SUBROUTINE surface_output_collect
3151   
3152   
3153   
3154!    SUBROUTINE data_output_surf_driver (time)
3155!       IMPLICIT NONE
3156!
3157!       REAL(kind=wp),INTENT(IN)    :: time
3158!       CHARACTER(LEN=64)           :: filename
3159!       INTEGER,SAVE                :: file_number=1
3160!
3161!       if(output_mode == "legacy_VTK")   then
3162!           CALL prepare_vtk_data
3163!
3164!           write(filename,'(a,i5.5,a)') 'surface_output_',file_number,'.vtk'
3165!           if(myid == surf_io_pe)   write(9,*) 'Legacy VTK surface output to file: ',trim(filename)
3166!           file_number = file_number+1
3167!
3168!           CALL write_ascii_vtk (filename, time)
3169!       END IF
3170!
3171!       RETURN
3172!    END SUBROUTINE data_output_surf_driver
3173   
3174   
3175!------------------------------------------------------------------------------!
3176! Description:
3177! ------------
3178!> Parin for output of surface parameters
3179!------------------------------------------------------------------------------!
3180    SUBROUTINE surface_output_parin
3181
3182       USE control_parameters,                                                 &
3183           ONLY:  message_string
3184
3185       IMPLICIT NONE
3186
3187       CHARACTER (LEN=80) ::  line  !< dummy string that contains the current line of the parameter file
3188       
3189
3190       NAMELIST /surface_output_parameters/                                    &
3191                                  averaging_interval_surf, data_output_surf,   &
3192                                  dt_dosurf, dt_dosurf_av,                     &
3193                                  skip_time_dosurf, skip_time_dosurf_av
3194                                 
3195       line = ' '
3196 
3197!
3198!--    Try to find the namelist
3199       REWIND ( 11 )
3200       line = ' '
3201       DO WHILE ( INDEX( line, '&surface_output_parameters' ) == 0 )
3202          READ ( 11, '(A)', END=14 )  line
3203       ENDDO
3204       BACKSPACE ( 11 )
3205
3206!
3207!--    Read namelist
3208       READ ( 11, surface_output_parameters, ERR = 10 )
3209!
3210!--    Set flag that indicates that surface data output is switched on
3211       surface_data_output = .TRUE.
3212       GOTO 14
3213
3214 10    BACKSPACE( 11 )
3215       READ( 11 , '(A)') line
3216       CALL parin_fail_message( 'surface_output_parameters', line )     
3217       
3218 14    CONTINUE
3219       
3220
3221    END SUBROUTINE surface_output_parin
3222   
3223   
3224!------------------------------------------------------------------------------!
3225! Description:
3226! ------------
3227!> Check the input parameters for consistency. Further pre-process the given
3228!> output variables, i.e. separate them into average and non-average output
3229!> variables and map them onto internal output array.
3230!------------------------------------------------------------------------------!
3231    SUBROUTINE surface_output_check_parameters
3232
3233       USE control_parameters,                                                 &
3234           ONLY:  averaging_interval, dt_data_output, dt_data_output_av,       &
3235                  message_string, skip_time_data_output,                       &
3236                  skip_time_data_output_av
3237
3238       IMPLICIT NONE
3239
3240       CHARACTER(LEN=100) ::  trimvar !< dummy for single output variable
3241
3242       INTEGER(iwp) ::  av    !< id indicating average or non-average data output
3243       INTEGER(iwp) ::  ilen  !< string length
3244       INTEGER(iwp) ::  n_out !< running index for number of output variables
3245
3246!
3247!--    Check the average interval
3248       IF ( averaging_interval_surf == 9999999.9_wp )  THEN
3249          averaging_interval_surf = averaging_interval
3250       ENDIF
3251!
3252!--    Set the default data-output interval dt_data_output if necessary
3253       IF ( dt_dosurf    == 9999999.9_wp )  dt_dosurf    = dt_data_output
3254       IF ( dt_dosurf_av == 9999999.9_wp )  dt_dosurf_av = dt_data_output
3255
3256       IF ( averaging_interval_surf > dt_dosurf_av )  THEN
3257          WRITE( message_string, * )  'averaging_interval_surf = ',            &
3258                averaging_interval_surf, ' must be <= dt_dosurf_av = ',        &
3259                dt_dosurf_av
3260          CALL message( 'surface_output_check_parameters',                     &
3261                        'PA0087', 1, 2, 0, 6, 0 )
3262       ENDIF
3263
3264!
3265!--   Count number of output variables and separate output strings for
3266!--   average and non-average output variables.
3267      n_out = 0
3268      DO WHILE ( data_output_surf(n_out+1)(1:1) /= ' ' )
3269
3270         n_out = n_out + 1
3271         ilen = LEN_TRIM( data_output_surf(n_out) )
3272         trimvar = TRIM( data_output_surf(n_out) )
3273
3274!
3275!--      Check for data averaging
3276         av = 0
3277         IF ( ilen > 3 )  THEN
3278            IF ( data_output_surf(n_out)(ilen-2:ilen) == '_av' )  THEN
3279               trimvar = data_output_surf(n_out)(1:ilen-3)
3280               av      = 1
3281            ENDIF
3282         ENDIF
3283
3284         dosurf_no(av) = dosurf_no(av) + 1
3285         dosurf(av,dosurf_no(av)) = TRIM( trimvar )
3286
3287!
3288!--      Check if all output variables are known
3289         SELECT CASE ( TRIM( trimvar ) )
3290
3291            CASE ( 'css', 'cssws', 'qsws_liq', 'qsws_soil', 'qsws_veg' )
3292               message_string = TRIM( trimvar ) //                             &
3293                             ' is not yet implemented in the surface output'
3294               CALL message( 'surface_output_check_parameters',                &
3295                             'PA0087', 1, 2, 0, 6, 0 )
3296
3297            CASE ( 'us', 'ts', 'qs', 'ss', 'qcs', 'ncs', 'qrs', 'nrs', 'ol' )
3298
3299            CASE ( 'z0', 'z0h', 'z0q' )
3300
3301            CASE ( 'theta1', 'qv1', 'thetav1' )
3302
3303            CASE ( 'usws', 'vsws' )
3304
3305            CASE ( 'qcsws', 'ncsws', 'qrsws', 'nrsws', 'sasws' )
3306
3307            CASE ( 'shf', 'qsws', 'ssws' )
3308
3309            CASE ( 'q_surface', 'theta_surface', 'thetav_surface' )
3310
3311            CASE ( 'rad_net' )
3312
3313            CASE ( 'rad_lw_in', 'rad_lw_out', 'rad_sw_in', 'rad_sw_out' )
3314
3315            CASE ( 'ghf' )
3316
3317            CASE ( 'r_a', 'r_canopy', 'r_soil', 'r_s' )
3318
3319            CASE DEFAULT
3320               message_string = TRIM( trimvar ) //                             &
3321                             ' is not part of the surface output'
3322               CALL message( 'surface_output_check_parameters',                &
3323                             'PA0087', 1, 2, 0, 6, 0 )
3324         END SELECT
3325
3326       ENDDO
3327
3328    END SUBROUTINE surface_output_check_parameters
3329   
3330   
3331!------------------------------------------------------------------------------!
3332! Description:
3333! ------------
3334!> Last action.
3335!------------------------------------------------------------------------------!
3336    SUBROUTINE surface_output_last_action( av )
3337 
3338      USE control_parameters,                                                  &
3339          ONLY:  io_blocks, io_group
3340
3341      USE pegrid,                                                              &
3342          ONLY:  comm2d, ierr
3343
3344      IMPLICIT NONE
3345     
3346      INTEGER(iwp) ::  av     !< id indicating average or non-average data output
3347      INTEGER(iwp) ::  i      !< loop index
3348
3349!
3350!--   Return, if nothing to output
3351      IF ( dosurf_no(av) == 0 )  RETURN
3352!
3353!--   Open files
3354      CALL check_open( 25+av )
3355!
3356!--   Write time coordinate
3357      DO  i = 0, io_blocks-1
3358         IF ( i == io_group )  THEN
3359            WRITE ( 25+av )  LEN_TRIM( 'END' )
3360            WRITE ( 25+av )  'END'
3361         ENDIF
3362#if defined( __parallel )
3363         CALL MPI_BARRIER( comm2d, ierr )
3364#endif
3365      ENDDO   
3366
3367    END SUBROUTINE surface_output_last_action   
3368
3369!--Private Subroutines (Only called from tinside his module)
3370
3371
3372
3373!    SUBROUTINE write_ascii_vtk (filename,time)
3374!       IMPLICIT NONE
3375!       CHARACTER(LEN=*),INTENT(IN)   :: filename
3376!       REAL(kind=wp),INTENT(IN)      :: time
3377!
3378!       INTEGER                       :: file_id=251
3379!       INTEGER                       :: n
3380!
3381!       IF(myid == surf_io_pe)  then
3382!          open(file_id,file=TRIM(filename))                                 ! Please use check_open in the final version
3383!
3384!          write(file_id,'(a)') "# vtk DataFile Version 3.0"
3385!          write(file_id,'(a,f8.2,a)') "legacy vtk File generated by PALM,  simulation time = ",time," sec"
3386!          write(file_id,'(a)') "ASCII"
3387!
3388!          write(file_id,'(a)') "DATASET POLYDATA"
3389!          write(file_id,'(a,i5,a)') "POINTS ",npo," float"
3390!          do n=1,npo
3391!             write (file_id,'(8f10.1)') points(1:3,n)
3392!          end do
3393!
3394!          write(file_id,'(a,8i10)') "POLYGONS ",npg, 5*npg
3395!          do n=1,npg
3396!             write(file_id,'(8i10)') polygons(:,n)
3397!          end do
3398!
3399!          write(file_id,'(a,i10)') "CELL_DATA ",npg
3400!          write(file_id,'(a,i10)') "SCALARS cell_scalars float 1 "
3401!          write(file_id,'(a,i10)') "LOOKUP_TABLE default "
3402!
3403!          do n=1,npg
3404!             write(file_id,'(f10.2)') cell_scalars(n)
3405!          end do
3406!
3407!          close (file_id)
3408!
3409!          DEALLOCATE (points,polygons,cell_scalars)
3410!       ENDIF
3411!
3412!       RETURN
3413!    END SUBROUTINE write_ascii_vtk
3414
3415END MODULE surface_output_mod
Note: See TracBrowser for help on using the repository browser.