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

Last change on this file since 4148 was 4148, checked in by suehring, 22 months ago

indoor_model_mod: Revision of indoor model; urban_surface_mod: parameters for waste heat from cooling and heating are introduced to building data base; initialization of building data base moved to an earlier point of time before indoor model will be initialized; radiation_model_mod: minor improvement in some comments; synthetic_turbulence_generator: unused variable removed

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