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

Last change on this file since 3710 was 3710, checked in by suehring, 2 years ago

Check if building-, water-, pavement-, vegetation- and soil types are within a valid range

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