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

Last change on this file since 4392 was 4392, checked in by pavelkrc, 20 months ago

Merge branch resler (updated documentation, optional flux tracing, bugfix)

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