Home

Context Navigation

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

Last change on this file since 4460 was 4444, checked in by raasch, 5 years ago
bugfix: cpp-directives for serial mode added
Property svn:keywords set to `Id`
File size: 245.0 KB

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