Home

Context Navigation

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

Last change on this file since 3765 was 3762, checked in by suehring, 5 years ago
remove unused variables in surface_data_output_mod and use compiler diretives to declare variables only if netcdf4-output is used
Property svn:keywords set to `Id`
File size: 243.3 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |