source: palm/trunk/SOURCE/urban_surface_mod.f90 @ 3637

Last change on this file since 3637 was 3637, checked in by knoop, 2 years ago

M Makefile

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