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

Last change on this file since 3832 was 3832, checked in by raasch, 2 years ago

some routines instrumented with openmp directives, loop reordering for performance optimization

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