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

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

Replace get_topography_top_index functions by pre-calculated arrays in order to save computational resources

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