Home

Context Navigation

urban_surface_mod.f90 @ 3223

Last change on this file since 3223 was 3223, checked in by suehring, 4 years ago
Bugfix for last commit
Property svn:keywords set to `Id`
File size: 434.1 KB

Line
1	!> @file urban_surface_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 2015-2018 Czech Technical University in Prague
18	! Copyright 2015-2018 Institute of Computer Science of the
19	! Czech Academy of Sciences, Prague
20	! Copyright 1997-2018 Leibniz Universitaet Hannover
21	!------------------------------------------------------------------------------!
22	!
23	! Current revisions:
24	! ------------------
25	!
26	!
27	! Former revisions:
28	! -----------------
29	! $Id: urban_surface_mod.f90 3223 2018-08-30 13:48:17Z suehring $
30	! Bugfix for commit 3222
31	!
32	! 3222 2018-08-30 13:35:35Z suehring
33	! Introduction of surface array for type and its name
34	!
35	! 3203 2018-08-23 10:48:36Z suehring
36	! Revise bulk parameter for emissivity at ground-floor level
37	!
38	! 3196 2018-08-13 12:26:14Z maronga
39	! Added maximum aerodynamic resistance of 300 for horiztonal surfaces.
40	!
41	! 3176 2018-07-26 17:12:48Z suehring
42	! Bugfix, update virtual potential surface temparture, else heat fluxes on
43	! roofs might become unphysical
44	!
45	! 3152 2018-07-19 13:26:52Z suehring
46	! Initialize q_surface, which might be used in surface_layer_fluxes
47	!
48	! 3151 2018-07-19 08:45:38Z raasch
49	! remaining preprocessor define strings __check removed
50	!
51	! 3136 2018-07-16 14:48:21Z suehring
52	! Limit also roughness length for heat and moisture where necessary
53	!
54	! 3123 2018-07-12 16:21:53Z suehring
55	! Correct working precision for INTEGER number
56	!
57	! 3115 2018-07-10 12:49:26Z suehring
58	! Additional building type to represent bridges
59	!
60	! 3091 2018-06-28 16:20:35Z suehring
61	! - Limit aerodynamic resistance at vertical walls.
62	! - Add check for local roughness length not exceeding surface-layer height and
63	! limit roughness length where necessary.
64	!
65	! 3065 2018-06-12 07:03:02Z Giersch
66	! Unused array dxdir was removed, dz was replaced by dzu to consider vertical
67	! grid stretching
68	!
69	! 3049 2018-05-29 13:52:36Z Giersch
70	! Error messages revised
71	!
72	! 3045 2018-05-28 07:55:41Z Giersch
73	! Error message added
74	!
75	! 3029 2018-05-23 12:19:17Z raasch
76	! bugfix: close unit 151 instead of 90
77	!
78	! 3014 2018-05-09 08:42:38Z maronga
79	! Added pc_transpiration_rate
80	!
81	! 2977 2018-04-17 10:27:57Z kanani
82	! Implement changes from branch radiation (r2948-2971) with minor modifications.
83	! (moh.hefny):
84	! Extended exn for all model domain height to avoid the need to get nzut.
85	!
86	! 2963 2018-04-12 14:47:44Z suehring
87	! Introduce index for vegetation/wall, pavement/green-wall and water/window
88	! surfaces, for clearer access of surface fraction, albedo, emissivity, etc. .
89	!
90	! 2943 2018-04-03 16:17:10Z suehring
91	! Calculate exner function at all height levels and remove some un-used
92	! variables.
93	!
94	! 2932 2018-03-26 09:39:22Z maronga
95	! renamed urban_surface_par to urban_surface_parameters
96	!
97	! 2921 2018-03-22 15:05:23Z Giersch
98	! The activation of spinup has been moved to parin
99	!
100	! 2920 2018-03-22 11:22:01Z kanani
101	! Remove unused pcbl, npcbl from ONLY list
102	! moh.hefny:
103	! Fixed bugs introduced by new structures and by moving radiation interaction
104	! into radiation_model_mod.f90.
105	! Bugfix: usm data output 3D didn't respect directions
106	!
107	! 2906 2018-03-19 08:56:40Z Giersch
108	! Local variable ids has to be initialized with a value of -1 in
109	! usm_average_3d_data
110	!
111	! 2894 2018-03-15 09:17:58Z Giersch
112	! Calculations of the index range of the subdomain on file which overlaps with
113	! the current subdomain are already done in read_restart_data_mod,
114	! usm_read/write_restart_data have been renamed to usm_r/wrd_local, variable
115	! named found has been introduced for checking if restart data was found,
116	! reading of restart strings has been moved completely to
117	! read_restart_data_mod, usm_rrd_local is already inside the overlap loop
118	! programmed in read_restart_data_mod, SAVE attribute added where necessary,
119	! deallocation and allocation of some arrays have been changed to take care of
120	! different restart files that can be opened (index i), the marker *** end usm
121	! *** is not necessary anymore, strings and their respective lengths are
122	! written out and read now in case of restart runs to get rid of prescribed
123	! character lengths
124	!
125	! 2805 2018-02-14 17:00:09Z suehring
126	! Initialization of resistances.
127	!
128	! 2797 2018-02-08 13:24:35Z suehring
129	! Comment concerning output of ground-heat flux added.
130	!
131	! 2766 2018-01-22 17:17:47Z kanani
132	! Removed redundant commas, added some blanks
133	!
134	! 2765 2018-01-22 11:34:58Z maronga
135	! Major bugfix in calculation of f_shf. Adjustment of roughness lengths in
136	! building_pars
137	!
138	! 2750 2018-01-15 16:26:51Z knoop
139	! Move flag plant canopy to modules
140	!
141	! 2737 2018-01-11 14:58:11Z kanani
142	! Removed unused variables t_surf_whole...
143	!
144	! 2735 2018-01-11 12:01:27Z suehring
145	! resistances are saved in surface attributes
146	!
147	! 2723 2018-01-05 09:27:03Z maronga
148	! Bugfix for spinups (end_time was increased twice in case of LSM + USM runs)
149	!
150	! 2720 2018-01-02 16:27:15Z kanani
151	! Correction of comment
152	!
153	! 2718 2018-01-02 08:49:38Z maronga
154	! Corrected "Former revisions" section
155	!
156	! 2705 2017-12-18 11:26:23Z maronga
157	! Changes from last commit documented
158	!
159	! 2703 2017-12-15 20:12:38Z maronga
160	! Workaround for calculation of r_a
161	!
162	! 2696 2017-12-14 17:12:51Z kanani
163	! - Change in file header (GPL part)
164	! - Bugfix in calculation of pt_surface and related fluxes. (BM)
165	! - Do not write surface temperatures onto pt array as this might cause
166	! problems with nesting. (MS)
167	! - Revised calculation of pt1 (now done in surface_layer_fluxes).
168	! Bugfix, f_shf_window and f_shf_green were not set at vertical surface
169	! elements. (MS)
170	! - merged with branch ebsolver
171	! green building surfaces do not evaporate yet
172	! properties of green wall layers and window layers are taken from wall layers
173	! this input data is missing. (RvT)
174	! - Merged with branch radiation (developed by Mohamed Salim)
175	! - Revised initialization. (MS)
176	! - Rename emiss_surf into emissivity, roughness_wall into z0, albedo_surf into
177	! albedo. (MS)
178	! - Move first call of usm_radiatin from usm_init to init_3d_model
179	! - fixed problem with near surface temperature
180	! - added near surface temperature t_surf_10cm_h(m), t_surf_10cm_v(l)%t(m)
181	! - does not work with temp profile including stability, ol
182	! t_surf_10cm = pt1 now
183	! - merged with 2357 bugfix, error message for nopointer version
184	! - added indoor model coupling with wall heat flux
185	! - added green substrate/ dry vegetation layer for buildings
186	! - merged with 2232 new surface-type structure
187	! - added transmissivity of window tiles
188	! - added MOSAIK tile approach for 3 different surfaces (RvT)
189	!
190	! 2583 2017-10-26 13:58:38Z knoop
191	! Bugfix: reverted MPI_Win_allocate_cptr introduction in last commit
192	!
193	! 2582 2017-10-26 13:19:46Z hellstea
194	! Workaround for gnufortran compiler added in usm_calc_svf. CALL MPI_Win_allocate is
195	! replaced by CALL MPI_Win_allocate_cptr if defined ( __gnufortran ).
196	!
197	! 2544 2017-10-13 18:09:32Z maronga
198	! Date and time quantities are now read from date_and_time_mod. Solar constant is
199	! read from radiation_model_mod
200	!
201	! 2516 2017-10-04 11:03:04Z suehring
202	! Remove tabs
203	!
204	! 2514 2017-10-04 09:52:37Z suehring
205	! upper bounds of 3d output changed from nx+1,ny+1 to nx,ny
206	! no output of ghost layer data
207	!
208	! 2350 2017-08-15 11:48:26Z kanani
209	! Bugfix and error message for nopointer version.
210	! Additional "! defined(__nopointer)" as workaround to enable compilation of
211	! nopointer version.
212	!
213	! 2318 2017-07-20 17:27:44Z suehring
214	! Get topography top index via Function call
215	!
216	! 2317 2017-07-20 17:27:19Z suehring
217	! Bugfix: adjust output of shf. Added support for spinups
218	!
219	! 2287 2017-06-15 16:46:30Z suehring
220	! Bugfix in determination topography-top index
221	!
222	! 2269 2017-06-09 11:57:32Z suehring
223	! Enable restart runs with different number of PEs
224	! Bugfixes nopointer branch
225	!
226	! 2258 2017-06-08 07:55:13Z suehring
227	! Bugfix, add pre-preprocessor directives to enable non-parrallel mode
228	!
229	! 2233 2017-05-30 18:08:54Z suehring
230	!
231	! 2232 2017-05-30 17:47:52Z suehring
232	! Adjustments according to new surface-type structure. Remove usm_wall_heat_flux;
233	! insteat, heat fluxes are directly applied in diffusion_s.
234	!
235	! 2213 2017-04-24 15:10:35Z kanani
236	! Removal of output quantities usm_lad and usm_canopy_hr
237	!
238	! 2209 2017-04-19 09:34:46Z kanani
239	! cpp switch __mpi3 removed,
240	! minor formatting,
241	! small bugfix for division by zero (Krc)
242	!
243	! 2113 2017-01-12 13:40:46Z kanani
244	! cpp switch __mpi3 added for MPI-3 standard code (Ketelsen)
245	!
246	! 2071 2016-11-17 11:22:14Z maronga
247	! Small bugfix (Resler)
248	!
249	! 2031 2016-10-21 15:11:58Z knoop
250	! renamed variable rho to rho_ocean
251	!
252	! 2024 2016-10-12 16:42:37Z kanani
253	! Bugfixes in deallocation of array plantt and reading of csf/csfsurf,
254	! optimization of MPI-RMA operations,
255	! declaration of pcbl as integer,
256	! renamed usm_radnet -> usm_rad_net, usm_canopy_khf -> usm_canopy_hr,
257	! splitted arrays svf -> svf & csf, svfsurf -> svfsurf & csfsurf,
258	! use of new control parameter varnamelength,
259	! added output variables usm_rad_ressw, usm_rad_reslw,
260	! minor formatting changes,
261	! minor optimizations.
262	!
263	! 2011 2016-09-19 17:29:57Z kanani
264	! Major reformatting according to PALM coding standard (comments, blanks,
265	! alphabetical ordering, etc.),
266	! removed debug_prints,
267	! removed auxiliary SUBROUTINE get_usm_info, instead, USM flag urban_surface is
268	! defined in MODULE control_parameters (modules.f90) to avoid circular
269	! dependencies,
270	! renamed canopy_heat_flux to pc_heating_rate, as meaning of quantity changed.
271	!
272	! 2007 2016-08-24 15:47:17Z kanani
273	! Initial revision
274	!
275	!
276	! Description:
277	! ------------
278	! 2016/6/9 - Initial version of the USM (Urban Surface Model)
279	! authors: Jaroslav Resler, Pavel Krc
280	! (Czech Technical University in Prague and Institute of
281	! Computer Science of the Czech Academy of Sciences, Prague)
282	! with contributions: Michal Belda, Nina Benesova, Ondrej Vlcek
283	! partly inspired by PALM LSM (B. Maronga)
284	! parameterizations of Ra checked with TUF3D (E. S. Krayenhoff)
285	!> Module for Urban Surface Model (USM)
286	!> The module includes:
287	!> 1. radiation model with direct/diffuse radiation, shading, reflections
288	!> and integration with plant canopy
289	!> 2. wall and wall surface model
290	!> 3. surface layer energy balance
291	!> 4. anthropogenic heat (only from transportation so far)
292	!> 5. necessary auxiliary subroutines (reading inputs, writing outputs,
293	!> restart simulations, ...)
294	!> It also make use of standard radiation and integrates it into
295	!> urban surface model.
296	!>
297	!> Further work:
298	!> -------------
299	!> 1. Remove global arrays surfouts, surfoutl and only keep track of radiosity
300	!> from surfaces that are visible from local surfaces (i.e. there is a SVF
301	!> where target is local). To do that, radiosity will be exchanged after each
302	!> reflection step using MPI_Alltoall instead of current MPI_Allgather.
303	!>
304	!> 2. Temporarily large values of surface heat flux can be observed, up to
305	!> 1.2 Km/s, which seem to be not realistic.
306	!>
307	!> @todo Output of _av variables in case of restarts
308	!> @todo Revise flux conversion in energy-balance solver
309	!> @todo Bugfixing in nopointer branch
310	!> @todo Check optimizations for RMA operations
311	!> @todo Alternatives for MPI_WIN_ALLOCATE? (causes problems with openmpi)
312	!> @todo Check for load imbalances in CPU measures, e.g. for exchange_horiz_prog
313	!> factor 3 between min and max time
314	!> @todo Move setting of flag indoor_model to indoor_model_mod once available
315	!> @todo Check divisions in wtend (etc.) calculations for possible division
316	!> by zero, e.g. in case fraq(0,m) + fraq(1,m) = 0?!
317	!> @todo Use unit 90 for OPEN/CLOSE of input files (FK)
318	!> @todo Move plant canopy stuff into plant canopy code
319	!------------------------------------------------------------------------------!
320	MODULE urban_surface_mod
321
322	#if ! defined( __nopointer )
323	USE arrays_3d, &
324	ONLY: dzu, hyp, zu, pt, pt_1, pt_2, p, u, v, w, hyp, tend
325	#endif
326
327	USE cloud_parameters, &
328	ONLY: cp, r_d
329
330	USE constants, &
331	ONLY: pi
332
333	USE control_parameters, &
334	ONLY: coupling_start_time, topography, dt_3d, humidity, &
335	intermediate_timestep_count, initializing_actions, &
336	intermediate_timestep_count_max, simulated_time, end_time, &
337	timestep_scheme, tsc, coupling_char, io_blocks, io_group, &
338	message_string, time_since_reference_point, surface_pressure, &
339	g, pt_surface, large_scale_forcing, lsf_surf, spinup, &
340	spinup_pt_mean, spinup_time, time_do3d, dt_do3d, &
341	average_count_3d, varnamelength, urban_surface, kappa, &
342	plant_canopy
343
344	USE cpulog, &
345	ONLY: cpu_log, log_point, log_point_s
346
347	USE date_and_time_mod, &
348	ONLY: time_utc_init
349
350	USE grid_variables, &
351	ONLY: dx, dy, ddx, ddy, ddx2, ddy2
352
353	USE indices, &
354	ONLY: nx, ny, nnx, nny, nnz, nxl, nxlg, nxr, nxrg, nyn, nyng, nys, &
355	nysg, nzb, nzt, nbgp, wall_flags_0
356
357	USE, INTRINSIC :: iso_c_binding
358
359	USE kinds
360
361	USE pegrid
362
363	USE plant_canopy_model_mod, &
364	ONLY: pc_heating_rate, pc_transpiration_rate
365
366	USE radiation_model_mod, &
367	ONLY: albedo_type, radiation_interaction, calc_zenith, zenith, &
368	radiation, rad_sw_in, rad_lw_in, rad_sw_out, rad_lw_out, &
369	sigma_sb, solar_constant, sun_direction, sun_dir_lat, &
370	sun_dir_lon, &
371	force_radiation_call, surfinsw, surfinlw, surfinswdir, &
372	surfinswdif, surfoutsw, surfoutlw, surfins,nsvfl, svf, svfsurf, &
373	surfinl, surfinlwdif, rad_sw_in_dir, rad_sw_in_diff, &
374	rad_lw_in_diff, surfouts, surfoutl, surfoutsl, surfoutll, surf, &
375	surfl, nsurfl, nsurfs, surfstart, pcbinsw, pcbinlw, &
376	iup_u, inorth_u, isouth_u, ieast_u, iwest_u, iup_l, &
377	inorth_l, isouth_l, ieast_l, iwest_l, id, &
378	iz, iy, ix, idir, jdir, kdir, nsurf_type, nsurf, idsvf, ndsvf, &
379	iup_a, idown_a, inorth_a, isouth_a, ieast_a, iwest_a, &
380	idcsf, ndcsf, kdcsf, pct, &
381	startland, endland, startwall, endwall, skyvf, skyvft
382
383	USE statistics, &
384	ONLY: hom, statistic_regions
385
386	USE surface_mod, &
387	ONLY: get_topography_top_index_ji, get_topography_top_index, &
388	ind_pav_green, ind_veg_wall, ind_wat_win, surf_usm_h, &
389	surf_usm_v, surface_restore_elements
390
391
392	IMPLICIT NONE
393
394
395	!-- configuration parameters (they can be setup in PALM config)
396	LOGICAL :: usm_material_model = .TRUE. !< flag parameter indicating wheather the model of heat in materials is used
397	LOGICAL :: usm_anthropogenic_heat = .FALSE. !< flag parameter indicating wheather the anthropogenic heat sources (e.g.transportation) are used
398	LOGICAL :: force_radiation_call_l = .FALSE. !< flag parameter for unscheduled radiation model calls
399	LOGICAL :: indoor_model = .FALSE. !< whether to use the indoor model
400	LOGICAL :: read_wall_temp_3d = .FALSE.
401
402
403	INTEGER(iwp) :: building_type = 1 !< default building type (preleminary setting)
404	INTEGER(iwp) :: land_category = 2 !< default category for land surface
405	INTEGER(iwp) :: wall_category = 2 !< default category for wall surface over pedestrian zone
406	INTEGER(iwp) :: pedestrian_category = 2 !< default category for wall surface in pedestrian zone
407	INTEGER(iwp) :: roof_category = 2 !< default category for root surface
408	REAL(wp) :: roughness_concrete = 0.001_wp !< roughness length of average concrete surface
409	!
410	!-- Indices of input attributes for (above) ground floor level
411	INTEGER(iwp) :: ind_alb_wall = 38 !< index in input list for albedo_type of wall fraction
412	INTEGER(iwp) :: ind_alb_green = 39 !< index in input list for albedo_type of green fraction
413	INTEGER(iwp) :: ind_alb_win = 40 !< index in input list for albedo_type of window fraction
414	INTEGER(iwp) :: ind_emis_wall_agfl = 14 !< index in input list for wall emissivity, above ground floor level
415	INTEGER(iwp) :: ind_emis_wall_gfl = 32 !< index in input list for wall emissivity, ground floor level
416	INTEGER(iwp) :: ind_emis_green_agfl = 15 !< index in input list for green emissivity, above ground floor level
417	INTEGER(iwp) :: ind_emis_green_gfl = 33 !< index in input list for green emissivity, ground floor level
418	INTEGER(iwp) :: ind_emis_win_agfl = 16 !< index in input list for window emissivity, above ground floor level
419	INTEGER(iwp) :: ind_emis_win_gfl = 34 !< index in input list for window emissivity, ground floor level
420	INTEGER(iwp) :: ind_green_frac_w_agfl = 2 !< index in input list for green fraction on wall, above ground floor level
421	INTEGER(iwp) :: ind_green_frac_w_gfl = 23 !< index in input list for green fraction on wall, ground floor level
422	INTEGER(iwp) :: ind_green_frac_r_agfl = 3 !< index in input list for green fraction on roof, above ground floor level
423	INTEGER(iwp) :: ind_green_frac_r_gfl = 24 !< index in input list for green fraction on roof, ground floor level
424	INTEGER(iwp) :: ind_hc1_agfl = 6 !< index in input list for heat capacity at first wall layer, above ground floor level
425	INTEGER(iwp) :: ind_hc1_gfl = 26 !< index in input list for heat capacity at first wall layer, ground floor level
426	INTEGER(iwp) :: ind_hc2_agfl = 7 !< index in input list for heat capacity at second wall layer, above ground floor level
427	INTEGER(iwp) :: ind_hc2_gfl = 27 !< index in input list for heat capacity at second wall layer, ground floor level
428	INTEGER(iwp) :: ind_hc3_agfl = 8 !< index in input list for heat capacity at third wall layer, above ground floor level
429	INTEGER(iwp) :: ind_hc3_gfl = 28 !< index in input list for heat capacity at third wall layer, ground floor level
430	INTEGER(iwp) :: ind_gflh = 20 !< index in input list for ground floor level height
431	INTEGER(iwp) :: ind_lai_r_agfl = 4 !< index in input list for LAI on roof, above ground floor level
432	INTEGER(iwp) :: ind_lai_r_gfl = 4 !< index in input list for LAI on roof, ground floor level
433	INTEGER(iwp) :: ind_lai_w_agfl = 5 !< index in input list for LAI on wall, above ground floor level
434	INTEGER(iwp) :: ind_lai_w_gfl = 25 !< index in input list for LAI on wall, ground floor level
435	INTEGER(iwp) :: ind_tc1_agfl = 9 !< index in input list for thermal conductivity at first wall layer, above ground floor level
436	INTEGER(iwp) :: ind_tc1_gfl = 29 !< index in input list for thermal conductivity at first wall layer, ground floor level
437	INTEGER(iwp) :: ind_tc2_agfl = 10 !< index in input list for thermal conductivity at second wall layer, above ground floor level
438	INTEGER(iwp) :: ind_tc2_gfl = 30 !< index in input list for thermal conductivity at second wall layer, ground floor level
439	INTEGER(iwp) :: ind_tc3_agfl = 11 !< index in input list for thermal conductivity at third wall layer, above ground floor level
440	INTEGER(iwp) :: ind_tc3_gfl = 31 !< index in input list for thermal conductivity at third wall layer, ground floor level
441	INTEGER(iwp) :: ind_thick_1 = 41 !< index for wall layer thickness - 1st layer
442	INTEGER(iwp) :: ind_thick_2 = 42 !< index for wall layer thickness - 2nd layer
443	INTEGER(iwp) :: ind_thick_3 = 43 !< index for wall layer thickness - 3rd layer
444	INTEGER(iwp) :: ind_thick_4 = 44 !< index for wall layer thickness - 4th layer
445	INTEGER(iwp) :: ind_trans_agfl = 17 !< index in input list for window transmissivity, above ground floor level
446	INTEGER(iwp) :: ind_trans_gfl = 35 !< index in input list for window transmissivity, ground floor level
447	INTEGER(iwp) :: ind_wall_frac_agfl = 0 !< index in input list for wall fraction, above ground floor level
448	INTEGER(iwp) :: ind_wall_frac_gfl = 21 !< index in input list for wall fraction, ground floor level
449	INTEGER(iwp) :: ind_win_frac_agfl = 1 !< index in input list for window fraction, above ground floor level
450	INTEGER(iwp) :: ind_win_frac_gfl = 22 !< index in input list for window fraction, ground floor level
451	INTEGER(iwp) :: ind_z0_agfl = 18 !< index in input list for z0, above ground floor level
452	INTEGER(iwp) :: ind_z0_gfl = 36 !< index in input list for z0, ground floor level
453	INTEGER(iwp) :: ind_z0qh_agfl = 19 !< index in input list for z0h / z0q, above ground floor level
454	INTEGER(iwp) :: ind_z0qh_gfl = 37 !< index in input list for z0h / z0q, ground floor level
455
456
457	REAL(wp) :: roof_height_limit = 4._wp !< height for distinguish between land surfaces and roofs
458	REAL(wp) :: ground_floor_level = 4.0_wp !< default ground floor level
459
460
461	CHARACTER(37), DIMENSION(0:7), PARAMETER :: building_type_name = (/ &
462	'user-defined ', & ! 0
463	'residential - 1950 ', & ! 1
464	'residential 1951 - 2000 ', & ! 2
465	'residential 2001 - ', & ! 3
466	'office - 1950 ', & ! 4
467	'office 1951 - 2000 ', & ! 5
468	'office 2001 - ', & ! 6
469	'bridges ' & ! 7
470	/)
471	!
472	!-- building parameters, 4 different types
473	!-- 0 - wall fraction, 1- window fraction, 2 - green fraction on wall, 3- green fraction
474	!-- at roof, 4 - lai of green fraction at roof, 5 - lai of green fraction at wall,
475	!-- 6 - heat capacity of wall layer 1, 7 - heat capacity of wall layer 2,
476	!-- 8 - heat capacity of wall layer 3, 9 - thermal conductivity of wall layer 1,
477	!-- 10 - thermal conductivity of wall layer 2, 11 - thermal conductivity of wall layer 3,
478	!-- 12 - indoor target summer temperature ( K ), 13 - indoor target winter temperature (K),
479	!-- 14 - emissivity of wall fraction, 15 - emissivity of green fraction, 16 - emissivity of window fraction,
480	!-- 17 - transmissivity of window fraction, 18 - z0, 19 - z0h/z0q, 20 - ground floor height,
481	!-- 21 - ground floor wall fraction, 22 - ground floor window fraction, 23 ground floor green fraction,
482	!-- 24 - ground floor green fraction on roof, 25 - ground floor lai of green fraction,
483	!-- 26 - ground floor heat capacity of wall layer 1, 27 - ground floor heat capacity of wall layer 1,
484	!-- 28 - ground floor heat capacity of wall layer 3, 29 - ground floor thermal conductivity of wall layer 1,
485	!-- 30 - ground floor thermal conductivity of wall layer 2, 31 - ground floor thermal conductivity of wall layer 3,
486	!-- 32 - ground floor emissivity of wall fraction, 33 - ground floor emissivity of green fraction,
487	!-- 34 - ground floor emissivity of window fraction, 35 - ground floor transmissivity of window fraction,
488	!-- 36 - ground floor z0, 37 - ground floor z0h/z0q, 38 - albedo type wall fraction
489	!-- 39 - albedo type green fraction, 40 - albedo type window fraction
490	!-- 41 - wall layer thickness - 1st layer, 42 - wall layer thickness - 2nd layer,
491	!-- 43 - wall layer thickness - 3rd layer, 44 - wall layer thickness - 4th layer,
492	!-- 45 - heat capacity of the wall surface, 46 - heat conductivity
493	!-- Please note, only preleminary dummy values so far!
494	REAL(wp), DIMENSION(0:46,1:7), PARAMETER :: building_pars = RESHAPE( (/ &
495	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
496	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
497	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
498	0.001_wp, 0.0001_wp, 4.0_wp, & !parameter 18-20
499	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
500	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
501	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
502	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
503	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
504	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
505	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 1
506	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
507	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
508	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
509	0.001_wp, 0.0001_wp, 4.0_wp, & !parameter 18-20
510	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
511	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
512	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
513	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
514	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
515	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
516	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 2
517	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
518	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
519	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
520	0.001_wp, 0.0001_wp, 4.0_wp, & !parameter 18-20
521	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
522	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
523	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
524	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
525	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
526	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
527	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 3
528	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
529	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
530	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
531	0.01_wp, 0.001_wp, 4.0_wp, & !parameter 18-20
532	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
533	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
534	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
535	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
536	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
537	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
538	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 4
539	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
540	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
541	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
542	0.001_wp, 0.0001_wp, 4.0_wp, & !parameter 18-20
543	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
544	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
545	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
546	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
547	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
548	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
549	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 5
550	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
551	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
552	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
553	0.001_wp, 0.0001_wp, 4.0_wp, & !parameter 18-20
554	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
555	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
556	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
557	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
558	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
559	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
560	20000.0_wp, 10.0_wp, & !parameter 45-46 - end of type 6
561	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp, & !parameter 0-5
562	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp, & !parameter 6-11
563	296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp, & !parameter 12-17
564	0.001_wp, 0.0001_wp, 0.0_wp, & !parameter 18-20
565	1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp, & !parameter 21-25
566	1000000.0_wp, 1000000.0_wp, 1000000.0_wp, & !parameter 26-28
567	0.3_wp, 0.3_wp, 0.3_wp, & !parameter 29-31
568	0.9_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp, & !parameter 32-37
569	24.0_wp, 24.0_wp, 24.0_wp, & !parameter 38-40
570	0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp, & !parameter 41-44
571	20000.0_wp, 10.0_wp & !parameter 45-46 - end of type 7 (bridges)
572	/), &
573	(/47, 7/) )
574
575	!
576	!-- Type for surface temperatures at vertical walls. Is not necessary for horizontal walls.
577	TYPE t_surf_vertical
578	REAL(wp), DIMENSION(:), ALLOCATABLE :: t
579	END TYPE t_surf_vertical
580	!
581	!-- Type for wall temperatures at vertical walls. Is not necessary for horizontal walls.
582	TYPE t_wall_vertical
583	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t
584	END TYPE t_wall_vertical
585
586
587	!-- arrays for time averages
588	!-- Attention: the variable rad_net_av is also used in the 3d field variable in radiation_model_mod.f90. It may be better to rename it
589	REAL(wp), DIMENSION(:), ALLOCATABLE :: rad_net_av !< average of rad_net_l
590	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinsw_av !< average of sw radiation falling to local surface including radiation from reflections
591	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlw_av !< average of lw radiation falling to local surface including radiation from reflections
592	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdir_av !< average of direct sw radiation falling to local surface
593	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswdif_av !< average of diffuse sw radiation from sky and model boundary falling to local surface
594	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwdif_av !< average of diffuse lw radiation from sky and model boundary falling to local surface
595	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinswref_av !< average of sw radiation falling to surface from reflections
596	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinlwref_av !< average of lw radiation falling to surface from reflections
597	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutsw_av !< average of total sw radiation outgoing from nonvirtual surfaces surfaces after all reflection
598	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw_av !< average of total lw radiation outgoing from nonvirtual surfaces surfaces after all reflection
599	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfins_av !< average of array of residua of sw radiation absorbed in surface after last reflection
600	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfinl_av !< average of array of residua of lw radiation absorbed in surface after last reflection
601	REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf_av !< average of total radiation flux incoming to minus outgoing from local surface
602	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_av !< average of wghf_eb
603	REAL(wp), DIMENSION(:), ALLOCATABLE :: wshf_eb_av !< average of wshf_eb
604	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_wall_av !< Average of t_wall
605	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_green_av !< average of wghf_eb_green
606	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_green_av !< Average of t_green
607	REAL(wp), DIMENSION(:), ALLOCATABLE :: wghf_eb_window_av !< average of wghf_eb_window
608	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: t_window_av !< Average of t_window
609
610
611	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
612	!-- anthropogenic heat sources
613	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
614	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: aheat !< daily average of anthropogenic heat (W/m2)
615	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: aheatprof !< diurnal profiles of anthropogenic heat for particular layers
616	INTEGER(iwp) :: naheatlayers = 1 !< number of layers of anthropogenic heat
617
618	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
619	!-- wall surface model
620	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
621	!-- wall surface model constants
622	INTEGER(iwp), PARAMETER :: nzb_wall = 0 !< inner side of the wall model (to be switched)
623	INTEGER(iwp), PARAMETER :: nzt_wall = 3 !< outer side of the wall model (to be switched)
624	INTEGER(iwp), PARAMETER :: nzw = 4 !< number of wall layers (fixed for now)
625
626	REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
627	!< normalized soil, wall and roof layer depths (m/m)
628	! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default = (/0.33_wp, 0.66_wp, 1.0_wp /)
629	REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /)
630	! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.33_wp, 0.66_wp, 1.0_wp /)
631	! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_window = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
632	!< normalized window layer depths (m/m)
633	! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
634	!< normalized green layer depths (m/m)
635	REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /)
636	! REAL(wp), DIMENSION(nzb_wall:nzt_wall) :: zwn_default_green = (/0.33_wp, 0.66_wp, 1.0_wp /)
637
638
639	REAL(wp) :: wall_inner_temperature = 295.0_wp !< temperature of the inner wall surface (~22 degrees C) (K)
640	REAL(wp) :: roof_inner_temperature = 295.0_wp !< temperature of the inner roof surface (~22 degrees C) (K)
641	REAL(wp) :: soil_inner_temperature = 288.0_wp !< temperature of the deep soil (~15 degrees C) (K)
642	REAL(wp) :: window_inner_temperature = 295.0_wp !< temperature of the inner window surface (~22 degrees C) (K)
643
644	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
645	!-- surface and material model variables for walls, ground, roofs
646	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
647	REAL(wp), DIMENSION(:), ALLOCATABLE :: zwn !< normalized wall layer depths (m)
648	REAL(wp), DIMENSION(:), ALLOCATABLE :: zwn_window !< normalized window layer depths (m)
649	REAL(wp), DIMENSION(:), ALLOCATABLE :: zwn_green !< normalized green layer depths (m)
650
651	#if defined( __nopointer )
652	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_h !< wall surface temperature (K) at horizontal walls
653	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_h_p !< progn. wall surface temperature (K) at horizontal walls
654	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_window_h !< window surface temperature (K) at horizontal walls
655	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_window_h_p !< progn. window surface temperature (K) at horizontal walls
656	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_h !< green surface temperature (K) at horizontal walls
657	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_h_p !< progn. green surface temperature (K) at horizontal walls
658	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_10cm_h !< near surface temperature (10cm) (K) at horizontal walls
659	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_10cm_h_p !< progn. near surface temperature (10cm) (K) at horizontal walls
660	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_v
661	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_v_p
662	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_window_v
663	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_window_v_p
664	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_green_v
665	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_green_v_p
666	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_10cm_v
667	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_10cm_v_p
668	#else
669	REAL(wp), DIMENSION(:), POINTER :: t_surf_h
670	REAL(wp), DIMENSION(:), POINTER :: t_surf_h_p
671	REAL(wp), DIMENSION(:), POINTER :: t_surf_window_h
672	REAL(wp), DIMENSION(:), POINTER :: t_surf_window_h_p
673	REAL(wp), DIMENSION(:), POINTER :: t_surf_green_h
674	REAL(wp), DIMENSION(:), POINTER :: t_surf_green_h_p
675	REAL(wp), DIMENSION(:), POINTER :: t_surf_10cm_h
676	REAL(wp), DIMENSION(:), POINTER :: t_surf_10cm_h_p
677
678	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_h_1
679	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_h_2
680	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_window_h_1
681	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_window_h_2
682	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_h_1
683	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_h_2
684	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_10cm_h_1
685	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_10cm_h_2
686
687	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_v
688	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_v_p
689	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v
690	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_window_v_p
691	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v
692	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_green_v_p
693	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_10cm_v
694	TYPE(t_surf_vertical), DIMENSION(:), POINTER :: t_surf_10cm_v_p
695
696	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_v_1
697	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_v_2
698	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_window_v_1
699	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_window_v_2
700	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_green_v_1
701	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_green_v_2
702	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_10cm_v_1
703	TYPE(t_surf_vertical), DIMENSION(0:3), TARGET :: t_surf_10cm_v_2
704
705	#endif
706	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_av !< average of wall surface temperature (K)
707	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_window_av !< average of window surface temperature (K)
708	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_green_av !< average of green wall surface temperature (K)
709	REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET :: t_surf_10cm_av !< average of whole wall surface temperature (K)
710
711	!-- Temporal tendencies for time stepping
712	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_m !< surface temperature tendency of wall (K)
713	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_window_m !< surface temperature tendency of window (K)
714	REAL(wp), DIMENSION(:), ALLOCATABLE :: tt_surface_green_m !< surface temperature tendency of green wall (K)
715
716	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
717	!-- Energy balance variables
718	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
719	!-- parameters of the land, roof and wall surfaces
720
721	#if defined( __nopointer )
722	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_wall_h !< Wall temperature (K)
723	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_wall_h_av !< Average of t_wall
724	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_wall_h_p !< Prog. wall temperature (K)
725	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_window_h !< Window temperature (K)
726	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_window_h_av !< Average of t_window
727	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_window_h_p !< Prog. window temperature (K)
728	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_green_h !< Green temperature (K)
729	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_green_h_av !< Average of t_green
730	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_green_h_p !< Prog. green temperature (K)
731
732	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_wall_v !< Wall temperature (K)
733	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_wall_v_av !< Average of t_wall
734	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_wall_v_p !< Prog. wall temperature (K)
735	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_window_v !< Window temperature (K)
736	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_window_v_av !< Average of t_window
737	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_window_v_p !< Prog. window temperature (K)
738	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_green_v !< Green temperature (K)
739	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_green_v_av !< Average of t_green
740	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_green_v_p !< Prog. green temperature (K)
741	#else
742	REAL(wp), DIMENSION(:,:), POINTER :: t_wall_h, t_wall_h_p
743	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_wall_h_av, t_wall_h_1, t_wall_h_2
744	REAL(wp), DIMENSION(:,:), POINTER :: t_window_h, t_window_h_p
745	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_window_h_av, t_window_h_1, t_window_h_2
746	REAL(wp), DIMENSION(:,:), POINTER :: t_green_h, t_green_h_p
747	REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_green_h_av, t_green_h_1, t_green_h_2
748
749	TYPE(t_wall_vertical), DIMENSION(:), POINTER :: t_wall_v, t_wall_v_p
750	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_wall_v_av, t_wall_v_1, t_wall_v_2
751	TYPE(t_wall_vertical), DIMENSION(:), POINTER :: t_window_v, t_window_v_p
752	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_window_v_av, t_window_v_1, t_window_v_2
753	TYPE(t_wall_vertical), DIMENSION(:), POINTER :: t_green_v, t_green_v_p
754	TYPE(t_wall_vertical), DIMENSION(0:3), TARGET :: t_green_v_av, t_green_v_1, t_green_v_2
755	#endif
756
757	!-- Wall temporal tendencies for time stepping
758	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_wall_m !< t_wall prognostic array
759	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_window_m !< t_window prognostic array
760	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_green_m !< t_green prognostic array
761
762	!-- Surface and material parameters classes (surface_type)
763	!-- albedo, emissivity, lambda_surf, roughness, thickness, volumetric heat capacity, thermal conductivity
764	INTEGER(iwp) :: n_surface_types !< number of the wall type categories
765	INTEGER(iwp), PARAMETER :: n_surface_params = 9 !< number of parameters for each type of the wall
766	INTEGER(iwp), PARAMETER :: ialbedo = 1 !< albedo of the surface
767	INTEGER(iwp), PARAMETER :: iemiss = 2 !< emissivity of the surface
768	INTEGER(iwp), PARAMETER :: ilambdas = 3 !< heat conductivity lambda S between surface and material ( W m-2 K-1 )
769	INTEGER(iwp), PARAMETER :: irough = 4 !< roughness length z0 for movements
770	INTEGER(iwp), PARAMETER :: iroughh = 5 !< roughness length z0h for scalars (heat, humidity,...)
771	INTEGER(iwp), PARAMETER :: icsurf = 6 !< Surface skin layer heat capacity (J m-2 K-1 )
772	INTEGER(iwp), PARAMETER :: ithick = 7 !< thickness of the surface (wall, roof, land) ( m )
773	INTEGER(iwp), PARAMETER :: irhoC = 8 !< volumetric heat capacity rho*C of the material ( J m-3 K-1 )
774	INTEGER(iwp), PARAMETER :: ilambdah = 9 !< thermal conductivity lambda H of the wall (W m-1 K-1 )
775	CHARACTER(12), DIMENSION(:), ALLOCATABLE :: surface_type_names !< names of wall types (used only for reports)
776	INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: surface_type_codes !< codes of wall types
777	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: surface_params !< parameters of wall types
778
779
780	!-- interfaces of subroutines accessed from outside of this module
781	INTERFACE usm_boundary_condition
782	MODULE PROCEDURE usm_boundary_condition
783	END INTERFACE usm_boundary_condition
784
785	INTERFACE usm_check_data_output
786	MODULE PROCEDURE usm_check_data_output
787	END INTERFACE usm_check_data_output
788
789	INTERFACE usm_check_parameters
790	MODULE PROCEDURE usm_check_parameters
791	END INTERFACE usm_check_parameters
792
793	INTERFACE usm_data_output_3d
794	MODULE PROCEDURE usm_data_output_3d
795	END INTERFACE usm_data_output_3d
796
797	INTERFACE usm_define_netcdf_grid
798	MODULE PROCEDURE usm_define_netcdf_grid
799	END INTERFACE usm_define_netcdf_grid
800
801	INTERFACE usm_init_urban_surface
802	MODULE PROCEDURE usm_init_urban_surface
803	END INTERFACE usm_init_urban_surface
804
805	INTERFACE usm_material_heat_model
806	MODULE PROCEDURE usm_material_heat_model
807	END INTERFACE usm_material_heat_model
808
809	INTERFACE usm_green_heat_model
810	MODULE PROCEDURE usm_green_heat_model
811	END INTERFACE usm_green_heat_model
812
813	INTERFACE usm_parin
814	MODULE PROCEDURE usm_parin
815	END INTERFACE usm_parin
816
817	INTERFACE usm_temperature_near_surface
818	MODULE PROCEDURE usm_temperature_near_surface
819	END INTERFACE usm_temperature_near_surface
820
821	INTERFACE usm_rrd_local
822	MODULE PROCEDURE usm_rrd_local
823	END INTERFACE usm_rrd_local
824
825	INTERFACE usm_surface_energy_balance
826	MODULE PROCEDURE usm_surface_energy_balance
827	END INTERFACE usm_surface_energy_balance
828
829	INTERFACE usm_swap_timelevel
830	MODULE PROCEDURE usm_swap_timelevel
831	END INTERFACE usm_swap_timelevel
832
833	INTERFACE usm_wrd_local
834	MODULE PROCEDURE usm_wrd_local
835	END INTERFACE usm_wrd_local
836
837	INTERFACE usm_allocate_surface
838	MODULE PROCEDURE usm_allocate_surface
839	END INTERFACE usm_allocate_surface
840
841	INTERFACE usm_average_3d_data
842	MODULE PROCEDURE usm_average_3d_data
843	END INTERFACE usm_average_3d_data
844
845
846	SAVE
847
848	PRIVATE
849
850	!-- Public functions
851	PUBLIC usm_boundary_condition, usm_check_parameters, usm_init_urban_surface,&
852	usm_rrd_local, &
853	usm_surface_energy_balance, usm_material_heat_model, &
854	usm_swap_timelevel, usm_check_data_output, usm_average_3d_data, &
855	usm_data_output_3d, usm_define_netcdf_grid, usm_parin, &
856	usm_wrd_local, usm_allocate_surface
857
858	!-- Public parameters, constants and initial values
859	PUBLIC usm_anthropogenic_heat, usm_material_model, &
860	usm_green_heat_model, usm_temperature_near_surface
861
862
863
864	CONTAINS
865
866	!------------------------------------------------------------------------------!
867	! Description:
868	! ------------
869	!> This subroutine creates the necessary indices of the urban surfaces
870	!> and plant canopy and it allocates the needed arrays for USM
871	!------------------------------------------------------------------------------!
872	SUBROUTINE usm_allocate_surface
873
874	IMPLICIT NONE
875
876	INTEGER(iwp) :: l
877
878	!
879	!-- Allocate radiation arrays which are part of the new data type.
880	!-- For horizontal surfaces.
881	ALLOCATE( surf_usm_h%surfhf(1:surf_usm_h%ns) )
882	ALLOCATE( surf_usm_h%rad_net_l(1:surf_usm_h%ns) )
883	!
884	!-- For vertical surfaces
885	DO l = 0, 3
886	ALLOCATE( surf_usm_v(l)%surfhf(1:surf_usm_v(l)%ns) )
887	ALLOCATE( surf_usm_v(l)%rad_net_l(1:surf_usm_v(l)%ns) )
888	ENDDO
889
890	!-- Wall surface model
891	!-- allocate arrays for wall surface model and define pointers
892
893	!-- allocate array of wall types and wall parameters
894	ALLOCATE ( surf_usm_h%surface_types(1:surf_usm_h%ns) )
895	ALLOCATE ( surf_usm_h%building_type(1:surf_usm_h%ns) )
896	ALLOCATE ( surf_usm_h%building_type_name(1:surf_usm_h%ns) )
897	surf_usm_h%building_type = 0
898	surf_usm_h%building_type_name = 'none'
899	DO l = 0, 3
900	ALLOCATE( surf_usm_v(l)%surface_types(1:surf_usm_v(l)%ns) )
901	ALLOCATE ( surf_usm_v(l)%building_type(1:surf_usm_v(l)%ns) )
902	ALLOCATE ( surf_usm_v(l)%building_type_name(1:surf_usm_v(l)%ns) )
903	surf_usm_v(l)%building_type = 0
904	surf_usm_v(l)%building_type_name = 'none'
905	ENDDO
906	!
907	!-- Allocate albedo_type and albedo. Each surface element
908	!-- has 3 values, 0: wall fraction, 1: green fraction, 2: window fraction.
909	ALLOCATE( surf_usm_h%albedo_type(0:2,1:surf_usm_h%ns) )
910	ALLOCATE( surf_usm_h%albedo(0:2,1:surf_usm_h%ns) )
911	surf_usm_h%albedo_type = albedo_type
912	DO l = 0, 3
913	ALLOCATE( surf_usm_v(l)%albedo_type(0:2,1:surf_usm_v(l)%ns) )
914	ALLOCATE( surf_usm_v(l)%albedo(0:2,1:surf_usm_v(l)%ns) )
915	surf_usm_v(l)%albedo_type = albedo_type
916	ENDDO
917
918
919	!
920	!-- Allocate indoor target temperature for summer and winter
921	ALLOCATE( surf_usm_h%target_temp_summer(1:surf_usm_h%ns) )
922	ALLOCATE( surf_usm_h%target_temp_winter(1:surf_usm_h%ns) )
923	DO l = 0, 3
924	ALLOCATE( surf_usm_v(l)%target_temp_summer(1:surf_usm_v(l)%ns) )
925	ALLOCATE( surf_usm_v(l)%target_temp_winter(1:surf_usm_v(l)%ns) )
926	ENDDO
927	!
928	!-- Allocate flag indicating ground floor level surface elements
929	ALLOCATE ( surf_usm_h%ground_level(1:surf_usm_h%ns) )
930	DO l = 0, 3
931	ALLOCATE( surf_usm_v(l)%ground_level(1:surf_usm_v(l)%ns) )
932	ENDDO
933	!
934	!-- Allocate arrays for relative surface fraction.
935	!-- 0 - wall fraction, 1 - green fraction, 2 - window fraction
936	ALLOCATE( surf_usm_h%frac(0:2,1:surf_usm_h%ns) )
937	surf_usm_h%frac = 0.0_wp
938	DO l = 0, 3
939	ALLOCATE( surf_usm_v(l)%frac(0:2,1:surf_usm_v(l)%ns) )
940	surf_usm_v(l)%frac = 0.0_wp
941	ENDDO
942
943	!-- wall and roof surface parameters. First for horizontal surfaces
944	ALLOCATE ( surf_usm_h%isroof_surf(1:surf_usm_h%ns) )
945	ALLOCATE ( surf_usm_h%lambda_surf(1:surf_usm_h%ns) )
946	ALLOCATE ( surf_usm_h%lambda_surf_window(1:surf_usm_h%ns) )
947	ALLOCATE ( surf_usm_h%lambda_surf_green(1:surf_usm_h%ns) )
948	ALLOCATE ( surf_usm_h%c_surface(1:surf_usm_h%ns) )
949	ALLOCATE ( surf_usm_h%c_surface_window(1:surf_usm_h%ns) )
950	ALLOCATE ( surf_usm_h%c_surface_green(1:surf_usm_h%ns) )
951	ALLOCATE ( surf_usm_h%transmissivity(1:surf_usm_h%ns) )
952	ALLOCATE ( surf_usm_h%lai(1:surf_usm_h%ns) )
953	ALLOCATE ( surf_usm_h%emissivity(0:2,1:surf_usm_h%ns) )
954	ALLOCATE ( surf_usm_h%r_a(1:surf_usm_h%ns) )
955	ALLOCATE ( surf_usm_h%r_a_green(1:surf_usm_h%ns) )
956	ALLOCATE ( surf_usm_h%r_a_window(1:surf_usm_h%ns) )
957
958	!
959	!-- For vertical surfaces.
960	DO l = 0, 3
961	ALLOCATE ( surf_usm_v(l)%lambda_surf(1:surf_usm_v(l)%ns) )
962	ALLOCATE ( surf_usm_v(l)%c_surface(1:surf_usm_v(l)%ns) )
963	ALLOCATE ( surf_usm_v(l)%lambda_surf_window(1:surf_usm_v(l)%ns) )
964	ALLOCATE ( surf_usm_v(l)%c_surface_window(1:surf_usm_v(l)%ns) )
965	ALLOCATE ( surf_usm_v(l)%lambda_surf_green(1:surf_usm_v(l)%ns) )
966	ALLOCATE ( surf_usm_v(l)%c_surface_green(1:surf_usm_v(l)%ns) )
967	ALLOCATE ( surf_usm_v(l)%transmissivity(1:surf_usm_v(l)%ns) )
968	ALLOCATE ( surf_usm_v(l)%lai(1:surf_usm_v(l)%ns) )
969	ALLOCATE ( surf_usm_v(l)%emissivity(0:2,1:surf_usm_v(l)%ns) )
970	ALLOCATE ( surf_usm_v(l)%r_a(1:surf_usm_v(l)%ns) )
971	ALLOCATE ( surf_usm_v(l)%r_a_green(1:surf_usm_v(l)%ns) )
972	ALLOCATE ( surf_usm_v(l)%r_a_window(1:surf_usm_v(l)%ns) )
973	ENDDO
974
975	!
976	!-- allocate wall and roof material parameters. First for horizontal surfaces
977	ALLOCATE ( surf_usm_h%thickness_wall(1:surf_usm_h%ns) )
978	ALLOCATE ( surf_usm_h%thickness_window(1:surf_usm_h%ns) )
979	ALLOCATE ( surf_usm_h%thickness_green(1:surf_usm_h%ns) )
980	ALLOCATE ( surf_usm_h%lambda_h(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
981	ALLOCATE ( surf_usm_h%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
982	ALLOCATE ( surf_usm_h%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
983	ALLOCATE ( surf_usm_h%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
984	ALLOCATE ( surf_usm_h%lambda_h_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
985	ALLOCATE ( surf_usm_h%rho_c_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
986
987	!
988	!-- For vertical surfaces.
989	DO l = 0, 3
990	ALLOCATE ( surf_usm_v(l)%thickness_wall(1:surf_usm_v(l)%ns) )
991	ALLOCATE ( surf_usm_v(l)%thickness_window(1:surf_usm_v(l)%ns) )
992	ALLOCATE ( surf_usm_v(l)%thickness_green(1:surf_usm_v(l)%ns) )
993	ALLOCATE ( surf_usm_v(l)%lambda_h(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
994	ALLOCATE ( surf_usm_v(l)%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
995	ALLOCATE ( surf_usm_v(l)%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
996	ALLOCATE ( surf_usm_v(l)%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
997	ALLOCATE ( surf_usm_v(l)%lambda_h_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
998	ALLOCATE ( surf_usm_v(l)%rho_c_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
999	ENDDO
1000
1001	!-- allocate wall and roof layers sizes. For horizontal surfaces.
1002	ALLOCATE ( zwn(nzb_wall:nzt_wall) )
1003	ALLOCATE ( surf_usm_h%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1004	ALLOCATE ( zwn_window(nzb_wall:nzt_wall) )
1005	ALLOCATE ( surf_usm_h%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1006	ALLOCATE ( zwn_green(nzb_wall:nzt_wall) )
1007	ALLOCATE ( surf_usm_h%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1008	ALLOCATE ( surf_usm_h%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1009	ALLOCATE ( surf_usm_h%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1010	ALLOCATE ( surf_usm_h%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1011	ALLOCATE ( surf_usm_h%zw(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1012	ALLOCATE ( surf_usm_h%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1013	ALLOCATE ( surf_usm_h%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1014	ALLOCATE ( surf_usm_h%ddz_window_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1015	ALLOCATE ( surf_usm_h%zw_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1016	ALLOCATE ( surf_usm_h%ddz_green(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1017	ALLOCATE ( surf_usm_h%dz_green_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1018	ALLOCATE ( surf_usm_h%ddz_green_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1019	ALLOCATE ( surf_usm_h%zw_green(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1020	!
1021	!-- For vertical surfaces.
1022	DO l = 0, 3
1023	ALLOCATE ( surf_usm_v(l)%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1024	ALLOCATE ( surf_usm_v(l)%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1025	ALLOCATE ( surf_usm_v(l)%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1026	ALLOCATE ( surf_usm_v(l)%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1027	ALLOCATE ( surf_usm_v(l)%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1028	ALLOCATE ( surf_usm_v(l)%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1029	ALLOCATE ( surf_usm_v(l)%zw(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1030	ALLOCATE ( surf_usm_v(l)%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1031	ALLOCATE ( surf_usm_v(l)%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1032	ALLOCATE ( surf_usm_v(l)%ddz_window_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1033	ALLOCATE ( surf_usm_v(l)%zw_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1034	ALLOCATE ( surf_usm_v(l)%ddz_green(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1035	ALLOCATE ( surf_usm_v(l)%dz_green_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1036	ALLOCATE ( surf_usm_v(l)%ddz_green_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1037	ALLOCATE ( surf_usm_v(l)%zw_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1038	ENDDO
1039
1040	!-- allocate wall and roof temperature arrays, for horizontal walls
1041	#if defined( __nopointer )
1042	IF ( .NOT. ALLOCATED( t_surf_h ) ) &
1043	ALLOCATE ( t_surf_h(1:surf_usm_h%ns) )
1044	IF ( .NOT. ALLOCATED( t_surf_h_p ) ) &
1045	ALLOCATE ( t_surf_h_p(1:surf_usm_h%ns) )
1046	IF ( .NOT. ALLOCATED( t_wall_h ) ) &
1047	ALLOCATE ( t_wall_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1048	IF ( .NOT. ALLOCATED( t_wall_h_p ) ) &
1049	ALLOCATE ( t_wall_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1050	IF ( .NOT. ALLOCATED( t_surf_window_h ) ) &
1051	ALLOCATE ( t_surf_window_h(1:surf_usm_h%ns) )
1052	IF ( .NOT. ALLOCATED( t_surf_window_h_p ) ) &
1053	ALLOCATE ( t_surf_window_h_p(1:surf_usm_h%ns) )
1054	IF ( .NOT. ALLOCATED( t_window_h ) ) &
1055	ALLOCATE ( t_window_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1056	IF ( .NOT. ALLOCATED( t_window_h_p ) ) &
1057	ALLOCATE ( t_window_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1058	IF ( .NOT. ALLOCATED( t_surf_green_h ) ) &
1059	ALLOCATE ( t_surf_green_h(1:surf_usm_h%ns) )
1060	IF ( .NOT. ALLOCATED( t_surf_green_h_p ) ) &
1061	ALLOCATE ( t_surf_green_h_p(1:surf_usm_h%ns) )
1062	IF ( .NOT. ALLOCATED( t_green_h ) ) &
1063	ALLOCATE ( t_green_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1064	IF ( .NOT. ALLOCATED( t_green_h_p ) ) &
1065	ALLOCATE ( t_green_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1066	IF ( .NOT. ALLOCATED( t_surf_10cm_h ) ) &
1067	ALLOCATE ( t_surf_10cm_h(1:surf_usm_h%ns) )
1068	IF ( .NOT. ALLOCATED( t_surf_10cm_h_p ) ) &
1069	ALLOCATE ( t_surf_10cm_h_p(1:surf_usm_h%ns) )
1070	#else
1071	!
1072	!-- Allocate if required. Note, in case of restarts, some of these arrays
1073	!-- might be already allocated.
1074	IF ( .NOT. ALLOCATED( t_surf_h_1 ) ) &
1075	ALLOCATE ( t_surf_h_1(1:surf_usm_h%ns) )
1076	IF ( .NOT. ALLOCATED( t_surf_h_2 ) ) &
1077	ALLOCATE ( t_surf_h_2(1:surf_usm_h%ns) )
1078	IF ( .NOT. ALLOCATED( t_wall_h_1 ) ) &
1079	ALLOCATE ( t_wall_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1080	IF ( .NOT. ALLOCATED( t_wall_h_2 ) ) &
1081	ALLOCATE ( t_wall_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1082	IF ( .NOT. ALLOCATED( t_surf_window_h_1 ) ) &
1083	ALLOCATE ( t_surf_window_h_1(1:surf_usm_h%ns) )
1084	IF ( .NOT. ALLOCATED( t_surf_window_h_2 ) ) &
1085	ALLOCATE ( t_surf_window_h_2(1:surf_usm_h%ns) )
1086	IF ( .NOT. ALLOCATED( t_window_h_1 ) ) &
1087	ALLOCATE ( t_window_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1088	IF ( .NOT. ALLOCATED( t_window_h_2 ) ) &
1089	ALLOCATE ( t_window_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1090	IF ( .NOT. ALLOCATED( t_surf_green_h_1 ) ) &
1091	ALLOCATE ( t_surf_green_h_1(1:surf_usm_h%ns) )
1092	IF ( .NOT. ALLOCATED( t_surf_green_h_2 ) ) &
1093	ALLOCATE ( t_surf_green_h_2(1:surf_usm_h%ns) )
1094	IF ( .NOT. ALLOCATED( t_green_h_1 ) ) &
1095	ALLOCATE ( t_green_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1096	IF ( .NOT. ALLOCATED( t_green_h_2 ) ) &
1097	ALLOCATE ( t_green_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1098	IF ( .NOT. ALLOCATED( t_surf_10cm_h_1 ) ) &
1099	ALLOCATE ( t_surf_10cm_h_1(1:surf_usm_h%ns) )
1100	IF ( .NOT. ALLOCATED( t_surf_10cm_h_2 ) ) &
1101	ALLOCATE ( t_surf_10cm_h_2(1:surf_usm_h%ns) )
1102	!
1103	!-- initial assignment of the pointers
1104	t_wall_h => t_wall_h_1; t_wall_h_p => t_wall_h_2
1105	t_window_h => t_window_h_1; t_window_h_p => t_window_h_2
1106	t_green_h => t_green_h_1; t_green_h_p => t_green_h_2
1107	t_surf_h => t_surf_h_1; t_surf_h_p => t_surf_h_2
1108	t_surf_window_h => t_surf_window_h_1; t_surf_window_h_p => t_surf_window_h_2
1109	t_surf_green_h => t_surf_green_h_1; t_surf_green_h_p => t_surf_green_h_2
1110	t_surf_10cm_h => t_surf_10cm_h_1; t_surf_10cm_h_p => t_surf_10cm_h_2
1111
1112	#endif
1113
1114	!-- allocate wall and roof temperature arrays, for vertical walls if required
1115	#if defined( __nopointer )
1116	DO l = 0, 3
1117	IF ( .NOT. ALLOCATED( t_surf_v(l)%t ) ) &
1118	ALLOCATE ( t_surf_v(l)%t(1:surf_usm_v(l)%ns) )
1119	IF ( .NOT. ALLOCATED( t_surf_v_p(l)%t ) ) &
1120	ALLOCATE ( t_surf_v_p(l)%t(1:surf_usm_v(l)%ns) )
1121	IF ( .NOT. ALLOCATED( t_wall_v(l)%t ) ) &
1122	ALLOCATE ( t_wall_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1123	IF ( .NOT. ALLOCATED( t_wall_v_p(l)%t ) ) &
1124	ALLOCATE ( t_wall_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1125	IF ( .NOT. ALLOCATED( t_surf_window_v(l)%t ) ) &
1126	ALLOCATE ( t_surf_window_v(l)%t(1:surf_usm_v(l)%ns) )
1127	IF ( .NOT. ALLOCATED( t_surf_window_v_p(l)%t ) ) &
1128	ALLOCATE ( t_surf_window_v_p(l)%t(1:surf_usm_v(l)%ns) )
1129	IF ( .NOT. ALLOCATED( t_window_v(l)%t ) ) &
1130	ALLOCATE ( t_window_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1131	IF ( .NOT. ALLOCATED( t_window_v_p(l)%t ) ) &
1132	ALLOCATE ( t_window_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1133	IF ( .NOT. ALLOCATED( t_green_v(l)%t ) ) &
1134	ALLOCATE ( t_green_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1135	IF ( .NOT. ALLOCATED( t_green_v_p(l)%t ) ) &
1136	ALLOCATE ( t_green_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1137	IF ( .NOT. ALLOCATED( t_surf_green_v(l)%t ) ) &
1138	ALLOCATE ( t_surf_green_v(l)%t(1:surf_usm_v(l)%ns) )
1139	IF ( .NOT. ALLOCATED( t_surf_green_v_p(l)%t ) ) &
1140	ALLOCATE ( t_surf_green_v_p(l)%t(1:surf_usm_v(l)%ns) )
1141	IF ( .NOT. ALLOCATED( t_surf_10cm_v(l)%t ) ) &
1142	ALLOCATE ( t_surf_10cm_v(l)%t(1:surf_usm_v(l)%ns) )
1143	IF ( .NOT. ALLOCATED( t_surf_10cm_v_p(l)%t ) ) &
1144	ALLOCATE ( t_surf_10cm_v_p(l)%t(1:surf_usm_v(l)%ns) )
1145	ENDDO
1146	#else
1147	!
1148	!-- Allocate if required. Note, in case of restarts, some of these arrays
1149	!-- might be already allocated.
1150	DO l = 0, 3
1151	IF ( .NOT. ALLOCATED( t_surf_v_1(l)%t ) ) &
1152	ALLOCATE ( t_surf_v_1(l)%t(1:surf_usm_v(l)%ns) )
1153	IF ( .NOT. ALLOCATED( t_surf_v_2(l)%t ) ) &
1154	ALLOCATE ( t_surf_v_2(l)%t(1:surf_usm_v(l)%ns) )
1155	IF ( .NOT. ALLOCATED( t_wall_v_1(l)%t ) ) &
1156	ALLOCATE ( t_wall_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1157	IF ( .NOT. ALLOCATED( t_wall_v_2(l)%t ) ) &
1158	ALLOCATE ( t_wall_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1159	IF ( .NOT. ALLOCATED( t_surf_window_v_1(l)%t ) ) &
1160	ALLOCATE ( t_surf_window_v_1(l)%t(1:surf_usm_v(l)%ns) )
1161	IF ( .NOT. ALLOCATED( t_surf_window_v_2(l)%t ) ) &
1162	ALLOCATE ( t_surf_window_v_2(l)%t(1:surf_usm_v(l)%ns) )
1163	IF ( .NOT. ALLOCATED( t_window_v_1(l)%t ) ) &
1164	ALLOCATE ( t_window_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1165	IF ( .NOT. ALLOCATED( t_window_v_2(l)%t ) ) &
1166	ALLOCATE ( t_window_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1167	IF ( .NOT. ALLOCATED( t_surf_green_v_1(l)%t ) ) &
1168	ALLOCATE ( t_surf_green_v_1(l)%t(1:surf_usm_v(l)%ns) )
1169	IF ( .NOT. ALLOCATED( t_surf_green_v_2(l)%t ) ) &
1170	ALLOCATE ( t_surf_green_v_2(l)%t(1:surf_usm_v(l)%ns) )
1171	IF ( .NOT. ALLOCATED( t_green_v_1(l)%t ) ) &
1172	ALLOCATE ( t_green_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1173	IF ( .NOT. ALLOCATED( t_green_v_2(l)%t ) ) &
1174	ALLOCATE ( t_green_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1175	IF ( .NOT. ALLOCATED( t_surf_10cm_v_1(l)%t ) ) &
1176	ALLOCATE ( t_surf_10cm_v_1(l)%t(1:surf_usm_v(l)%ns) )
1177	IF ( .NOT. ALLOCATED( t_surf_10cm_v_2(l)%t ) ) &
1178	ALLOCATE ( t_surf_10cm_v_2(l)%t(1:surf_usm_v(l)%ns) )
1179	ENDDO
1180	!
1181	!-- initial assignment of the pointers
1182	t_wall_v => t_wall_v_1; t_wall_v_p => t_wall_v_2
1183	t_surf_v => t_surf_v_1; t_surf_v_p => t_surf_v_2
1184	t_window_v => t_window_v_1; t_window_v_p => t_window_v_2
1185	t_green_v => t_green_v_1; t_green_v_p => t_green_v_2
1186	t_surf_window_v => t_surf_window_v_1; t_surf_window_v_p => t_surf_window_v_2
1187	t_surf_green_v => t_surf_green_v_1; t_surf_green_v_p => t_surf_green_v_2
1188	t_surf_10cm_v => t_surf_10cm_v_1; t_surf_10cm_v_p => t_surf_10cm_v_2
1189
1190	#endif
1191	!
1192	!-- Allocate intermediate timestep arrays. For horizontal surfaces.
1193	ALLOCATE ( surf_usm_h%tt_surface_m(1:surf_usm_h%ns) )
1194	ALLOCATE ( surf_usm_h%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1195	ALLOCATE ( surf_usm_h%tt_surface_window_m(1:surf_usm_h%ns) )
1196	ALLOCATE ( surf_usm_h%tt_window_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1197	ALLOCATE ( surf_usm_h%tt_green_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1198	ALLOCATE ( surf_usm_h%tt_surface_green_m(1:surf_usm_h%ns) )
1199
1200	!
1201	!-- Set inital values for prognostic quantities
1202	IF ( ALLOCATED( surf_usm_h%tt_surface_m ) ) surf_usm_h%tt_surface_m = 0.0_wp
1203	IF ( ALLOCATED( surf_usm_h%tt_wall_m ) ) surf_usm_h%tt_wall_m = 0.0_wp
1204	IF ( ALLOCATED( surf_usm_h%tt_surface_window_m ) ) surf_usm_h%tt_surface_window_m = 0.0_wp
1205	IF ( ALLOCATED( surf_usm_h%tt_window_m ) ) surf_usm_h%tt_window_m = 0.0_wp
1206	IF ( ALLOCATED( surf_usm_h%tt_green_m ) ) surf_usm_h%tt_green_m = 0.0_wp
1207	IF ( ALLOCATED( surf_usm_h%tt_surface_green_m ) ) surf_usm_h%tt_surface_green_m = 0.0_wp
1208	!
1209	!-- Now, for vertical surfaces
1210	DO l = 0, 3
1211	ALLOCATE ( surf_usm_v(l)%tt_surface_m(1:surf_usm_v(l)%ns) )
1212	ALLOCATE ( surf_usm_v(l)%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1213	IF ( ALLOCATED( surf_usm_v(l)%tt_surface_m ) ) surf_usm_v(l)%tt_surface_m = 0.0_wp
1214	IF ( ALLOCATED( surf_usm_v(l)%tt_wall_m ) ) surf_usm_v(l)%tt_wall_m = 0.0_wp
1215	ALLOCATE ( surf_usm_v(l)%tt_surface_window_m(1:surf_usm_v(l)%ns) )
1216	ALLOCATE ( surf_usm_v(l)%tt_window_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1217	IF ( ALLOCATED( surf_usm_v(l)%tt_surface_window_m ) ) surf_usm_v(l)%tt_surface_window_m = 0.0_wp
1218	IF ( ALLOCATED( surf_usm_v(l)%tt_window_m ) ) surf_usm_v(l)%tt_window_m = 0.0_wp
1219	ALLOCATE ( surf_usm_v(l)%tt_surface_green_m(1:surf_usm_v(l)%ns) )
1220	IF ( ALLOCATED( surf_usm_v(l)%tt_surface_green_m ) ) surf_usm_v(l)%tt_surface_green_m = 0.0_wp
1221	ALLOCATE ( surf_usm_v(l)%tt_green_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1222	IF ( ALLOCATED( surf_usm_v(l)%tt_green_m ) ) surf_usm_v(l)%tt_green_m = 0.0_wp
1223	ENDDO
1224
1225	!-- allocate wall heat flux output array and set initial values. For horizontal surfaces
1226	! ALLOCATE ( surf_usm_h%wshf(1:surf_usm_h%ns) ) !can be removed
1227	ALLOCATE ( surf_usm_h%wshf_eb(1:surf_usm_h%ns) )
1228	ALLOCATE ( surf_usm_h%wghf_eb(1:surf_usm_h%ns) )
1229	ALLOCATE ( surf_usm_h%wghf_eb_window(1:surf_usm_h%ns) )
1230	ALLOCATE ( surf_usm_h%wghf_eb_green(1:surf_usm_h%ns) )
1231	ALLOCATE ( surf_usm_h%iwghf_eb(1:surf_usm_h%ns) )
1232	ALLOCATE ( surf_usm_h%iwghf_eb_window(1:surf_usm_h%ns) )
1233	IF ( ALLOCATED( surf_usm_h%wshf ) ) surf_usm_h%wshf = 0.0_wp
1234	IF ( ALLOCATED( surf_usm_h%wshf_eb ) ) surf_usm_h%wshf_eb = 0.0_wp
1235	IF ( ALLOCATED( surf_usm_h%wghf_eb ) ) surf_usm_h%wghf_eb = 0.0_wp
1236	IF ( ALLOCATED( surf_usm_h%wghf_eb_window ) ) surf_usm_h%wghf_eb_window = 0.0_wp
1237	IF ( ALLOCATED( surf_usm_h%wghf_eb_green ) ) surf_usm_h%wghf_eb_green = 0.0_wp
1238	IF ( ALLOCATED( surf_usm_h%iwghf_eb ) ) surf_usm_h%iwghf_eb = 0.0_wp
1239	IF ( ALLOCATED( surf_usm_h%iwghf_eb_window ) ) surf_usm_h%iwghf_eb_window = 0.0_wp
1240	!
1241	!-- Now, for vertical surfaces
1242	DO l = 0, 3
1243	! ALLOCATE ( surf_usm_v(l)%wshf(1:surf_usm_v(l)%ns) ) ! can be removed
1244	ALLOCATE ( surf_usm_v(l)%wshf_eb(1:surf_usm_v(l)%ns) )
1245	ALLOCATE ( surf_usm_v(l)%wghf_eb(1:surf_usm_v(l)%ns) )
1246	ALLOCATE ( surf_usm_v(l)%wghf_eb_window(1:surf_usm_v(l)%ns) )
1247	ALLOCATE ( surf_usm_v(l)%wghf_eb_green(1:surf_usm_v(l)%ns) )
1248	ALLOCATE ( surf_usm_v(l)%iwghf_eb(1:surf_usm_v(l)%ns) )
1249	ALLOCATE ( surf_usm_v(l)%iwghf_eb_window(1:surf_usm_v(l)%ns) )
1250	IF ( ALLOCATED( surf_usm_v(l)%wshf ) ) surf_usm_v(l)%wshf = 0.0_wp
1251	IF ( ALLOCATED( surf_usm_v(l)%wshf_eb ) ) surf_usm_v(l)%wshf_eb = 0.0_wp
1252	IF ( ALLOCATED( surf_usm_v(l)%wghf_eb ) ) surf_usm_v(l)%wghf_eb = 0.0_wp
1253	IF ( ALLOCATED( surf_usm_v(l)%wghf_eb_window ) ) surf_usm_v(l)%wghf_eb_window = 0.0_wp
1254	IF ( ALLOCATED( surf_usm_v(l)%wghf_eb_green ) ) surf_usm_v(l)%wghf_eb_green = 0.0_wp
1255	IF ( ALLOCATED( surf_usm_v(l)%iwghf_eb ) ) surf_usm_v(l)%iwghf_eb = 0.0_wp
1256	IF ( ALLOCATED( surf_usm_v(l)%iwghf_eb_window ) ) surf_usm_v(l)%iwghf_eb_window = 0.0_wp
1257	ENDDO
1258
1259	END SUBROUTINE usm_allocate_surface
1260
1261
1262	!------------------------------------------------------------------------------!
1263	! Description:
1264	! ------------
1265	!> Sum up and time-average urban surface output quantities as well as allocate
1266	!> the array necessary for storing the average.
1267	!------------------------------------------------------------------------------!
1268	SUBROUTINE usm_average_3d_data( mode, variable )
1269
1270	IMPLICIT NONE
1271
1272	CHARACTER (len=*), INTENT(IN) :: mode
1273	CHARACTER (len=*), INTENT(IN) :: variable
1274
1275	INTEGER(iwp) :: i, j, k, l, m, ids, idsint, iwl, istat
1276	CHARACTER (len=varnamelength) :: var, surfid
1277	INTEGER(iwp), PARAMETER :: nd = 5
1278	CHARACTER(len=6), DIMENSION(0:nd-1), PARAMETER :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /)
1279	INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER :: dirint = (/ iup_u, isouth_u, inorth_u, iwest_u, ieast_u /)
1280
1281	!-- find the real name of the variable
1282	ids = -1
1283	var = TRIM(variable)
1284	DO i = 0, nd-1
1285	k = len(TRIM(var))
1286	j = len(TRIM(dirname(i)))
1287	IF ( var(k-j+1:k) == dirname(i) ) THEN
1288	ids = i
1289	idsint = dirint(ids)
1290	var = var(:k-j)
1291	EXIT
1292	ENDIF
1293	ENDDO
1294	IF ( ids == -1 ) THEN
1295	var = TRIM(variable)
1296	ENDIF
1297	IF ( var(1:11) == 'usm_t_wall_' .AND. len(TRIM(var)) >= 12 ) THEN
1298	!-- wall layers
1299	READ(var(12:12), '(I1)', iostat=istat ) iwl
1300	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
1301	var = var(1:10)
1302	ELSE
1303	!-- wrong wall layer index
1304	RETURN
1305	ENDIF
1306	ENDIF
1307	IF ( var(1:13) == 'usm_t_window_' .AND. len(TRIM(var)) >= 14 ) THEN
1308	!-- wall layers
1309	READ(var(14:14), '(I1)', iostat=istat ) iwl
1310	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
1311	var = var(1:12)
1312	ELSE
1313	!-- wrong window layer index
1314	RETURN
1315	ENDIF
1316	ENDIF
1317	IF ( var(1:12) == 'usm_t_green_' .AND. len(TRIM(var)) >= 13 ) THEN
1318	!-- wall layers
1319	READ(var(13:13), '(I1)', iostat=istat ) iwl
1320	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
1321	var = var(1:11)
1322	ELSE
1323	!-- wrong green layer index
1324	RETURN
1325	ENDIF
1326	ENDIF
1327
1328	IF ( mode == 'allocate' ) THEN
1329
1330	SELECT CASE ( TRIM( var ) )
1331
1332	CASE ( 'usm_rad_net' )
1333	!-- array of complete radiation balance
1334	IF ( .NOT. ALLOCATED(surf_usm_h%rad_net_av) ) THEN
1335	ALLOCATE( surf_usm_h%rad_net_av(1:surf_usm_h%ns) )
1336	surf_usm_h%rad_net_av = 0.0_wp
1337	ENDIF
1338	DO l = 0, 3
1339	IF ( .NOT. ALLOCATED(surf_usm_v(l)%rad_net_av) ) THEN
1340	ALLOCATE( surf_usm_v(l)%rad_net_av(1:surf_usm_v(l)%ns) )
1341	surf_usm_v(l)%rad_net_av = 0.0_wp
1342	ENDIF
1343	ENDDO
1344
1345	CASE ( 'usm_rad_insw' )
1346	!-- array of sw radiation falling to surface after i-th reflection
1347	IF ( .NOT. ALLOCATED(surfinsw_av) ) THEN
1348	ALLOCATE( surfinsw_av(nsurfl) )
1349	surfinsw_av = 0.0_wp
1350	ENDIF
1351
1352	CASE ( 'usm_rad_inlw' )
1353	!-- array of lw radiation falling to surface after i-th reflection
1354	IF ( .NOT. ALLOCATED(surfinlw_av) ) THEN
1355	ALLOCATE( surfinlw_av(nsurfl) )
1356	surfinlw_av = 0.0_wp
1357	ENDIF
1358
1359	CASE ( 'usm_rad_inswdir' )
1360	!-- array of direct sw radiation falling to surface from sun
1361	IF ( .NOT. ALLOCATED(surfinswdir_av) ) THEN
1362	ALLOCATE( surfinswdir_av(nsurfl) )
1363	surfinswdir_av = 0.0_wp
1364	ENDIF
1365
1366	CASE ( 'usm_rad_inswdif' )
1367	!-- array of difusion sw radiation falling to surface from sky and borders of the domain
1368	IF ( .NOT. ALLOCATED(surfinswdif_av) ) THEN
1369	ALLOCATE( surfinswdif_av(nsurfl) )
1370	surfinswdif_av = 0.0_wp
1371	ENDIF
1372
1373	CASE ( 'usm_rad_inswref' )
1374	!-- array of sw radiation falling to surface from reflections
1375	IF ( .NOT. ALLOCATED(surfinswref_av) ) THEN
1376	ALLOCATE( surfinswref_av(nsurfl) )
1377	surfinswref_av = 0.0_wp
1378	ENDIF
1379
1380	CASE ( 'usm_rad_inlwdif' )
1381	!-- array of sw radiation falling to surface after i-th reflection
1382	IF ( .NOT. ALLOCATED(surfinlwdif_av) ) THEN
1383	ALLOCATE( surfinlwdif_av(nsurfl) )
1384	surfinlwdif_av = 0.0_wp
1385	ENDIF
1386
1387	CASE ( 'usm_rad_inlwref' )
1388	!-- array of lw radiation falling to surface from reflections
1389	IF ( .NOT. ALLOCATED(surfinlwref_av) ) THEN
1390	ALLOCATE( surfinlwref_av(nsurfl) )
1391	surfinlwref_av = 0.0_wp
1392	ENDIF
1393
1394	CASE ( 'usm_rad_outsw' )
1395	!-- array of sw radiation emitted from surface after i-th reflection
1396	IF ( .NOT. ALLOCATED(surfoutsw_av) ) THEN
1397	ALLOCATE( surfoutsw_av(nsurfl) )
1398	surfoutsw_av = 0.0_wp
1399	ENDIF
1400
1401	CASE ( 'usm_rad_outlw' )
1402	!-- array of lw radiation emitted from surface after i-th reflection
1403	IF ( .NOT. ALLOCATED(surfoutlw_av) ) THEN
1404	ALLOCATE( surfoutlw_av(nsurfl) )
1405	surfoutlw_av = 0.0_wp
1406	ENDIF
1407	CASE ( 'usm_rad_ressw' )
1408	!-- array of residua of sw radiation absorbed in surface after last reflection
1409	IF ( .NOT. ALLOCATED(surfins_av) ) THEN
1410	ALLOCATE( surfins_av(nsurfl) )
1411	surfins_av = 0.0_wp
1412	ENDIF
1413
1414	CASE ( 'usm_rad_reslw' )
1415	!-- array of residua of lw radiation absorbed in surface after last reflection
1416	IF ( .NOT. ALLOCATED(surfinl_av) ) THEN
1417	ALLOCATE( surfinl_av(nsurfl) )
1418	surfinl_av = 0.0_wp
1419	ENDIF
1420
1421	CASE ( 'usm_rad_hf' )
1422	!-- array of heat flux from radiation for surfaces after i-th reflection
1423	IF ( .NOT. ALLOCATED(surf_usm_h%surfhf_av) ) THEN
1424	ALLOCATE( surf_usm_h%surfhf_av(1:surf_usm_h%ns) )
1425	surf_usm_h%surfhf_av = 0.0_wp
1426	ENDIF
1427	DO l = 0, 3
1428	IF ( .NOT. ALLOCATED(surf_usm_v(l)%surfhf_av) ) THEN
1429	ALLOCATE( surf_usm_v(l)%surfhf_av(1:surf_usm_v(l)%ns) )
1430	surf_usm_v(l)%surfhf_av = 0.0_wp
1431	ENDIF
1432	ENDDO
1433
1434	CASE ( 'usm_wshf' )
1435	!-- array of sensible heat flux from surfaces
1436	!-- land surfaces
1437	IF ( .NOT. ALLOCATED(surf_usm_h%wshf_eb_av) ) THEN
1438	ALLOCATE( surf_usm_h%wshf_eb_av(1:surf_usm_h%ns) )
1439	surf_usm_h%wshf_eb_av = 0.0_wp
1440	ENDIF
1441	DO l = 0, 3
1442	IF ( .NOT. ALLOCATED(surf_usm_v(l)%wshf_eb_av) ) THEN
1443	ALLOCATE( surf_usm_v(l)%wshf_eb_av(1:surf_usm_v(l)%ns) )
1444	surf_usm_v(l)%wshf_eb_av = 0.0_wp
1445	ENDIF
1446	ENDDO
1447	!
1448	!-- Please note, the following output quantities belongs to the
1449	!-- individual tile fractions - ground heat flux at wall-, window-,
1450	!-- and green fraction. Aggregated ground-heat flux is treated
1451	!-- accordingly in average_3d_data, sum_up_3d_data, etc..
1452	CASE ( 'usm_wghf' )
1453	!-- array of heat flux from ground (wall, roof, land)
1454	IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_av) ) THEN
1455	ALLOCATE( surf_usm_h%wghf_eb_av(1:surf_usm_h%ns) )
1456	surf_usm_h%wghf_eb_av = 0.0_wp
1457	ENDIF
1458	DO l = 0, 3
1459	IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_av) ) THEN
1460	ALLOCATE( surf_usm_v(l)%wghf_eb_av(1:surf_usm_v(l)%ns) )
1461	surf_usm_v(l)%wghf_eb_av = 0.0_wp
1462	ENDIF
1463	ENDDO
1464
1465	CASE ( 'usm_wghf_window' )
1466	!-- array of heat flux from window ground (wall, roof, land)
1467	IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_window_av) ) THEN
1468	ALLOCATE( surf_usm_h%wghf_eb_window_av(1:surf_usm_h%ns) )
1469	surf_usm_h%wghf_eb_window_av = 0.0_wp
1470	ENDIF
1471	DO l = 0, 3
1472	IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_window_av) ) THEN
1473	ALLOCATE( surf_usm_v(l)%wghf_eb_window_av(1:surf_usm_v(l)%ns) )
1474	surf_usm_v(l)%wghf_eb_window_av = 0.0_wp
1475	ENDIF
1476	ENDDO
1477
1478	CASE ( 'usm_wghf_green' )
1479	!-- array of heat flux from green ground (wall, roof, land)
1480	IF ( .NOT. ALLOCATED(surf_usm_h%wghf_eb_green_av) ) THEN
1481	ALLOCATE( surf_usm_h%wghf_eb_green_av(1:surf_usm_h%ns) )
1482	surf_usm_h%wghf_eb_green_av = 0.0_wp
1483	ENDIF
1484	DO l = 0, 3
1485	IF ( .NOT. ALLOCATED(surf_usm_v(l)%wghf_eb_green_av) ) THEN
1486	ALLOCATE( surf_usm_v(l)%wghf_eb_green_av(1:surf_usm_v(l)%ns) )
1487	surf_usm_v(l)%wghf_eb_green_av = 0.0_wp
1488	ENDIF
1489	ENDDO
1490
1491	CASE ( 'usm_iwghf' )
1492	!-- array of heat flux from indoor ground (wall, roof, land)
1493	IF ( .NOT. ALLOCATED(surf_usm_h%iwghf_eb_av) ) THEN
1494	ALLOCATE( surf_usm_h%iwghf_eb_av(1:surf_usm_h%ns) )
1495	surf_usm_h%iwghf_eb_av = 0.0_wp
1496	ENDIF
1497	DO l = 0, 3
1498	IF ( .NOT. ALLOCATED(surf_usm_v(l)%iwghf_eb_av) ) THEN
1499	ALLOCATE( surf_usm_v(l)%iwghf_eb_av(1:surf_usm_v(l)%ns) )
1500	surf_usm_v(l)%iwghf_eb_av = 0.0_wp
1501	ENDIF
1502	ENDDO
1503
1504	CASE ( 'usm_iwghf_window' )
1505	!-- array of heat flux from indoor window ground (wall, roof, land)
1506	IF ( .NOT. ALLOCATED(surf_usm_h%iwghf_eb_window_av) ) THEN
1507	ALLOCATE( surf_usm_h%iwghf_eb_window_av(1:surf_usm_h%ns) )
1508	surf_usm_h%iwghf_eb_window_av = 0.0_wp
1509	ENDIF
1510	DO l = 0, 3
1511	IF ( .NOT. ALLOCATED(surf_usm_v(l)%iwghf_eb_window_av) ) THEN
1512	ALLOCATE( surf_usm_v(l)%iwghf_eb_window_av(1:surf_usm_v(l)%ns) )
1513	surf_usm_v(l)%iwghf_eb_window_av = 0.0_wp
1514	ENDIF
1515	ENDDO
1516
1517	CASE ( 'usm_t_surf' )
1518	!-- surface temperature for surfaces
1519	IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_av) ) THEN
1520	ALLOCATE( surf_usm_h%t_surf_av(1:surf_usm_h%ns) )
1521	surf_usm_h%t_surf_av = 0.0_wp
1522	ENDIF
1523	DO l = 0, 3
1524	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_av) ) THEN
1525	ALLOCATE( surf_usm_v(l)%t_surf_av(1:surf_usm_v(l)%ns) )
1526	surf_usm_v(l)%t_surf_av = 0.0_wp
1527	ENDIF
1528	ENDDO
1529
1530	CASE ( 'usm_t_surf_window' )
1531	!-- surface temperature for window surfaces
1532	IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_window_av) ) THEN
1533	ALLOCATE( surf_usm_h%t_surf_window_av(1:surf_usm_h%ns) )
1534	surf_usm_h%t_surf_window_av = 0.0_wp
1535	ENDIF
1536	DO l = 0, 3
1537	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_window_av) ) THEN
1538	ALLOCATE( surf_usm_v(l)%t_surf_window_av(1:surf_usm_v(l)%ns) )
1539	surf_usm_v(l)%t_surf_window_av = 0.0_wp
1540	ENDIF
1541	ENDDO
1542
1543	CASE ( 'usm_t_surf_green' )
1544	!-- surface temperature for green surfaces
1545	IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_green_av) ) THEN
1546	ALLOCATE( surf_usm_h%t_surf_green_av(1:surf_usm_h%ns) )
1547	surf_usm_h%t_surf_green_av = 0.0_wp
1548	ENDIF
1549	DO l = 0, 3
1550	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_green_av) ) THEN
1551	ALLOCATE( surf_usm_v(l)%t_surf_green_av(1:surf_usm_v(l)%ns) )
1552	surf_usm_v(l)%t_surf_green_av = 0.0_wp
1553	ENDIF
1554	ENDDO
1555
1556	CASE ( 'usm_t_surf_10cm' )
1557	!-- near surface temperature for whole surfaces
1558	IF ( .NOT. ALLOCATED(surf_usm_h%t_surf_10cm_av) ) THEN
1559	ALLOCATE( surf_usm_h%t_surf_10cm_av(1:surf_usm_h%ns) )
1560	surf_usm_h%t_surf_10cm_av = 0.0_wp
1561	ENDIF
1562	DO l = 0, 3
1563	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_surf_10cm_av) ) THEN
1564	ALLOCATE( surf_usm_v(l)%t_surf_10cm_av(1:surf_usm_v(l)%ns) )
1565	surf_usm_v(l)%t_surf_10cm_av = 0.0_wp
1566	ENDIF
1567	ENDDO
1568
1569	CASE ( 'usm_t_wall' )
1570	!-- wall temperature for iwl layer of walls and land
1571	IF ( .NOT. ALLOCATED(surf_usm_h%t_wall_av) ) THEN
1572	ALLOCATE( surf_usm_h%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1573	surf_usm_h%t_wall_av = 0.0_wp
1574	ENDIF
1575	DO l = 0, 3
1576	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_wall_av) ) THEN
1577	ALLOCATE( surf_usm_v(l)%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1578	surf_usm_v(l)%t_wall_av = 0.0_wp
1579	ENDIF
1580	ENDDO
1581
1582	CASE ( 'usm_t_window' )
1583	!-- window temperature for iwl layer of walls and land
1584	IF ( .NOT. ALLOCATED(surf_usm_h%t_window_av) ) THEN
1585	ALLOCATE( surf_usm_h%t_window_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1586	surf_usm_h%t_window_av = 0.0_wp
1587	ENDIF
1588	DO l = 0, 3
1589	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_window_av) ) THEN
1590	ALLOCATE( surf_usm_v(l)%t_window_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1591	surf_usm_v(l)%t_window_av = 0.0_wp
1592	ENDIF
1593	ENDDO
1594
1595	CASE ( 'usm_t_green' )
1596	!-- green temperature for iwl layer of walls and land
1597	IF ( .NOT. ALLOCATED(surf_usm_h%t_green_av) ) THEN
1598	ALLOCATE( surf_usm_h%t_green_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1599	surf_usm_h%t_green_av = 0.0_wp
1600	ENDIF
1601	DO l = 0, 3
1602	IF ( .NOT. ALLOCATED(surf_usm_v(l)%t_green_av) ) THEN
1603	ALLOCATE( surf_usm_v(l)%t_green_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1604	surf_usm_v(l)%t_green_av = 0.0_wp
1605	ENDIF
1606	ENDDO
1607
1608	CASE DEFAULT
1609	CONTINUE
1610
1611	END SELECT
1612
1613	ELSEIF ( mode == 'sum' ) THEN
1614
1615	SELECT CASE ( TRIM( var ) )
1616
1617	CASE ( 'usm_rad_net' )
1618	!-- array of complete radiation balance
1619	DO m = 1, surf_usm_h%ns
1620	surf_usm_h%rad_net_av(m) = &
1621	surf_usm_h%rad_net_av(m) + &
1622	surf_usm_h%rad_net_l(m)
1623	ENDDO
1624	DO l = 0, 3
1625	DO m = 1, surf_usm_v(l)%ns
1626	surf_usm_v(l)%rad_net_av(m) = &
1627	surf_usm_v(l)%rad_net_av(m) + &
1628	surf_usm_v(l)%rad_net_l(m)
1629	ENDDO
1630	ENDDO
1631
1632	CASE ( 'usm_rad_insw' )
1633	!-- array of sw radiation falling to surface after i-th reflection
1634	DO l = 1, nsurfl
1635	IF ( surfl(id,l) == idsint ) THEN
1636	surfinsw_av(l) = surfinsw_av(l) + surfinsw(l)
1637	ENDIF
1638	ENDDO
1639
1640	CASE ( 'usm_rad_inlw' )
1641	!-- array of lw radiation falling to surface after i-th reflection
1642	DO l = 1, nsurfl
1643	IF ( surfl(id,l) == idsint ) THEN
1644	surfinlw_av(l) = surfinlw_av(l) + surfinlw(l)
1645	ENDIF
1646	ENDDO
1647
1648	CASE ( 'usm_rad_inswdir' )
1649	!-- array of direct sw radiation falling to surface from sun
1650	DO l = 1, nsurfl
1651	IF ( surfl(id,l) == idsint ) THEN
1652	surfinswdir_av(l) = surfinswdir_av(l) + surfinswdir(l)
1653	ENDIF
1654	ENDDO
1655
1656	CASE ( 'usm_rad_inswdif' )
1657	!-- array of difusion sw radiation falling to surface from sky and borders of the domain
1658	DO l = 1, nsurfl
1659	IF ( surfl(id,l) == idsint ) THEN
1660	surfinswdif_av(l) = surfinswdif_av(l) + surfinswdif(l)
1661	ENDIF
1662	ENDDO
1663
1664	CASE ( 'usm_rad_inswref' )
1665	!-- array of sw radiation falling to surface from reflections
1666	DO l = 1, nsurfl
1667	IF ( surfl(id,l) == idsint ) THEN
1668	surfinswref_av(l) = surfinswref_av(l) + surfinsw(l) - &
1669	surfinswdir(l) - surfinswdif(l)
1670	ENDIF
1671	ENDDO
1672
1673
1674	CASE ( 'usm_rad_inlwdif' )
1675	!-- array of sw radiation falling to surface after i-th reflection
1676	DO l = 1, nsurfl
1677	IF ( surfl(id,l) == idsint ) THEN
1678	surfinlwdif_av(l) = surfinlwdif_av(l) + surfinlwdif(l)
1679	ENDIF
1680	ENDDO
1681	!
1682	CASE ( 'usm_rad_inlwref' )
1683	!-- array of lw radiation falling to surface from reflections
1684	DO l = 1, nsurfl
1685	IF ( surfl(id,l) == idsint ) THEN
1686	surfinlwref_av(l) = surfinlwref_av(l) + &
1687	surfinlw(l) - surfinlwdif(l)
1688	ENDIF
1689	ENDDO
1690
1691	CASE ( 'usm_rad_outsw' )
1692	!-- array of sw radiation emitted from surface after i-th reflection
1693	DO l = 1, nsurfl
1694	IF ( surfl(id,l) == idsint ) THEN
1695	surfoutsw_av(l) = surfoutsw_av(l) + surfoutsw(l)
1696	ENDIF
1697	ENDDO
1698
1699	CASE ( 'usm_rad_outlw' )
1700	!-- array of lw radiation emitted from surface after i-th reflection
1701	DO l = 1, nsurfl
1702	IF ( surfl(id,l) == idsint ) THEN
1703	surfoutlw_av(l) = surfoutlw_av(l) + surfoutlw(l)
1704	ENDIF
1705	ENDDO
1706
1707	CASE ( 'usm_rad_ressw' )
1708	!-- array of residua of sw radiation absorbed in surface after last reflection
1709	DO l = 1, nsurfl
1710	IF ( surfl(id,l) == idsint ) THEN
1711	surfins_av(l) = surfins_av(l) + surfins(l)
1712	ENDIF
1713	ENDDO
1714
1715	CASE ( 'usm_rad_reslw' )
1716	!-- array of residua of lw radiation absorbed in surface after last reflection
1717	DO l = 1, nsurfl
1718	IF ( surfl(id,l) == idsint ) THEN
1719	surfinl_av(l) = surfinl_av(l) + surfinl(l)
1720	ENDIF
1721	ENDDO
1722
1723	CASE ( 'usm_rad_hf' )
1724	!-- array of heat flux from radiation for surfaces after i-th reflection
1725	DO m = 1, surf_usm_h%ns
1726	surf_usm_h%surfhf_av(m) = &
1727	surf_usm_h%surfhf_av(m) + &
1728	surf_usm_h%surfhf(m)
1729	ENDDO
1730	DO l = 0, 3
1731	DO m = 1, surf_usm_v(l)%ns
1732	surf_usm_v(l)%surfhf_av(m) = &
1733	surf_usm_v(l)%surfhf_av(m) + &
1734	surf_usm_v(l)%surfhf(m)
1735	ENDDO
1736	ENDDO
1737
1738	CASE ( 'usm_wshf' )
1739	!-- array of sensible heat flux from surfaces (land, roof, wall)
1740	DO m = 1, surf_usm_h%ns
1741	surf_usm_h%wshf_eb_av(m) = &
1742	surf_usm_h%wshf_eb_av(m) + &
1743	surf_usm_h%wshf_eb(m)
1744	ENDDO
1745	DO l = 0, 3
1746	DO m = 1, surf_usm_v(l)%ns
1747	surf_usm_v(l)%wshf_eb_av(m) = &
1748	surf_usm_v(l)%wshf_eb_av(m) + &
1749	surf_usm_v(l)%wshf_eb(m)
1750	ENDDO
1751	ENDDO
1752
1753	CASE ( 'usm_wghf' )
1754	!-- array of heat flux from ground (wall, roof, land)
1755	DO m = 1, surf_usm_h%ns
1756	surf_usm_h%wghf_eb_av(m) = &
1757	surf_usm_h%wghf_eb_av(m) + &
1758	surf_usm_h%wghf_eb(m)
1759	ENDDO
1760	DO l = 0, 3
1761	DO m = 1, surf_usm_v(l)%ns
1762	surf_usm_v(l)%wghf_eb_av(m) = &
1763	surf_usm_v(l)%wghf_eb_av(m) + &
1764	surf_usm_v(l)%wghf_eb(m)
1765	ENDDO
1766	ENDDO
1767
1768	CASE ( 'usm_wghf_window' )
1769	!-- array of heat flux from window ground (wall, roof, land)
1770	DO m = 1, surf_usm_h%ns
1771	surf_usm_h%wghf_eb_window_av(m) = &
1772	surf_usm_h%wghf_eb_window_av(m) + &
1773	surf_usm_h%wghf_eb_window(m)
1774	ENDDO
1775	DO l = 0, 3
1776	DO m = 1, surf_usm_v(l)%ns
1777	surf_usm_v(l)%wghf_eb_window_av(m) = &
1778	surf_usm_v(l)%wghf_eb_window_av(m) + &
1779	surf_usm_v(l)%wghf_eb_window(m)
1780	ENDDO
1781	ENDDO
1782
1783	CASE ( 'usm_wghf_green' )
1784	!-- array of heat flux from green ground (wall, roof, land)
1785	DO m = 1, surf_usm_h%ns
1786	surf_usm_h%wghf_eb_green_av(m) = &
1787	surf_usm_h%wghf_eb_green_av(m) + &
1788	surf_usm_h%wghf_eb_green(m)
1789	ENDDO
1790	DO l = 0, 3
1791	DO m = 1, surf_usm_v(l)%ns
1792	surf_usm_v(l)%wghf_eb_green_av(m) = &
1793	surf_usm_v(l)%wghf_eb_green_av(m) + &
1794	surf_usm_v(l)%wghf_eb_green(m)
1795	ENDDO
1796	ENDDO
1797
1798	CASE ( 'usm_iwghf' )
1799	!-- array of heat flux from indoor ground (wall, roof, land)
1800	DO m = 1, surf_usm_h%ns
1801	surf_usm_h%iwghf_eb_av(m) = &
1802	surf_usm_h%iwghf_eb_av(m) + &
1803	surf_usm_h%iwghf_eb(m)
1804	ENDDO
1805	DO l = 0, 3
1806	DO m = 1, surf_usm_v(l)%ns
1807	surf_usm_v(l)%iwghf_eb_av(m) = &
1808	surf_usm_v(l)%iwghf_eb_av(m) + &
1809	surf_usm_v(l)%iwghf_eb(m)
1810	ENDDO
1811	ENDDO
1812
1813	CASE ( 'usm_iwghf_window' )
1814	!-- array of heat flux from indoor window ground (wall, roof, land)
1815	DO m = 1, surf_usm_h%ns
1816	surf_usm_h%iwghf_eb_window_av(m) = &
1817	surf_usm_h%iwghf_eb_window_av(m) + &
1818	surf_usm_h%iwghf_eb_window(m)
1819	ENDDO
1820	DO l = 0, 3
1821	DO m = 1, surf_usm_v(l)%ns
1822	surf_usm_v(l)%iwghf_eb_window_av(m) = &
1823	surf_usm_v(l)%iwghf_eb_window_av(m) + &
1824	surf_usm_v(l)%iwghf_eb_window(m)
1825	ENDDO
1826	ENDDO
1827
1828	CASE ( 'usm_t_surf' )
1829	!-- surface temperature for surfaces
1830	DO m = 1, surf_usm_h%ns
1831	surf_usm_h%t_surf_av(m) = &
1832	surf_usm_h%t_surf_av(m) + &
1833	t_surf_h(m)
1834	ENDDO
1835	DO l = 0, 3
1836	DO m = 1, surf_usm_v(l)%ns
1837	surf_usm_v(l)%t_surf_av(m) = &
1838	surf_usm_v(l)%t_surf_av(m) + &
1839	t_surf_v(l)%t(m)
1840	ENDDO
1841	ENDDO
1842
1843	CASE ( 'usm_t_surf_window' )
1844	!-- surface temperature for window surfaces
1845	DO m = 1, surf_usm_h%ns
1846	surf_usm_h%t_surf_window_av(m) = &
1847	surf_usm_h%t_surf_window_av(m) + &
1848	t_surf_window_h(m)
1849	ENDDO
1850	DO l = 0, 3
1851	DO m = 1, surf_usm_v(l)%ns
1852	surf_usm_v(l)%t_surf_window_av(m) = &
1853	surf_usm_v(l)%t_surf_window_av(m) + &
1854	t_surf_window_v(l)%t(m)
1855	ENDDO
1856	ENDDO
1857
1858	CASE ( 'usm_t_surf_green' )
1859	!-- surface temperature for green surfaces
1860	DO m = 1, surf_usm_h%ns
1861	surf_usm_h%t_surf_green_av(m) = &
1862	surf_usm_h%t_surf_green_av(m) + &
1863	t_surf_green_h(m)
1864	ENDDO
1865	DO l = 0, 3
1866	DO m = 1, surf_usm_v(l)%ns
1867	surf_usm_v(l)%t_surf_green_av(m) = &
1868	surf_usm_v(l)%t_surf_green_av(m) + &
1869	t_surf_green_v(l)%t(m)
1870	ENDDO
1871	ENDDO
1872
1873	CASE ( 'usm_t_surf_10cm' )
1874	!-- near surface temperature for whole surfaces
1875	DO m = 1, surf_usm_h%ns
1876	surf_usm_h%t_surf_10cm_av(m) = &
1877	surf_usm_h%t_surf_10cm_av(m) + &
1878	t_surf_10cm_h(m)
1879	ENDDO
1880	DO l = 0, 3
1881	DO m = 1, surf_usm_v(l)%ns
1882	surf_usm_v(l)%t_surf_10cm_av(m) = &
1883	surf_usm_v(l)%t_surf_10cm_av(m) + &
1884	t_surf_10cm_v(l)%t(m)
1885	ENDDO
1886	ENDDO
1887
1888
1889	CASE ( 'usm_t_wall' )
1890	!-- wall temperature for iwl layer of walls and land
1891	DO m = 1, surf_usm_h%ns
1892	surf_usm_h%t_wall_av(iwl,m) = &
1893	surf_usm_h%t_wall_av(iwl,m) + &
1894	t_wall_h(iwl,m)
1895	ENDDO
1896	DO l = 0, 3
1897	DO m = 1, surf_usm_v(l)%ns
1898	surf_usm_v(l)%t_wall_av(iwl,m) = &
1899	surf_usm_v(l)%t_wall_av(iwl,m) + &
1900	t_wall_v(l)%t(iwl,m)
1901	ENDDO
1902	ENDDO
1903
1904	CASE ( 'usm_t_window' )
1905	!-- window temperature for iwl layer of walls and land
1906	DO m = 1, surf_usm_h%ns
1907	surf_usm_h%t_window_av(iwl,m) = &
1908	surf_usm_h%t_window_av(iwl,m) + &
1909	t_window_h(iwl,m)
1910	ENDDO
1911	DO l = 0, 3
1912	DO m = 1, surf_usm_v(l)%ns
1913	surf_usm_v(l)%t_window_av(iwl,m) = &
1914	surf_usm_v(l)%t_window_av(iwl,m) + &
1915	t_window_v(l)%t(iwl,m)
1916	ENDDO
1917	ENDDO
1918
1919	CASE ( 'usm_t_green' )
1920	!-- green temperature for iwl layer of walls and land
1921	DO m = 1, surf_usm_h%ns
1922	surf_usm_h%t_green_av(iwl,m) = &
1923	surf_usm_h%t_green_av(iwl,m) + &
1924	t_green_h(iwl,m)
1925	ENDDO
1926	DO l = 0, 3
1927	DO m = 1, surf_usm_v(l)%ns
1928	surf_usm_v(l)%t_green_av(iwl,m) = &
1929	surf_usm_v(l)%t_green_av(iwl,m) + &
1930	t_green_v(l)%t(iwl,m)
1931	ENDDO
1932	ENDDO
1933
1934	CASE DEFAULT
1935	CONTINUE
1936
1937	END SELECT
1938
1939	ELSEIF ( mode == 'average' ) THEN
1940
1941	SELECT CASE ( TRIM( var ) )
1942
1943	CASE ( 'usm_rad_net' )
1944	!-- array of complete radiation balance
1945	DO m = 1, surf_usm_h%ns
1946	surf_usm_h%rad_net_av(m) = &
1947	surf_usm_h%rad_net_av(m) / &
1948	REAL( average_count_3d, kind=wp )
1949	ENDDO
1950	DO l = 0, 3
1951	DO m = 1, surf_usm_v(l)%ns
1952	surf_usm_v(l)%rad_net_av(m) = &
1953	surf_usm_v(l)%rad_net_av(m) / &
1954	REAL( average_count_3d, kind=wp )
1955	ENDDO
1956	ENDDO
1957
1958	CASE ( 'usm_rad_insw' )
1959	!-- array of sw radiation falling to surface after i-th reflection
1960	DO l = 1, nsurfl
1961	IF ( surfl(id,l) == idsint ) THEN
1962	surfinsw_av(l) = surfinsw_av(l) / REAL( average_count_3d, kind=wp )
1963	ENDIF
1964	ENDDO
1965
1966	CASE ( 'usm_rad_inlw' )
1967	!-- array of lw radiation falling to surface after i-th reflection
1968	DO l = 1, nsurfl
1969	IF ( surfl(id,l) == idsint ) THEN
1970	surfinlw_av(l) = surfinlw_av(l) / REAL( average_count_3d, kind=wp )
1971	ENDIF
1972	ENDDO
1973
1974	CASE ( 'usm_rad_inswdir' )
1975	!-- array of direct sw radiation falling to surface from sun
1976	DO l = 1, nsurfl
1977	IF ( surfl(id,l) == idsint ) THEN
1978	surfinswdir_av(l) = surfinswdir_av(l) / REAL( average_count_3d, kind=wp )
1979	ENDIF
1980	ENDDO
1981
1982	CASE ( 'usm_rad_inswdif' )
1983	!-- array of difusion sw radiation falling to surface from sky and borders of the domain
1984	DO l = 1, nsurfl
1985	IF ( surfl(id,l) == idsint ) THEN
1986	surfinswdif_av(l) = surfinswdif_av(l) / REAL( average_count_3d, kind=wp )
1987	ENDIF
1988	ENDDO
1989
1990	CASE ( 'usm_rad_inswref' )
1991	!-- array of sw radiation falling to surface from reflections
1992	DO l = 1, nsurfl
1993	IF ( surfl(id,l) == idsint ) THEN
1994	surfinswref_av(l) = surfinswref_av(l) / REAL( average_count_3d, kind=wp )
1995	ENDIF
1996	ENDDO
1997
1998	CASE ( 'usm_rad_inlwdif' )
1999	!-- array of sw radiation falling to surface after i-th reflection
2000	DO l = 1, nsurfl
2001	IF ( surfl(id,l) == idsint ) THEN
2002	surfinlwdif_av(l) = surfinlwdif_av(l) / REAL( average_count_3d, kind=wp )
2003	ENDIF
2004	ENDDO
2005
2006	CASE ( 'usm_rad_inlwref' )
2007	!-- array of lw radiation falling to surface from reflections
2008	DO l = 1, nsurfl
2009	IF ( surfl(id,l) == idsint ) THEN
2010	surfinlwref_av(l) = surfinlwref_av(l) / REAL( average_count_3d, kind=wp )
2011	ENDIF
2012	ENDDO
2013
2014	CASE ( 'usm_rad_outsw' )
2015	!-- array of sw radiation emitted from surface after i-th reflection
2016	DO l = 1, nsurfl
2017	IF ( surfl(id,l) == idsint ) THEN
2018	surfoutsw_av(l) = surfoutsw_av(l) / REAL( average_count_3d, kind=wp )
2019	ENDIF
2020	ENDDO
2021
2022	CASE ( 'usm_rad_outlw' )
2023	!-- array of lw radiation emitted from surface after i-th reflection
2024	DO l = 1, nsurfl
2025	IF ( surfl(id,l) == idsint ) THEN
2026	surfoutlw_av(l) = surfoutlw_av(l) / REAL( average_count_3d, kind=wp )
2027	ENDIF
2028	ENDDO
2029
2030	CASE ( 'usm_rad_ressw' )
2031	!-- array of residua of sw radiation absorbed in surface after last reflection
2032	DO l = 1, nsurfl
2033	IF ( surfl(id,l) == idsint ) THEN
2034	surfins_av(l) = surfins_av(l) / REAL( average_count_3d, kind=wp )
2035	ENDIF
2036	ENDDO
2037
2038	CASE ( 'usm_rad_reslw' )
2039	!-- array of residua of lw radiation absorbed in surface after last reflection
2040	DO l = 1, nsurfl
2041	IF ( surfl(id,l) == idsint ) THEN
2042	surfinl_av(l) = surfinl_av(l) / REAL( average_count_3d, kind=wp )
2043	ENDIF
2044	ENDDO
2045
2046	CASE ( 'usm_rad_hf' )
2047	!-- array of heat flux from radiation for surfaces after i-th reflection
2048	DO m = 1, surf_usm_h%ns
2049	surf_usm_h%surfhf_av(m) = &
2050	surf_usm_h%surfhf_av(m) / &
2051	REAL( average_count_3d, kind=wp )
2052	ENDDO
2053	DO l = 0, 3
2054	DO m = 1, surf_usm_v(l)%ns
2055	surf_usm_v(l)%surfhf_av(m) = &
2056	surf_usm_v(l)%surfhf_av(m) / &
2057	REAL( average_count_3d, kind=wp )
2058	ENDDO
2059	ENDDO
2060
2061	CASE ( 'usm_wshf' )
2062	!-- array of sensible heat flux from surfaces (land, roof, wall)
2063	DO m = 1, surf_usm_h%ns
2064	surf_usm_h%wshf_eb_av(m) = &
2065	surf_usm_h%wshf_eb_av(m) / &
2066	REAL( average_count_3d, kind=wp )
2067	ENDDO
2068	DO l = 0, 3
2069	DO m = 1, surf_usm_v(l)%ns
2070	surf_usm_v(l)%wshf_eb_av(m) = &
2071	surf_usm_v(l)%wshf_eb_av(m) / &
2072	REAL( average_count_3d, kind=wp )
2073	ENDDO
2074	ENDDO
2075
2076	CASE ( 'usm_wghf' )
2077	!-- array of heat flux from ground (wall, roof, land)
2078	DO m = 1, surf_usm_h%ns
2079	surf_usm_h%wghf_eb_av(m) = &
2080	surf_usm_h%wghf_eb_av(m) / &
2081	REAL( average_count_3d, kind=wp )
2082	ENDDO
2083	DO l = 0, 3
2084	DO m = 1, surf_usm_v(l)%ns
2085	surf_usm_v(l)%wghf_eb_av(m) = &
2086	surf_usm_v(l)%wghf_eb_av(m) / &
2087	REAL( average_count_3d, kind=wp )
2088	ENDDO
2089	ENDDO
2090
2091	CASE ( 'usm_wghf_window' )
2092	!-- array of heat flux from window ground (wall, roof, land)
2093	DO m = 1, surf_usm_h%ns
2094	surf_usm_h%wghf_eb_window_av(m) = &
2095	surf_usm_h%wghf_eb_window_av(m) / &
2096	REAL( average_count_3d, kind=wp )
2097	ENDDO
2098	DO l = 0, 3
2099	DO m = 1, surf_usm_v(l)%ns
2100	surf_usm_v(l)%wghf_eb_window_av(m) = &
2101	surf_usm_v(l)%wghf_eb_window_av(m) / &
2102	REAL( average_count_3d, kind=wp )
2103	ENDDO
2104	ENDDO
2105
2106	CASE ( 'usm_wghf_green' )
2107	!-- array of heat flux from green ground (wall, roof, land)
2108	DO m = 1, surf_usm_h%ns
2109	surf_usm_h%wghf_eb_green_av(m) = &
2110	surf_usm_h%wghf_eb_green_av(m) / &
2111	REAL( average_count_3d, kind=wp )
2112	ENDDO
2113	DO l = 0, 3
2114	DO m = 1, surf_usm_v(l)%ns
2115	surf_usm_v(l)%wghf_eb_green_av(m) = &
2116	surf_usm_v(l)%wghf_eb_green_av(m) / &
2117	REAL( average_count_3d, kind=wp )
2118	ENDDO
2119	ENDDO
2120
2121	CASE ( 'usm_iwghf' )
2122	!-- array of heat flux from indoor ground (wall, roof, land)
2123	DO m = 1, surf_usm_h%ns
2124	surf_usm_h%iwghf_eb_av(m) = &
2125	surf_usm_h%iwghf_eb_av(m) / &
2126	REAL( average_count_3d, kind=wp )
2127	ENDDO
2128	DO l = 0, 3
2129	DO m = 1, surf_usm_v(l)%ns
2130	surf_usm_v(l)%iwghf_eb_av(m) = &
2131	surf_usm_v(l)%iwghf_eb_av(m) / &
2132	REAL( average_count_3d, kind=wp )
2133	ENDDO
2134	ENDDO
2135
2136	CASE ( 'usm_iwghf_window' )
2137	!-- array of heat flux from indoor window ground (wall, roof, land)
2138	DO m = 1, surf_usm_h%ns
2139	surf_usm_h%iwghf_eb_window_av(m) = &
2140	surf_usm_h%iwghf_eb_window_av(m) / &
2141	REAL( average_count_3d, kind=wp )
2142	ENDDO
2143	DO l = 0, 3
2144	DO m = 1, surf_usm_v(l)%ns
2145	surf_usm_v(l)%iwghf_eb_window_av(m) = &
2146	surf_usm_v(l)%iwghf_eb_window_av(m) / &
2147	REAL( average_count_3d, kind=wp )
2148	ENDDO
2149	ENDDO
2150
2151	CASE ( 'usm_t_surf' )
2152	!-- surface temperature for surfaces
2153	DO m = 1, surf_usm_h%ns
2154	surf_usm_h%t_surf_av(m) = &
2155	surf_usm_h%t_surf_av(m) / &
2156	REAL( average_count_3d, kind=wp )
2157	ENDDO
2158	DO l = 0, 3
2159	DO m = 1, surf_usm_v(l)%ns
2160	surf_usm_v(l)%t_surf_av(m) = &
2161	surf_usm_v(l)%t_surf_av(m) / &
2162	REAL( average_count_3d, kind=wp )
2163	ENDDO
2164	ENDDO
2165
2166	CASE ( 'usm_t_surf_window' )
2167	!-- surface temperature for window surfaces
2168	DO m = 1, surf_usm_h%ns
2169	surf_usm_h%t_surf_window_av(m) = &
2170	surf_usm_h%t_surf_window_av(m) / &
2171	REAL( average_count_3d, kind=wp )
2172	ENDDO
2173	DO l = 0, 3
2174	DO m = 1, surf_usm_v(l)%ns
2175	surf_usm_v(l)%t_surf_window_av(m) = &
2176	surf_usm_v(l)%t_surf_window_av(m) / &
2177	REAL( average_count_3d, kind=wp )
2178	ENDDO
2179	ENDDO
2180
2181	CASE ( 'usm_t_surf_green' )
2182	!-- surface temperature for green surfaces
2183	DO m = 1, surf_usm_h%ns
2184	surf_usm_h%t_surf_green_av(m) = &
2185	surf_usm_h%t_surf_green_av(m) / &
2186	REAL( average_count_3d, kind=wp )
2187	ENDDO
2188	DO l = 0, 3
2189	DO m = 1, surf_usm_v(l)%ns
2190	surf_usm_v(l)%t_surf_green_av(m) = &
2191	surf_usm_v(l)%t_surf_green_av(m) / &
2192	REAL( average_count_3d, kind=wp )
2193	ENDDO
2194	ENDDO
2195
2196	CASE ( 'usm_t_surf_10cm' )
2197	!-- near surface temperature for whole surfaces
2198	DO m = 1, surf_usm_h%ns
2199	surf_usm_h%t_surf_10cm_av(m) = &
2200	surf_usm_h%t_surf_10cm_av(m) / &
2201	REAL( average_count_3d, kind=wp )
2202	ENDDO
2203	DO l = 0, 3
2204	DO m = 1, surf_usm_v(l)%ns
2205	surf_usm_v(l)%t_surf_10cm_av(m) = &
2206	surf_usm_v(l)%t_surf_10cm_av(m) / &
2207	REAL( average_count_3d, kind=wp )
2208	ENDDO
2209	ENDDO
2210
2211	CASE ( 'usm_t_wall' )
2212	!-- wall temperature for iwl layer of walls and land
2213	DO m = 1, surf_usm_h%ns
2214	surf_usm_h%t_wall_av(iwl,m) = &
2215	surf_usm_h%t_wall_av(iwl,m) / &
2216	REAL( average_count_3d, kind=wp )
2217	ENDDO
2218	DO l = 0, 3
2219	DO m = 1, surf_usm_v(l)%ns
2220	surf_usm_v(l)%t_wall_av(iwl,m) = &
2221	surf_usm_v(l)%t_wall_av(iwl,m) / &
2222	REAL( average_count_3d, kind=wp )
2223	ENDDO
2224	ENDDO
2225
2226	CASE ( 'usm_t_window' )
2227	!-- window temperature for iwl layer of walls and land
2228	DO m = 1, surf_usm_h%ns
2229	surf_usm_h%t_window_av(iwl,m) = &
2230	surf_usm_h%t_window_av(iwl,m) / &
2231	REAL( average_count_3d, kind=wp )
2232	ENDDO
2233	DO l = 0, 3
2234	DO m = 1, surf_usm_v(l)%ns
2235	surf_usm_v(l)%t_window_av(iwl,m) = &
2236	surf_usm_v(l)%t_window_av(iwl,m) / &
2237	REAL( average_count_3d, kind=wp )
2238	ENDDO
2239	ENDDO
2240
2241	CASE ( 'usm_t_green' )
2242	!-- green temperature for iwl layer of walls and land
2243	DO m = 1, surf_usm_h%ns
2244	surf_usm_h%t_green_av(iwl,m) = &
2245	surf_usm_h%t_green_av(iwl,m) / &
2246	REAL( average_count_3d, kind=wp )
2247	ENDDO
2248	DO l = 0, 3
2249	DO m = 1, surf_usm_v(l)%ns
2250	surf_usm_v(l)%t_green_av(iwl,m) = &
2251	surf_usm_v(l)%t_green_av(iwl,m) / &
2252	REAL( average_count_3d, kind=wp )
2253	ENDDO
2254	ENDDO
2255
2256
2257	END SELECT
2258
2259	ENDIF
2260
2261	END SUBROUTINE usm_average_3d_data
2262
2263
2264
2265	!------------------------------------------------------------------------------!
2266	! Description:
2267	! ------------
2268	!> Set internal Neumann boundary condition at outer soil grid points
2269	!> for temperature and humidity.
2270	!------------------------------------------------------------------------------!
2271	SUBROUTINE usm_boundary_condition
2272
2273	IMPLICIT NONE
2274
2275	INTEGER(iwp) :: i !< grid index x-direction
2276	INTEGER(iwp) :: ioff !< offset index x-direction indicating location of soil grid point
2277	INTEGER(iwp) :: j !< grid index y-direction
2278	INTEGER(iwp) :: joff !< offset index x-direction indicating location of soil grid point
2279	INTEGER(iwp) :: k !< grid index z-direction
2280	INTEGER(iwp) :: koff !< offset index x-direction indicating location of soil grid point
2281	INTEGER(iwp) :: l !< running index surface-orientation
2282	INTEGER(iwp) :: m !< running index surface elements
2283
2284	koff = surf_usm_h%koff
2285	DO m = 1, surf_usm_h%ns
2286	i = surf_usm_h%i(m)
2287	j = surf_usm_h%j(m)
2288	k = surf_usm_h%k(m)
2289	pt(k+koff,j,i) = pt(k,j,i)
2290	ENDDO
2291
2292	DO l = 0, 3
2293	ioff = surf_usm_v(l)%ioff
2294	joff = surf_usm_v(l)%joff
2295	DO m = 1, surf_usm_v(l)%ns
2296	i = surf_usm_v(l)%i(m)
2297	j = surf_usm_v(l)%j(m)
2298	k = surf_usm_v(l)%k(m)
2299	pt(k,j+joff,i+ioff) = pt(k,j,i)
2300	ENDDO
2301	ENDDO
2302
2303	END SUBROUTINE usm_boundary_condition
2304
2305
2306	!------------------------------------------------------------------------------!
2307	!
2308	! Description:
2309	! ------------
2310	!> Subroutine checks variables and assigns units.
2311	!> It is called out from subroutine check_parameters.
2312	!------------------------------------------------------------------------------!
2313	SUBROUTINE usm_check_data_output( variable, unit )
2314
2315	IMPLICIT NONE
2316
2317	CHARACTER (len=*),INTENT(IN) :: variable !:
2318	CHARACTER (len=*),INTENT(OUT) :: unit !:
2319
2320	CHARACTER (len=varnamelength) :: var
2321
2322	var = TRIM(variable)
2323	IF ( var(1:12) == 'usm_rad_net_' .OR. var(1:13) == 'usm_rad_insw_' .OR. &
2324	var(1:13) == 'usm_rad_inlw_' .OR. var(1:16) == 'usm_rad_inswdir_' .OR. &
2325	var(1:16) == 'usm_rad_inswdif_' .OR. var(1:16) == 'usm_rad_inswref_' .OR. &
2326	var(1:16) == 'usm_rad_inlwdif_' .OR. var(1:16) == 'usm_rad_inlwref_' .OR. &
2327	var(1:14) == 'usm_rad_outsw_' .OR. var(1:14) == 'usm_rad_outlw_' .OR. &
2328	var(1:14) == 'usm_rad_ressw_' .OR. var(1:14) == 'usm_rad_reslw_' .OR. &
2329	var(1:11) == 'usm_rad_hf_' .OR. &
2330	var(1:9) == 'usm_wshf_' .OR. var(1:9) == 'usm_wghf_' .OR. &
2331	var(1:16) == 'usm_wghf_window_' .OR. var(1:15) == 'usm_wghf_green_' .OR. &
2332	var(1:10) == 'usm_iwghf_' .OR. var(1:17) == 'usm_iwghf_window_' ) THEN
2333	unit = 'W/m2'
2334	ELSE IF ( var(1:10) == 'usm_t_surf' .OR. var(1:10) == 'usm_t_wall' .OR. &
2335	var(1:12) == 'usm_t_window' .OR. var(1:17) == 'usm_t_surf_window' .OR. &
2336	var(1:16) == 'usm_t_surf_green' .OR. &
2337	var(1:11) == 'usm_t_green' .OR. &
2338	var(1:15) == 'usm_t_surf_10cm') THEN
2339	unit = 'K'
2340	ELSE IF ( var(1:9) == 'usm_surfz' .OR. var(1:7) == 'usm_svf' .OR. &
2341	var(1:7) == 'usm_dif' .OR. var(1:11) == 'usm_surfcat' .OR. &
2342	var(1:11) == 'usm_surfalb' .OR. var(1:12) == 'usm_surfemis' .OR. &
2343	var(1:9) == 'usm_skyvf' .OR. var(1:9) == 'usm_skyvft' ) THEN
2344	unit = '1'
2345	ELSE
2346	unit = 'illegal'
2347	ENDIF
2348
2349	END SUBROUTINE usm_check_data_output
2350
2351
2352	!------------------------------------------------------------------------------!
2353	! Description:
2354	! ------------
2355	!> Check parameters routine for urban surface model
2356	!------------------------------------------------------------------------------!
2357	SUBROUTINE usm_check_parameters
2358
2359	USE control_parameters, &
2360	ONLY: bc_pt_b, bc_q_b, constant_flux_layer, large_scale_forcing, &
2361	lsf_surf, topography
2362
2363	!
2364	!-- Dirichlet boundary conditions are required as the surface fluxes are
2365	!-- calculated from the temperature/humidity gradients in the urban surface
2366	!-- model
2367	IF ( bc_pt_b == 'neumann' .OR. bc_q_b == 'neumann' ) THEN
2368	message_string = 'urban surface model requires setting of '// &
2369	'bc_pt_b = "dirichlet" and '// &
2370	'bc_q_b = "dirichlet"'
2371	CALL message( 'usm_check_parameters', 'PA0590', 1, 2, 0, 6, 0 )
2372	ENDIF
2373
2374	IF ( .NOT. constant_flux_layer ) THEN
2375	message_string = 'urban surface model requires '// &
2376	'constant_flux_layer = .T.'
2377	CALL message( 'usm_check_parameters', 'PA0084', 1, 2, 0, 6, 0 )
2378	ENDIF
2379
2380	IF ( .NOT. radiation ) THEN
2381	message_string = 'urban surface model requires '// &
2382	'the radiation model to be switched on'
2383	CALL message( 'usm_check_parameters', 'PA0084', 1, 2, 0, 6, 0 )
2384	ENDIF
2385	!
2386	!-- Surface forcing has to be disabled for LSF in case of enabled
2387	!-- urban surface module
2388	IF ( large_scale_forcing ) THEN
2389	lsf_surf = .FALSE.
2390	ENDIF
2391	!
2392	!-- Topography
2393	IF ( topography == 'flat' ) THEN
2394	message_string = 'topography /= "flat" is required '// &
2395	'when using the urban surface model'
2396	CALL message( 'check_parameters', 'PA0592', 1, 2, 0, 6, 0 )
2397	ENDIF
2398	!
2399	!-- naheatlayers
2400	IF ( naheatlayers > nzt ) THEN
2401	message_string = 'number of anthropogenic heat layers '// &
2402	'"naheatlayers" can not be larger than'// &
2403	' number of domain layers "nzt"'
2404	CALL message( 'check_parameters', 'PA0593', 1, 2, 0, 6, 0 )
2405	ENDIF
2406
2407	END SUBROUTINE usm_check_parameters
2408
2409
2410	!------------------------------------------------------------------------------!
2411	!
2412	! Description:
2413	! ------------
2414	!> Output of the 3D-arrays in netCDF and/or AVS format
2415	!> for variables of urban_surface model.
2416	!> It resorts the urban surface module output quantities from surf style
2417	!> indexing into temporary 3D array with indices (i,j,k).
2418	!> It is called from subroutine data_output_3d.
2419	!------------------------------------------------------------------------------!
2420	SUBROUTINE usm_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
2421
2422	IMPLICIT NONE
2423
2424	INTEGER(iwp), INTENT(IN) :: av !<
2425	CHARACTER (len=*), INTENT(IN) :: variable !<
2426	INTEGER(iwp), INTENT(IN) :: nzb_do !< lower limit of the data output (usually 0)
2427	INTEGER(iwp), INTENT(IN) :: nzt_do !< vertical upper limit of the data output (usually nz_do3d)
2428	LOGICAL, INTENT(OUT) :: found !<
2429	REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) :: local_pf !< sp - it has to correspond to module data_output_3d
2430	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: temp_pf !< temp array for urban surface output procedure
2431
2432	CHARACTER (len=varnamelength) :: var, surfid
2433	INTEGER(iwp), PARAMETER :: nd = 5
2434	CHARACTER(len=6), DIMENSION(0:nd-1), PARAMETER :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /)
2435	INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER :: dirint = (/ iup_u, isouth_u, inorth_u, iwest_u, ieast_u /)
2436	INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER :: diridx = (/ -1, 1, 0, 3, 2 /)
2437	!< index for surf_*_v: 0:3 = (North, South, East, West)
2438	INTEGER(iwp), DIMENSION(0:nd-1) :: dirstart
2439	INTEGER(iwp), DIMENSION(0:nd-1) :: dirend
2440	INTEGER(iwp) :: ids,idsint,idsidx,isurf,isvf,isurfs,isurflt
2441	INTEGER(iwp) :: is,js,ks,i,j,k,iwl,istat, l, m
2442	INTEGER(iwp) :: k_topo !< topography top index
2443
2444	dirstart = (/ startland, startwall, startwall, startwall, startwall /)
2445	dirend = (/ endland, endwall, endwall, endwall, endwall /)
2446
2447	found = .TRUE.
2448	temp_pf = -1._wp
2449
2450	ids = -1
2451	var = TRIM(variable)
2452	DO i = 0, nd-1
2453	k = len(TRIM(var))
2454	j = len(TRIM(dirname(i)))
2455	IF ( var(k-j+1:k) == dirname(i) ) THEN
2456	ids = i
2457	idsint = dirint(ids)
2458	idsidx = diridx(ids)
2459	var = var(:k-j)
2460	EXIT
2461	ENDIF
2462	ENDDO
2463	IF ( ids == -1 ) THEN
2464	var = TRIM(variable)
2465	ENDIF
2466	IF ( var(1:11) == 'usm_t_wall_' .AND. len(TRIM(var)) >= 12 ) THEN
2467	!-- wall layers
2468	READ(var(12:12), '(I1)', iostat=istat ) iwl
2469	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
2470	var = var(1:10)
2471	ENDIF
2472	ENDIF
2473	IF ( var(1:13) == 'usm_t_window_' .AND. len(TRIM(var)) >= 14 ) THEN
2474	!-- window layers
2475	READ(var(14:14), '(I1)', iostat=istat ) iwl
2476	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
2477	var = var(1:12)
2478	ENDIF
2479	ENDIF
2480	IF ( var(1:12) == 'usm_t_green_' .AND. len(TRIM(var)) >= 13 ) THEN
2481	!-- green layers
2482	READ(var(13:13), '(I1)', iostat=istat ) iwl
2483	IF ( istat == 0 .AND. iwl >= nzb_wall .AND. iwl <= nzt_wall ) THEN
2484	var = var(1:11)
2485	ENDIF
2486	ENDIF
2487	IF ( (var(1:8) == 'usm_svf_' .OR. var(1:8) == 'usm_dif_') .AND. len(TRIM(var)) >= 13 ) THEN
2488	!-- svf values to particular surface
2489	surfid = var(9:)
2490	i = index(surfid,'_')
2491	j = index(surfid(i+1:),'_')
2492	READ(surfid(1:i-1),*, iostat=istat ) is
2493	IF ( istat == 0 ) THEN
2494	READ(surfid(i+1:i+j-1),*, iostat=istat ) js
2495	ENDIF
2496	IF ( istat == 0 ) THEN
2497	READ(surfid(i+j+1:),*, iostat=istat ) ks
2498	ENDIF
2499	IF ( istat == 0 ) THEN
2500	var = var(1:7)
2501	ENDIF
2502	ENDIF
2503
2504	SELECT CASE ( TRIM(var) )
2505
2506	CASE ( 'usm_surfz' )
2507	!-- array of lw radiation falling to local surface after i-th reflection
2508	IF ( idsint == iup_u ) THEN
2509	DO m = 1, surf_usm_h%ns
2510	i = surf_usm_h%i(m)
2511	j = surf_usm_h%j(m)
2512	k = surf_usm_h%k(m)
2513	temp_pf(0,j,i) = MAX( temp_pf(0,j,i), REAL( k, kind=wp) )
2514	ENDDO
2515	ELSE
2516	l = idsidx
2517	DO m = 1, surf_usm_v(l)%ns
2518	i = surf_usm_v(l)%i(m)
2519	j = surf_usm_v(l)%j(m)
2520	k = surf_usm_v(l)%k(m)
2521	temp_pf(0,j,i) = MAX( temp_pf(0,j,i), REAL( k, kind=wp) + 1.0_wp )
2522	ENDDO
2523	ENDIF
2524
2525	CASE ( 'usm_surfcat' )
2526	!-- surface category
2527	IF ( idsint == iup_u ) THEN
2528	DO m = 1, surf_usm_h%ns
2529	i = surf_usm_h%i(m)
2530	j = surf_usm_h%j(m)
2531	k = surf_usm_h%k(m)
2532	temp_pf(k,j,i) = surf_usm_h%surface_types(m)
2533	ENDDO
2534	ELSE
2535	l = idsidx
2536	DO m = 1, surf_usm_v(l)%ns
2537	i = surf_usm_v(l)%i(m)
2538	j = surf_usm_v(l)%j(m)
2539	k = surf_usm_v(l)%k(m)
2540	temp_pf(k,j,i) = surf_usm_v(l)%surface_types(m)
2541	ENDDO
2542	ENDIF
2543
2544	CASE ( 'usm_surfalb' )
2545	!-- surface albedo, weighted average
2546	IF ( idsint == iup_u ) THEN
2547	DO m = 1, surf_usm_h%ns
2548	i = surf_usm_h%i(m)
2549	j = surf_usm_h%j(m)
2550	k = surf_usm_h%k(m)
2551	temp_pf(k,j,i) = surf_usm_h%frac(ind_veg_wall,m) * &
2552	surf_usm_h%albedo(ind_veg_wall,m) + &
2553	surf_usm_h%frac(ind_pav_green,m) * &
2554	surf_usm_h%albedo(ind_pav_green,m) + &
2555	surf_usm_h%frac(ind_wat_win,m) * &
2556	surf_usm_h%albedo(ind_wat_win,m)
2557	ENDDO
2558	ELSE
2559	l = idsidx
2560	DO m = 1, surf_usm_v(l)%ns
2561	i = surf_usm_v(l)%i(m)
2562	j = surf_usm_v(l)%j(m)
2563	k = surf_usm_v(l)%k(m)
2564	temp_pf(k,j,i) = surf_usm_v(l)%frac(ind_veg_wall,m) * &
2565	surf_usm_v(l)%albedo(ind_veg_wall,m) + &
2566	surf_usm_v(l)%frac(ind_pav_green,m) * &
2567	surf_usm_v(l)%albedo(ind_pav_green,m) + &
2568	surf_usm_v(l)%frac(ind_wat_win,m) * &
2569	surf_usm_v(l)%albedo(ind_wat_win,m)
2570	ENDDO
2571	ENDIF
2572
2573	CASE ( 'usm_surfemis' )
2574	!-- surface emissivity, weighted average
2575	IF ( idsint == iup_u ) THEN
2576	DO m = 1, surf_usm_h%ns
2577	i = surf_usm_h%i(m)
2578	j = surf_usm_h%j(m)
2579	k = surf_usm_h%k(m)
2580	temp_pf(k,j,i) = surf_usm_h%frac(ind_veg_wall,m) * &
2581	surf_usm_h%emissivity(ind_veg_wall,m) + &
2582	surf_usm_h%frac(ind_pav_green,m) * &
2583	surf_usm_h%emissivity(ind_pav_green,m) + &
2584	surf_usm_h%frac(ind_wat_win,m) * &
2585	surf_usm_h%emissivity(ind_wat_win,m)
2586	ENDDO
2587	ELSE
2588	l = idsidx
2589	DO m = 1, surf_usm_v(l)%ns
2590	i = surf_usm_v(l)%i(m)
2591	j = surf_usm_v(l)%j(m)
2592	k = surf_usm_v(l)%k(m)
2593	temp_pf(k,j,i) = surf_usm_v(l)%frac(ind_veg_wall,m) *&
2594	surf_usm_v(l)%emissivity(ind_veg_wall,m) +&
2595	surf_usm_v(l)%frac(ind_pav_green,m) *&
2596	surf_usm_v(l)%emissivity(ind_pav_green,m)+&
2597	surf_usm_v(l)%frac(ind_wat_win,m) *&
2598	surf_usm_v(l)%emissivity(ind_wat_win,m)
2599	ENDDO
2600	ENDIF
2601
2602	CASE ( 'usm_surfwintrans' )
2603	!-- transmissivity window tiles
2604	IF ( idsint == iup_u ) THEN
2605	DO m = 1, surf_usm_h%ns
2606	i = surf_usm_h%i(m)
2607	j = surf_usm_h%j(m)
2608	k = surf_usm_h%k(m)
2609	temp_pf(k,j,i) = surf_usm_h%transmissivity(m)
2610	ENDDO
2611	ELSE
2612	l = idsidx
2613	DO m = 1, surf_usm_v(l)%ns
2614	i = surf_usm_v(l)%i(m)
2615	j = surf_usm_v(l)%j(m)
2616	k = surf_usm_v(l)%k(m)
2617	temp_pf(k,j,i) = surf_usm_v(l)%transmissivity(m)
2618	ENDDO
2619	ENDIF
2620
2621	CASE ( 'usm_skyvf' )
2622	!-- sky view factor
2623	DO isurf = dirstart(ids), dirend(ids)
2624	IF ( surfl(id,isurf) == idsint ) THEN
2625	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = skyvf(isurf)
2626	ENDIF
2627	ENDDO
2628
2629	CASE ( 'usm_skyvft' )
2630	!-- sky view factor
2631	DO isurf = dirstart(ids), dirend(ids)
2632	IF ( surfl(id,isurf) == ids ) THEN
2633	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = skyvft(isurf)
2634	ENDIF
2635	ENDDO
2636
2637	!
2638	!-- Not adjusted so far
2639	CASE ( 'usm_svf', 'usm_dif' )
2640	!-- shape view factors or iradiance factors to selected surface
2641	IF ( TRIM(var)=='usm_svf' ) THEN
2642	k = 1
2643	ELSE
2644	k = 2
2645	ENDIF
2646	DO isvf = 1, nsvfl
2647	isurflt = svfsurf(1, isvf)
2648	isurfs = svfsurf(2, isvf)
2649
2650	IF ( surf(ix,isurfs) == is .AND. surf(iy,isurfs) == js .AND. &
2651	surf(iz,isurfs) == ks .AND. surf(id,isurfs) == idsint ) THEN
2652	!-- correct source surface
2653	temp_pf(surfl(iz,isurflt),surfl(iy,isurflt),surfl(ix,isurflt)) = svf(k,isvf)
2654	ENDIF
2655	ENDDO
2656
2657	CASE ( 'usm_rad_net' )
2658	!-- array of complete radiation balance
2659	IF ( av == 0 ) THEN
2660	IF ( idsint == iup_u ) THEN
2661	DO m = 1, surf_usm_h%ns
2662	i = surf_usm_h%i(m)
2663	j = surf_usm_h%j(m)
2664	k = surf_usm_h%k(m)
2665	temp_pf(k,j,i) = surf_usm_h%rad_net_l(m)
2666	ENDDO
2667	ELSE
2668	l = idsidx
2669	DO m = 1, surf_usm_v(l)%ns
2670	i = surf_usm_v(l)%i(m)
2671	j = surf_usm_v(l)%j(m)
2672	k = surf_usm_v(l)%k(m)
2673	temp_pf(k,j,i) = surf_usm_v(l)%rad_net_l(m)
2674	ENDDO
2675	ENDIF
2676	ELSE
2677	IF ( idsint == iup_u ) THEN
2678	DO m = 1, surf_usm_h%ns
2679	i = surf_usm_h%i(m)
2680	j = surf_usm_h%j(m)
2681	k = surf_usm_h%k(m)
2682	temp_pf(k,j,i) = surf_usm_h%rad_net_av(m)
2683	ENDDO
2684	ELSE
2685	l = idsidx
2686	DO m = 1, surf_usm_v(l)%ns
2687	i = surf_usm_v(l)%i(m)
2688	j = surf_usm_v(l)%j(m)
2689	k = surf_usm_v(l)%k(m)
2690	temp_pf(k,j,i) = surf_usm_v(l)%rad_net_av(m)
2691	ENDDO
2692	ENDIF
2693	ENDIF
2694
2695	CASE ( 'usm_rad_insw' )
2696	!-- array of sw radiation falling to surface after i-th reflection
2697	DO isurf = dirstart(ids), dirend(ids)
2698	IF ( surfl(id,isurf) == idsint ) THEN
2699	IF ( av == 0 ) THEN
2700	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinsw(isurf)
2701	ELSE
2702	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinsw_av(isurf)
2703	ENDIF
2704	ENDIF
2705	ENDDO
2706
2707	CASE ( 'usm_rad_inlw' )
2708	!-- array of lw radiation falling to surface after i-th reflection
2709	DO isurf = dirstart(ids), dirend(ids)
2710	IF ( surfl(id,isurf) == idsint ) THEN
2711	IF ( av == 0 ) THEN
2712	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw(isurf)
2713	ELSE
2714	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw_av(isurf)
2715	ENDIF
2716	ENDIF
2717	ENDDO
2718
2719	CASE ( 'usm_rad_inswdir' )
2720	!-- array of direct sw radiation falling to surface from sun
2721	DO isurf = dirstart(ids), dirend(ids)
2722	IF ( surfl(id,isurf) == idsint ) THEN
2723	IF ( av == 0 ) THEN
2724	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdir(isurf)
2725	ELSE
2726	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdir_av(isurf)
2727	ENDIF
2728	ENDIF
2729	ENDDO
2730
2731	CASE ( 'usm_rad_inswdif' )
2732	!-- array of difusion sw radiation falling to surface from sky and borders of the domain
2733	DO isurf = dirstart(ids), dirend(ids)
2734	IF ( surfl(id,isurf) == idsint ) THEN
2735	IF ( av == 0 ) THEN
2736	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdif(isurf)
2737	ELSE
2738	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdif_av(isurf)
2739	ENDIF
2740	ENDIF
2741	ENDDO
2742
2743	CASE ( 'usm_rad_inswref' )
2744	!-- array of sw radiation falling to surface from reflections
2745	DO isurf = dirstart(ids), dirend(ids)
2746	IF ( surfl(id,isurf) == idsint ) THEN
2747	IF ( av == 0 ) THEN
2748	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = &
2749	surfinsw(isurf) - surfinswdir(isurf) - surfinswdif(isurf)
2750	ELSE
2751	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswref_av(isurf)
2752	ENDIF
2753	ENDIF
2754	ENDDO
2755
2756	CASE ( 'usm_rad_inlwdif' )
2757	!-- array of difusion lw radiation falling to surface from sky and borders of the domain
2758	DO isurf = dirstart(ids), dirend(ids)
2759	IF ( surfl(id,isurf) == idsint ) THEN
2760	IF ( av == 0 ) THEN
2761	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwdif(isurf)
2762	ELSE
2763	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwdif_av(isurf)
2764	ENDIF
2765	ENDIF
2766	ENDDO
2767
2768	CASE ( 'usm_rad_inlwref' )
2769	!-- array of lw radiation falling to surface from reflections
2770	DO isurf = dirstart(ids), dirend(ids)
2771	IF ( surfl(id,isurf) == idsint ) THEN
2772	IF ( av == 0 ) THEN
2773	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw(isurf) - surfinlwdif(isurf)
2774	ELSE
2775	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwref_av(isurf)
2776	ENDIF
2777	ENDIF
2778	ENDDO
2779
2780	CASE ( 'usm_rad_outsw' )
2781	!-- array of sw radiation emitted from surface after i-th reflection
2782	DO isurf = dirstart(ids), dirend(ids)
2783	IF ( surfl(id,isurf) == idsint ) THEN
2784	IF ( av == 0 ) THEN
2785	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutsw(isurf)
2786	ELSE
2787	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutsw_av(isurf)
2788	ENDIF
2789	ENDIF
2790	ENDDO
2791
2792	CASE ( 'usm_rad_outlw' )
2793	!-- array of lw radiation emitted from surface after i-th reflection
2794	DO isurf = dirstart(ids), dirend(ids)
2795	IF ( surfl(id,isurf) == idsint ) THEN
2796	IF ( av == 0 ) THEN
2797	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutlw(isurf)
2798	ELSE
2799	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutlw_av(isurf)
2800	ENDIF
2801	ENDIF
2802	ENDDO
2803
2804	CASE ( 'usm_rad_ressw' )
2805	!-- average of array of residua of sw radiation absorbed in surface after last reflection
2806	DO isurf = dirstart(ids), dirend(ids)
2807	IF ( surfl(id,isurf) == idsint ) THEN
2808	IF ( av == 0 ) THEN
2809	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfins(isurf)
2810	ELSE
2811	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfins_av(isurf)
2812	ENDIF
2813	ENDIF
2814	ENDDO
2815
2816	CASE ( 'usm_rad_reslw' )
2817	!-- average of array of residua of lw radiation absorbed in surface after last reflection
2818	DO isurf = dirstart(ids), dirend(ids)
2819	IF ( surfl(id,isurf) == idsint ) THEN
2820	IF ( av == 0 ) THEN
2821	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinl(isurf)
2822	ELSE
2823	temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinl_av(isurf)
2824	ENDIF
2825	ENDIF
2826	ENDDO
2827
2828	CASE ( 'usm_rad_hf' )
2829	!-- array of heat flux from radiation for surfaces after all reflections
2830	IF ( av == 0 ) THEN
2831	IF ( idsint == iup_u ) THEN
2832	DO m = 1, surf_usm_h%ns
2833	i = surf_usm_h%i(m)
2834	j = surf_usm_h%j(m)
2835	k = surf_usm_h%k(m)
2836	temp_pf(k,j,i) = surf_usm_h%surfhf(m)
2837	ENDDO
2838	ELSE
2839	l = idsidx
2840	DO m = 1, surf_usm_v(l)%ns
2841	i = surf_usm_v(l)%i(m)
2842	j = surf_usm_v(l)%j(m)
2843	k = surf_usm_v(l)%k(m)
2844	temp_pf(k,j,i) = surf_usm_v(l)%surfhf(m)
2845	ENDDO
2846	ENDIF
2847	ELSE
2848	IF ( idsint == iup_u ) THEN
2849	DO m = 1, surf_usm_h%ns
2850	i = surf_usm_h%i(m)
2851	j = surf_usm_h%j(m)
2852	k = surf_usm_h%k(m)
2853	temp_pf(k,j,i) = surf_usm_h%surfhf_av(m)
2854	ENDDO
2855	ELSE
2856	l = idsidx
2857	DO m = 1, surf_usm_v(l)%ns
2858	i = surf_usm_v(l)%i(m)
2859	j = surf_usm_v(l)%j(m)
2860	k = surf_usm_v(l)%k(m)
2861	temp_pf(k,j,i) = surf_usm_v(l)%surfhf_av(m)
2862	ENDDO
2863	ENDIF
2864	ENDIF
2865
2866	CASE ( 'usm_wshf' )
2867	!-- array of sensible heat flux from surfaces
2868	IF ( av == 0 ) THEN
2869	IF ( idsint == iup_u ) THEN
2870	DO m = 1, surf_usm_h%ns
2871	i = surf_usm_h%i(m)
2872	j = surf_usm_h%j(m)
2873	k = surf_usm_h%k(m)
2874	temp_pf(k,j,i) = surf_usm_h%wshf_eb(m)
2875	ENDDO
2876	ELSE
2877	l = idsidx
2878	DO m = 1, surf_usm_v(l)%ns
2879	i = surf_usm_v(l)%i(m)
2880	j = surf_usm_v(l)%j(m)
2881	k = surf_usm_v(l)%k(m)
2882	temp_pf(k,j,i) = surf_usm_v(l)%wshf_eb(m)
2883	ENDDO
2884	ENDIF
2885	ELSE
2886	IF ( idsint == iup_u ) THEN
2887	DO m = 1, surf_usm_h%ns
2888	i = surf_usm_h%i(m)
2889	j = surf_usm_h%j(m)
2890	k = surf_usm_h%k(m)
2891	temp_pf(k,j,i) = surf_usm_h%wshf_eb_av(m)
2892	ENDDO
2893	ELSE
2894	l = idsidx
2895	DO m = 1, surf_usm_v(l)%ns
2896	i = surf_usm_v(l)%i(m)
2897	j = surf_usm_v(l)%j(m)
2898	k = surf_usm_v(l)%k(m)
2899	temp_pf(k,j,i) = surf_usm_v(l)%wshf_eb_av(m)
2900	ENDDO
2901	ENDIF
2902	ENDIF
2903
2904
2905	CASE ( 'usm_wghf' )
2906	!-- array of heat flux from ground (land, wall, roof)
2907	IF ( av == 0 ) THEN
2908	IF ( idsint == iup_u ) THEN
2909	DO m = 1, surf_usm_h%ns
2910	i = surf_usm_h%i(m)
2911	j = surf_usm_h%j(m)
2912	k = surf_usm_h%k(m)
2913	temp_pf(k,j,i) = surf_usm_h%wghf_eb(m)
2914	ENDDO
2915	ELSE
2916	l = idsidx
2917	DO m = 1, surf_usm_v(l)%ns
2918	i = surf_usm_v(l)%i(m)
2919	j = surf_usm_v(l)%j(m)
2920	k = surf_usm_v(l)%k(m)
2921	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb(m)
2922	ENDDO
2923	ENDIF
2924	ELSE
2925	IF ( idsint == iup_u ) THEN
2926	DO m = 1, surf_usm_h%ns
2927	i = surf_usm_h%i(m)
2928	j = surf_usm_h%j(m)
2929	k = surf_usm_h%k(m)
2930	temp_pf(k,j,i) = surf_usm_h%wghf_eb_av(m)
2931	ENDDO
2932	ELSE
2933	l = idsidx
2934	DO m = 1, surf_usm_v(l)%ns
2935	i = surf_usm_v(l)%i(m)
2936	j = surf_usm_v(l)%j(m)
2937	k = surf_usm_v(l)%k(m)
2938	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_av(m)
2939	ENDDO
2940	ENDIF
2941	ENDIF
2942
2943	CASE ( 'usm_wghf_window' )
2944	!-- array of heat flux from window ground (land, wall, roof)
2945
2946	IF ( av == 0 ) THEN
2947	IF ( idsint == iup_u ) THEN
2948	DO m = 1, surf_usm_h%ns
2949	i = surf_usm_h%i(m)
2950	j = surf_usm_h%j(m)
2951	k = surf_usm_h%k(m)
2952	temp_pf(k,j,i) = surf_usm_h%wghf_eb_window(m)
2953	ENDDO
2954	ELSE
2955	l = idsidx
2956	DO m = 1, surf_usm_v(l)%ns
2957	i = surf_usm_v(l)%i(m)
2958	j = surf_usm_v(l)%j(m)
2959	k = surf_usm_v(l)%k(m)
2960	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_window(m)
2961	ENDDO
2962	ENDIF
2963	ELSE
2964	IF ( idsint == iup_u ) THEN
2965	DO m = 1, surf_usm_h%ns
2966	i = surf_usm_h%i(m)
2967	j = surf_usm_h%j(m)
2968	k = surf_usm_h%k(m)
2969	temp_pf(k,j,i) = surf_usm_h%wghf_eb_window_av(m)
2970	ENDDO
2971	ELSE
2972	l = idsidx
2973	DO m = 1, surf_usm_v(l)%ns
2974	i = surf_usm_v(l)%i(m)
2975	j = surf_usm_v(l)%j(m)
2976	k = surf_usm_v(l)%k(m)
2977	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_window_av(m)
2978	ENDDO
2979	ENDIF
2980	ENDIF
2981
2982	CASE ( 'usm_wghf_green' )
2983	!-- array of heat flux from green ground (land, wall, roof)
2984
2985	IF ( av == 0 ) THEN
2986	IF ( idsint == iup_u ) THEN
2987	DO m = 1, surf_usm_h%ns
2988	i = surf_usm_h%i(m)
2989	j = surf_usm_h%j(m)
2990	k = surf_usm_h%k(m)
2991	temp_pf(k,j,i) = surf_usm_h%wghf_eb_green(m)
2992	ENDDO
2993	ELSE
2994	l = idsidx
2995	DO m = 1, surf_usm_v(l)%ns
2996	i = surf_usm_v(l)%i(m)
2997	j = surf_usm_v(l)%j(m)
2998	k = surf_usm_v(l)%k(m)
2999	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_green(m)
3000	ENDDO
3001	ENDIF
3002	ELSE
3003	IF ( idsint == iup_u ) THEN
3004	DO m = 1, surf_usm_h%ns
3005	i = surf_usm_h%i(m)
3006	j = surf_usm_h%j(m)
3007	k = surf_usm_h%k(m)
3008	temp_pf(k,j,i) = surf_usm_h%wghf_eb_green_av(m)
3009	ENDDO
3010	ELSE
3011	l = idsidx
3012	DO m = 1, surf_usm_v(l)%ns
3013	i = surf_usm_v(l)%i(m)
3014	j = surf_usm_v(l)%j(m)
3015	k = surf_usm_v(l)%k(m)
3016	temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_green_av(m)
3017	ENDDO
3018	ENDIF
3019	ENDIF
3020
3021	CASE ( 'usm_iwghf' )
3022	!-- array of heat flux from indoor ground (land, wall, roof)
3023	IF ( av == 0 ) THEN
3024	IF ( idsint == iup_u ) THEN
3025	DO m = 1, surf_usm_h%ns
3026	i = surf_usm_h%i(m)
3027	j = surf_usm_h%j(m)
3028	k = surf_usm_h%k(m)
3029	temp_pf(k,j,i) = surf_usm_h%iwghf_eb(m)
3030	ENDDO
3031	ELSE
3032	l = idsidx
3033	DO m = 1, surf_usm_v(l)%ns
3034	i = surf_usm_v(l)%i(m)
3035	j = surf_usm_v(l)%j(m)
3036	k = surf_usm_v(l)%k(m)
3037	temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb(m)
3038	ENDDO
3039	ENDIF
3040	ELSE
3041	IF ( idsint == iup_u ) THEN
3042	DO m = 1, surf_usm_h%ns
3043	i = surf_usm_h%i(m)
3044	j = surf_usm_h%j(m)
3045	k = surf_usm_h%k(m)
3046	temp_pf(k,j,i) = surf_usm_h%iwghf_eb_av(m)
3047	ENDDO
3048	ELSE
3049	l = idsidx
3050	DO m = 1, surf_usm_v(l)%ns
3051	i = surf_usm_v(l)%i(m)
3052	j = surf_usm_v(l)%j(m)
3053	k = surf_usm_v(l)%k(m)
3054	temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_av(m)
3055	ENDDO
3056	ENDIF
3057	ENDIF
3058
3059	CASE ( 'usm_iwghf_window' )
3060	!-- array of heat flux from indoor window ground (land, wall, roof)
3061
3062	IF ( av == 0 ) THEN
3063	IF ( idsint == iup_u ) THEN
3064	DO m = 1, surf_usm_h%ns
3065	i = surf_usm_h%i(m)
3066	j = surf_usm_h%j(m)
3067	k = surf_usm_h%k(m)
3068	temp_pf(k,j,i) = surf_usm_h%iwghf_eb_window(m)
3069	ENDDO
3070	ELSE
3071	l = idsidx
3072	DO m = 1, surf_usm_v(l)%ns
3073	i = surf_usm_v(l)%i(m)
3074	j = surf_usm_v(l)%j(m)
3075	k = surf_usm_v(l)%k(m)
3076	temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_window(m)
3077	ENDDO
3078	ENDIF
3079	ELSE
3080	IF ( idsint == iup_u ) THEN
3081	DO m = 1, surf_usm_h%ns
3082	i = surf_usm_h%i(m)
3083	j = surf_usm_h%j(m)
3084	k = surf_usm_h%k(m)
3085	temp_pf(k,j,i) = surf_usm_h%iwghf_eb_window_av(m)
3086	ENDDO
3087	ELSE
3088	l = idsidx
3089	DO m = 1, surf_usm_v(l)%ns
3090	i = surf_usm_v(l)%i(m)
3091	j = surf_usm_v(l)%j(m)
3092	k = surf_usm_v(l)%k(m)
3093	temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_window_av(m)
3094	ENDDO
3095	ENDIF
3096	ENDIF
3097
3098	CASE ( 'usm_t_surf' )
3099	!-- surface temperature for surfaces
3100	IF ( av == 0 ) THEN
3101	IF ( idsint == iup_u ) THEN
3102	DO m = 1, surf_usm_h%ns
3103	i = surf_usm_h%i(m)
3104	j = surf_usm_h%j(m)
3105	k = surf_usm_h%k(m)
3106	temp_pf(k,j,i) = t_surf_h(m)
3107	ENDDO
3108	ELSE
3109	l = idsidx
3110	DO m = 1, surf_usm_v(l)%ns
3111	i = surf_usm_v(l)%i(m)
3112	j = surf_usm_v(l)%j(m)
3113	k = surf_usm_v(l)%k(m)
3114	temp_pf(k,j,i) = t_surf_v(l)%t(m)
3115	ENDDO
3116	ENDIF
3117	ELSE
3118	IF ( idsint == iup_u ) THEN
3119	DO m = 1, surf_usm_h%ns
3120	i = surf_usm_h%i(m)
3121	j = surf_usm_h%j(m)
3122	k = surf_usm_h%k(m)
3123	temp_pf(k,j,i) = surf_usm_h%t_surf_av(m)
3124	ENDDO
3125	ELSE
3126	l = idsidx
3127	DO m = 1, surf_usm_v(l)%ns
3128	i = surf_usm_v(l)%i(m)
3129	j = surf_usm_v(l)%j(m)
3130	k = surf_usm_v(l)%k(m)
3131	temp_pf(k,j,i) = surf_usm_v(l)%t_surf_av(m)
3132	ENDDO
3133	ENDIF
3134	ENDIF
3135
3136	CASE ( 'usm_t_surf_window' )
3137	!-- surface temperature for window surfaces
3138
3139	IF ( av == 0 ) THEN
3140	IF ( idsint == iup_u ) THEN
3141	DO m = 1, surf_usm_h%ns
3142	i = surf_usm_h%i(m)
3143	j = surf_usm_h%j(m)
3144	k = surf_usm_h%k(m)
3145	temp_pf(k,j,i) = t_surf_window_h(m)
3146	ENDDO
3147	ELSE
3148	l = idsidx
3149	DO m = 1, surf_usm_v(l)%ns
3150	i = surf_usm_v(l)%i(m)
3151	j = surf_usm_v(l)%j(m)
3152	k = surf_usm_v(l)%k(m)
3153	temp_pf(k,j,i) = t_surf_window_v(l)%t(m)
3154	ENDDO
3155	ENDIF
3156
3157	ELSE
3158	IF ( idsint == iup_u ) THEN
3159	DO m = 1, surf_usm_h%ns
3160	i = surf_usm_h%i(m)
3161	j = surf_usm_h%j(m)
3162	k = surf_usm_h%k(m)
3163	temp_pf(k,j,i) = surf_usm_h%t_surf_window_av(m)
3164	ENDDO
3165	ELSE
3166	l = idsidx
3167	DO m = 1, surf_usm_v(l)%ns
3168	i = surf_usm_v(l)%i(m)
3169	j = surf_usm_v(l)%j(m)
3170	k = surf_usm_v(l)%k(m)
3171	temp_pf(k,j,i) = surf_usm_v(l)%t_surf_window_av(m)
3172	ENDDO
3173
3174	ENDIF
3175
3176	ENDIF
3177
3178	CASE ( 'usm_t_surf_green' )
3179	!-- surface temperature for green surfaces
3180
3181	IF ( av == 0 ) THEN
3182	IF ( idsint == iup_u ) THEN
3183	DO m = 1, surf_usm_h%ns
3184	i = surf_usm_h%i(m)
3185	j = surf_usm_h%j(m)
3186	k = surf_usm_h%k(m)
3187	temp_pf(k,j,i) = t_surf_green_h(m)
3188	ENDDO
3189	ELSE
3190	l = idsidx
3191	DO m = 1, surf_usm_v(l)%ns
3192	i = surf_usm_v(l)%i(m)
3193	j = surf_usm_v(l)%j(m)
3194	k = surf_usm_v(l)%k(m)
3195	temp_pf(k,j,i) = t_surf_green_v(l)%t(m)
3196	ENDDO
3197	ENDIF
3198
3199	ELSE
3200	IF ( idsint == iup_u ) THEN
3201	DO m = 1, surf_usm_h%ns
3202	i = surf_usm_h%i(m)
3203	j = surf_usm_h%j(m)
3204	k = surf_usm_h%k(m)
3205	temp_pf(k,j,i) = surf_usm_h%t_surf_green_av(m)
3206	ENDDO
3207	ELSE
3208	l = idsidx
3209	DO m = 1, surf_usm_v(l)%ns
3210	i = surf_usm_v(l)%i(m)
3211	j = surf_usm_v(l)%j(m)
3212	k = surf_usm_v(l)%k(m)
3213	temp_pf(k,j,i) = surf_usm_v(l)%t_surf_green_av(m)
3214	ENDDO
3215
3216	ENDIF
3217
3218	ENDIF
3219
3220	CASE ( 'usm_t_surf_10cm' )
3221	!-- near surface temperature for whole surfaces
3222
3223	IF ( av == 0 ) THEN
3224	IF ( idsint == iup_u ) THEN
3225	DO m = 1, surf_usm_h%ns
3226	i = surf_usm_h%i(m)
3227	j = surf_usm_h%j(m)
3228	k = surf_usm_h%k(m)
3229	temp_pf(k,j,i) = t_surf_10cm_h(m)
3230	ENDDO
3231	ELSE
3232	l = idsidx
3233	DO m = 1, surf_usm_v(l)%ns
3234	i = surf_usm_v(l)%i(m)
3235	j = surf_usm_v(l)%j(m)
3236	k = surf_usm_v(l)%k(m)
3237	temp_pf(k,j,i) = t_surf_10cm_v(l)%t(m)
3238	ENDDO
3239	ENDIF
3240
3241	ELSE
3242	IF ( idsint == iup_u ) THEN
3243	DO m = 1, surf_usm_h%ns
3244	i = surf_usm_h%i(m)
3245	j = surf_usm_h%j(m)
3246	k = surf_usm_h%k(m)
3247	temp_pf(k,j,i) = surf_usm_h%t_surf_10cm_av(m)
3248	ENDDO
3249	ELSE
3250	l = idsidx
3251	DO m = 1, surf_usm_v(l)%ns
3252	i = surf_usm_v(l)%i(m)
3253	j = surf_usm_v(l)%j(m)
3254	k = surf_usm_v(l)%k(m)
3255	temp_pf(k,j,i) = surf_usm_v(l)%t_surf_10cm_av(m)
3256	ENDDO
3257
3258	ENDIF
3259
3260	ENDIF
3261
3262
3263	CASE ( 'usm_t_wall' )
3264	!-- wall temperature for iwl layer of walls and land
3265	IF ( av == 0 ) THEN
3266	IF ( idsint == iup_u ) THEN
3267	DO m = 1, surf_usm_h%ns
3268	i = surf_usm_h%i(m)
3269	j = surf_usm_h%j(m)
3270	k = surf_usm_h%k(m)
3271	temp_pf(k,j,i) = t_wall_h(iwl,m)
3272	ENDDO
3273	ELSE
3274	l = idsidx
3275	DO m = 1, surf_usm_v(l)%ns
3276	i = surf_usm_v(l)%i(m)
3277	j = surf_usm_v(l)%j(m)
3278	k = surf_usm_v(l)%k(m)
3279	temp_pf(k,j,i) = t_wall_v(l)%t(iwl,m)
3280	ENDDO
3281	ENDIF
3282	ELSE
3283	IF ( idsint == iup_u ) THEN
3284	DO m = 1, surf_usm_h%ns
3285	i = surf_usm_h%i(m)
3286	j = surf_usm_h%j(m)
3287	k = surf_usm_h%k(m)
3288	temp_pf(k,j,i) = surf_usm_h%t_wall_av(iwl,m)
3289	ENDDO
3290	ELSE
3291	l = idsidx
3292	DO m = 1, surf_usm_v(l)%ns
3293	i = surf_usm_v(l)%i(m)
3294	j = surf_usm_v(l)%j(m)
3295	k = surf_usm_v(l)%k(m)
3296	temp_pf(k,j,i) = surf_usm_v(l)%t_wall_av(iwl,m)
3297	ENDDO
3298	ENDIF
3299	ENDIF
3300
3301	CASE ( 'usm_t_window' )
3302	!-- window temperature for iwl layer of walls and land
3303	IF ( av == 0 ) THEN
3304	IF ( idsint == iup_u ) THEN
3305	DO m = 1, surf_usm_h%ns
3306	i = surf_usm_h%i(m)
3307	j = surf_usm_h%j(m)
3308	k = surf_usm_h%k(m)
3309	temp_pf(k,j,i) = t_window_h(iwl,m)
3310	ENDDO
3311	ELSE
3312	l = idsidx
3313	DO m = 1, surf_usm_v(l)%ns
3314	i = surf_usm_v(l)%i(m)
3315	j = surf_usm_v(l)%j(m)
3316	k = surf_usm_v(l)%k(m)
3317	temp_pf(k,j,i) = t_window_v(l)%t(iwl,m)
3318	ENDDO
3319	ENDIF
3320	ELSE
3321	IF ( idsint == iup_u ) THEN
3322	DO m = 1, surf_usm_h%ns
3323	i = surf_usm_h%i(m)
3324	j = surf_usm_h%j(m)
3325	k = surf_usm_h%k(m)
3326	temp_pf(k,j,i) = surf_usm_h%t_window_av(iwl,m)
3327	ENDDO
3328	ELSE
3329	l = idsidx
3330	DO m = 1, surf_usm_v(l)%ns
3331	i = surf_usm_v(l)%i(m)
3332	j = surf_usm_v(l)%j(m)
3333	k = surf_usm_v(l)%k(m)
3334	temp_pf(k,j,i) = surf_usm_v(l)%t_window_av(iwl,m)
3335	ENDDO
3336	ENDIF
3337	ENDIF
3338
3339	CASE ( 'usm_t_green' )
3340	!-- green temperature for iwl layer of walls and land
3341	IF ( av == 0 ) THEN
3342	IF ( idsint == iup_u ) THEN
3343	DO m = 1, surf_usm_h%ns
3344	i = surf_usm_h%i(m)
3345	j = surf_usm_h%j(m)
3346	k = surf_usm_h%k(m)
3347	temp_pf(k,j,i) = t_green_h(iwl,m)
3348	ENDDO
3349	ELSE
3350	l = idsidx
3351	DO m = 1, surf_usm_v(l)%ns
3352	i = surf_usm_v(l)%i(m)
3353	j = surf_usm_v(l)%j(m)
3354	k = surf_usm_v(l)%k(m)
3355	temp_pf(k,j,i) = t_green_v(l)%t(iwl,m)
3356	ENDDO
3357	ENDIF
3358	ELSE
3359	IF ( idsint == iup_u ) THEN
3360	DO m = 1, surf_usm_h%ns
3361	i = surf_usm_h%i(m)
3362	j = surf_usm_h%j(m)
3363	k = surf_usm_h%k(m)
3364	temp_pf(k,j,i) = surf_usm_h%t_green_av(iwl,m)
3365	ENDDO
3366	ELSE
3367	l = idsidx
3368	DO m = 1, surf_usm_v(l)%ns
3369	i = surf_usm_v(l)%i(m)
3370	j = surf_usm_v(l)%j(m)
3371	k = surf_usm_v(l)%k(m)
3372	temp_pf(k,j,i) = surf_usm_v(l)%t_green_av(iwl,m)
3373	ENDDO
3374	ENDIF
3375	ENDIF
3376
3377
3378	CASE DEFAULT
3379	found = .FALSE.
3380
3381	END SELECT
3382
3383	!
3384	!-- Rearrange dimensions for NetCDF output
3385	!-- FIXME: this may generate FPE overflow upon conversion from DP to SP
3386	DO j = nys, nyn
3387	DO i = nxl, nxr
3388	DO k = nzb_do, nzt_do
3389	local_pf(i,j,k) = temp_pf(k,j,i)
3390	ENDDO
3391	ENDDO
3392	ENDDO
3393
3394	END SUBROUTINE usm_data_output_3d
3395
3396
3397	!------------------------------------------------------------------------------!
3398	!
3399	! Description:
3400	! ------------
3401	!> Soubroutine defines appropriate grid for netcdf variables.
3402	!> It is called out from subroutine netcdf.
3403	!------------------------------------------------------------------------------!
3404	SUBROUTINE usm_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
3405
3406	IMPLICIT NONE
3407
3408	CHARACTER (len=*), INTENT(IN) :: variable !<
3409	LOGICAL, INTENT(OUT) :: found !<
3410	CHARACTER (len=*), INTENT(OUT) :: grid_x !<
3411	CHARACTER (len=*), INTENT(OUT) :: grid_y !<
3412	CHARACTER (len=*), INTENT(OUT) :: grid_z !<
3413
3414	CHARACTER (len=varnamelength) :: var
3415
3416	var = TRIM(variable)
3417	IF ( var(1:12) == 'usm_rad_net_' .OR. var(1:13) == 'usm_rad_insw_' .OR. &
3418	var(1:13) == 'usm_rad_inlw_' .OR. var(1:16) == 'usm_rad_inswdir_' .OR. &
3419	var(1:16) == 'usm_rad_inswdif_' .OR. var(1:16) == 'usm_rad_inswref_' .OR. &
3420	var(1:16) == 'usm_rad_inlwdif_' .OR. var(1:16) == 'usm_rad_inlwref_' .OR. &
3421	var(1:14) == 'usm_rad_outsw_' .OR. var(1:14) == 'usm_rad_outlw_' .OR. &
3422	var(1:14) == 'usm_rad_ressw_' .OR. var(1:14) == 'usm_rad_reslw_' .OR. &
3423	var(1:11) == 'usm_rad_hf_' .OR. &
3424	var(1:9) == 'usm_wshf_' .OR. var(1:9) == 'usm_wghf_' .OR. &
3425	var(1:16) == 'usm_wghf_window_' .OR. var(1:15) == 'usm_wghf_green_' .OR. &
3426	var(1:10) == 'usm_iwghf_' .OR. var(1:17) == 'usm_iwghf_window_' .OR. &
3427	var(1:10) == 'usm_t_surf' .OR. var(1:10) == 'usm_t_wall' .OR. &
3428	var(1:17) == 'usm_t_surf_window' .OR. var(1:12) == 'usm_t_window' .OR. &
3429	var(1:16) == 'usm_t_surf_green' .OR. &
3430	var(1:15) == 'usm_t_surf_10cm' .OR. &
3431	var(1:9) == 'usm_surfz' .OR. var(1:7) == 'usm_svf' .OR. &
3432	var(1:7) == 'usm_dif' .OR. var(1:11) == 'usm_surfcat' .OR. &
3433	var(1:11) == 'usm_surfalb' .OR. var(1:12) == 'usm_surfemis' .OR. &
3434	var(1:16) == 'usm_surfwintrans' .OR. &
3435	var(1:9) == 'usm_skyvf' .OR. var(1:9) == 'usm_skyvft' ) THEN
3436
3437	found = .TRUE.
3438	grid_x = 'x'
3439	grid_y = 'y'
3440	grid_z = 'zu'
3441	ELSE
3442	found = .FALSE.
3443	grid_x = 'none'
3444	grid_y = 'none'
3445	grid_z = 'none'
3446	ENDIF
3447
3448	END SUBROUTINE usm_define_netcdf_grid
3449
3450
3451	!------------------------------------------------------------------------------!
3452	! Description:
3453	! ------------
3454	!> Initialization of the wall surface model
3455	!------------------------------------------------------------------------------!
3456	SUBROUTINE usm_init_material_model
3457
3458	IMPLICIT NONE
3459
3460	INTEGER(iwp) :: k, l, m !< running indices
3461
3462	CALL location_message( ' initialization of wall surface model', .TRUE. )
3463
3464	!-- Calculate wall grid spacings.
3465	!-- Temperature is defined at the center of the wall layers,
3466	!-- whereas gradients/fluxes are defined at the edges (_stag)
3467	!-- apply for all particular surface grids. First for horizontal surfaces
3468	DO m = 1, surf_usm_h%ns
3469
3470	surf_usm_h%dz_wall(nzb_wall,m) = surf_usm_h%zw(nzb_wall,m)
3471	DO k = nzb_wall+1, nzt_wall
3472	surf_usm_h%dz_wall(k,m) = surf_usm_h%zw(k,m) - &
3473	surf_usm_h%zw(k-1,m)
3474	ENDDO
3475	surf_usm_h%dz_window(nzb_wall,m) = surf_usm_h%zw_window(nzb_wall,m)
3476	DO k = nzb_wall+1, nzt_wall
3477	surf_usm_h%dz_window(k,m) = surf_usm_h%zw_window(k,m) - &
3478	surf_usm_h%zw_window(k-1,m)
3479	ENDDO
3480	surf_usm_h%dz_green(nzb_wall,m) = surf_usm_h%zw_green(nzb_wall,m)
3481	DO k = nzb_wall+1, nzt_wall
3482	surf_usm_h%dz_green(k,m) = surf_usm_h%zw_green(k,m) - &
3483	surf_usm_h%zw_green(k-1,m)
3484	ENDDO
3485
3486	surf_usm_h%dz_wall(nzt_wall+1,m) = surf_usm_h%dz_wall(nzt_wall,m)
3487
3488	DO k = nzb_wall, nzt_wall-1
3489	surf_usm_h%dz_wall_stag(k,m) = 0.5 * ( &
3490	surf_usm_h%dz_wall(k+1,m) + surf_usm_h%dz_wall(k,m) )
3491	ENDDO
3492	surf_usm_h%dz_wall_stag(nzt_wall,m) = surf_usm_h%dz_wall(nzt_wall,m)
3493
3494	surf_usm_h%dz_window(nzt_wall+1,m) = surf_usm_h%dz_window(nzt_wall,m)
3495
3496	DO k = nzb_wall, nzt_wall-1
3497	surf_usm_h%dz_window_stag(k,m) = 0.5 * ( &
3498	surf_usm_h%dz_window(k+1,m) + surf_usm_h%dz_window(k,m) )
3499	ENDDO
3500	surf_usm_h%dz_window_stag(nzt_wall,m) = surf_usm_h%dz_window(nzt_wall,m)
3501
3502	surf_usm_h%dz_green(nzt_wall+1,m) = surf_usm_h%dz_green(nzt_wall,m)
3503
3504	DO k = nzb_wall, nzt_wall-1
3505	surf_usm_h%dz_green_stag(k,m) = 0.5 * ( &
3506	surf_usm_h%dz_green(k+1,m) + surf_usm_h%dz_green(k,m) )
3507	ENDDO
3508	surf_usm_h%dz_green_stag(nzt_wall,m) = surf_usm_h%dz_green(nzt_wall,m)
3509	ENDDO
3510	surf_usm_h%ddz_wall = 1.0_wp / surf_usm_h%dz_wall
3511	surf_usm_h%ddz_wall_stag = 1.0_wp / surf_usm_h%dz_wall_stag
3512	surf_usm_h%ddz_window = 1.0_wp / surf_usm_h%dz_window
3513	surf_usm_h%ddz_window_stag = 1.0_wp / surf_usm_h%dz_window_stag
3514	surf_usm_h%ddz_green = 1.0_wp / surf_usm_h%dz_green
3515	surf_usm_h%ddz_green_stag = 1.0_wp / surf_usm_h%dz_green_stag
3516	!
3517	!-- For vertical surfaces
3518	DO l = 0, 3
3519	DO m = 1, surf_usm_v(l)%ns
3520	surf_usm_v(l)%dz_wall(nzb_wall,m) = surf_usm_v(l)%zw(nzb_wall,m)
3521	DO k = nzb_wall+1, nzt_wall
3522	surf_usm_v(l)%dz_wall(k,m) = surf_usm_v(l)%zw(k,m) - &
3523	surf_usm_v(l)%zw(k-1,m)
3524	ENDDO
3525	surf_usm_v(l)%dz_window(nzb_wall,m) = surf_usm_v(l)%zw_window(nzb_wall,m)
3526	DO k = nzb_wall+1, nzt_wall
3527	surf_usm_v(l)%dz_window(k,m) = surf_usm_v(l)%zw_window(k,m) - &
3528	surf_usm_v(l)%zw_window(k-1,m)
3529	ENDDO
3530	surf_usm_v(l)%dz_green(nzb_wall,m) = surf_usm_v(l)%zw_green(nzb_wall,m)
3531	DO k = nzb_wall+1, nzt_wall
3532	surf_usm_v(l)%dz_green(k,m) = surf_usm_v(l)%zw_green(k,m) - &
3533	surf_usm_v(l)%zw_green(k-1,m)
3534	ENDDO
3535
3536	surf_usm_v(l)%dz_wall(nzt_wall+1,m) = &
3537	surf_usm_v(l)%dz_wall(nzt_wall,m)
3538
3539	DO k = nzb_wall, nzt_wall-1
3540	surf_usm_v(l)%dz_wall_stag(k,m) = 0.5 * ( &
3541	surf_usm_v(l)%dz_wall(k+1,m) + &
3542	surf_usm_v(l)%dz_wall(k,m) )
3543	ENDDO
3544	surf_usm_v(l)%dz_wall_stag(nzt_wall,m) = &
3545	surf_usm_v(l)%dz_wall(nzt_wall,m)
3546	surf_usm_v(l)%dz_window(nzt_wall+1,m) = &
3547	surf_usm_v(l)%dz_window(nzt_wall,m)
3548
3549	DO k = nzb_wall, nzt_wall-1
3550	surf_usm_v(l)%dz_window_stag(k,m) = 0.5 * ( &
3551	surf_usm_v(l)%dz_window(k+1,m) + &
3552	surf_usm_v(l)%dz_window(k,m) )
3553	ENDDO
3554	surf_usm_v(l)%dz_window_stag(nzt_wall,m) = &
3555	surf_usm_v(l)%dz_window(nzt_wall,m)
3556	surf_usm_v(l)%dz_green(nzt_wall+1,m) = &
3557	surf_usm_v(l)%dz_green(nzt_wall,m)
3558
3559	DO k = nzb_wall, nzt_wall-1
3560	surf_usm_v(l)%dz_green_stag(k,m) = 0.5 * ( &
3561	surf_usm_v(l)%dz_green(k+1,m) + &
3562	surf_usm_v(l)%dz_green(k,m) )
3563	ENDDO
3564	surf_usm_v(l)%dz_green_stag(nzt_wall,m) = &
3565	surf_usm_v(l)%dz_green(nzt_wall,m)
3566	ENDDO
3567	surf_usm_v(l)%ddz_wall = 1.0_wp / surf_usm_v(l)%dz_wall
3568	surf_usm_v(l)%ddz_wall_stag = 1.0_wp / surf_usm_v(l)%dz_wall_stag
3569	surf_usm_v(l)%ddz_window = 1.0_wp / surf_usm_v(l)%dz_window
3570	surf_usm_v(l)%ddz_window_stag = 1.0_wp / surf_usm_v(l)%dz_window_stag
3571	surf_usm_v(l)%ddz_green = 1.0_wp / surf_usm_v(l)%dz_green
3572	surf_usm_v(l)%ddz_green_stag = 1.0_wp / surf_usm_v(l)%dz_green_stag
3573	ENDDO
3574
3575
3576	CALL location_message( ' wall structures filed out', .TRUE. )
3577
3578	CALL location_message( ' initialization of wall surface model finished', .TRUE. )
3579
3580	END SUBROUTINE usm_init_material_model
3581
3582
3583	!------------------------------------------------------------------------------!
3584	! Description:
3585	! ------------
3586	!> Initialization of the urban surface model
3587	!------------------------------------------------------------------------------!
3588	SUBROUTINE usm_init_urban_surface
3589
3590	USE arrays_3d, &
3591	ONLY: zw
3592
3593	USE netcdf_data_input_mod, &
3594	ONLY: building_pars_f, building_type_f, terrain_height_f
3595
3596	IMPLICIT NONE
3597
3598	INTEGER(iwp) :: i !< loop index x-dirction
3599	INTEGER(iwp) :: ind_emis_wall !< index in input list for wall emissivity
3600	INTEGER(iwp) :: ind_emis_green !< index in input list for green emissivity
3601	INTEGER(iwp) :: ind_emis_win !< index in input list for window emissivity
3602	INTEGER(iwp) :: ind_green_frac_w !< index in input list for green fraction on wall
3603	INTEGER(iwp) :: ind_green_frac_r !< index in input list for green fraction on roof
3604	INTEGER(iwp) :: ind_hc1 !< index in input list for heat capacity at first wall layer
3605	INTEGER(iwp) :: ind_hc2 !< index in input list for heat capacity at second wall layer
3606	INTEGER(iwp) :: ind_hc3 !< index in input list for heat capacity at third wall layer
3607	INTEGER(iwp) :: ind_lai_r !< index in input list for LAI on roof
3608	INTEGER(iwp) :: ind_lai_w !< index in input list for LAI on wall
3609	INTEGER(iwp) :: ind_tc1 !< index in input list for thermal conductivity at first wall layer
3610	INTEGER(iwp) :: ind_tc2 !< index in input list for thermal conductivity at second wall layer
3611	INTEGER(iwp) :: ind_tc3 !< index in input list for thermal conductivity at third wall layer
3612	INTEGER(iwp) :: ind_trans !< index in input list for window transmissivity
3613	INTEGER(iwp) :: ind_wall_frac !< index in input list for wall fraction
3614	INTEGER(iwp) :: ind_win_frac !< index in input list for window fraction
3615	INTEGER(iwp) :: ind_z0 !< index in input list for z0
3616	INTEGER(iwp) :: ind_z0qh !< index in input list for z0h / z0q
3617	INTEGER(iwp) :: j !< loop index y-dirction
3618	INTEGER(iwp) :: k !< loop index z-dirction
3619	INTEGER(iwp) :: l !< loop index surface orientation
3620	INTEGER(iwp) :: m !< loop index surface element
3621	INTEGER(iwp) :: st !< dummy
3622
3623	REAL(wp) :: c, d, tin, twin
3624	REAL(wp) :: ground_floor_level_l !< local height of ground floor level
3625	REAL(wp) :: z_agl !< height above ground
3626	REAL(wp), DIMENSION(nzb:nzt) :: exn !< value of the Exner function in layers
3627
3628	!
3629	!-- NOPOINTER version not implemented yet
3630	#if defined( __nopointer )
3631	message_string = 'The urban surface module only runs with POINTER version'
3632	CALL message( 'urban_surface_mod', 'PA0452', 1, 2, 0, 6, 0 )
3633	#endif
3634
3635	CALL cpu_log( log_point_s(78), 'usm_init', 'start' )
3636	!-- surface forcing have to be disabled for LSF
3637	!-- in case of enabled urban surface module
3638	IF ( large_scale_forcing ) THEN
3639	lsf_surf = .FALSE.
3640	ENDIF
3641
3642	!
3643	!-- Flag surface elements belonging to the ground floor level. Therefore,
3644	!-- use terrain height array from file, if available. This flag is later used
3645	!-- to control initialization of surface attributes.
3646	surf_usm_h%ground_level = .FALSE.
3647	DO m = 1, surf_usm_h%ns
3648	i = surf_usm_h%i(m)
3649	j = surf_usm_h%j(m)
3650	k = surf_usm_h%k(m)
3651	!
3652	!-- Get local ground level. If no ground level is given in input file,
3653	!-- use default value.
3654	ground_floor_level_l = ground_floor_level
3655	IF ( building_pars_f%from_file ) THEN
3656	IF ( building_pars_f%pars_xy(ind_gflh,j,i) /= &
3657	building_pars_f%fill ) &
3658	ground_floor_level_l = building_pars_f%pars_xy(ind_gflh,j,i)
3659	ENDIF
3660	!
3661	!-- Determine height of surface element above ground level
3662	IF ( terrain_height_f%from_file ) THEN
3663	z_agl = zw(k) - terrain_height_f%var(j,i)
3664	ELSE
3665	z_agl = zw(k)
3666	ENDIF
3667	!
3668	!-- Set flag for ground level
3669	IF ( z_agl <= ground_floor_level_l ) &
3670	surf_usm_h%ground_level(m) = .TRUE.
3671	ENDDO
3672
3673	DO l = 0, 3
3674	surf_usm_v(l)%ground_level = .FALSE.
3675	DO m = 1, surf_usm_v(l)%ns
3676	i = surf_usm_v(l)%i(m) + surf_usm_v(l)%ioff
3677	j = surf_usm_v(l)%j(m) + surf_usm_v(l)%joff
3678	k = surf_usm_v(l)%k(m)
3679	!
3680	!-- Get local ground level. If no ground level is given in input file,
3681	!-- use default value.
3682	ground_floor_level_l = ground_floor_level
3683	IF ( building_pars_f%from_file ) THEN
3684	IF ( building_pars_f%pars_xy(ind_gflh,j,i) /= &
3685	building_pars_f%fill ) &
3686	ground_floor_level_l = building_pars_f%pars_xy(ind_gflh,j,i)
3687	ENDIF
3688	!
3689	!-- Determine height of surface element above ground level. Please
3690	!-- note, height of surface element is determined with respect to
3691	!-- its height of the adjoing atmospheric grid point.
3692	IF ( terrain_height_f%from_file ) THEN
3693	z_agl = zw(k) - terrain_height_f%var(j-surf_usm_v(l)%joff, &
3694	i-surf_usm_v(l)%ioff)
3695	ELSE
3696	z_agl = zw(k)
3697	ENDIF
3698	!
3699	!-- Set flag for ground level
3700	IF ( z_agl <= ground_floor_level_l ) &
3701	surf_usm_v(l)%ground_level(m) = .TRUE.
3702
3703	ENDDO
3704	ENDDO
3705	!
3706	!-- Initialization of resistances.
3707	DO m = 1, surf_usm_h%ns
3708	surf_usm_h%r_a(m) = 50.0_wp
3709	surf_usm_h%r_a_green(m) = 50.0_wp
3710	surf_usm_h%r_a_window(m) = 50.0_wp
3711	ENDDO
3712	DO l = 0, 3
3713	DO m = 1, surf_usm_v(l)%ns
3714	surf_usm_v(l)%r_a(m) = 50.0_wp
3715	surf_usm_v(l)%r_a_green(m) = 50.0_wp
3716	surf_usm_v(l)%r_a_window(m) = 50.0_wp
3717	ENDDO
3718	ENDDO
3719	!
3720	!-- Initialize urban-type surface attribute. According to initialization in
3721	!-- land-surface model, follow a 3-level approach.
3722	!-- Level 1 - initialization via default attributes
3723	DO m = 1, surf_usm_h%ns
3724	!
3725	!-- Now, all horizontal surfaces are roof surfaces (?)
3726	surf_usm_h%isroof_surf(m) = .TRUE.
3727	surf_usm_h%surface_types(m) = roof_category !< default category for root surface
3728	!
3729	!-- In order to distinguish between ground floor level and
3730	!-- above-ground-floor level surfaces, set input indices.
3731	ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, &
3732	surf_usm_h%ground_level(m) )
3733	ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, &
3734	surf_usm_h%ground_level(m) )
3735	ind_green_frac_w = MERGE( ind_green_frac_w_gfl, ind_green_frac_w_agfl, &
3736	surf_usm_h%ground_level(m) )
3737	ind_green_frac_r = MERGE( ind_green_frac_r_gfl, ind_green_frac_r_agfl, &
3738	surf_usm_h%ground_level(m) )
3739	ind_lai_r = MERGE( ind_lai_r_gfl, ind_lai_r_agfl, &
3740	surf_usm_h%ground_level(m) )
3741	ind_lai_w = MERGE( ind_lai_w_gfl, ind_lai_w_agfl, &
3742	surf_usm_h%ground_level(m) )
3743	ind_hc1 = MERGE( ind_hc1_gfl, ind_hc1_agfl, &
3744	surf_usm_h%ground_level(m) )
3745	ind_hc2 = MERGE( ind_hc2_gfl, ind_hc2_agfl, &
3746	surf_usm_h%ground_level(m) )
3747	ind_hc3 = MERGE( ind_hc3_gfl, ind_hc3_agfl, &
3748	surf_usm_h%ground_level(m) )
3749	ind_tc1 = MERGE( ind_tc1_gfl, ind_tc1_agfl, &
3750	surf_usm_h%ground_level(m) )
3751	ind_tc2 = MERGE( ind_tc2_gfl, ind_tc2_agfl, &
3752	surf_usm_h%ground_level(m) )
3753	ind_tc3 = MERGE( ind_tc3_gfl, ind_tc3_agfl, &
3754	surf_usm_h%ground_level(m) )
3755	ind_emis_wall = MERGE( ind_emis_wall_gfl, ind_emis_wall_agfl, &
3756	surf_usm_h%ground_level(m) )
3757	ind_emis_green = MERGE( ind_emis_green_gfl, ind_emis_green_agfl, &
3758	surf_usm_h%ground_level(m) )
3759	ind_emis_win = MERGE( ind_emis_win_gfl, ind_emis_win_agfl, &
3760	surf_usm_h%ground_level(m) )
3761	ind_trans = MERGE( ind_trans_gfl, ind_trans_agfl, &
3762	surf_usm_h%ground_level(m) )
3763	ind_z0 = MERGE( ind_z0_gfl, ind_z0_agfl, &
3764	surf_usm_h%ground_level(m) )
3765	ind_z0qh = MERGE( ind_z0qh_gfl, ind_z0qh_agfl, &
3766	surf_usm_h%ground_level(m) )
3767	!
3768	!-- Store building type and its name on each surface element
3769	surf_usm_h%building_type(m) = building_type
3770	surf_usm_h%building_type_name(m) = building_type_name(building_type)
3771	!
3772	!-- Initialize relatvie wall- (0), green- (1) and window (2) fractions
3773	surf_usm_h%frac(ind_veg_wall,m) = building_pars(ind_wall_frac,building_type)
3774	surf_usm_h%frac(ind_pav_green,m) = building_pars(ind_green_frac_r,building_type)
3775	surf_usm_h%frac(ind_wat_win,m) = building_pars(ind_win_frac,building_type)
3776	surf_usm_h%lai(m) = building_pars(ind_green_frac_r,building_type)
3777
3778	surf_usm_h%rho_c_wall(nzb_wall,m) = building_pars(ind_hc1,building_type)
3779	surf_usm_h%rho_c_wall(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3780	surf_usm_h%rho_c_wall(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3781	surf_usm_h%rho_c_wall(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3782	surf_usm_h%lambda_h(nzb_wall,m) = building_pars(ind_tc1,building_type)
3783	surf_usm_h%lambda_h(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3784	surf_usm_h%lambda_h(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3785	surf_usm_h%lambda_h(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3786	surf_usm_h%rho_c_green(nzb_wall,m) = building_pars(ind_hc1,building_type)
3787	surf_usm_h%rho_c_green(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3788	surf_usm_h%rho_c_green(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3789	surf_usm_h%rho_c_green(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3790	surf_usm_h%lambda_h_green(nzb_wall,m) = building_pars(ind_tc1,building_type)
3791	surf_usm_h%lambda_h_green(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3792	surf_usm_h%lambda_h_green(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3793	surf_usm_h%lambda_h_green(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3794	surf_usm_h%rho_c_window(nzb_wall,m) = building_pars(ind_hc1,building_type)
3795	surf_usm_h%rho_c_window(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3796	surf_usm_h%rho_c_window(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3797	surf_usm_h%rho_c_window(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3798	surf_usm_h%lambda_h_window(nzb_wall,m) = building_pars(ind_tc1,building_type)
3799	surf_usm_h%lambda_h_window(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3800	surf_usm_h%lambda_h_window(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3801	surf_usm_h%lambda_h_window(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3802
3803	surf_usm_h%target_temp_summer(m) = building_pars(12,building_type)
3804	surf_usm_h%target_temp_winter(m) = building_pars(13,building_type)
3805	!
3806	!-- emissivity of wall-, green- and window fraction
3807	surf_usm_h%emissivity(ind_veg_wall,m) = building_pars(ind_emis_wall,building_type)
3808	surf_usm_h%emissivity(ind_pav_green,m) = building_pars(ind_emis_green,building_type)
3809	surf_usm_h%emissivity(ind_wat_win,m) = building_pars(ind_emis_win,building_type)
3810
3811	surf_usm_h%transmissivity(m) = building_pars(ind_trans,building_type)
3812
3813	surf_usm_h%z0(m) = building_pars(ind_z0,building_type)
3814	surf_usm_h%z0h(m) = building_pars(ind_z0qh,building_type)
3815	surf_usm_h%z0q(m) = building_pars(ind_z0qh,building_type)
3816	!
3817	!-- albedo type for wall fraction, green fraction, window fraction
3818	surf_usm_h%albedo_type(ind_veg_wall,m) = INT( building_pars(ind_alb_wall,building_type) )
3819	surf_usm_h%albedo_type(ind_pav_green,m) = INT( building_pars(ind_alb_green,building_type) )
3820	surf_usm_h%albedo_type(ind_wat_win,m) = INT( building_pars(ind_alb_win,building_type) )
3821
3822	surf_usm_h%zw(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3823	surf_usm_h%zw(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3824	surf_usm_h%zw(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3825	surf_usm_h%zw(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3826
3827	surf_usm_h%zw_green(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3828	surf_usm_h%zw_green(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3829	surf_usm_h%zw_green(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3830	surf_usm_h%zw_green(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3831
3832	surf_usm_h%zw_window(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3833	surf_usm_h%zw_window(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3834	surf_usm_h%zw_window(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3835	surf_usm_h%zw_window(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3836
3837	surf_usm_h%c_surface(m) = building_pars(45,building_type)
3838	surf_usm_h%lambda_surf(m) = building_pars(46,building_type)
3839	surf_usm_h%c_surface_green(m) = building_pars(45,building_type)
3840	surf_usm_h%lambda_surf_green(m) = building_pars(46,building_type)
3841	surf_usm_h%c_surface_window(m) = building_pars(45,building_type)
3842	surf_usm_h%lambda_surf_window(m) = building_pars(46,building_type)
3843
3844	ENDDO
3845
3846	DO l = 0, 3
3847	DO m = 1, surf_usm_v(l)%ns
3848
3849	surf_usm_v(l)%surface_types(m) = wall_category !< default category for root surface
3850	!
3851	!-- In order to distinguish between ground floor level and
3852	!-- above-ground-floor level surfaces, set input indices.
3853	ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, &
3854	surf_usm_v(l)%ground_level(m) )
3855	ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, &
3856	surf_usm_v(l)%ground_level(m) )
3857	ind_green_frac_w = MERGE( ind_green_frac_w_gfl, ind_green_frac_w_agfl, &
3858	surf_usm_v(l)%ground_level(m) )
3859	ind_green_frac_r = MERGE( ind_green_frac_r_gfl, ind_green_frac_r_agfl, &
3860	surf_usm_v(l)%ground_level(m) )
3861	ind_lai_r = MERGE( ind_lai_r_gfl, ind_lai_r_agfl, &
3862	surf_usm_v(l)%ground_level(m) )
3863	ind_lai_w = MERGE( ind_lai_w_gfl, ind_lai_w_agfl, &
3864	surf_usm_v(l)%ground_level(m) )
3865	ind_hc1 = MERGE( ind_hc1_gfl, ind_hc1_agfl, &
3866	surf_usm_v(l)%ground_level(m) )
3867	ind_hc2 = MERGE( ind_hc2_gfl, ind_hc2_agfl, &
3868	surf_usm_v(l)%ground_level(m) )
3869	ind_hc3 = MERGE( ind_hc3_gfl, ind_hc3_agfl, &
3870	surf_usm_v(l)%ground_level(m) )
3871	ind_tc1 = MERGE( ind_tc1_gfl, ind_tc1_agfl, &
3872	surf_usm_v(l)%ground_level(m) )
3873	ind_tc2 = MERGE( ind_tc2_gfl, ind_tc2_agfl, &
3874	surf_usm_v(l)%ground_level(m) )
3875	ind_tc3 = MERGE( ind_tc3_gfl, ind_tc3_agfl, &
3876	surf_usm_v(l)%ground_level(m) )
3877	ind_emis_wall = MERGE( ind_emis_wall_gfl, ind_emis_wall_agfl, &
3878	surf_usm_v(l)%ground_level(m) )
3879	ind_emis_green = MERGE( ind_emis_green_gfl, ind_emis_green_agfl, &
3880	surf_usm_v(l)%ground_level(m) )
3881	ind_emis_win = MERGE( ind_emis_win_gfl, ind_emis_win_agfl, &
3882	surf_usm_v(l)%ground_level(m) )
3883	ind_trans = MERGE( ind_trans_gfl, ind_trans_agfl, &
3884	surf_usm_v(l)%ground_level(m) )
3885	ind_z0 = MERGE( ind_z0_gfl, ind_z0_agfl, &
3886	surf_usm_v(l)%ground_level(m) )
3887	ind_z0qh = MERGE( ind_z0qh_gfl, ind_z0qh_agfl, &
3888	surf_usm_v(l)%ground_level(m) )
3889	!
3890	!-- Store building type and its name on each surface element
3891	surf_usm_v(l)%building_type(m) = building_type
3892	surf_usm_v(l)%building_type_name(m) = building_type_name(building_type)
3893	!
3894	!-- Initialize relatvie wall- (0), green- (1) and window (2) fractions
3895	surf_usm_v(l)%frac(ind_veg_wall,m) = building_pars(ind_wall_frac,building_type)
3896	surf_usm_v(l)%frac(ind_pav_green,m) = building_pars(ind_green_frac_w,building_type)
3897	surf_usm_v(l)%frac(ind_wat_win,m) = building_pars(ind_win_frac,building_type)
3898	surf_usm_v(l)%lai(m) = building_pars(ind_lai_w,building_type)
3899
3900	surf_usm_v(l)%rho_c_wall(nzb_wall,m) = building_pars(ind_hc1,building_type)
3901	surf_usm_v(l)%rho_c_wall(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3902	surf_usm_v(l)%rho_c_wall(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3903	surf_usm_v(l)%rho_c_wall(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3904
3905	surf_usm_v(l)%rho_c_green(nzb_wall,m) = building_pars(ind_hc1,building_type)
3906	surf_usm_v(l)%rho_c_green(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3907	surf_usm_v(l)%rho_c_green(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3908	surf_usm_v(l)%rho_c_green(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3909
3910	surf_usm_v(l)%rho_c_window(nzb_wall,m) = building_pars(ind_hc1,building_type)
3911	surf_usm_v(l)%rho_c_window(nzb_wall+1,m) = building_pars(ind_hc1,building_type)
3912	surf_usm_v(l)%rho_c_window(nzb_wall+2,m) = building_pars(ind_hc2,building_type)
3913	surf_usm_v(l)%rho_c_window(nzb_wall+3,m) = building_pars(ind_hc3,building_type)
3914
3915	surf_usm_v(l)%lambda_h(nzb_wall,m) = building_pars(ind_tc1,building_type)
3916	surf_usm_v(l)%lambda_h(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3917	surf_usm_v(l)%lambda_h(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3918	surf_usm_v(l)%lambda_h(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3919
3920	surf_usm_v(l)%lambda_h_green(nzb_wall,m) = building_pars(ind_tc1,building_type)
3921	surf_usm_v(l)%lambda_h_green(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3922	surf_usm_v(l)%lambda_h_green(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3923	surf_usm_v(l)%lambda_h_green(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3924
3925	surf_usm_v(l)%lambda_h_window(nzb_wall,m) = building_pars(ind_tc1,building_type)
3926	surf_usm_v(l)%lambda_h_window(nzb_wall+1,m) = building_pars(ind_tc1,building_type)
3927	surf_usm_v(l)%lambda_h_window(nzb_wall+2,m) = building_pars(ind_tc2,building_type)
3928	surf_usm_v(l)%lambda_h_window(nzb_wall+3,m) = building_pars(ind_tc3,building_type)
3929
3930	surf_usm_v(l)%target_temp_summer(m) = building_pars(12,building_type)
3931	surf_usm_v(l)%target_temp_winter(m) = building_pars(13,building_type)
3932	!
3933	!-- emissivity of wall-, green- and window fraction
3934	surf_usm_v(l)%emissivity(ind_veg_wall,m) = building_pars(ind_emis_wall,building_type)
3935	surf_usm_v(l)%emissivity(ind_pav_green,m) = building_pars(ind_emis_green,building_type)
3936	surf_usm_v(l)%emissivity(ind_wat_win,m) = building_pars(ind_emis_win,building_type)
3937
3938	surf_usm_v(l)%transmissivity(m) = building_pars(ind_trans,building_type)
3939
3940	surf_usm_v(l)%z0(m) = building_pars(ind_z0,building_type)
3941	surf_usm_v(l)%z0h(m) = building_pars(ind_z0qh,building_type)
3942	surf_usm_v(l)%z0q(m) = building_pars(ind_z0qh,building_type)
3943
3944	surf_usm_v(l)%albedo_type(ind_veg_wall,m) = INT( building_pars(ind_alb_wall,building_type) )
3945	surf_usm_v(l)%albedo_type(ind_pav_green,m) = INT( building_pars(ind_alb_green,building_type) )
3946	surf_usm_v(l)%albedo_type(ind_wat_win,m) = INT( building_pars(ind_alb_win,building_type) )
3947
3948	surf_usm_v(l)%zw(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3949	surf_usm_v(l)%zw(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3950	surf_usm_v(l)%zw(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3951	surf_usm_v(l)%zw(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3952
3953	surf_usm_v(l)%zw_green(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3954	surf_usm_v(l)%zw_green(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3955	surf_usm_v(l)%zw_green(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3956	surf_usm_v(l)%zw_green(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3957
3958	surf_usm_v(l)%zw_window(nzb_wall,m) = building_pars(ind_thick_1,building_type)
3959	surf_usm_v(l)%zw_window(nzb_wall+1,m) = building_pars(ind_thick_2,building_type)
3960	surf_usm_v(l)%zw_window(nzb_wall+2,m) = building_pars(ind_thick_3,building_type)
3961	surf_usm_v(l)%zw_window(nzb_wall+3,m) = building_pars(ind_thick_4,building_type)
3962
3963	surf_usm_v(l)%c_surface(m) = building_pars(45,building_type)
3964	surf_usm_v(l)%lambda_surf(m) = building_pars(46,building_type)
3965	surf_usm_v(l)%c_surface_green(m) = building_pars(45,building_type)
3966	surf_usm_v(l)%lambda_surf_green(m) = building_pars(46,building_type)
3967	surf_usm_v(l)%c_surface_window(m) = building_pars(45,building_type)
3968	surf_usm_v(l)%lambda_surf_window(m) = building_pars(46,building_type)
3969
3970	ENDDO
3971	ENDDO
3972	!
3973	!-- Level 2 - initialization via building type read from file
3974	IF ( building_type_f%from_file ) THEN
3975	DO m = 1, surf_usm_h%ns
3976	i = surf_usm_h%i(m)
3977	j = surf_usm_h%j(m)
3978	!
3979	!-- For the moment, limit building type to 6 (to overcome errors in input file).
3980	st = building_type_f%var(j,i)
3981	IF ( st /= building_type_f%fill ) THEN
3982
3983	!
3984	!-- In order to distinguish between ground floor level and
3985	!-- above-ground-floor level surfaces, set input indices.
3986	ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, &
3987	surf_usm_h%ground_level(m) )
3988	ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, &
3989	surf_usm_h%ground_level(m) )
3990	ind_green_frac_w = MERGE( ind_green_frac_w_gfl, ind_green_frac_w_agfl, &
3991	surf_usm_h%ground_level(m) )
3992	ind_green_frac_r = MERGE( ind_green_frac_r_gfl, ind_green_frac_r_agfl, &
3993	surf_usm_h%ground_level(m) )
3994	ind_lai_r = MERGE( ind_lai_r_gfl, ind_lai_r_agfl, &
3995	surf_usm_h%ground_level(m) )
3996	ind_lai_w = MERGE( ind_lai_w_gfl, ind_lai_w_agfl, &
3997	surf_usm_h%ground_level(m) )
3998	ind_hc1 = MERGE( ind_hc1_gfl, ind_hc1_agfl, &
3999	surf_usm_h%ground_level(m) )
4000	ind_hc2 = MERGE( ind_hc2_gfl, ind_hc2_agfl, &
4001	surf_usm_h%ground_level(m) )
4002	ind_hc3 = MERGE( ind_hc3_gfl, ind_hc3_agfl, &
4003	surf_usm_h%ground_level(m) )
4004	ind_tc1 = MERGE( ind_tc1_gfl, ind_tc1_agfl, &
4005	surf_usm_h%ground_level(m) )
4006	ind_tc2 = MERGE( ind_tc2_gfl, ind_tc2_agfl, &
4007	surf_usm_h%ground_level(m) )
4008	ind_tc3 = MERGE( ind_tc3_gfl, ind_tc3_agfl, &
4009	surf_usm_h%ground_level(m) )
4010	ind_emis_wall = MERGE( ind_emis_wall_gfl, ind_emis_wall_agfl, &
4011	surf_usm_h%ground_level(m) )
4012	ind_emis_green = MERGE( ind_emis_green_gfl, ind_emis_green_agfl, &
4013	surf_usm_h%ground_level(m) )
4014	ind_emis_win = MERGE( ind_emis_win_gfl, ind_emis_win_agfl, &
4015	surf_usm_h%ground_level(m) )
4016	ind_trans = MERGE( ind_trans_gfl, ind_trans_agfl, &
4017	surf_usm_h%ground_level(m) )
4018	ind_z0 = MERGE( ind_z0_gfl, ind_z0_agfl, &
4019	surf_usm_h%ground_level(m) )
4020	ind_z0qh = MERGE( ind_z0qh_gfl, ind_z0qh_agfl, &
4021	surf_usm_h%ground_level(m) )
4022	!
4023	!-- Store building type and its name on each surface element
4024	surf_usm_h%building_type(m) = st
4025	surf_usm_h%building_type_name(m) = building_type_name(st)
4026	!
4027	!-- Initialize relatvie wall- (0), green- (1) and window (2) fractions
4028	surf_usm_h%frac(ind_veg_wall,m) = building_pars(ind_wall_frac,st)
4029	surf_usm_h%frac(ind_pav_green,m) = building_pars(ind_green_frac_r,st)
4030	surf_usm_h%frac(ind_wat_win,m) = building_pars(ind_win_frac,st)
4031	surf_usm_h%lai(m) = building_pars(ind_green_frac_r,st)
4032
4033	surf_usm_h%rho_c_wall(nzb_wall,m) = building_pars(ind_hc1,st)
4034	surf_usm_h%rho_c_wall(nzb_wall+1,m) = building_pars(ind_hc1,st)
4035	surf_usm_h%rho_c_wall(nzb_wall+2,m) = building_pars(ind_hc2,st)
4036	surf_usm_h%rho_c_wall(nzb_wall+3,m) = building_pars(ind_hc3,st)
4037	surf_usm_h%lambda_h(nzb_wall,m) = building_pars(ind_tc1,st)
4038	surf_usm_h%lambda_h(nzb_wall+1,m) = building_pars(ind_tc1,st)
4039	surf_usm_h%lambda_h(nzb_wall+2,m) = building_pars(ind_tc2,st)
4040	surf_usm_h%lambda_h(nzb_wall+3,m) = building_pars(ind_tc3,st)
4041
4042	surf_usm_h%rho_c_green(nzb_wall,m) = building_pars(ind_hc1,st)
4043	surf_usm_h%rho_c_green(nzb_wall+1,m) = building_pars(ind_hc1,st)
4044	surf_usm_h%rho_c_green(nzb_wall+2,m) = building_pars(ind_hc2,st)
4045	surf_usm_h%rho_c_green(nzb_wall+3,m) = building_pars(ind_hc3,st)
4046	surf_usm_h%lambda_h_green(nzb_wall,m) = building_pars(ind_tc1,st)
4047	surf_usm_h%lambda_h_green(nzb_wall+1,m) = building_pars(ind_tc1,st)
4048	surf_usm_h%lambda_h_green(nzb_wall+2,m) = building_pars(ind_tc2,st)
4049	surf_usm_h%lambda_h_green(nzb_wall+3,m) = building_pars(ind_tc3,st)
4050
4051	surf_usm_h%rho_c_window(nzb_wall,m) = building_pars(ind_hc1,st)
4052	surf_usm_h%rho_c_window(nzb_wall+1,m) = building_pars(ind_hc1,st)
4053	surf_usm_h%rho_c_window(nzb_wall+2,m) = building_pars(ind_hc2,st)
4054	surf_usm_h%rho_c_window(nzb_wall+3,m) = building_pars(ind_hc3,st)
4055	surf_usm_h%lambda_h_window(nzb_wall,m) = building_pars(ind_tc1,st)
4056	surf_usm_h%lambda_h_window(nzb_wall+1,m) = building_pars(ind_tc1,st)
4057	surf_usm_h%lambda_h_window(nzb_wall+2,m) = building_pars(ind_tc2,st)
4058	surf_usm_h%lambda_h_window(nzb_wall+3,m) = building_pars(ind_tc3,st)
4059
4060	surf_usm_h%target_temp_summer(m) = building_pars(12,st)
4061	surf_usm_h%target_temp_winter(m) = building_pars(13,st)
4062	!
4063	!-- emissivity of wall-, green- and window fraction
4064	surf_usm_h%emissivity(ind_veg_wall,m) = building_pars(ind_emis_wall,st)
4065	surf_usm_h%emissivity(ind_pav_green,m) = building_pars(ind_emis_green,st)
4066	surf_usm_h%emissivity(ind_wat_win,m) = building_pars(ind_emis_win,st)
4067
4068	surf_usm_h%transmissivity(m) = building_pars(ind_trans,st)
4069
4070	surf_usm_h%z0(m) = building_pars(ind_z0,st)
4071	surf_usm_h%z0h(m) = building_pars(ind_z0qh,st)
4072	surf_usm_h%z0q(m) = building_pars(ind_z0qh,st)
4073	!
4074	!-- albedo type for wall fraction, green fraction, window fraction
4075	surf_usm_h%albedo_type(ind_veg_wall,m) = INT( building_pars(ind_alb_wall,st) )
4076	surf_usm_h%albedo_type(ind_pav_green,m) = INT( building_pars(ind_alb_green,st) )
4077	surf_usm_h%albedo_type(ind_wat_win,m) = INT( building_pars(ind_alb_win,st) )
4078
4079	surf_usm_h%zw(nzb_wall,m) = building_pars(ind_thick_1,st)
4080	surf_usm_h%zw(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4081	surf_usm_h%zw(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4082	surf_usm_h%zw(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4083
4084	surf_usm_h%zw_green(nzb_wall,m) = building_pars(ind_thick_1,st)
4085	surf_usm_h%zw_green(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4086	surf_usm_h%zw_green(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4087	surf_usm_h%zw_green(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4088
4089	surf_usm_h%zw_window(nzb_wall,m) = building_pars(ind_thick_1,st)
4090	surf_usm_h%zw_window(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4091	surf_usm_h%zw_window(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4092	surf_usm_h%zw_window(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4093
4094	surf_usm_h%c_surface(m) = building_pars(45,st)
4095	surf_usm_h%lambda_surf(m) = building_pars(46,st)
4096	surf_usm_h%c_surface_green(m) = building_pars(45,st)
4097	surf_usm_h%lambda_surf_green(m) = building_pars(46,st)
4098	surf_usm_h%c_surface_window(m) = building_pars(45,st)
4099	surf_usm_h%lambda_surf_window(m) = building_pars(46,st)
4100
4101	ENDIF
4102	ENDDO
4103
4104	DO l = 0, 3
4105	DO m = 1, surf_usm_v(l)%ns
4106	i = surf_usm_v(l)%i(m) + surf_usm_v(l)%ioff
4107	j = surf_usm_v(l)%j(m) + surf_usm_v(l)%joff
4108	!
4109	!-- For the moment, limit building type to 6 (to overcome errors in input file).
4110
4111	st = building_type_f%var(j,i)
4112	IF ( st /= building_type_f%fill ) THEN
4113
4114	!
4115	!-- In order to distinguish between ground floor level and
4116	!-- above-ground-floor level surfaces, set input indices.
4117	ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, &
4118	surf_usm_v(l)%ground_level(m) )
4119	ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, &
4120	surf_usm_v(l)%ground_level(m) )
4121	ind_green_frac_w = MERGE( ind_green_frac_w_gfl, ind_green_frac_w_agfl, &
4122	surf_usm_v(l)%ground_level(m) )
4123	ind_green_frac_r = MERGE( ind_green_frac_r_gfl, ind_green_frac_r_agfl, &
4124	surf_usm_v(l)%ground_level(m) )
4125	ind_lai_r = MERGE( ind_lai_r_gfl, ind_lai_r_agfl, &
4126	surf_usm_v(l)%ground_level(m) )
4127	ind_lai_w = MERGE( ind_lai_w_gfl, ind_lai_w_agfl, &
4128	surf_usm_v(l)%ground_level(m) )
4129	ind_hc1 = MERGE( ind_hc1_gfl, ind_hc1_agfl, &
4130	surf_usm_v(l)%ground_level(m) )
4131	ind_hc2 = MERGE( ind_hc2_gfl, ind_hc2_agfl, &
4132	surf_usm_v(l)%ground_level(m) )
4133	ind_hc3 = MERGE( ind_hc3_gfl, ind_hc3_agfl, &
4134	surf_usm_v(l)%ground_level(m) )
4135	ind_tc1 = MERGE( ind_tc1_gfl, ind_tc1_agfl, &
4136	surf_usm_v(l)%ground_level(m) )
4137	ind_tc2 = MERGE( ind_tc2_gfl, ind_tc2_agfl, &
4138	surf_usm_v(l)%ground_level(m) )
4139	ind_tc3 = MERGE( ind_tc3_gfl, ind_tc3_agfl, &
4140	surf_usm_v(l)%ground_level(m) )
4141	ind_emis_wall = MERGE( ind_emis_wall_gfl, ind_emis_wall_agfl, &
4142	surf_usm_v(l)%ground_level(m) )
4143	ind_emis_green = MERGE( ind_emis_green_gfl, ind_emis_green_agfl, &
4144	surf_usm_v(l)%ground_level(m) )
4145	ind_emis_win = MERGE( ind_emis_win_gfl, ind_emis_win_agfl, &
4146	surf_usm_v(l)%ground_level(m) )
4147	ind_trans = MERGE( ind_trans_gfl, ind_trans_agfl, &
4148	surf_usm_v(l)%ground_level(m) )
4149	ind_z0 = MERGE( ind_z0_gfl, ind_z0_agfl, &
4150	surf_usm_v(l)%ground_level(m) )
4151	ind_z0qh = MERGE( ind_z0qh_gfl, ind_z0qh_agfl, &
4152	surf_usm_v(l)%ground_level(m) )
4153	!
4154	!-- Store building type and its name on each surface element
4155	surf_usm_v(l)%building_type(m) = st
4156	surf_usm_v(l)%building_type_name(m) = building_type_name(st)
4157	!
4158	!-- Initialize relatvie wall- (0), green- (1) and window (2) fractions
4159	surf_usm_v(l)%frac(ind_veg_wall,m) = building_pars(ind_wall_frac,st)
4160	surf_usm_v(l)%frac(ind_pav_green,m) = building_pars(ind_green_frac_w,st)
4161	surf_usm_v(l)%frac(ind_wat_win,m) = building_pars(ind_win_frac,st)
4162	surf_usm_v(l)%lai(m) = building_pars(ind_lai_w,st)
4163
4164	surf_usm_v(l)%rho_c_wall(nzb_wall,m) = building_pars(ind_hc1,st)
4165	surf_usm_v(l)%rho_c_wall(nzb_wall+1,m) = building_pars(ind_hc1,st)
4166	surf_usm_v(l)%rho_c_wall(nzb_wall+2,m) = building_pars(ind_hc2,st)
4167	surf_usm_v(l)%rho_c_wall(nzb_wall+3,m) = building_pars(ind_hc3,st)
4168
4169	surf_usm_v(l)%rho_c_green(nzb_wall,m) = building_pars(ind_hc1,st)
4170	surf_usm_v(l)%rho_c_green(nzb_wall+1,m) = building_pars(ind_hc1,st)
4171	surf_usm_v(l)%rho_c_green(nzb_wall+2,m) = building_pars(ind_hc2,st)
4172	surf_usm_v(l)%rho_c_green(nzb_wall+3,m) = building_pars(ind_hc3,st)
4173
4174	surf_usm_v(l)%rho_c_window(nzb_wall,m) = building_pars(ind_hc1,st)
4175	surf_usm_v(l)%rho_c_window(nzb_wall+1,m) = building_pars(ind_hc1,st)
4176	surf_usm_v(l)%rho_c_window(nzb_wall+2,m) = building_pars(ind_hc2,st)
4177	surf_usm_v(l)%rho_c_window(nzb_wall+3,m) = building_pars(ind_hc3,st)
4178
4179	surf_usm_v(l)%lambda_h(nzb_wall,m) = building_pars(ind_tc1,st)
4180	surf_usm_v(l)%lambda_h(nzb_wall+1,m) = building_pars(ind_tc1,st)
4181	surf_usm_v(l)%lambda_h(nzb_wall+2,m) = building_pars(ind_tc2,st)
4182	surf_usm_v(l)%lambda_h(nzb_wall+3,m) = building_pars(ind_tc3,st)
4183
4184	surf_usm_v(l)%lambda_h_green(nzb_wall,m) = building_pars(ind_tc1,st)
4185	surf_usm_v(l)%lambda_h_green(nzb_wall+1,m) = building_pars(ind_tc1,st)
4186	surf_usm_v(l)%lambda_h_green(nzb_wall+2,m) = building_pars(ind_tc2,st)
4187	surf_usm_v(l)%lambda_h_green(nzb_wall+3,m) = building_pars(ind_tc3,st)
4188
4189	surf_usm_v(l)%lambda_h_window(nzb_wall,m) = building_pars(ind_tc1,st)
4190	surf_usm_v(l)%lambda_h_window(nzb_wall+1,m) = building_pars(ind_tc1,st)
4191	surf_usm_v(l)%lambda_h_window(nzb_wall+2,m) = building_pars(ind_tc2,st)
4192	surf_usm_v(l)%lambda_h_window(nzb_wall+3,m) = building_pars(ind_tc3,st)
4193
4194	surf_usm_v(l)%target_temp_summer(m) = building_pars(12,st)
4195	surf_usm_v(l)%target_temp_winter(m) = building_pars(13,st)
4196	!
4197	!-- emissivity of wall-, green- and window fraction
4198	surf_usm_v(l)%emissivity(ind_veg_wall,m) = building_pars(ind_emis_wall,st)
4199	surf_usm_v(l)%emissivity(ind_pav_green,m) = building_pars(ind_emis_green,st)
4200	surf_usm_v(l)%emissivity(ind_wat_win,m) = building_pars(ind_emis_win,st)
4201
4202	surf_usm_v(l)%transmissivity(m) = building_pars(ind_trans,st)
4203
4204	surf_usm_v(l)%z0(m) = building_pars(ind_z0,st)
4205	surf_usm_v(l)%z0h(m) = building_pars(ind_z0qh,st)
4206	surf_usm_v(l)%z0q(m) = building_pars(ind_z0qh,st)
4207
4208	surf_usm_v(l)%albedo_type(ind_veg_wall,m) = INT( building_pars(ind_alb_wall,st) )
4209	surf_usm_v(l)%albedo_type(ind_pav_green,m) = INT( building_pars(ind_alb_green,st) )
4210	surf_usm_v(l)%albedo_type(ind_wat_win,m) = INT( building_pars(ind_alb_win,st) )
4211
4212	surf_usm_v(l)%zw(nzb_wall,m) = building_pars(ind_thick_1,st)
4213	surf_usm_v(l)%zw(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4214	surf_usm_v(l)%zw(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4215	surf_usm_v(l)%zw(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4216
4217	surf_usm_v(l)%zw_green(nzb_wall,m) = building_pars(ind_thick_1,st)
4218	surf_usm_v(l)%zw_green(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4219	surf_usm_v(l)%zw_green(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4220	surf_usm_v(l)%zw_green(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4221
4222	surf_usm_v(l)%zw_window(nzb_wall,m) = building_pars(ind_thick_1,st)
4223	surf_usm_v(l)%zw_window(nzb_wall+1,m) = building_pars(ind_thick_2,st)
4224	surf_usm_v(l)%zw_window(nzb_wall+2,m) = building_pars(ind_thick_3,st)
4225	surf_usm_v(l)%zw_window(nzb_wall+3,m) = building_pars(ind_thick_4,st)
4226
4227	surf_usm_v(l)%c_surface(m) = building_pars(45,st)
4228	surf_usm_v(l)%lambda_surf(m) = building_pars(46,st)
4229	surf_usm_v(l)%c_surface_green(m) = building_pars(45,st)
4230	surf_usm_v(l)%lambda_surf_green(m) = building_pars(46,st)
4231	surf_usm_v(l)%c_surface_window(m) = building_pars(45,st)
4232	surf_usm_v(l)%lambda_surf_window(m) = building_pars(46,st)
4233
4234
4235	ENDIF
4236	ENDDO
4237	ENDDO
4238	ENDIF
4239
4240	!
4241	!-- Level 3 - initialization via building_pars read from file
4242	IF ( building_pars_f%from_file ) THEN
4243	DO m = 1, surf_usm_h%ns
4244	i = surf_usm_h%i(m)
4245	j = surf_usm_h%j(m)
4246
4247	!
4248	!-- In order to distinguish between ground floor level and
4249	!-- above-ground-floor level surfaces, set input indices.
4250	ind_wall_frac = MERGE( ind_wall_frac_gfl, ind_wall_frac_agfl, &
4251	surf_usm_h%ground_level(m) )
4252	ind_win_frac = MERGE( ind_win_frac_gfl, ind_win_frac_agfl, &
4253	surf_usm_h%ground_level(m) )
4254	ind_green_frac_w = MERGE( ind_green_frac_w_gfl, ind_green_frac_w_agfl, &
4255	surf_usm_h%ground_level(m) )
4256	ind_green_frac_r = MERGE( ind_green_frac_r_gfl, ind_green_frac_r_agfl, &
4257	surf_usm_h%ground_level(m) )
4258	ind_lai_r = MERGE( ind_lai_r_gfl, ind_lai_r_agfl, &
4259	surf_usm_h%ground_level(m) )
4260	ind_lai_w = MERGE( ind_lai_w_gfl, ind_lai_w_agfl, &
4261	surf_usm_h%ground_level(m) )
4262	ind_hc1 = MERGE( ind_hc1_gfl, ind_hc1_agfl, &
4263	surf_usm_h%ground_level(m) )
4264	ind_hc2 = MERGE( ind_hc2_gfl, ind_hc2_agfl, &
4265	surf_usm_h%ground_level(m) )
4266	ind_hc3 = MERGE( ind_hc3_gfl, ind_hc3_agfl, &
4267	surf_usm_h%ground_level(m) )
4268	ind_tc1 = MERGE( ind_tc1_gfl, ind_tc1_agfl, &
4269	surf_usm_h%ground_level(m) )
4270	ind_tc2 = MERGE( ind_tc2_gfl, ind_tc2_agfl, &
4271	surf_usm_h%ground_level(m) )
4272	ind_tc3 = MERGE( ind_tc3_gfl, ind_tc3_agfl, &
4273	surf_usm_h%ground_level(m) )
4274	ind_emis_wall = MERGE( ind_emis_wall_gfl, ind_emis_wall_agfl, &
4275	surf_usm_h%ground_level(m) )
4276	ind_emis_green = MERGE( ind_emis_green_gfl, ind_emis_green_agfl, &
4277	surf_usm_h%ground_level(m) )
4278	ind_emis_win = MERGE( ind_emis_win_gfl, ind_emis_win_agfl, &
4279	surf_usm_h%ground_level(m) )
4280	ind_trans = MERGE( ind_trans_gfl, ind_trans_agfl, &
4281	surf_usm_h%ground_level(m) )
4282	ind_z0 = MERGE( ind_z0_gfl, ind_z0_agfl, &
4283	surf_usm_h%ground_level(m) )
4284	ind_z0qh = MERGE( ind_z0qh_gfl, ind_z0qh_agfl, &
4285	surf_usm_h%ground_level(m) )
4286
4287	!
4288	!-- Initialize relatvie wall- (0), green- (1) and window (2) fractions
4289	IF ( building_pars_f%pars_xy(ind_wall_frac,j,i) /= building_pars_f%fill ) &
4290	surf_usm_h%frac(ind_veg_wall,m) = building_pars_f%pars_xy(ind_wall_frac,j,i)
4291	IF ( building_pars_f%pars_xy(ind_green_frac_r,j,i) /= building_pars_f%fill ) &
4292	surf_usm_h%frac(ind_pav_green,m) = building_pars_f%pars_xy(ind_green_frac_r,j,i)
4293	IF ( building_pars_f%pars_xy(ind_win_frac,j,i) /= building_pars_f%fill ) &
4294	surf_usm_h%frac(ind_wat_win,m) = building_pars_f%pars_xy(ind_win_frac,j,i)
4295
4296
4297	IF ( building_pars_f%pars_xy(ind_lai_r,j,i) /= building_pars_f%fill ) &
4298	surf_usm_h%lai(m) = building_pars_f%pars_xy(ind_lai_r,j,i)
4299
4300	IF ( building_pars_f%pars_xy(ind_hc1,j,i) /= building_pars_f%fill ) THEN
4301	surf_usm_h%rho_c_wall(nzb_wall,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4302	surf_usm_h%rho_c_wall(nzb_wall+1,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4303	ENDIF
4304	IF ( building_pars_f%pars_xy(ind_hc2,j,i) /= building_pars_f%fill ) &
4305	surf_usm_h%rho_c_wall(nzb_wall+2,m) = building_pars_f%pars_xy(ind_hc2,j,i)
4306	IF ( building_pars_f%pars_xy(ind_hc3,j,i) /= building_pars_f%fill ) &
4307	surf_usm_h%rho_c_wall(nzb_wall+3,m) = building_pars_f%pars_xy(ind_hc3,j,i)
4308	IF ( building_pars_f%pars_xy(ind_hc1,j,i) /= building_pars_f%fill ) THEN
4309	surf_usm_h%rho_c_green(nzb_wall,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4310	surf_usm_h%rho_c_green(nzb_wall+1,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4311	ENDIF
4312	IF ( building_pars_f%pars_xy(ind_hc2,j,i) /= building_pars_f%fill ) &
4313	surf_usm_h%rho_c_green(nzb_wall+2,m) = building_pars_f%pars_xy(ind_hc2,j,i)
4314	IF ( building_pars_f%pars_xy(ind_hc3,j,i) /= building_pars_f%fill ) &
4315	surf_usm_h%rho_c_green(nzb_wall+3,m) = building_pars_f%pars_xy(ind_hc3,j,i)
4316	IF ( building_pars_f%pars_xy(ind_hc1,j,i) /= building_pars_f%fill ) THEN
4317	surf_usm_h%rho_c_window(nzb_wall,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4318	surf_usm_h%rho_c_window(nzb_wall+1,m) = building_pars_f%pars_xy(ind_hc1,j,i)
4319	ENDIF
4320	IF ( building_pars_f%pars_xy(ind_hc2,j,i) /= building_pars_f%fill ) &
4321	surf_usm_h%rho_c_window(nzb_wall+2,m) = building_pars_f%pars_xy(ind_hc2,j,i)
4322	IF ( building_pars_f%pars_xy(ind_hc3,j,i) /= building_pars_f%fill ) &
4323	surf_usm_h%rho_c_window(nzb_wall+3,m) = building_pars_f%pars_xy(ind_hc3,j,i)
4324
4325	IF ( building_pars_f%pars_xy(ind_tc1,j,i) /= building_pars_f%fill ) THEN
4326	surf_usm_h%lambda_h(nzb_wall,m) = building_pars_f%pars_xy(ind_tc1,j,i)
4327	surf_usm_h%lambda_h(nzb_wall+1,m) = building_pars_f%pars_xy(ind_tc1,j,i)
4328	ENDIF
4329	IF ( building_pars_f%pars_xy(ind_tc2,j,i) /= building_pars_f%fill ) &
4330	surf_usm_h%lambda_h(nzb_wall+2,m) = building_pars_f%pars_xy(ind_tc2,j,i)
4331	IF ( building_pars_f%pars_xy(ind_tc3,j,i) /= building_pars_f%fill ) &
4332	surf_usm_h%lambda_h(nzb_wall+3,m) = building_pars_f%pars_xy(ind_tc3,j,i)
4333	IF ( building_pars_f%pars_xy(ind_tc1,j,i) /= building_pars_f%fill ) THEN
4334	surf_usm_h%lambda_h_green(nzb_wall,m) = building_pars_f%pars_xy(ind_tc1,j,i)
4335	surf_usm_h%lambda_h_green(nzb_wall+1,m) = building_pars_f%pars_xy(ind_tc1,j,i)
4336	ENDIF
4337	IF ( building_pars_f%pars_xy(ind_tc2,j,i) /= building_pars_f%fill ) &
4338	surf_usm_h%lambda_h_green(nzb_wall+2,m) = building_pars_f%pars_xy(ind_tc2,j,i)
4339	IF ( building_pars_f%pars_xy(ind_tc3,j,i) /= building_pars_f%fill ) &
4340	surf_usm_h%lambda_h_green(nzb_wall+3,m) = building_pars_f%pars_xy(ind_tc3,j,i)
4341	IF ( building_pars_f%pars_xy(ind_tc1,j,i)