Home

Context Navigation

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

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