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

Last change on this file since 3712 was 3712, checked in by kanani, 2 years ago

Formatting and clean-up (urban_surface_mod)

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