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

Last change on this file since 2977 was 2977, checked in by kanani, 3 years ago

Fixes for radiative transfer model

  • Property svn:keywords set to Id
File size: 425.0 KB
RevLine 
[2737]1!> @file urban_surface_mod.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 2015-2018 Czech Technical University in Prague
[2920]18! Copyright 2015-2018 Institute of Computer Science of the
19!                     Czech Academy of Sciences, Prague
[2737]20! Copyright 1997-2018 Leibniz Universitaet Hannover
21!------------------------------------------------------------------------------!
22!
23! Current revisions:
24! ------------------
[2705]25!
26!
[2737]27! Former revisions:
28! -----------------
29! $Id: urban_surface_mod.f90 2977 2018-04-17 10:27:57Z kanani $
[2977]30! Implement changes from branch radiation (r2948-2971) with minor modifications.
31! (moh.hefny):
32! Extended exn for all model domain height to avoid the need to get nzut.
33!
34! 2963 2018-04-12 14:47:44Z suehring
[2963]35! Introduce index for vegetation/wall, pavement/green-wall and water/window
36! surfaces, for clearer access of surface fraction, albedo, emissivity, etc. .
37!
38! 2943 2018-04-03 16:17:10Z suehring
[2943]39! Calculate exner function at all height levels and remove some un-used
40! variables.
41!
42! 2932 2018-03-26 09:39:22Z maronga
[2932]43! renamed urban_surface_par to urban_surface_parameters
44!
45! 2921 2018-03-22 15:05:23Z Giersch
[2921]46! The activation of spinup has been moved to parin
47!
48! 2920 2018-03-22 11:22:01Z kanani
[2920]49! Remove unused pcbl, npcbl from ONLY list
50! moh.hefny:
51! Fixed bugs introduced by new structures and by moving radiation interaction
52! into radiation_model_mod.f90.
53! Bugfix: usm data output 3D didn't respect directions
54!
55! 2906 2018-03-19 08:56:40Z Giersch
[2906]56! Local variable ids has to be initialized with a value of -1 in
57! usm_average_3d_data
58!
59! 2894 2018-03-15 09:17:58Z Giersch
[2894]60! Calculations of the index range of the subdomain on file which overlaps with
61! the current subdomain are already done in read_restart_data_mod,
62! usm_read/write_restart_data have been renamed to usm_r/wrd_local, variable
63! named found has been introduced for checking if restart data was found,
64! reading of restart strings has been moved completely to
65! read_restart_data_mod, usm_rrd_local is already inside the overlap loop
66! programmed in read_restart_data_mod, SAVE attribute added where necessary,
67! deallocation and allocation of some arrays have been changed to take care of
68! different restart files that can be opened (index i), the marker *** end usm
69! *** is not necessary anymore, strings and their respective lengths are
70! written out and read now in case of restart runs to get rid of prescribed
71! character lengths
72!
73! 2805 2018-02-14 17:00:09Z suehring
[2805]74! Initialization of resistances.
75!
76! 2797 2018-02-08 13:24:35Z suehring
[2797]77! Comment concerning output of ground-heat flux added.
78!
79! 2766 2018-01-22 17:17:47Z kanani
[2766]80! Removed redundant commas, added some blanks
81!
82! 2765 2018-01-22 11:34:58Z maronga
[2765]83! Major bugfix in calculation of f_shf. Adjustment of roughness lengths in
84! building_pars
85!
86! 2750 2018-01-15 16:26:51Z knoop
[2746]87! Move flag plant canopy to modules
88!
89! 2737 2018-01-11 14:58:11Z kanani
[2737]90! Removed unused variables t_surf_whole...
91!
92! 2735 2018-01-11 12:01:27Z suehring
93! resistances are saved in surface attributes
94!
[2735]95! 2723 2018-01-05 09:27:03Z maronga
[2737]96! Bugfix for spinups (end_time was increased twice in case of LSM + USM runs)
97!
98! 2720 2018-01-02 16:27:15Z kanani
99! Correction of comment
100!
101! 2718 2018-01-02 08:49:38Z maronga
102! Corrected "Former revisions" section
103!
104! 2705 2017-12-18 11:26:23Z maronga
105! Changes from last commit documented
106!
107! 2703 2017-12-15 20:12:38Z maronga
108! Workaround for calculation of r_a
109!
110! 2696 2017-12-14 17:12:51Z kanani
111! - Change in file header (GPL part)
112! - Bugfix in calculation of pt_surface and related fluxes. (BM)
113! - Do not write surface temperatures onto pt array as this might cause
114!   problems with nesting. (MS)
115! - Revised calculation of pt1 (now done in surface_layer_fluxes).
116!   Bugfix, f_shf_window and f_shf_green were not set at vertical surface
117!   elements. (MS)
118! - merged with branch ebsolver
119!   green building surfaces do not evaporate yet
120!   properties of green wall layers and window layers are taken from wall layers
121!   this input data is missing. (RvT)
122! - Merged with branch radiation (developed by Mohamed Salim)
123! - Revised initialization. (MS)
124! - Rename emiss_surf into emissivity, roughness_wall into z0, albedo_surf into
125!   albedo. (MS)
126! - Move first call of usm_radiatin from usm_init to init_3d_model
127! - fixed problem with near surface temperature
128! - added near surface temperature t_surf_10cm_h(m), t_surf_10cm_v(l)%t(m)
129! - does not work with temp profile including stability, ol
130!   t_surf_10cm = pt1 now
131! - merged with 2357 bugfix, error message for nopointer version
132! - added indoor model coupling with wall heat flux
133! - added green substrate/ dry vegetation layer for buildings
134! - merged with 2232 new surface-type structure
135! - added transmissivity of window tiles
136! - added MOSAIK tile approach for 3 different surfaces (RvT)
137!
138! 2583 2017-10-26 13:58:38Z knoop
139! Bugfix: reverted MPI_Win_allocate_cptr introduction in last commit
140!
141! 2582 2017-10-26 13:19:46Z hellstea
142! Workaround for gnufortran compiler added in usm_calc_svf. CALL MPI_Win_allocate is
143! replaced by CALL MPI_Win_allocate_cptr if defined ( __gnufortran ).
144!
145! 2544 2017-10-13 18:09:32Z maronga
146! Date and time quantities are now read from date_and_time_mod. Solar constant is
147! read from radiation_model_mod
148!
149! 2516 2017-10-04 11:03:04Z suehring
150! Remove tabs
151!
152! 2514 2017-10-04 09:52:37Z suehring
153! upper bounds of 3d output changed from nx+1,ny+1 to nx,ny
154! no output of ghost layer data
155!
156! 2350 2017-08-15 11:48:26Z kanani
157! Bugfix and error message for nopointer version.
158! Additional "! defined(__nopointer)" as workaround to enable compilation of
159! nopointer version.
160!
161! 2318 2017-07-20 17:27:44Z suehring
162! Get topography top index via Function call
163!
164! 2317 2017-07-20 17:27:19Z suehring
165! Bugfix: adjust output of shf. Added support for spinups
166!
167! 2287 2017-06-15 16:46:30Z suehring
168! Bugfix in determination topography-top index
169!
170! 2269 2017-06-09 11:57:32Z suehring
171! Enable restart runs with different number of PEs
172! Bugfixes nopointer branch
173!
174! 2258 2017-06-08 07:55:13Z suehring
175! Bugfix, add pre-preprocessor directives to enable non-parrallel mode
176!
177! 2233 2017-05-30 18:08:54Z suehring
178!
179! 2232 2017-05-30 17:47:52Z suehring
180! Adjustments according to new surface-type structure. Remove usm_wall_heat_flux;
181! insteat, heat fluxes are directly applied in diffusion_s.
182!
183! 2213 2017-04-24 15:10:35Z kanani
184! Removal of output quantities usm_lad and usm_canopy_hr
185!
186! 2209 2017-04-19 09:34:46Z kanani
187! cpp switch __mpi3 removed,
188! minor formatting,
189! small bugfix for division by zero (Krc)
190!
191! 2113 2017-01-12 13:40:46Z kanani
192! cpp switch __mpi3 added for MPI-3 standard code (Ketelsen)
193!
194! 2071 2016-11-17 11:22:14Z maronga
195! Small bugfix (Resler)
196!
197! 2031 2016-10-21 15:11:58Z knoop
198! renamed variable rho to rho_ocean
199!
200! 2024 2016-10-12 16:42:37Z kanani
201! Bugfixes in deallocation of array plantt and reading of csf/csfsurf,
202! optimization of MPI-RMA operations,
203! declaration of pcbl as integer,
204! renamed usm_radnet -> usm_rad_net, usm_canopy_khf -> usm_canopy_hr,
205! splitted arrays svf -> svf & csf, svfsurf -> svfsurf & csfsurf,
206! use of new control parameter varnamelength,
207! added output variables usm_rad_ressw, usm_rad_reslw,
208! minor formatting changes,
209! minor optimizations.
210!
211! 2011 2016-09-19 17:29:57Z kanani
212! Major reformatting according to PALM coding standard (comments, blanks,
213! alphabetical ordering, etc.),
214! removed debug_prints,
215! removed auxiliary SUBROUTINE get_usm_info, instead, USM flag urban_surface is
216! defined in MODULE control_parameters (modules.f90) to avoid circular
217! dependencies,
218! renamed canopy_heat_flux to pc_heating_rate, as meaning of quantity changed.
219!
220! 2007 2016-08-24 15:47:17Z kanani
221! Initial revision
222!
223!
224! Description:
225! ------------
226! 2016/6/9 - Initial version of the USM (Urban Surface Model)
227!            authors: Jaroslav Resler, Pavel Krc
228!                     (Czech Technical University in Prague and Institute of
229!                      Computer Science of the Czech Academy of Sciences, Prague)
230!            with contributions: Michal Belda, Nina Benesova, Ondrej Vlcek
231!            partly inspired by PALM LSM (B. Maronga)
232!            parameterizations of Ra checked with TUF3D (E. S. Krayenhoff)
233!> Module for Urban Surface Model (USM)
234!> The module includes:
235!>    1. radiation model with direct/diffuse radiation, shading, reflections
236!>       and integration with plant canopy
237!>    2. wall and wall surface model
238!>    3. surface layer energy balance
239!>    4. anthropogenic heat (only from transportation so far)
240!>    5. necessary auxiliary subroutines (reading inputs, writing outputs,
241!>       restart simulations, ...)
242!> It also make use of standard radiation and integrates it into
243!> urban surface model.
244!>
245!> Further work:
246!> -------------
[2920]247!> 1. Remove global arrays surfouts, surfoutl and only keep track of radiosity
[2737]248!>    from surfaces that are visible from local surfaces (i.e. there is a SVF
249!>    where target is local). To do that, radiosity will be exchanged after each
250!>    reflection step using MPI_Alltoall instead of current MPI_Allgather.
251!>
[2920]252!> 2. Temporarily large values of surface heat flux can be observed, up to
[2737]253!>    1.2 Km/s, which seem to be not realistic.
254!>
255!> @todo Output of _av variables in case of restarts
256!> @todo Revise flux conversion in energy-balance solver
257!> @todo Bugfixing in nopointer branch
258!> @todo Check optimizations for RMA operations
259!> @todo Alternatives for MPI_WIN_ALLOCATE? (causes problems with openmpi)
260!> @todo Check for load imbalances in CPU measures, e.g. for exchange_horiz_prog
261!>       factor 3 between min and max time
262!> @todo Move setting of flag indoor_model to indoor_model_mod once available
263!> @todo Check divisions in wtend (etc.) calculations for possible division
264!>       by zero, e.g. in case fraq(0,m) + fraq(1,m) = 0?!
265!> @todo Use unit 90 for OPEN/CLOSE of input files (FK)
266!------------------------------------------------------------------------------!
267 MODULE urban_surface_mod
268
269#if ! defined( __nopointer )
270    USE arrays_3d,                                                             &
[2920]271        ONLY:  hyp, zu, pt, pt_1, pt_2, p, u, v, w, hyp, tend
[2737]272#endif
273
274    USE cloud_parameters,                                                      &
275        ONLY:  cp, r_d
276
277    USE constants,                                                             &
278        ONLY:  pi
279   
280    USE control_parameters,                                                    &
281        ONLY:  coupling_start_time, dz, topography, dt_3d,                     &
282               intermediate_timestep_count, initializing_actions,              &
283               intermediate_timestep_count_max, simulated_time, end_time,      &
284               timestep_scheme, tsc, coupling_char, io_blocks, io_group,       &
285               message_string, time_since_reference_point, surface_pressure,   &
286               g, pt_surface, large_scale_forcing, lsf_surf, spinup,           &
287               spinup_pt_mean, spinup_time, time_do3d, dt_do3d,                &
[2746]288               average_count_3d, varnamelength, urban_surface, kappa,          &
289               plant_canopy
[2737]290
291    USE cpulog,                                                                &
292        ONLY:  cpu_log, log_point, log_point_s
293
294    USE date_and_time_mod,                                                     &
[2920]295        ONLY:  time_utc_init
[2737]296
297    USE grid_variables,                                                        &
298        ONLY:  dx, dy, ddx, ddy, ddx2, ddy2
299
300    USE indices,                                                               &
301        ONLY:  nx, ny, nnx, nny, nnz, nxl, nxlg, nxr, nxrg, nyn, nyng, nys,    &
302               nysg, nzb, nzt, nbgp, wall_flags_0
303
304    USE, INTRINSIC :: iso_c_binding
305
306    USE kinds
307             
308    USE pegrid
309   
310    USE plant_canopy_model_mod,                                                &
[2920]311        ONLY:  pc_heating_rate
[2737]312   
313    USE radiation_model_mod,                                                   &
[2943]314        ONLY:  albedo_type, radiation_interaction, calc_zenith, zenith,        &
[2737]315               rad_sw_in, rad_lw_in, rad_sw_out, rad_lw_out,                   &
316               sigma_sb, solar_constant, sun_direction, sun_dir_lat,           &
317               sun_dir_lon,                                                    &
318               force_radiation_call, surfinsw, surfinlw, surfinswdir,          &
319               surfinswdif, surfoutsw, surfoutlw, surfins,nsvfl, svf, svfsurf, &
[2920]320               surfinl, surfinlwdif, rad_sw_in_dir, rad_sw_in_diff,            &
[2737]321               rad_lw_in_diff, surfouts, surfoutl, surfoutsl, surfoutll, surf, &
[2920]322               surfl, nsurfl, nsurfs, surfstart, pcbinsw, pcbinlw,             &
323               iup_u, inorth_u, isouth_u, ieast_u, iwest_u, iup_l,             &
324               inorth_l, isouth_l, ieast_l, iwest_l, id,                       &
[2943]325               iz, iy, ix, idir, jdir, kdir,  nsurf_type, nsurf, idsvf, ndsvf, &
[2920]326               iup_a, idown_a, inorth_a, isouth_a, ieast_a, iwest_a,           &
327               idcsf, ndcsf, kdcsf, pct,                                       &
328               startland, endland, startwall, endwall, skyvf, skyvft
[2737]329
330    USE statistics,                                                            &
331        ONLY:  hom, statistic_regions
332
[2963]333    USE surface_mod,                                                           &
334        ONLY:  get_topography_top_index_ji, get_topography_top_index,          &
335               ind_pav_green, ind_veg_wall, ind_wat_win, surf_usm_h,           &
336               surf_usm_v, surface_restore_elements
[2737]337
338
339    IMPLICIT NONE
340
341
342!-- configuration parameters (they can be setup in PALM config)
343    LOGICAL                                        ::  usm_material_model = .TRUE.        !< flag parameter indicating wheather the  model of heat in materials is used
344    LOGICAL                                        ::  usm_anthropogenic_heat = .FALSE.   !< flag parameter indicating wheather the anthropogenic heat sources (e.g.transportation) are used
345    LOGICAL                                        ::  force_radiation_call_l = .FALSE.   !< flag parameter for unscheduled radiation model calls
346    LOGICAL                                        ::  indoor_model = .FALSE.              !< whether to use the indoor model
[2920]347    LOGICAL                                        ::  read_wall_temp_3d = .FALSE.
[2737]348
[2920]349
[2737]350    INTEGER(iwp)                                   ::  building_type = 1                  !< default building type (preleminary setting)
351    INTEGER(iwp)                                   ::  land_category = 2                  !< default category for land surface
352    INTEGER(iwp)                                   ::  wall_category = 2                  !< default category for wall surface over pedestrian zone
[2920]353    INTEGER(iwp)                                   ::  pedestrian_category = 2            !< default category for wall surface in pedestrian zone
[2737]354    INTEGER(iwp)                                   ::  roof_category = 2                  !< default category for root surface
[2920]355    REAL(wp)                                       ::  roughness_concrete = 0.001_wp      !< roughness length of average concrete surface
[2737]356!
357!-- Indices of input attributes for (above) ground floor level
358    INTEGER(iwp) ::  ind_alb_wall          = 38 !< index in input list for albedo_type of wall fraction
359    INTEGER(iwp) ::  ind_alb_green         = 39 !< index in input list for albedo_type of green fraction
360    INTEGER(iwp) ::  ind_alb_win           = 40 !< index in input list for albedo_type of window fraction
361    INTEGER(iwp) ::  ind_emis_wall_agfl    = 14 !< index in input list for wall emissivity, above ground floor level
362    INTEGER(iwp) ::  ind_emis_wall_gfl     = 32 !< index in input list for wall emissivity, ground floor level
363    INTEGER(iwp) ::  ind_emis_green_agfl   = 15 !< index in input list for green emissivity, above ground floor level
364    INTEGER(iwp) ::  ind_emis_green_gfl    = 33 !< index in input list for green emissivity, ground floor level
365    INTEGER(iwp) ::  ind_emis_win_agfl     = 16 !< index in input list for window emissivity, above ground floor level
366    INTEGER(iwp) ::  ind_emis_win_gfl      = 34 !< index in input list for window emissivity, ground floor level
367    INTEGER(iwp) ::  ind_green_frac_w_agfl = 2  !< index in input list for green fraction on wall, above ground floor level
368    INTEGER(iwp) ::  ind_green_frac_w_gfl  = 23 !< index in input list for green fraction on wall, ground floor level
369    INTEGER(iwp) ::  ind_green_frac_r_agfl = 3  !< index in input list for green fraction on roof, above ground floor level
370    INTEGER(iwp) ::  ind_green_frac_r_gfl  = 24 !< index in input list for green fraction on roof, ground floor level
371    INTEGER(iwp) ::  ind_hc1_agfl          =  6 !< index in input list for heat capacity at first wall layer, above ground floor level
372    INTEGER(iwp) ::  ind_hc1_gfl           = 26 !< index in input list for heat capacity at first wall layer, ground floor level
373    INTEGER(iwp) ::  ind_hc2_agfl          = 7  !< index in input list for heat capacity at second wall layer, above ground floor level
374    INTEGER(iwp) ::  ind_hc2_gfl           = 27 !< index in input list for heat capacity at second wall layer, ground floor level
375    INTEGER(iwp) ::  ind_hc3_agfl          = 8  !< index in input list for heat capacity at third wall layer, above ground floor level
376    INTEGER(iwp) ::  ind_hc3_gfl           = 28 !< index in input list for heat capacity at third wall layer, ground floor level
377    INTEGER(iwp) ::  ind_gflh              = 20 !< index in input list for ground floor level height
378    INTEGER(iwp) ::  ind_lai_r_agfl        = 4  !< index in input list for LAI on roof, above ground floor level
379    INTEGER(iwp) ::  ind_lai_r_gfl         = 4  !< index in input list for LAI on roof, ground floor level
380    INTEGER(iwp) ::  ind_lai_w_agfl        = 5  !< index in input list for LAI on wall, above ground floor level
381    INTEGER(iwp) ::  ind_lai_w_gfl         = 25 !< index in input list for LAI on wall, ground floor level
382    INTEGER(iwp) ::  ind_tc1_agfl          = 9  !< index in input list for thermal conductivity at first wall layer, above ground floor level
383    INTEGER(iwp) ::  ind_tc1_gfl           = 29 !< index in input list for thermal conductivity at first wall layer, ground floor level
384    INTEGER(iwp) ::  ind_tc2_agfl          = 10 !< index in input list for thermal conductivity at second wall layer, above ground floor level
385    INTEGER(iwp) ::  ind_tc2_gfl           = 30 !< index in input list for thermal conductivity at second wall layer, ground floor level
386    INTEGER(iwp) ::  ind_tc3_agfl          = 11 !< index in input list for thermal conductivity at third wall layer, above ground floor level
387    INTEGER(iwp) ::  ind_tc3_gfl           = 31 !< index in input list for thermal conductivity at third wall layer, ground floor level
388    INTEGER(iwp) ::  ind_thick_1           = 41 !< index for wall layer thickness - 1st layer
389    INTEGER(iwp) ::  ind_thick_2           = 42 !< index for wall layer thickness - 2nd layer
390    INTEGER(iwp) ::  ind_thick_3           = 43 !< index for wall layer thickness - 3rd layer
391    INTEGER(iwp) ::  ind_thick_4           = 44 !< index for wall layer thickness - 4th layer
392    INTEGER(iwp) ::  ind_trans_agfl        = 17 !< index in input list for window transmissivity, above ground floor level
393    INTEGER(iwp) ::  ind_trans_gfl         = 35 !< index in input list for window transmissivity, ground floor level
394    INTEGER(iwp) ::  ind_wall_frac_agfl    = 0  !< index in input list for wall fraction, above ground floor level
395    INTEGER(iwp) ::  ind_wall_frac_gfl     = 21 !< index in input list for wall fraction, ground floor level
396    INTEGER(iwp) ::  ind_win_frac_agfl     = 1  !< index in input list for window fraction, above ground floor level
397    INTEGER(iwp) ::  ind_win_frac_gfl      = 22 !< index in input list for window fraction, ground floor level
398    INTEGER(iwp) ::  ind_z0_agfl           = 18 !< index in input list for z0, above ground floor level
399    INTEGER(iwp) ::  ind_z0_gfl            = 36 !< index in input list for z0, ground floor level
400    INTEGER(iwp) ::  ind_z0qh_agfl         = 19 !< index in input list for z0h / z0q, above ground floor level
401    INTEGER(iwp) ::  ind_z0qh_gfl          = 37 !< index in input list for z0h / z0q, ground floor level
402
403
404    REAL(wp)  ::  roof_height_limit = 4._wp          !< height for distinguish between land surfaces and roofs
405    REAL(wp)  ::  ground_floor_level = 4.0_wp        !< default ground floor level
406
[2920]407
[2737]408    CHARACTER(37), DIMENSION(0:6), PARAMETER :: building_type_name = (/     &
409                                   'user-defined                         ', & !  0
410                                   'residential - 1950                   ', & !  1
411                                   'residential 1951 - 2000              ', & !  2
412                                   'residential 2001 -                   ', & !  3
413                                   'office - 1950                        ', & !  4
414                                   'office 1951 - 2000                   ', & !  5
415                                   'office 2001 -                        '  & !  6
416                                                                     /)
417!
418!-- building parameters, 4 different types
419!-- 0 - wall fraction, 1- window fraction, 2 - green fraction on wall, 3- green fraction
420!-- at roof, 4 - lai of green fraction at roof,  5 - lai of green fraction at wall,
421!-- 6 - heat capacity of wall layer 1, 7 - heat capacity of wall layer 2,
422!-- 8 - heat capacity of wall layer 3, 9 - thermal conductivity of wall layer 1,
423!-- 10 - thermal conductivity of wall layer 2, 11 - thermal conductivity of wall layer 3, 
424!-- 12 - indoor target summer temperature ( K ), 13 - indoor target winter temperature (K),
425!-- 14 - emissivity of wall fraction, 15 - emissivity of green fraction, 16 - emissivity of window fraction,
426!-- 17 - transmissivity of window fraction, 18 - z0, 19 - z0h/z0q, 20 - ground floor height,
427!-- 21 - ground floor wall fraction, 22 - ground floor window fraction, 23 ground floor green fraction,
428!-- 24 - ground floor green fraction on roof, 25 - ground floor lai of green fraction,
429!-- 26 - ground floor heat capacity of wall layer 1, 27 - ground floor heat capacity of wall layer 1,
430!-- 28 - ground floor heat capacity of wall layer 3, 29 - ground floor thermal conductivity of wall layer 1,
431!-- 30 - ground floor thermal conductivity of wall layer 2, 31 - ground floor thermal conductivity of wall layer 3,
432!-- 32 - ground floor emissivity of wall fraction, 33 - ground floor emissivity of green fraction,
433!-- 34 - ground floor emissivity of window fraction, 35 - ground floor transmissivity of window fraction,
434!-- 36 - ground floor z0, 37 - ground floor z0h/z0q, 38 - albedo type wall fraction
435!-- 39 - albedo type green fraction, 40 - albedo type window fraction
436!-- 41 - wall layer thickness - 1st layer, 42 - wall layer thickness - 2nd layer,
437!-- 43 - wall layer thickness - 3rd layer, 44 - wall layer thickness - 4th layer,
438!-- 45 - heat capacity of the wall surface, 46 - heat conductivity
439!-- Please note, only preleminary dummy values so far!
440    REAL(wp), DIMENSION(0:46,1:6), PARAMETER :: building_pars = RESHAPE( (/    &
441        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
442        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
443        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
[2765]444        0.001_wp, 0.0001_wp, 4.0_wp,                                             & !parameter 18-20
[2737]445        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
446        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
447        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
448        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
449        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
450        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
451        20000.0_wp, 10.0_wp,                                                   & !parameter 45-46 - end of type 1
452        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
453        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
454        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
[2765]455        0.001_wp, 0.0001_wp, 4.0_wp,                                             & !parameter 18-20
[2737]456        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
457        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
458        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
459        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
460        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
461        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
462        20000.0_wp, 10.0_wp,                                                   & !parameter 45-46 - end of type 2
463        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
464        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
465        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
[2765]466        0.001_wp, 0.0001_wp, 4.0_wp,                                             & !parameter 18-20
[2737]467        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
468        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
469        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
470        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
471        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
472        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
473        20000.0_wp, 10.0_wp,                                                   & !parameter 45-46 - end of type 3
474        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
475        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
476        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
477        0.01_wp, 0.001_wp, 4.0_wp,                                             & !parameter 18-20
478        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
479        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
480        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
481        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
482        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
483        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
484        20000.0_wp, 10.0_wp,                                                   & !parameter 45-46 - end of type 4
485        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
486        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
487        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
[2765]488        0.001_wp, 0.0001_wp, 4.0_wp,                                             & !parameter 18-20
[2737]489        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
490        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
491        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
492        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
493        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
494        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
495        20000.0_wp, 10.0_wp,                                                   & !parameter 45-46 - end of type 5
496        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 1.0_wp, 1.0_wp,                        & !parameter 0-5
497        1000000.0_wp, 1000000.0_wp, 1000000.0_wp, 0.3_wp, 0.3_wp, 0.3_wp,      & !parameter 6-11
498        296.15_wp, 293.15_wp, 0.9_wp, 0.9_wp, 0.01_wp, 0.99_wp,                & !parameter 12-17
[2765]499        0.001_wp, 0.0001_wp, 4.0_wp,                                             & !parameter 18-20
[2737]500        1.0_wp, 0.0_wp, 0.0_wp, 0.0_wp, 3.0_wp,                                & !parameter 21-25
501        1000000.0_wp, 1000000.0_wp, 1000000.0_wp,                              & !parameter 26-28                     
502        0.3_wp, 0.3_wp, 0.3_wp,                                                & !parameter 29-31       
503        0.4_wp, 0.4_wp, 0.4_wp, 0.4_wp, 0.01_wp, 0.001_wp,                     & !parameter 32-37
504        24.0_wp, 24.0_wp, 24.0_wp,                                             & !parameter 38-40
505        0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp,                                & !parameter 41-44
506        20000.0_wp, 10.0_wp                                                    & !parameter 45-46 - end of type 6
507                                                                          /),  &
508                                                               (/47, 6/) )
509
510!
511!-- Type for surface temperatures at vertical walls. Is not necessary for horizontal walls.
512    TYPE t_surf_vertical
513       REAL(wp), DIMENSION(:), ALLOCATABLE         :: t
514    END TYPE t_surf_vertical
515!
516!-- Type for wall temperatures at vertical walls. Is not necessary for horizontal walls.
517    TYPE t_wall_vertical
518       REAL(wp), DIMENSION(:,:), ALLOCATABLE       :: t
519    END TYPE t_wall_vertical
520
521
522!-- arrays for time averages
523!-- Attention: the variable rad_net_av is also used in the 3d field variable in radiation_model_mod.f90. It may be better to rename it
524    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  rad_net_av       !< average of rad_net_l
525    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinsw_av      !< average of sw radiation falling to local surface including radiation from reflections
526    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinlw_av      !< average of lw radiation falling to local surface including radiation from reflections
527    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinswdir_av   !< average of direct sw radiation falling to local surface
528    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinswdif_av   !< average of diffuse sw radiation from sky and model boundary falling to local surface
529    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinlwdif_av   !< average of diffuse lw radiation from sky and model boundary falling to local surface
530    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinswref_av   !< average of sw radiation falling to surface from reflections
531    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinlwref_av   !< average of lw radiation falling to surface from reflections
532    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfoutsw_av     !< average of total sw radiation outgoing from nonvirtual surfaces surfaces after all reflection
533    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfoutlw_av     !< average of total lw radiation outgoing from nonvirtual surfaces surfaces after all reflection
534    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfins_av       !< average of array of residua of sw radiation absorbed in surface after last reflection
535    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinl_av       !< average of array of residua of lw radiation absorbed in surface after last reflection
536    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfhf_av        !< average of total radiation flux incoming to minus outgoing from local surface 
537    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  wghf_eb_av       !< average of wghf_eb
538    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  wshf_eb_av       !< average of wshf_eb
539    REAL(wp), DIMENSION(:,:), ALLOCATABLE          ::  t_wall_av        !< Average of t_wall
540    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  wghf_eb_green_av !< average of wghf_eb_green
541    REAL(wp), DIMENSION(:,:), ALLOCATABLE          ::  t_green_av       !< Average of t_green
542    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  wghf_eb_window_av !< average of wghf_eb_window
543    REAL(wp), DIMENSION(:,:), ALLOCATABLE          ::  t_window_av      !< Average of t_window   
544   
545
546!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
547!-- anthropogenic heat sources
548!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
[2920]549    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE        ::  aheat             !< daily average of anthropogenic heat (W/m2)
550    REAL(wp), DIMENSION(:,:), ALLOCATABLE          ::  aheatprof         !< diurnal profiles of anthropogenic heat for particular layers
551    INTEGER(wp)                                    ::  naheatlayers = 1  !< number of layers of anthropogenic heat
[2737]552
553!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
554!-- wall surface model
555!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
556!-- wall surface model constants
557    INTEGER(iwp), PARAMETER                        :: nzb_wall = 0       !< inner side of the wall model (to be switched)
558    INTEGER(iwp), PARAMETER                        :: nzt_wall = 3       !< outer side of the wall model (to be switched)
559    INTEGER(iwp), PARAMETER                        :: nzw = 4            !< number of wall layers (fixed for now)
560
561    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
562                                                                         !< normalized soil, wall and roof layer depths (m/m)
563!    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default = (/0.33_wp, 0.66_wp, 1.0_wp /)
564    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_window = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /)
565!    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_window = (/0.33_wp, 0.66_wp, 1.0_wp /)
566!    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_window = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
567                                                                         !< normalized window layer depths (m/m)
568!    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_green = (/0.0242_wp, 0.0969_wp, 0.346_wp, 1.0_wp /)
569                                                                         !< normalized green layer depths (m/m)
570    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_green = (/0.25_wp, 0.5_wp, 0.75_wp, 1.0_wp /)
571!    REAL(wp), DIMENSION(nzb_wall:nzt_wall)         :: zwn_default_green = (/0.33_wp, 0.66_wp, 1.0_wp /)
572
573
[2920]574    REAL(wp)                                       :: wall_inner_temperature = 295.0_wp    !< temperature of the inner wall surface (~22 degrees C) (K)
575    REAL(wp)                                       :: roof_inner_temperature = 295.0_wp    !< temperature of the inner roof surface (~22 degrees C) (K)
576    REAL(wp)                                       :: soil_inner_temperature = 288.0_wp    !< temperature of the deep soil (~15 degrees C) (K)
577    REAL(wp)                                       :: window_inner_temperature = 295.0_wp  !< temperature of the inner window surface (~22 degrees C) (K)
578
[2737]579!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
580!-- surface and material model variables for walls, ground, roofs
581!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
582    REAL(wp), DIMENSION(:), ALLOCATABLE            :: zwn                !< normalized wall layer depths (m)
583    REAL(wp), DIMENSION(:), ALLOCATABLE            :: zwn_window         !< normalized window layer depths (m)
584    REAL(wp), DIMENSION(:), ALLOCATABLE            :: zwn_green          !< normalized green layer depths (m)
585
586#if defined( __nopointer )
587    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_h           !< wall surface temperature (K) at horizontal walls
588    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_h_p         !< progn. wall surface temperature (K) at horizontal walls
589    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_window_h    !< window surface temperature (K) at horizontal walls
590    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_window_h_p  !< progn. window surface temperature (K) at horizontal walls
591    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_green_h     !< green surface temperature (K) at horizontal walls
592    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_green_h_p   !< progn. green surface temperature (K) at horizontal walls
593    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_10cm_h      !< near surface temperature (10cm) (K) at horizontal walls
594    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_10cm_h_p    !< progn. near surface temperature (10cm) (K) at horizontal walls
595    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_v
596    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_v_p
597    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_window_v
598    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_window_v_p
599    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_green_v
600    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_green_v_p
601    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_10cm_v
602    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  ::  t_surf_10cm_v_p
603#else
604    REAL(wp), DIMENSION(:), POINTER                :: t_surf_h
605    REAL(wp), DIMENSION(:), POINTER                :: t_surf_h_p
606    REAL(wp), DIMENSION(:), POINTER                :: t_surf_window_h
607    REAL(wp), DIMENSION(:), POINTER                :: t_surf_window_h_p
608    REAL(wp), DIMENSION(:), POINTER                :: t_surf_green_h
609    REAL(wp), DIMENSION(:), POINTER                :: t_surf_green_h_p
610    REAL(wp), DIMENSION(:), POINTER                :: t_surf_10cm_h
611    REAL(wp), DIMENSION(:), POINTER                :: t_surf_10cm_h_p
612
613    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_h_1
614    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_h_2
615    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_window_h_1
616    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_window_h_2
617    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_green_h_1
618    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_green_h_2
619    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_10cm_h_1
620    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_10cm_h_2
621
622    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_v
623    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_v_p
624    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_window_v
625    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_window_v_p
626    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_green_v
627    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_green_v_p
628    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_10cm_v
629    TYPE(t_surf_vertical), DIMENSION(:), POINTER ::  t_surf_10cm_v_p
630
631    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_v_1
632    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_v_2
633    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_window_v_1
634    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_window_v_2
635    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_green_v_1
636    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_green_v_2
637    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_10cm_v_1
638    TYPE(t_surf_vertical), DIMENSION(0:3), TARGET  :: t_surf_10cm_v_2
639   
640#endif
641    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_av          !< average of wall surface temperature (K)
642    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_window_av   !< average of window surface temperature (K)
643    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_green_av    !< average of green wall surface temperature (K)
644    REAL(wp), DIMENSION(:), ALLOCATABLE, TARGET    :: t_surf_10cm_av    !< average of whole wall surface temperature (K)
645
646!-- Temporal tendencies for time stepping           
647    REAL(wp), DIMENSION(:), ALLOCATABLE            :: tt_surface_m       !< surface temperature tendency of wall (K)
648    REAL(wp), DIMENSION(:), ALLOCATABLE            :: tt_surface_window_m !< surface temperature tendency of window (K)
649    REAL(wp), DIMENSION(:), ALLOCATABLE            :: tt_surface_green_m !< surface temperature tendency of green wall (K)
650
651!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
652!-- Energy balance variables
653!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
654!-- parameters of the land, roof and wall surfaces
655
656#if defined( __nopointer )
657    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_wall_h             !< Wall temperature (K)
658    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_wall_h_av          !< Average of t_wall
659    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_wall_h_p           !< Prog. wall temperature (K)
660    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_window_h           !< Window temperature (K)
661    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_window_h_av        !< Average of t_window
662    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_window_h_p         !< Prog. window temperature (K)
663    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_green_h            !< Green temperature (K)
664    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_green_h_av         !< Average of t_green
665    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_green_h_p          !< Prog. green temperature (K)
666
667    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_wall_v             !< Wall temperature (K)
668    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_wall_v_av          !< Average of t_wall
669    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_wall_v_p           !< Prog. wall temperature (K)
670    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_window_v           !< Window temperature (K)
671    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_window_v_av        !< Average of t_window
672    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_window_v_p         !< Prog. window temperature (K)
673    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_green_v            !< Green temperature (K)
674    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_green_v_av         !< Average of t_green
675    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_green_v_p          !< Prog. green temperature (K)
676#else
677    REAL(wp), DIMENSION(:,:), POINTER                :: t_wall_h, t_wall_h_p
678    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_wall_h_av, t_wall_h_1, t_wall_h_2
679    REAL(wp), DIMENSION(:,:), POINTER                :: t_window_h, t_window_h_p
680    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_window_h_av, t_window_h_1, t_window_h_2
681    REAL(wp), DIMENSION(:,:), POINTER                :: t_green_h, t_green_h_p
682    REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET    :: t_green_h_av, t_green_h_1, t_green_h_2
683
684    TYPE(t_wall_vertical), DIMENSION(:), POINTER   :: t_wall_v, t_wall_v_p
685    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_wall_v_av, t_wall_v_1, t_wall_v_2
686    TYPE(t_wall_vertical), DIMENSION(:), POINTER   :: t_window_v, t_window_v_p
687    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_window_v_av, t_window_v_1, t_window_v_2
688    TYPE(t_wall_vertical), DIMENSION(:), POINTER   :: t_green_v, t_green_v_p
689    TYPE(t_wall_vertical), DIMENSION(0:3), TARGET  :: t_green_v_av, t_green_v_1, t_green_v_2
690#endif
691
692!-- Wall temporal tendencies for time stepping
693    REAL(wp), DIMENSION(:,:), ALLOCATABLE          :: tt_wall_m          !< t_wall prognostic array
694    REAL(wp), DIMENSION(:,:), ALLOCATABLE          :: tt_window_m        !< t_window prognostic array
695    REAL(wp), DIMENSION(:,:), ALLOCATABLE          :: tt_green_m         !< t_green prognostic array
696
697!-- Surface and material parameters classes (surface_type)
698!-- albedo, emissivity, lambda_surf, roughness, thickness, volumetric heat capacity, thermal conductivity
699    INTEGER(iwp)                                   :: n_surface_types      !< number of the wall type categories
[2920]700    INTEGER(iwp), PARAMETER                        :: n_surface_params = 9 !< number of parameters for each type of the wall
[2737]701    INTEGER(iwp), PARAMETER                        :: ialbedo  = 1         !< albedo of the surface
702    INTEGER(iwp), PARAMETER                        :: iemiss   = 2         !< emissivity of the surface
[2920]703    INTEGER(iwp), PARAMETER                        :: ilambdas = 3         !< heat conductivity lambda S between surface and material ( W m-2 K-1 )
704    INTEGER(iwp), PARAMETER                        :: irough   = 4         !< roughness length z0 for movements
705    INTEGER(iwp), PARAMETER                        :: iroughh  = 5         !< roughness length z0h for scalars (heat, humidity,...)
706    INTEGER(iwp), PARAMETER                        :: icsurf   = 6         !< Surface skin layer heat capacity (J m-2 K-1 )
707    INTEGER(iwp), PARAMETER                        :: ithick   = 7         !< thickness of the surface (wall, roof, land)  ( m )
708    INTEGER(iwp), PARAMETER                        :: irhoC    = 8         !< volumetric heat capacity rho*C of the material ( J m-3 K-1 )
709    INTEGER(iwp), PARAMETER                        :: ilambdah = 9         !< thermal conductivity lambda H of the wall (W m-1 K-1 )
[2737]710    CHARACTER(12), DIMENSION(:), ALLOCATABLE       :: surface_type_names   !< names of wall types (used only for reports)
711    INTEGER(iwp), DIMENSION(:), ALLOCATABLE        :: surface_type_codes   !< codes of wall types
712    REAL(wp), DIMENSION(:,:), ALLOCATABLE          :: surface_params       !< parameters of wall types
713
714   
715!-- interfaces of subroutines accessed from outside of this module
716    INTERFACE usm_boundary_condition
717       MODULE PROCEDURE usm_boundary_condition
718    END INTERFACE usm_boundary_condition
719
720    INTERFACE usm_check_data_output
721       MODULE PROCEDURE usm_check_data_output
722    END INTERFACE usm_check_data_output
723   
724    INTERFACE usm_check_parameters
725       MODULE PROCEDURE usm_check_parameters
726    END INTERFACE usm_check_parameters
727   
728    INTERFACE usm_data_output_3d
729       MODULE PROCEDURE usm_data_output_3d
730    END INTERFACE usm_data_output_3d
731   
732    INTERFACE usm_define_netcdf_grid
733       MODULE PROCEDURE usm_define_netcdf_grid
734    END INTERFACE usm_define_netcdf_grid
735
736    INTERFACE usm_init_urban_surface
737       MODULE PROCEDURE usm_init_urban_surface
738    END INTERFACE usm_init_urban_surface
739
740    INTERFACE usm_material_heat_model
741       MODULE PROCEDURE usm_material_heat_model
742    END INTERFACE usm_material_heat_model
743   
744    INTERFACE usm_green_heat_model
745       MODULE PROCEDURE usm_green_heat_model
746    END INTERFACE usm_green_heat_model
747   
748    INTERFACE usm_parin
749       MODULE PROCEDURE usm_parin
750    END INTERFACE usm_parin
751   
752    INTERFACE usm_temperature_near_surface
753       MODULE PROCEDURE usm_temperature_near_surface
754    END INTERFACE usm_temperature_near_surface
755
[2894]756    INTERFACE usm_rrd_local
757       MODULE PROCEDURE usm_rrd_local
758    END INTERFACE usm_rrd_local
[2737]759
760    INTERFACE usm_surface_energy_balance
761       MODULE PROCEDURE usm_surface_energy_balance
762    END INTERFACE usm_surface_energy_balance
763   
764    INTERFACE usm_swap_timelevel
765       MODULE PROCEDURE usm_swap_timelevel
766    END INTERFACE usm_swap_timelevel
767       
[2894]768    INTERFACE usm_wrd_local
769       MODULE PROCEDURE usm_wrd_local
770    END INTERFACE usm_wrd_local
[2737]771
772    INTERFACE usm_allocate_surface
773       MODULE PROCEDURE usm_allocate_surface
774    END INTERFACE usm_allocate_surface
775
776    INTERFACE usm_average_3d_data
777       MODULE PROCEDURE usm_average_3d_data
778    END INTERFACE usm_average_3d_data
779
780   
781    SAVE
782
783    PRIVATE 
784   
785!-- Public functions
786    PUBLIC usm_boundary_condition, usm_check_parameters, usm_init_urban_surface,&
[2920]787           usm_rrd_local,                                                      & 
788           usm_surface_energy_balance, usm_material_heat_model,                &
789           usm_swap_timelevel, usm_check_data_output, usm_average_3d_data,     &
790           usm_data_output_3d, usm_define_netcdf_grid, usm_parin,              &
[2894]791           usm_wrd_local, usm_allocate_surface
[2737]792
793!-- Public parameters, constants and initial values
[2920]794    PUBLIC usm_anthropogenic_heat, usm_material_model,                          &
[2737]795           usm_green_heat_model, usm_temperature_near_surface
796
797
798
799 CONTAINS
800
801!------------------------------------------------------------------------------!
802! Description:
803! ------------
804!> This subroutine creates the necessary indices of the urban surfaces
805!> and plant canopy and it allocates the needed arrays for USM
806!------------------------------------------------------------------------------!
807    SUBROUTINE usm_allocate_surface
808   
809        IMPLICIT NONE
810       
811        INTEGER(iwp) ::  l
812
813!
814!--     Allocate radiation arrays which are part of the new data type.
815!--     For horizontal surfaces.
816        ALLOCATE( surf_usm_h%surfhf(1:surf_usm_h%ns)    )
817        ALLOCATE( surf_usm_h%rad_net_l(1:surf_usm_h%ns) )
818!
819!--     For vertical surfaces
820        DO  l = 0, 3
821           ALLOCATE( surf_usm_v(l)%surfhf(1:surf_usm_v(l)%ns)    )
822           ALLOCATE( surf_usm_v(l)%rad_net_l(1:surf_usm_v(l)%ns) )
823        ENDDO
824
825!--     Wall surface model
826!--     allocate arrays for wall surface model and define pointers
827       
828!--     allocate array of wall types and wall parameters
829        ALLOCATE ( surf_usm_h%surface_types(1:surf_usm_h%ns) )
830        DO  l = 0, 3
831           ALLOCATE( surf_usm_v(l)%surface_types(1:surf_usm_v(l)%ns) )
832        ENDDO
833!
834!--     Allocate albedo_type and albedo. Each surface element
835!--     has 3 values, 0: wall fraction, 1: green fraction, 2: window fraction.
836        ALLOCATE( surf_usm_h%albedo_type(0:2,1:surf_usm_h%ns) )
837        ALLOCATE( surf_usm_h%albedo(0:2,1:surf_usm_h%ns)      )
838        surf_usm_h%albedo_type = albedo_type
839        DO  l = 0, 3
840           ALLOCATE( surf_usm_v(l)%albedo_type(0:2,1:surf_usm_v(l)%ns) )
841           ALLOCATE( surf_usm_v(l)%albedo(0:2,1:surf_usm_v(l)%ns)      )
842           surf_usm_v(l)%albedo_type = albedo_type
843        ENDDO       
844
845
846!
847!--     Allocate indoor target temperature for summer and winter
848        ALLOCATE( surf_usm_h%target_temp_summer(1:surf_usm_h%ns) )
849        ALLOCATE( surf_usm_h%target_temp_winter(1:surf_usm_h%ns) )
850        DO  l = 0, 3
851           ALLOCATE( surf_usm_v(l)%target_temp_summer(1:surf_usm_v(l)%ns) )
852           ALLOCATE( surf_usm_v(l)%target_temp_winter(1:surf_usm_v(l)%ns) )
853        ENDDO   
854!
855!--     Allocate flag indicating ground floor level surface elements
856        ALLOCATE ( surf_usm_h%ground_level(1:surf_usm_h%ns) ) 
857        DO  l = 0, 3
858           ALLOCATE( surf_usm_v(l)%ground_level(1:surf_usm_v(l)%ns) )
859        ENDDO   
860!
861!--      Allocate arrays for relative surface fraction.
862!--      0 - wall fraction, 1 - green fraction, 2 - window fraction
863         ALLOCATE( surf_usm_h%frac(0:2,1:surf_usm_h%ns) )
864         surf_usm_h%frac = 0.0_wp
865         DO  l = 0, 3
866            ALLOCATE( surf_usm_v(l)%frac(0:2,1:surf_usm_v(l)%ns) )
867            surf_usm_v(l)%frac = 0.0_wp
868         ENDDO
869       
870!--     wall and roof surface parameters. First for horizontal surfaces
871        ALLOCATE ( surf_usm_h%isroof_surf(1:surf_usm_h%ns)     )
872        ALLOCATE ( surf_usm_h%lambda_surf(1:surf_usm_h%ns)     )
873        ALLOCATE ( surf_usm_h%lambda_surf_window(1:surf_usm_h%ns) )
874        ALLOCATE ( surf_usm_h%lambda_surf_green(1:surf_usm_h%ns)  )
875        ALLOCATE ( surf_usm_h%c_surface(1:surf_usm_h%ns)       )
876        ALLOCATE ( surf_usm_h%c_surface_window(1:surf_usm_h%ns)   )
877        ALLOCATE ( surf_usm_h%c_surface_green(1:surf_usm_h%ns)    )
878        ALLOCATE ( surf_usm_h%transmissivity(1:surf_usm_h%ns)  )
879        ALLOCATE ( surf_usm_h%lai(1:surf_usm_h%ns)             )
880        ALLOCATE ( surf_usm_h%emissivity(0:2,1:surf_usm_h%ns)  )
881        ALLOCATE ( surf_usm_h%r_a(1:surf_usm_h%ns)             )
882        ALLOCATE ( surf_usm_h%r_a_green(1:surf_usm_h%ns)       )
883        ALLOCATE ( surf_usm_h%r_a_window(1:surf_usm_h%ns)      )
884
885!
886!--     For vertical surfaces.
887        DO  l = 0, 3
888           ALLOCATE ( surf_usm_v(l)%lambda_surf(1:surf_usm_v(l)%ns)     )
889           ALLOCATE ( surf_usm_v(l)%c_surface(1:surf_usm_v(l)%ns)       )
890           ALLOCATE ( surf_usm_v(l)%lambda_surf_window(1:surf_usm_v(l)%ns) )
891           ALLOCATE ( surf_usm_v(l)%c_surface_window(1:surf_usm_v(l)%ns)   )
892           ALLOCATE ( surf_usm_v(l)%lambda_surf_green(1:surf_usm_v(l)%ns)  )
893           ALLOCATE ( surf_usm_v(l)%c_surface_green(1:surf_usm_v(l)%ns)    )
894           ALLOCATE ( surf_usm_v(l)%transmissivity(1:surf_usm_v(l)%ns)  )
895           ALLOCATE ( surf_usm_v(l)%lai(1:surf_usm_v(l)%ns)             )
896           ALLOCATE ( surf_usm_v(l)%emissivity(0:2,1:surf_usm_v(l)%ns)  )
897           ALLOCATE ( surf_usm_v(l)%r_a(1:surf_usm_v(l)%ns)             )
898           ALLOCATE ( surf_usm_v(l)%r_a_green(1:surf_usm_v(l)%ns)       )
899           ALLOCATE ( surf_usm_v(l)%r_a_window(1:surf_usm_v(l)%ns)      )
900        ENDDO
901
902!       
903!--     allocate wall and roof material parameters. First for horizontal surfaces
904        ALLOCATE ( surf_usm_h%thickness_wall(1:surf_usm_h%ns)               )
905        ALLOCATE ( surf_usm_h%thickness_window(1:surf_usm_h%ns)                  )
906        ALLOCATE ( surf_usm_h%thickness_green(1:surf_usm_h%ns)                   )
907        ALLOCATE ( surf_usm_h%lambda_h(nzb_wall:nzt_wall,1:surf_usm_h%ns)   )
908        ALLOCATE ( surf_usm_h%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
909        ALLOCATE ( surf_usm_h%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
910        ALLOCATE ( surf_usm_h%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_h%ns)    )
911        ALLOCATE ( surf_usm_h%lambda_h_green(nzb_wall:nzt_wall,1:surf_usm_h%ns)  )
912        ALLOCATE ( surf_usm_h%rho_c_green(nzb_wall:nzt_wall,1:surf_usm_h%ns)     )
913
914!
915!--     For vertical surfaces.
916        DO  l = 0, 3
917           ALLOCATE ( surf_usm_v(l)%thickness_wall(1:surf_usm_v(l)%ns)               )
918           ALLOCATE ( surf_usm_v(l)%thickness_window(1:surf_usm_v(l)%ns)                  )
919           ALLOCATE ( surf_usm_v(l)%thickness_green(1:surf_usm_v(l)%ns)                   )
920           ALLOCATE ( surf_usm_v(l)%lambda_h(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)   )
921           ALLOCATE ( surf_usm_v(l)%rho_c_wall(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
922           ALLOCATE ( surf_usm_v(l)%lambda_h_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
923           ALLOCATE ( surf_usm_v(l)%rho_c_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)    )
924           ALLOCATE ( surf_usm_v(l)%lambda_h_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)  )
925           ALLOCATE ( surf_usm_v(l)%rho_c_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)     )
926        ENDDO
927
928!--     allocate wall and roof layers sizes. For horizontal surfaces.
929        ALLOCATE ( zwn(nzb_wall:nzt_wall) )
930        ALLOCATE ( surf_usm_h%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)     )
931        ALLOCATE ( zwn_window(nzb_wall:nzt_wall) )
932        ALLOCATE ( surf_usm_h%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)     )
933        ALLOCATE ( zwn_green(nzb_wall:nzt_wall) )
934        ALLOCATE ( surf_usm_h%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)      )
935        ALLOCATE ( surf_usm_h%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)    )
936        ALLOCATE ( surf_usm_h%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns)  )
937        ALLOCATE ( surf_usm_h%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
938        ALLOCATE ( surf_usm_h%zw(nzb_wall:nzt_wall,1:surf_usm_h%ns)            )
939        ALLOCATE ( surf_usm_h%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)    )
940        ALLOCATE ( surf_usm_h%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns)  )
941        ALLOCATE ( surf_usm_h%ddz_window_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
942        ALLOCATE ( surf_usm_h%zw_window(nzb_wall:nzt_wall,1:surf_usm_h%ns)       )
943        ALLOCATE ( surf_usm_h%ddz_green(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)     )
944        ALLOCATE ( surf_usm_h%dz_green_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns)   )
945        ALLOCATE ( surf_usm_h%ddz_green_stag(nzb_wall:nzt_wall,1:surf_usm_h%ns)  )
946        ALLOCATE ( surf_usm_h%zw_green(nzb_wall:nzt_wall,1:surf_usm_h%ns)        )
947!
948!--     For vertical surfaces.
949        DO  l = 0, 3
950           ALLOCATE ( surf_usm_v(l)%dz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)     )
951           ALLOCATE ( surf_usm_v(l)%dz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)     )
952           ALLOCATE ( surf_usm_v(l)%dz_green(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)      )
953           ALLOCATE ( surf_usm_v(l)%ddz_wall(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)    )
954           ALLOCATE ( surf_usm_v(l)%dz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)  )
955           ALLOCATE ( surf_usm_v(l)%ddz_wall_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
956           ALLOCATE ( surf_usm_v(l)%zw(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)            )
957           ALLOCATE ( surf_usm_v(l)%ddz_window(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)    )
958           ALLOCATE ( surf_usm_v(l)%dz_window_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)  )
959           ALLOCATE ( surf_usm_v(l)%ddz_window_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
960           ALLOCATE ( surf_usm_v(l)%zw_window(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)       )
961           ALLOCATE ( surf_usm_v(l)%ddz_green(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)     )
962           ALLOCATE ( surf_usm_v(l)%dz_green_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)   )
963           ALLOCATE ( surf_usm_v(l)%ddz_green_stag(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)  )
964           ALLOCATE ( surf_usm_v(l)%zw_green(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns)        )
965        ENDDO
966
967!--     allocate wall and roof temperature arrays, for horizontal walls
968#if defined( __nopointer )
969        IF ( .NOT. ALLOCATED( t_surf_h ) )                                     &
970           ALLOCATE ( t_surf_h(1:surf_usm_h%ns) )
971        IF ( .NOT. ALLOCATED( t_surf_h_p ) )                                   &
972           ALLOCATE ( t_surf_h_p(1:surf_usm_h%ns) )
973        IF ( .NOT. ALLOCATED( t_wall_h ) )                                     &           
974           ALLOCATE ( t_wall_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
975        IF ( .NOT. ALLOCATED( t_wall_h_p ) )                                   &           
976           ALLOCATE ( t_wall_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
977        IF ( .NOT. ALLOCATED( t_surf_window_h ) )                              &
978           ALLOCATE ( t_surf_window_h(1:surf_usm_h%ns) )
979        IF ( .NOT. ALLOCATED( t_surf_window_h_p ) )                            &
980           ALLOCATE ( t_surf_window_h_p(1:surf_usm_h%ns) )
981        IF ( .NOT. ALLOCATED( t_window_h ) )                                   &           
982           ALLOCATE ( t_window_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
983        IF ( .NOT. ALLOCATED( t_window_h_p ) )                                 &           
984           ALLOCATE ( t_window_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
985        IF ( .NOT. ALLOCATED( t_surf_green_h ) )                               &
986           ALLOCATE ( t_surf_green_h(1:surf_usm_h%ns) )
987        IF ( .NOT. ALLOCATED( t_surf_green_h_p ) )                             &
988           ALLOCATE ( t_surf_green_h_p(1:surf_usm_h%ns) )           
989        IF ( .NOT. ALLOCATED( t_green_h ) )                                    &           
990           ALLOCATE ( t_green_h(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
991        IF ( .NOT. ALLOCATED( t_green_h_p ) )                                  &           
992           ALLOCATE ( t_green_h_p(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
993        IF ( .NOT. ALLOCATED( t_surf_10cm_h ) )                                &
994           ALLOCATE ( t_surf_10cm_h(1:surf_usm_h%ns) )
995        IF ( .NOT. ALLOCATED( t_surf_10cm_h_p ) )                              &
996           ALLOCATE ( t_surf_10cm_h_p(1:surf_usm_h%ns) )
997#else
998!
999!--     Allocate if required. Note, in case of restarts, some of these arrays
1000!--     might be already allocated.
1001        IF ( .NOT. ALLOCATED( t_surf_h_1 ) )                                   &
1002           ALLOCATE ( t_surf_h_1(1:surf_usm_h%ns) )
1003        IF ( .NOT. ALLOCATED( t_surf_h_2 ) )                                   &
1004           ALLOCATE ( t_surf_h_2(1:surf_usm_h%ns) )
1005        IF ( .NOT. ALLOCATED( t_wall_h_1 ) )                                   &           
1006           ALLOCATE ( t_wall_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
1007        IF ( .NOT. ALLOCATED( t_wall_h_2 ) )                                   &           
1008           ALLOCATE ( t_wall_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )         
1009        IF ( .NOT. ALLOCATED( t_surf_window_h_1 ) )                            &
1010           ALLOCATE ( t_surf_window_h_1(1:surf_usm_h%ns) )
1011        IF ( .NOT. ALLOCATED( t_surf_window_h_2 ) )                            &
1012           ALLOCATE ( t_surf_window_h_2(1:surf_usm_h%ns) )
1013        IF ( .NOT. ALLOCATED( t_window_h_1 ) )                                 &           
1014           ALLOCATE ( t_window_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
1015        IF ( .NOT. ALLOCATED( t_window_h_2 ) )                                 &           
1016           ALLOCATE ( t_window_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )         
1017        IF ( .NOT. ALLOCATED( t_surf_green_h_1 ) )                             &
1018           ALLOCATE ( t_surf_green_h_1(1:surf_usm_h%ns) )
1019        IF ( .NOT. ALLOCATED( t_surf_green_h_2 ) )                             &
1020           ALLOCATE ( t_surf_green_h_2(1:surf_usm_h%ns) )
1021        IF ( .NOT. ALLOCATED( t_green_h_1 ) )                                  &           
1022           ALLOCATE ( t_green_h_1(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) ) 
1023        IF ( .NOT. ALLOCATED( t_green_h_2 ) )                                  &           
1024           ALLOCATE ( t_green_h_2(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )         
1025        IF ( .NOT. ALLOCATED( t_surf_10cm_h_1 ) )                              &
1026           ALLOCATE ( t_surf_10cm_h_1(1:surf_usm_h%ns) )
1027        IF ( .NOT. ALLOCATED( t_surf_10cm_h_2 ) )                              &
1028           ALLOCATE ( t_surf_10cm_h_2(1:surf_usm_h%ns) )
1029!           
1030!--     initial assignment of the pointers
1031        t_wall_h    => t_wall_h_1;    t_wall_h_p    => t_wall_h_2
1032        t_window_h    => t_window_h_1;    t_window_h_p    => t_window_h_2
1033        t_green_h    => t_green_h_1;    t_green_h_p    => t_green_h_2
1034        t_surf_h => t_surf_h_1; t_surf_h_p => t_surf_h_2           
1035        t_surf_window_h => t_surf_window_h_1; t_surf_window_h_p => t_surf_window_h_2 
1036        t_surf_green_h => t_surf_green_h_1; t_surf_green_h_p => t_surf_green_h_2           
1037        t_surf_10cm_h => t_surf_10cm_h_1; t_surf_10cm_h_p => t_surf_10cm_h_2 
1038 
1039#endif
1040
1041!--     allocate wall and roof temperature arrays, for vertical walls if required
1042#if defined( __nopointer )
1043        DO  l = 0, 3
1044           IF ( .NOT. ALLOCATED( t_surf_v(l)%t ) )                             &
1045              ALLOCATE ( t_surf_v(l)%t(1:surf_usm_v(l)%ns) )
1046           IF ( .NOT. ALLOCATED( t_surf_v_p(l)%t ) )                           &
1047              ALLOCATE ( t_surf_v_p(l)%t(1:surf_usm_v(l)%ns) )
1048           IF ( .NOT. ALLOCATED( t_wall_v(l)%t ) )                             &
1049              ALLOCATE ( t_wall_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1050           IF ( .NOT. ALLOCATED( t_wall_v_p(l)%t ) )                           &                 
1051              ALLOCATE ( t_wall_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1052           IF ( .NOT. ALLOCATED( t_surf_window_v(l)%t ) )                      &
1053              ALLOCATE ( t_surf_window_v(l)%t(1:surf_usm_v(l)%ns) )
1054           IF ( .NOT. ALLOCATED( t_surf_window_v_p(l)%t ) )                    &
1055              ALLOCATE ( t_surf_window_v_p(l)%t(1:surf_usm_v(l)%ns) )
1056           IF ( .NOT. ALLOCATED( t_window_v(l)%t ) )                           &
1057              ALLOCATE ( t_window_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1058           IF ( .NOT. ALLOCATED( t_window_v_p(l)%t ) )                         &                 
1059              ALLOCATE ( t_window_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1060           IF ( .NOT. ALLOCATED( t_green_v(l)%t ) )                            &
1061              ALLOCATE ( t_green_v(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1062           IF ( .NOT. ALLOCATED( t_green_v_p(l)%t ) )                          &                 
1063              ALLOCATE ( t_green_v_p(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1064           IF ( .NOT. ALLOCATED( t_surf_green_v(l)%t ) )                       &
1065              ALLOCATE ( t_surf_green_v(l)%t(1:surf_usm_v(l)%ns) )
1066           IF ( .NOT. ALLOCATED( t_surf_green_v_p(l)%t ) )                     &
1067              ALLOCATE ( t_surf_green_v_p(l)%t(1:surf_usm_v(l)%ns) )
1068           IF ( .NOT. ALLOCATED( t_surf_10cm_v(l)%t ) )                        &
1069              ALLOCATE ( t_surf_10cm_v(l)%t(1:surf_usm_v(l)%ns) )
1070           IF ( .NOT. ALLOCATED( t_surf_10cm_v_p(l)%t ) )                        &
1071              ALLOCATE ( t_surf_10cm_v_p(l)%t(1:surf_usm_v(l)%ns) )
1072        ENDDO
1073#else
1074!
1075!--     Allocate if required. Note, in case of restarts, some of these arrays
1076!--     might be already allocated.
1077        DO  l = 0, 3
1078           IF ( .NOT. ALLOCATED( t_surf_v_1(l)%t ) )                           &
1079              ALLOCATE ( t_surf_v_1(l)%t(1:surf_usm_v(l)%ns) )
1080           IF ( .NOT. ALLOCATED( t_surf_v_2(l)%t ) )                           &
1081              ALLOCATE ( t_surf_v_2(l)%t(1:surf_usm_v(l)%ns) )
1082           IF ( .NOT. ALLOCATED( t_wall_v_1(l)%t ) )                           &           
1083              ALLOCATE ( t_wall_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1084           IF ( .NOT. ALLOCATED( t_wall_v_2(l)%t ) )                           &           
1085              ALLOCATE ( t_wall_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1086           IF ( .NOT. ALLOCATED( t_surf_window_v_1(l)%t ) )                    &
1087              ALLOCATE ( t_surf_window_v_1(l)%t(1:surf_usm_v(l)%ns) )
1088           IF ( .NOT. ALLOCATED( t_surf_window_v_2(l)%t ) )                    &
1089              ALLOCATE ( t_surf_window_v_2(l)%t(1:surf_usm_v(l)%ns) )
1090           IF ( .NOT. ALLOCATED( t_window_v_1(l)%t ) )                         &           
1091              ALLOCATE ( t_window_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1092           IF ( .NOT. ALLOCATED( t_window_v_2(l)%t ) )                         &           
1093              ALLOCATE ( t_window_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1094           IF ( .NOT. ALLOCATED( t_surf_green_v_1(l)%t ) )                     &
1095              ALLOCATE ( t_surf_green_v_1(l)%t(1:surf_usm_v(l)%ns) )
1096           IF ( .NOT. ALLOCATED( t_surf_green_v_2(l)%t ) )                     &
1097              ALLOCATE ( t_surf_green_v_2(l)%t(1:surf_usm_v(l)%ns) )
1098           IF ( .NOT. ALLOCATED( t_green_v_1(l)%t ) )                          &           
1099              ALLOCATE ( t_green_v_1(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1100           IF ( .NOT. ALLOCATED( t_green_v_2(l)%t ) )                          &           
1101              ALLOCATE ( t_green_v_2(l)%t(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) ) 
1102           IF ( .NOT. ALLOCATED( t_surf_10cm_v_1(l)%t ) )                     &
1103              ALLOCATE ( t_surf_10cm_v_1(l)%t(1:surf_usm_v(l)%ns) )
1104           IF ( .NOT. ALLOCATED( t_surf_10cm_v_2(l)%t ) )                     &
1105              ALLOCATE ( t_surf_10cm_v_2(l)%t(1:surf_usm_v(l)%ns) )
1106        ENDDO
1107!
1108!--     initial assignment of the pointers
1109        t_wall_v    => t_wall_v_1;    t_wall_v_p    => t_wall_v_2
1110        t_surf_v => t_surf_v_1; t_surf_v_p => t_surf_v_2
1111        t_window_v    => t_window_v_1;    t_window_v_p    => t_window_v_2
1112        t_green_v    => t_green_v_1;    t_green_v_p    => t_green_v_2
1113        t_surf_window_v => t_surf_window_v_1; t_surf_window_v_p => t_surf_window_v_2
1114        t_surf_green_v => t_surf_green_v_1; t_surf_green_v_p => t_surf_green_v_2
1115        t_surf_10cm_v => t_surf_10cm_v_1; t_surf_10cm_v_p => t_surf_10cm_v_2
1116
1117#endif
1118!
1119!--     Allocate intermediate timestep arrays. For horizontal surfaces.
1120        ALLOCATE ( surf_usm_h%tt_surface_m(1:surf_usm_h%ns)                  )
1121        ALLOCATE ( surf_usm_h%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1122        ALLOCATE ( surf_usm_h%tt_surface_window_m(1:surf_usm_h%ns)             )
1123        ALLOCATE ( surf_usm_h%tt_window_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns) )
1124        ALLOCATE ( surf_usm_h%tt_green_m(nzb_wall:nzt_wall+1,1:surf_usm_h%ns)  )
1125        ALLOCATE ( surf_usm_h%tt_surface_green_m(1:surf_usm_h%ns)              )
1126
1127!
1128!--     Set inital values for prognostic quantities
1129        IF ( ALLOCATED( surf_usm_h%tt_surface_m ) )  surf_usm_h%tt_surface_m = 0.0_wp
1130        IF ( ALLOCATED( surf_usm_h%tt_wall_m    ) )  surf_usm_h%tt_wall_m    = 0.0_wp
1131        IF ( ALLOCATED( surf_usm_h%tt_surface_window_m ) )  surf_usm_h%tt_surface_window_m = 0.0_wp
1132        IF ( ALLOCATED( surf_usm_h%tt_window_m    )      )  surf_usm_h%tt_window_m         = 0.0_wp
1133        IF ( ALLOCATED( surf_usm_h%tt_green_m    )       )  surf_usm_h%tt_green_m          = 0.0_wp
1134        IF ( ALLOCATED( surf_usm_h%tt_surface_green_m )  )  surf_usm_h%tt_surface_green_m  = 0.0_wp
1135!
1136!--     Now, for vertical surfaces
1137        DO  l = 0, 3
1138           ALLOCATE ( surf_usm_v(l)%tt_surface_m(1:surf_usm_v(l)%ns)                  )
1139           ALLOCATE ( surf_usm_v(l)%tt_wall_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1140           IF ( ALLOCATED( surf_usm_v(l)%tt_surface_m ) )  surf_usm_v(l)%tt_surface_m = 0.0_wp
1141           IF ( ALLOCATED( surf_usm_v(l)%tt_wall_m    ) )  surf_usm_v(l)%tt_wall_m    = 0.0_wp
1142           ALLOCATE ( surf_usm_v(l)%tt_surface_window_m(1:surf_usm_v(l)%ns)             )
1143           ALLOCATE ( surf_usm_v(l)%tt_window_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns) )
1144           IF ( ALLOCATED( surf_usm_v(l)%tt_surface_window_m ) )  surf_usm_v(l)%tt_surface_window_m = 0.0_wp
1145           IF ( ALLOCATED( surf_usm_v(l)%tt_window_m  ) )  surf_usm_v(l)%tt_window_m    = 0.0_wp
1146           ALLOCATE ( surf_usm_v(l)%tt_surface_green_m(1:surf_usm_v(l)%ns)              )
1147           IF ( ALLOCATED( surf_usm_v(l)%tt_surface_green_m ) )  surf_usm_v(l)%tt_surface_green_m = 0.0_wp
1148           ALLOCATE ( surf_usm_v(l)%tt_green_m(nzb_wall:nzt_wall+1,1:surf_usm_v(l)%ns)  )
1149           IF ( ALLOCATED( surf_usm_v(l)%tt_green_m   ) )  surf_usm_v(l)%tt_green_m    = 0.0_wp
1150        ENDDO
1151
1152!--     allocate wall heat flux output array and set initial values. For horizontal surfaces
1153!         ALLOCATE ( surf_usm_h%wshf(1:surf_usm_h%ns)    )  !can be removed
1154        ALLOCATE ( surf_usm_h%wshf_eb(1:surf_usm_h%ns) )
1155        ALLOCATE ( surf_usm_h%wghf_eb(1:surf_usm_h%ns) )
1156        ALLOCATE ( surf_usm_h%wghf_eb_window(1:surf_usm_h%ns) )
1157        ALLOCATE ( surf_usm_h%wghf_eb_green(1:surf_usm_h%ns) )
1158        ALLOCATE ( surf_usm_h%iwghf_eb(1:surf_usm_h%ns) )
1159        ALLOCATE ( surf_usm_h%iwghf_eb_window(1:surf_usm_h%ns) )
1160        IF ( ALLOCATED( surf_usm_h%wshf    ) )  surf_usm_h%wshf    = 0.0_wp
1161        IF ( ALLOCATED( surf_usm_h%wshf_eb ) )  surf_usm_h%wshf_eb = 0.0_wp
1162        IF ( ALLOCATED( surf_usm_h%wghf_eb ) )  surf_usm_h%wghf_eb = 0.0_wp
1163        IF ( ALLOCATED( surf_usm_h%wghf_eb_window ) )  surf_usm_h%wghf_eb_window = 0.0_wp
1164        IF ( ALLOCATED( surf_usm_h%wghf_eb_green ) )  surf_usm_h%wghf_eb_green = 0.0_wp
1165        IF ( ALLOCATED( surf_usm_h%iwghf_eb ) )  surf_usm_h%iwghf_eb = 0.0_wp
1166        IF ( ALLOCATED( surf_usm_h%iwghf_eb_window ) )  surf_usm_h%iwghf_eb_window = 0.0_wp
1167!
1168!--     Now, for vertical surfaces
1169        DO  l = 0, 3
1170!            ALLOCATE ( surf_usm_v(l)%wshf(1:surf_usm_v(l)%ns)    )    ! can be removed
1171           ALLOCATE ( surf_usm_v(l)%wshf_eb(1:surf_usm_v(l)%ns) )
1172           ALLOCATE ( surf_usm_v(l)%wghf_eb(1:surf_usm_v(l)%ns) )
1173           ALLOCATE ( surf_usm_v(l)%wghf_eb_window(1:surf_usm_v(l)%ns) )
1174           ALLOCATE ( surf_usm_v(l)%wghf_eb_green(1:surf_usm_v(l)%ns) )
1175           ALLOCATE ( surf_usm_v(l)%iwghf_eb(1:surf_usm_v(l)%ns) )
1176           ALLOCATE ( surf_usm_v(l)%iwghf_eb_window(1:surf_usm_v(l)%ns) )
1177           IF ( ALLOCATED( surf_usm_v(l)%wshf    ) )  surf_usm_v(l)%wshf    = 0.0_wp
1178           IF ( ALLOCATED( surf_usm_v(l)%wshf_eb ) )  surf_usm_v(l)%wshf_eb = 0.0_wp
1179           IF ( ALLOCATED( surf_usm_v(l)%wghf_eb ) )  surf_usm_v(l)%wghf_eb = 0.0_wp
1180           IF ( ALLOCATED( surf_usm_v(l)%wghf_eb_window ) )  surf_usm_v(l)%wghf_eb_window = 0.0_wp
1181           IF ( ALLOCATED( surf_usm_v(l)%wghf_eb_green ) )  surf_usm_v(l)%wghf_eb_green = 0.0_wp
1182           IF ( ALLOCATED( surf_usm_v(l)%iwghf_eb ) )  surf_usm_v(l)%iwghf_eb = 0.0_wp
1183           IF ( ALLOCATED( surf_usm_v(l)%iwghf_eb_window ) )  surf_usm_v(l)%iwghf_eb_window = 0.0_wp
1184        ENDDO
1185       
1186    END SUBROUTINE usm_allocate_surface
1187
1188
1189!------------------------------------------------------------------------------!
1190! Description:
1191! ------------
1192!> Sum up and time-average urban surface output quantities as well as allocate
1193!> the array necessary for storing the average.
1194!------------------------------------------------------------------------------!
1195    SUBROUTINE usm_average_3d_data( mode, variable )
1196
1197        IMPLICIT NONE
1198
1199        CHARACTER (len=*), INTENT(IN) ::  mode
1200        CHARACTER (len=*), INTENT(IN) :: variable
1201 
[2920]1202        INTEGER(iwp)                                       :: i, j, k, l, m, ids, idsint, iwl, istat
[2737]1203        CHARACTER (len=varnamelength)                      :: var, surfid
1204        INTEGER(iwp), PARAMETER                            :: nd = 5
1205        CHARACTER(len=6), DIMENSION(0:nd-1), PARAMETER     :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /)
[2920]1206        INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER         :: dirint = (/ iup_u, isouth_u, inorth_u, iwest_u, ieast_u /)
[2737]1207
1208!--     find the real name of the variable
[2906]1209        ids = -1
[2737]1210        var = TRIM(variable)
1211        DO i = 0, nd-1
1212            k = len(TRIM(var))
1213            j = len(TRIM(dirname(i)))
1214            IF ( var(k-j+1:k) == dirname(i) )  THEN
1215                ids = i
[2920]1216                idsint = dirint(ids)
[2737]1217                var = var(:k-j)
1218                EXIT
1219            ENDIF
1220        ENDDO
1221        IF ( ids == -1 )  THEN
1222            var = TRIM(variable)
1223        ENDIF
1224        IF ( var(1:11) == 'usm_t_wall_'  .AND.  len(TRIM(var)) >= 12 )  THEN
1225!--          wall layers
1226            READ(var(12:12), '(I1)', iostat=istat ) iwl
1227            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
1228                var = var(1:10)
1229            ELSE
1230!--             wrong wall layer index
1231                RETURN
1232            ENDIF
1233        ENDIF
1234        IF ( var(1:13) == 'usm_t_window_'  .AND.  len(TRIM(var)) >= 14 )  THEN
1235!--          wall layers
1236            READ(var(14:14), '(I1)', iostat=istat ) iwl
1237            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
1238                var = var(1:12)
1239            ELSE
1240!--             wrong window layer index
1241                RETURN
1242            ENDIF
1243        ENDIF
1244        IF ( var(1:12) == 'usm_t_green_'  .AND.  len(TRIM(var)) >= 13 )  THEN
1245!--          wall layers
1246            READ(var(13:13), '(I1)', iostat=istat ) iwl
1247            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
1248                var = var(1:11)
1249            ELSE
1250!--             wrong green layer index
1251                RETURN
1252            ENDIF
1253        ENDIF
1254
1255        IF ( mode == 'allocate' )  THEN
1256           
1257           SELECT CASE ( TRIM( var ) )
1258               
1259                CASE ( 'usm_rad_net' )
1260!--                 array of complete radiation balance
1261                    IF ( .NOT.  ALLOCATED(surf_usm_h%rad_net_av) )  THEN
1262                        ALLOCATE( surf_usm_h%rad_net_av(1:surf_usm_h%ns) )
1263                        surf_usm_h%rad_net_av = 0.0_wp
1264                    ENDIF
1265                    DO  l = 0, 3
1266                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%rad_net_av) )  THEN
1267                           ALLOCATE( surf_usm_v(l)%rad_net_av(1:surf_usm_v(l)%ns) )
1268                           surf_usm_v(l)%rad_net_av = 0.0_wp
1269                       ENDIF
1270                    ENDDO
1271                   
1272                CASE ( 'usm_rad_insw' )
1273!--                 array of sw radiation falling to surface after i-th reflection
[2920]1274                    IF ( .NOT.  ALLOCATED(surfinsw_av) )  THEN
1275                        ALLOCATE( surfinsw_av(nsurfl) )
1276                        surfinsw_av = 0.0_wp
[2737]1277                    ENDIF
[2920]1278
[2737]1279                CASE ( 'usm_rad_inlw' )
1280!--                 array of lw radiation falling to surface after i-th reflection
[2920]1281                    IF ( .NOT.  ALLOCATED(surfinlw_av) )  THEN
1282                        ALLOCATE( surfinlw_av(nsurfl) )
1283                        surfinlw_av = 0.0_wp
[2737]1284                    ENDIF
1285
1286                CASE ( 'usm_rad_inswdir' )
1287!--                 array of direct sw radiation falling to surface from sun
1288                    IF ( .NOT.  ALLOCATED(surfinswdir_av) )  THEN
[2920]1289                        ALLOCATE( surfinswdir_av(nsurfl) )
[2737]1290                        surfinswdir_av = 0.0_wp
1291                    ENDIF
1292
1293                CASE ( 'usm_rad_inswdif' )
1294!--                 array of difusion sw radiation falling to surface from sky and borders of the domain
1295                    IF ( .NOT.  ALLOCATED(surfinswdif_av) )  THEN
[2920]1296                        ALLOCATE( surfinswdif_av(nsurfl) )
[2737]1297                        surfinswdif_av = 0.0_wp
1298                    ENDIF
1299
1300                CASE ( 'usm_rad_inswref' )
1301!--                 array of sw radiation falling to surface from reflections
1302                    IF ( .NOT.  ALLOCATED(surfinswref_av) )  THEN
[2920]1303                        ALLOCATE( surfinswref_av(nsurfl) )
[2737]1304                        surfinswref_av = 0.0_wp
1305                    ENDIF
1306
1307                CASE ( 'usm_rad_inlwdif' )
1308!--                 array of sw radiation falling to surface after i-th reflection
1309                   IF ( .NOT.  ALLOCATED(surfinlwdif_av) )  THEN
[2920]1310                        ALLOCATE( surfinlwdif_av(nsurfl) )
[2737]1311                        surfinlwdif_av = 0.0_wp
1312                    ENDIF
1313
1314                CASE ( 'usm_rad_inlwref' )
1315!--                 array of lw radiation falling to surface from reflections
1316                    IF ( .NOT.  ALLOCATED(surfinlwref_av) )  THEN
[2920]1317                        ALLOCATE( surfinlwref_av(nsurfl) )
[2737]1318                        surfinlwref_av = 0.0_wp
1319                    ENDIF
1320
1321                CASE ( 'usm_rad_outsw' )
1322!--                 array of sw radiation emitted from surface after i-th reflection
1323                    IF ( .NOT.  ALLOCATED(surfoutsw_av) )  THEN
[2920]1324                        ALLOCATE( surfoutsw_av(nsurfl) )
[2737]1325                        surfoutsw_av = 0.0_wp
1326                    ENDIF
1327
1328                CASE ( 'usm_rad_outlw' )
1329!--                 array of lw radiation emitted from surface after i-th reflection
1330                    IF ( .NOT.  ALLOCATED(surfoutlw_av) )  THEN
[2920]1331                        ALLOCATE( surfoutlw_av(nsurfl) )
[2737]1332                        surfoutlw_av = 0.0_wp
1333                    ENDIF
1334                CASE ( 'usm_rad_ressw' )
1335!--                 array of residua of sw radiation absorbed in surface after last reflection
1336                    IF ( .NOT.  ALLOCATED(surfins_av) )  THEN
[2920]1337                        ALLOCATE( surfins_av(nsurfl) )
[2737]1338                        surfins_av = 0.0_wp
1339                    ENDIF
1340                                   
1341                CASE ( 'usm_rad_reslw' )
1342!--                 array of residua of lw radiation absorbed in surface after last reflection
1343                    IF ( .NOT.  ALLOCATED(surfinl_av) )  THEN
[2920]1344                        ALLOCATE( surfinl_av(nsurfl) )
[2737]1345                        surfinl_av = 0.0_wp
1346                    ENDIF
1347                                   
1348                CASE ( 'usm_rad_hf' )
1349!--                 array of heat flux from radiation for surfaces after i-th reflection
1350                    IF ( .NOT.  ALLOCATED(surf_usm_h%surfhf_av) )  THEN
1351                        ALLOCATE( surf_usm_h%surfhf_av(1:surf_usm_h%ns) )
1352                        surf_usm_h%surfhf_av = 0.0_wp
1353                    ENDIF
1354                    DO  l = 0, 3
1355                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%surfhf_av) )  THEN
1356                           ALLOCATE( surf_usm_v(l)%surfhf_av(1:surf_usm_v(l)%ns) )
1357                           surf_usm_v(l)%surfhf_av = 0.0_wp
1358                       ENDIF
1359                    ENDDO
1360
1361                CASE ( 'usm_wshf' )
1362!--                 array of sensible heat flux from surfaces
1363!--                 land surfaces
1364                    IF ( .NOT.  ALLOCATED(surf_usm_h%wshf_eb_av) )  THEN
1365                        ALLOCATE( surf_usm_h%wshf_eb_av(1:surf_usm_h%ns) )
1366                        surf_usm_h%wshf_eb_av = 0.0_wp
1367                    ENDIF
1368                    DO  l = 0, 3
1369                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%wshf_eb_av) )  THEN
1370                           ALLOCATE( surf_usm_v(l)%wshf_eb_av(1:surf_usm_v(l)%ns) )
1371                           surf_usm_v(l)%wshf_eb_av = 0.0_wp
1372                       ENDIF
1373                    ENDDO
[2797]1374!
1375!--             Please note, the following output quantities belongs to the
1376!--             individual tile fractions - ground heat flux at wall-, window-,
1377!--             and green fraction. Aggregated ground-heat flux is treated
1378!--             accordingly in average_3d_data, sum_up_3d_data, etc..
[2737]1379                CASE ( 'usm_wghf' )
1380!--                 array of heat flux from ground (wall, roof, land)
1381                    IF ( .NOT.  ALLOCATED(surf_usm_h%wghf_eb_av) )  THEN
1382                        ALLOCATE( surf_usm_h%wghf_eb_av(1:surf_usm_h%ns) )
1383                        surf_usm_h%wghf_eb_av = 0.0_wp
1384                    ENDIF
1385                    DO  l = 0, 3
1386                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%wghf_eb_av) )  THEN
1387                           ALLOCATE( surf_usm_v(l)%wghf_eb_av(1:surf_usm_v(l)%ns) )
1388                           surf_usm_v(l)%wghf_eb_av = 0.0_wp
1389                       ENDIF
1390                    ENDDO
1391
1392                CASE ( 'usm_wghf_window' )
1393!--                 array of heat flux from window ground (wall, roof, land)
1394                    IF ( .NOT.  ALLOCATED(surf_usm_h%wghf_eb_window_av) )  THEN
1395                        ALLOCATE( surf_usm_h%wghf_eb_window_av(1:surf_usm_h%ns) )
1396                        surf_usm_h%wghf_eb_window_av = 0.0_wp
1397                    ENDIF
1398                    DO  l = 0, 3
1399                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%wghf_eb_window_av) )  THEN
1400                           ALLOCATE( surf_usm_v(l)%wghf_eb_window_av(1:surf_usm_v(l)%ns) )
1401                           surf_usm_v(l)%wghf_eb_window_av = 0.0_wp
1402                       ENDIF
1403                    ENDDO
1404
1405                CASE ( 'usm_wghf_green' )
1406!--                 array of heat flux from green ground (wall, roof, land)
1407                    IF ( .NOT.  ALLOCATED(surf_usm_h%wghf_eb_green_av) )  THEN
1408                        ALLOCATE( surf_usm_h%wghf_eb_green_av(1:surf_usm_h%ns) )
1409                        surf_usm_h%wghf_eb_green_av = 0.0_wp
1410                    ENDIF
1411                    DO  l = 0, 3
1412                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%wghf_eb_green_av) )  THEN
1413                           ALLOCATE( surf_usm_v(l)%wghf_eb_green_av(1:surf_usm_v(l)%ns) )
1414                           surf_usm_v(l)%wghf_eb_green_av = 0.0_wp
1415                       ENDIF
1416                    ENDDO
1417
1418                CASE ( 'usm_iwghf' )
1419!--                 array of heat flux from indoor ground (wall, roof, land)
1420                    IF ( .NOT.  ALLOCATED(surf_usm_h%iwghf_eb_av) )  THEN
1421                        ALLOCATE( surf_usm_h%iwghf_eb_av(1:surf_usm_h%ns) )
1422                        surf_usm_h%iwghf_eb_av = 0.0_wp
1423                    ENDIF
1424                    DO  l = 0, 3
1425                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%iwghf_eb_av) )  THEN
1426                           ALLOCATE( surf_usm_v(l)%iwghf_eb_av(1:surf_usm_v(l)%ns) )
1427                           surf_usm_v(l)%iwghf_eb_av = 0.0_wp
1428                       ENDIF
1429                    ENDDO
1430
1431                CASE ( 'usm_iwghf_window' )
1432!--                 array of heat flux from indoor window ground (wall, roof, land)
1433                    IF ( .NOT.  ALLOCATED(surf_usm_h%iwghf_eb_window_av) )  THEN
1434                        ALLOCATE( surf_usm_h%iwghf_eb_window_av(1:surf_usm_h%ns) )
1435                        surf_usm_h%iwghf_eb_window_av = 0.0_wp
1436                    ENDIF
1437                    DO  l = 0, 3
1438                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%iwghf_eb_window_av) )  THEN
1439                           ALLOCATE( surf_usm_v(l)%iwghf_eb_window_av(1:surf_usm_v(l)%ns) )
1440                           surf_usm_v(l)%iwghf_eb_window_av = 0.0_wp
1441                       ENDIF
1442                    ENDDO
1443                   
1444                CASE ( 'usm_t_surf' )
1445!--                 surface temperature for surfaces
1446                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_surf_av) )  THEN
1447                        ALLOCATE( surf_usm_h%t_surf_av(1:surf_usm_h%ns) )
1448                        surf_usm_h%t_surf_av = 0.0_wp
1449                    ENDIF
1450                    DO  l = 0, 3
1451                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_surf_av) )  THEN
1452                           ALLOCATE( surf_usm_v(l)%t_surf_av(1:surf_usm_v(l)%ns) )
1453                           surf_usm_v(l)%t_surf_av = 0.0_wp
1454                       ENDIF
1455                    ENDDO
1456
1457                CASE ( 'usm_t_surf_window' )
1458!--                 surface temperature for window surfaces
1459                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_surf_window_av) )  THEN
1460                        ALLOCATE( surf_usm_h%t_surf_window_av(1:surf_usm_h%ns) )
1461                        surf_usm_h%t_surf_window_av = 0.0_wp
1462                    ENDIF
1463                    DO  l = 0, 3
1464                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_surf_window_av) )  THEN
1465                           ALLOCATE( surf_usm_v(l)%t_surf_window_av(1:surf_usm_v(l)%ns) )
1466                           surf_usm_v(l)%t_surf_window_av = 0.0_wp
1467                       ENDIF
1468                    ENDDO
1469                   
1470                CASE ( 'usm_t_surf_green' )
1471!--                 surface temperature for green surfaces
1472                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_surf_green_av) )  THEN
1473                        ALLOCATE( surf_usm_h%t_surf_green_av(1:surf_usm_h%ns) )
1474                        surf_usm_h%t_surf_green_av = 0.0_wp
1475                    ENDIF
1476                    DO  l = 0, 3
1477                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_surf_green_av) )  THEN
1478                           ALLOCATE( surf_usm_v(l)%t_surf_green_av(1:surf_usm_v(l)%ns) )
1479                           surf_usm_v(l)%t_surf_green_av = 0.0_wp
1480                       ENDIF
1481                    ENDDO
1482               
1483                CASE ( 'usm_t_surf_10cm' )
1484!--                 near surface temperature for whole surfaces
1485                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_surf_10cm_av) )  THEN
1486                        ALLOCATE( surf_usm_h%t_surf_10cm_av(1:surf_usm_h%ns) )
1487                        surf_usm_h%t_surf_10cm_av = 0.0_wp
1488                    ENDIF
1489                    DO  l = 0, 3
1490                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_surf_10cm_av) )  THEN
1491                           ALLOCATE( surf_usm_v(l)%t_surf_10cm_av(1:surf_usm_v(l)%ns) )
1492                           surf_usm_v(l)%t_surf_10cm_av = 0.0_wp
1493                       ENDIF
1494                    ENDDO
1495
1496                CASE ( 'usm_t_wall' )
1497!--                 wall temperature for iwl layer of walls and land
1498                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_wall_av) )  THEN
1499                        ALLOCATE( surf_usm_h%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1500                        surf_usm_h%t_wall_av = 0.0_wp
1501                    ENDIF
1502                    DO  l = 0, 3
1503                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_wall_av) )  THEN
1504                           ALLOCATE( surf_usm_v(l)%t_wall_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1505                           surf_usm_v(l)%t_wall_av = 0.0_wp
1506                       ENDIF
1507                    ENDDO
1508
1509                CASE ( 'usm_t_window' )
1510!--                 window temperature for iwl layer of walls and land
1511                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_window_av) )  THEN
1512                        ALLOCATE( surf_usm_h%t_window_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1513                        surf_usm_h%t_window_av = 0.0_wp
1514                    ENDIF
1515                    DO  l = 0, 3
1516                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_window_av) )  THEN
1517                           ALLOCATE( surf_usm_v(l)%t_window_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1518                           surf_usm_v(l)%t_window_av = 0.0_wp
1519                       ENDIF
1520                    ENDDO
1521
1522                CASE ( 'usm_t_green' )
1523!--                 green temperature for iwl layer of walls and land
1524                    IF ( .NOT.  ALLOCATED(surf_usm_h%t_green_av) )  THEN
1525                        ALLOCATE( surf_usm_h%t_green_av(nzb_wall:nzt_wall,1:surf_usm_h%ns) )
1526                        surf_usm_h%t_green_av = 0.0_wp
1527                    ENDIF
1528                    DO  l = 0, 3
1529                       IF ( .NOT.  ALLOCATED(surf_usm_v(l)%t_green_av) )  THEN
1530                           ALLOCATE( surf_usm_v(l)%t_green_av(nzb_wall:nzt_wall,1:surf_usm_v(l)%ns) )
1531                           surf_usm_v(l)%t_green_av = 0.0_wp
1532                       ENDIF
1533                    ENDDO
1534
1535               CASE DEFAULT
1536                   CONTINUE
1537
1538           END SELECT
1539
1540        ELSEIF ( mode == 'sum' )  THEN
1541           
1542           SELECT CASE ( TRIM( var ) )
1543               
1544                CASE ( 'usm_rad_net' )
1545!--                 array of complete radiation balance
1546                    DO  m = 1, surf_usm_h%ns
1547                       surf_usm_h%rad_net_av(m) =                              &
1548                                          surf_usm_h%rad_net_av(m) +           &
1549                                          surf_usm_h%rad_net_l(m)
1550                    ENDDO
1551                    DO  l = 0, 3
1552                       DO  m = 1, surf_usm_v(l)%ns
1553                          surf_usm_v(l)%rad_net_av(m) =                        &
1554                                          surf_usm_v(l)%rad_net_av(m) +        &
1555                                          surf_usm_v(l)%rad_net_l(m)
1556                       ENDDO
1557                    ENDDO
1558                   
1559                CASE ( 'usm_rad_insw' )
1560!--                 array of sw radiation falling to surface after i-th reflection
[2920]1561                    DO l = 1, nsurfl
1562                        IF ( surfl(id,l) == idsint )  THEN
[2737]1563                            surfinsw_av(l) = surfinsw_av(l) + surfinsw(l)
1564                        ENDIF
1565                    ENDDO
1566                             
1567                CASE ( 'usm_rad_inlw' )
1568!--                 array of lw radiation falling to surface after i-th reflection
[2920]1569                    DO l = 1, nsurfl
1570                        IF ( surfl(id,l) == idsint )  THEN
[2737]1571                            surfinlw_av(l) = surfinlw_av(l) + surfinlw(l)
1572                        ENDIF
1573                    ENDDO
1574                   
1575                CASE ( 'usm_rad_inswdir' )
1576!--                 array of direct sw radiation falling to surface from sun
[2920]1577                    DO l = 1, nsurfl
1578                        IF ( surfl(id,l) == idsint )  THEN
[2737]1579                            surfinswdir_av(l) = surfinswdir_av(l) + surfinswdir(l)
1580                        ENDIF
1581                    ENDDO
1582                   
1583                CASE ( 'usm_rad_inswdif' )
1584!--                 array of difusion sw radiation falling to surface from sky and borders of the domain
[2920]1585                    DO l = 1, nsurfl
1586                        IF ( surfl(id,l) == idsint )  THEN
[2737]1587                            surfinswdif_av(l) = surfinswdif_av(l) + surfinswdif(l)
1588                        ENDIF
1589                    ENDDO
1590                   
1591                CASE ( 'usm_rad_inswref' )
1592!--                 array of sw radiation falling to surface from reflections
[2920]1593                    DO l = 1, nsurfl
1594                        IF ( surfl(id,l) == idsint )  THEN
[2737]1595                            surfinswref_av(l) = surfinswref_av(l) + surfinsw(l) - &
1596                                                surfinswdir(l) - surfinswdif(l)
1597                        ENDIF
1598                    ENDDO
1599
1600                   
1601                CASE ( 'usm_rad_inlwdif' )
1602!--                 array of sw radiation falling to surface after i-th reflection
[2920]1603                    DO l = 1, nsurfl
1604                        IF ( surfl(id,l) == idsint )  THEN
1605                            surfinlwdif_av(l) = surfinlwdif_av(l) + surfinlwdif(l)
[2737]1606                        ENDIF
1607                    ENDDO
1608!                     
1609                CASE ( 'usm_rad_inlwref' )
1610!--                 array of lw radiation falling to surface from reflections
[2920]1611                    DO l = 1, nsurfl
1612                        IF ( surfl(id,l) == idsint )  THEN
1613                            surfinlwref_av(l) = surfinlwref_av(l) + &
1614                                                surfinlw(l) - surfinlwdif(l)
[2737]1615                        ENDIF
1616                    ENDDO
1617                   
1618                CASE ( 'usm_rad_outsw' )
1619!--                 array of sw radiation emitted from surface after i-th reflection
[2920]1620                    DO l = 1, nsurfl
1621                        IF ( surfl(id,l) == idsint )  THEN
1622                            surfoutsw_av(l) = surfoutsw_av(l) + surfoutsw(l)
[2737]1623                        ENDIF
1624                    ENDDO
1625                   
1626                CASE ( 'usm_rad_outlw' )
1627!--                 array of lw radiation emitted from surface after i-th reflection
[2920]1628                    DO l = 1, nsurfl
1629                        IF ( surfl(id,l) == idsint )  THEN
1630                            surfoutlw_av(l) = surfoutlw_av(l) + surfoutlw(l)
[2737]1631                        ENDIF
1632                    ENDDO
1633                   
1634                CASE ( 'usm_rad_ressw' )
1635!--                 array of residua of sw radiation absorbed in surface after last reflection
[2920]1636                    DO l = 1, nsurfl
1637                        IF ( surfl(id,l) == idsint )  THEN
1638                            surfins_av(l) = surfins_av(l) + surfins(l)
[2737]1639                        ENDIF
1640                    ENDDO
1641                                   
1642                CASE ( 'usm_rad_reslw' )
1643!--                 array of residua of lw radiation absorbed in surface after last reflection
[2920]1644                    DO l = 1, nsurfl
1645                        IF ( surfl(id,l) == idsint )  THEN
1646                            surfinl_av(l) = surfinl_av(l) + surfinl(l)
[2737]1647                        ENDIF
1648                    ENDDO
1649                   
1650                CASE ( 'usm_rad_hf' )
1651!--                 array of heat flux from radiation for surfaces after i-th reflection
1652                    DO  m = 1, surf_usm_h%ns
1653                       surf_usm_h%surfhf_av(m) =                               &
1654                                          surf_usm_h%surfhf_av(m) +            &
1655                                          surf_usm_h%surfhf(m)
1656                    ENDDO
1657                    DO  l = 0, 3
1658                       DO  m = 1, surf_usm_v(l)%ns
1659                          surf_usm_v(l)%surfhf_av(m) =                         &
1660                                          surf_usm_v(l)%surfhf_av(m) +         &
1661                                          surf_usm_v(l)%surfhf(m)
1662                       ENDDO
1663                    ENDDO
1664                   
1665                CASE ( 'usm_wshf' )
1666!--                 array of sensible heat flux from surfaces (land, roof, wall)
1667                    DO  m = 1, surf_usm_h%ns
1668                       surf_usm_h%wshf_eb_av(m) =                              &
1669                                          surf_usm_h%wshf_eb_av(m) +           &
1670                                          surf_usm_h%wshf_eb(m)
1671                    ENDDO
1672                    DO  l = 0, 3
1673                       DO  m = 1, surf_usm_v(l)%ns
1674                          surf_usm_v(l)%wshf_eb_av(m) =                        &
1675                                          surf_usm_v(l)%wshf_eb_av(m) +        &
1676                                          surf_usm_v(l)%wshf_eb(m)
1677                       ENDDO
1678                    ENDDO
1679                   
1680                CASE ( 'usm_wghf' )
1681!--                 array of heat flux from ground (wall, roof, land)
1682                    DO  m = 1, surf_usm_h%ns
1683                       surf_usm_h%wghf_eb_av(m) =                              &
1684                                          surf_usm_h%wghf_eb_av(m) +           &
1685                                          surf_usm_h%wghf_eb(m)
1686                    ENDDO
1687                    DO  l = 0, 3
1688                       DO  m = 1, surf_usm_v(l)%ns
1689                          surf_usm_v(l)%wghf_eb_av(m) =                        &
1690                                          surf_usm_v(l)%wghf_eb_av(m) +        &
1691                                          surf_usm_v(l)%wghf_eb(m)
1692                       ENDDO
1693                    ENDDO
1694                   
1695                CASE ( 'usm_wghf_window' )
1696!--                 array of heat flux from window ground (wall, roof, land)
1697                    DO  m = 1, surf_usm_h%ns
1698                       surf_usm_h%wghf_eb_window_av(m) =                              &
1699                                          surf_usm_h%wghf_eb_window_av(m) +           &
1700                                          surf_usm_h%wghf_eb_window(m)
1701                    ENDDO
1702                    DO  l = 0, 3
1703                       DO  m = 1, surf_usm_v(l)%ns
1704                          surf_usm_v(l)%wghf_eb_window_av(m) =                        &
1705                                          surf_usm_v(l)%wghf_eb_window_av(m) +        &
1706                                          surf_usm_v(l)%wghf_eb_window(m)
1707                       ENDDO
1708                    ENDDO
1709
1710                CASE ( 'usm_wghf_green' )
1711!--                 array of heat flux from green ground (wall, roof, land)
1712                    DO  m = 1, surf_usm_h%ns
1713                       surf_usm_h%wghf_eb_green_av(m) =                              &
1714                                          surf_usm_h%wghf_eb_green_av(m) +           &
1715                                          surf_usm_h%wghf_eb_green(m)
1716                    ENDDO
1717                    DO  l = 0, 3
1718                       DO  m = 1, surf_usm_v(l)%ns
1719                          surf_usm_v(l)%wghf_eb_green_av(m) =                        &
1720                                          surf_usm_v(l)%wghf_eb_green_av(m) +        &
1721                                          surf_usm_v(l)%wghf_eb_green(m)
1722                       ENDDO
1723                    ENDDO
1724                   
1725                CASE ( 'usm_iwghf' )
1726!--                 array of heat flux from indoor ground (wall, roof, land)
1727                    DO  m = 1, surf_usm_h%ns
1728                       surf_usm_h%iwghf_eb_av(m) =                              &
1729                                          surf_usm_h%iwghf_eb_av(m) +           &
1730                                          surf_usm_h%iwghf_eb(m)
1731                    ENDDO
1732                    DO  l = 0, 3
1733                       DO  m = 1, surf_usm_v(l)%ns
1734                          surf_usm_v(l)%iwghf_eb_av(m) =                        &
1735                                          surf_usm_v(l)%iwghf_eb_av(m) +        &
1736                                          surf_usm_v(l)%iwghf_eb(m)
1737                       ENDDO
1738                    ENDDO
1739                   
1740                CASE ( 'usm_iwghf_window' )
1741!--                 array of heat flux from indoor window ground (wall, roof, land)
1742                    DO  m = 1, surf_usm_h%ns
1743                       surf_usm_h%iwghf_eb_window_av(m) =                              &
1744                                          surf_usm_h%iwghf_eb_window_av(m) +           &
1745                                          surf_usm_h%iwghf_eb_window(m)
1746                    ENDDO
1747                    DO  l = 0, 3
1748                       DO  m = 1, surf_usm_v(l)%ns
1749                          surf_usm_v(l)%iwghf_eb_window_av(m) =                        &
1750                                          surf_usm_v(l)%iwghf_eb_window_av(m) +        &
1751                                          surf_usm_v(l)%iwghf_eb_window(m)
1752                       ENDDO
1753                    ENDDO
1754                   
1755                CASE ( 'usm_t_surf' )
1756!--                 surface temperature for surfaces
1757                    DO  m = 1, surf_usm_h%ns
1758                       surf_usm_h%t_surf_av(m) =                               & 
1759                                          surf_usm_h%t_surf_av(m) +            &
1760                                          t_surf_h(m)
1761                    ENDDO
1762                    DO  l = 0, 3
1763                       DO  m = 1, surf_usm_v(l)%ns
1764                          surf_usm_v(l)%t_surf_av(m) =                         &
1765                                          surf_usm_v(l)%t_surf_av(m) +         &
1766                                          t_surf_v(l)%t(m)
1767                       ENDDO
1768                    ENDDO
1769                   
1770                CASE ( 'usm_t_surf_window' )
1771!--                 surface temperature for window surfaces
1772                    DO  m = 1, surf_usm_h%ns
1773                       surf_usm_h%t_surf_window_av(m) =                               & 
1774                                          surf_usm_h%t_surf_window_av(m) +            &
1775                                          t_surf_window_h(m)
1776                    ENDDO
1777                    DO  l = 0, 3
1778                       DO  m = 1, surf_usm_v(l)%ns
1779                          surf_usm_v(l)%t_surf_window_av(m) =                         &
1780                                          surf_usm_v(l)%t_surf_window_av(m) +         &
1781                                          t_surf_window_v(l)%t(m)
1782                       ENDDO
1783                    ENDDO
1784                   
1785                CASE ( 'usm_t_surf_green' )
1786!--                 surface temperature for green surfaces
1787                    DO  m = 1, surf_usm_h%ns
1788                       surf_usm_h%t_surf_green_av(m) =                               & 
1789                                          surf_usm_h%t_surf_green_av(m) +            &
1790                                          t_surf_green_h(m)
1791                    ENDDO
1792                    DO  l = 0, 3
1793                       DO  m = 1, surf_usm_v(l)%ns
1794                          surf_usm_v(l)%t_surf_green_av(m) =                         &
1795                                          surf_usm_v(l)%t_surf_green_av(m) +         &
1796                                          t_surf_green_v(l)%t(m)
1797                       ENDDO
1798                    ENDDO
1799               
1800                CASE ( 'usm_t_surf_10cm' )
1801!--                 near surface temperature for whole surfaces
1802                    DO  m = 1, surf_usm_h%ns
1803                       surf_usm_h%t_surf_10cm_av(m) =                               & 
1804                                          surf_usm_h%t_surf_10cm_av(m) +            &
1805                                          t_surf_10cm_h(m)
1806                    ENDDO
1807                    DO  l = 0, 3
1808                       DO  m = 1, surf_usm_v(l)%ns
1809                          surf_usm_v(l)%t_surf_10cm_av(m) =                         &
1810                                          surf_usm_v(l)%t_surf_10cm_av(m) +         &
1811                                          t_surf_10cm_v(l)%t(m)
1812                       ENDDO
1813                    ENDDO
1814
1815                   
1816                CASE ( 'usm_t_wall' )
1817!--                 wall temperature for  iwl layer of walls and land
1818                    DO  m = 1, surf_usm_h%ns
1819                       surf_usm_h%t_wall_av(iwl,m) =                           &
1820                                          surf_usm_h%t_wall_av(iwl,m) +        &
1821                                          t_wall_h(iwl,m)
1822                    ENDDO
1823                    DO  l = 0, 3
1824                       DO  m = 1, surf_usm_v(l)%ns
1825                          surf_usm_v(l)%t_wall_av(iwl,m) =                     &
1826                                          surf_usm_v(l)%t_wall_av(iwl,m) +     &
1827                                          t_wall_v(l)%t(iwl,m)
1828                       ENDDO
1829                    ENDDO
1830                   
1831                CASE ( 'usm_t_window' )
1832!--                 window temperature for  iwl layer of walls and land
1833                    DO  m = 1, surf_usm_h%ns
1834                       surf_usm_h%t_window_av(iwl,m) =                           &
1835                                          surf_usm_h%t_window_av(iwl,m) +        &
1836                                          t_window_h(iwl,m)
1837                    ENDDO
1838                    DO  l = 0, 3
1839                       DO  m = 1, surf_usm_v(l)%ns
1840                          surf_usm_v(l)%t_window_av(iwl,m) =                     &
1841                                          surf_usm_v(l)%t_window_av(iwl,m) +     &
1842                                          t_window_v(l)%t(iwl,m)
1843                       ENDDO
1844                    ENDDO
1845
1846                CASE ( 'usm_t_green' )
1847!--                 green temperature for  iwl layer of walls and land
1848                    DO  m = 1, surf_usm_h%ns
1849                       surf_usm_h%t_green_av(iwl,m) =                           &
1850                                          surf_usm_h%t_green_av(iwl,m) +        &
1851                                          t_green_h(iwl,m)
1852                    ENDDO
1853                    DO  l = 0, 3
1854                       DO  m = 1, surf_usm_v(l)%ns
1855                          surf_usm_v(l)%t_green_av(iwl,m) =                     &
1856                                          surf_usm_v(l)%t_green_av(iwl,m) +     &
1857                                          t_green_v(l)%t(iwl,m)
1858                       ENDDO
1859                    ENDDO
1860
1861                CASE DEFAULT
1862                    CONTINUE
1863
1864           END SELECT
1865
1866        ELSEIF ( mode == 'average' )  THEN
1867           
1868           SELECT CASE ( TRIM( var ) )
1869               
1870                CASE ( 'usm_rad_net' )
1871!--                 array of complete radiation balance
1872                    DO  m = 1, surf_usm_h%ns
1873                       surf_usm_h%rad_net_av(m) =                              &
1874                                          surf_usm_h%rad_net_av(m) /           &
1875                                          REAL( average_count_3d, kind=wp )
1876                    ENDDO
1877                    DO  l = 0, 3
1878                       DO  m = 1, surf_usm_v(l)%ns
1879                          surf_usm_v(l)%rad_net_av(m) =                        &
1880                                          surf_usm_v(l)%rad_net_av(m) /        &
1881                                          REAL( average_count_3d, kind=wp )
1882                       ENDDO
1883                    ENDDO
1884                   
1885                CASE ( 'usm_rad_insw' )
1886!--                 array of sw radiation falling to surface after i-th reflection
[2920]1887                    DO l = 1, nsurfl
1888                        IF ( surfl(id,l) == idsint )  THEN
[2737]1889                            surfinsw_av(l) = surfinsw_av(l) / REAL( average_count_3d, kind=wp )
1890                        ENDIF
1891                    ENDDO
1892                             
1893                CASE ( 'usm_rad_inlw' )
1894!--                 array of lw radiation falling to surface after i-th reflection
[2920]1895                    DO l = 1, nsurfl
1896                        IF ( surfl(id,l) == idsint )  THEN
[2737]1897                            surfinlw_av(l) = surfinlw_av(l) / REAL( average_count_3d, kind=wp )
1898                        ENDIF
1899                    ENDDO
1900                   
1901                CASE ( 'usm_rad_inswdir' )
1902!--                 array of direct sw radiation falling to surface from sun
[2920]1903                    DO l = 1, nsurfl
1904                        IF ( surfl(id,l) == idsint )  THEN
[2737]1905                            surfinswdir_av(l) = surfinswdir_av(l) / REAL( average_count_3d, kind=wp )
1906                        ENDIF
1907                    ENDDO
1908                   
1909                CASE ( 'usm_rad_inswdif' )
1910!--                 array of difusion sw radiation falling to surface from sky and borders of the domain
[2920]1911                    DO l = 1, nsurfl
1912                        IF ( surfl(id,l) == idsint )  THEN
[2737]1913                            surfinswdif_av(l) = surfinswdif_av(l) / REAL( average_count_3d, kind=wp )
1914                        ENDIF
1915                    ENDDO
1916                   
1917                CASE ( 'usm_rad_inswref' )
1918!--                 array of sw radiation falling to surface from reflections
[2920]1919                    DO l = 1, nsurfl
1920                        IF ( surfl(id,l) == idsint )  THEN
[2737]1921                            surfinswref_av(l) = surfinswref_av(l) / REAL( average_count_3d, kind=wp )
1922                        ENDIF
1923                    ENDDO
1924                   
1925                CASE ( 'usm_rad_inlwdif' )
1926!--                 array of sw radiation falling to surface after i-th reflection
[2920]1927                    DO l = 1, nsurfl
1928                        IF ( surfl(id,l) == idsint )  THEN
[2737]1929                            surfinlwdif_av(l) = surfinlwdif_av(l) / REAL( average_count_3d, kind=wp )
1930                        ENDIF
1931                    ENDDO
1932                   
1933                CASE ( 'usm_rad_inlwref' )
1934!--                 array of lw radiation falling to surface from reflections
[2920]1935                    DO l = 1, nsurfl
1936                        IF ( surfl(id,l) == idsint )  THEN
[2737]1937                            surfinlwref_av(l) = surfinlwref_av(l) / REAL( average_count_3d, kind=wp )
1938                        ENDIF
1939                    ENDDO
1940                   
1941                CASE ( 'usm_rad_outsw' )
1942!--                 array of sw radiation emitted from surface after i-th reflection
[2920]1943                    DO l = 1, nsurfl
1944                        IF ( surfl(id,l) == idsint )  THEN
[2737]1945                            surfoutsw_av(l) = surfoutsw_av(l) / REAL( average_count_3d, kind=wp )
1946                        ENDIF
1947                    ENDDO
1948                   
1949                CASE ( 'usm_rad_outlw' )
1950!--                 array of lw radiation emitted from surface after i-th reflection
[2920]1951                    DO l = 1, nsurfl
1952                        IF ( surfl(id,l) == idsint )  THEN
[2737]1953                            surfoutlw_av(l) = surfoutlw_av(l) / REAL( average_count_3d, kind=wp )
1954                        ENDIF
1955                    ENDDO
1956                   
1957                CASE ( 'usm_rad_ressw' )
1958!--                 array of residua of sw radiation absorbed in surface after last reflection
[2920]1959                    DO l = 1, nsurfl
1960                        IF ( surfl(id,l) == idsint )  THEN
[2737]1961                            surfins_av(l) = surfins_av(l) / REAL( average_count_3d, kind=wp )
1962                        ENDIF
1963                    ENDDO
1964                                   
1965                CASE ( 'usm_rad_reslw' )
1966!--                 array of residua of lw radiation absorbed in surface after last reflection
[2920]1967                    DO l = 1, nsurfl
1968                        IF ( surfl(id,l) == idsint )  THEN
[2737]1969                            surfinl_av(l) = surfinl_av(l) / REAL( average_count_3d, kind=wp )
1970                        ENDIF
1971                    ENDDO
1972                   
1973                CASE ( 'usm_rad_hf' )
1974!--                 array of heat flux from radiation for surfaces after i-th reflection
1975                    DO  m = 1, surf_usm_h%ns
1976                       surf_usm_h%surfhf_av(m) =                               &
1977                                          surf_usm_h%surfhf_av(m) /            &
1978                                          REAL( average_count_3d, kind=wp )
1979                    ENDDO
1980                    DO  l = 0, 3
1981                       DO  m = 1, surf_usm_v(l)%ns
1982                          surf_usm_v(l)%surfhf_av(m) =                         &
1983                                          surf_usm_v(l)%surfhf_av(m) /         &
1984                                          REAL( average_count_3d, kind=wp )
1985                       ENDDO
1986                    ENDDO
1987                   
1988                CASE ( 'usm_wshf' )
1989!--                 array of sensible heat flux from surfaces (land, roof, wall)
1990                    DO  m = 1, surf_usm_h%ns
1991                       surf_usm_h%wshf_eb_av(m) =                              &
1992                                          surf_usm_h%wshf_eb_av(m) /           &
1993                                          REAL( average_count_3d, kind=wp )
1994                    ENDDO
1995                    DO  l = 0, 3
1996                       DO  m = 1, surf_usm_v(l)%ns
1997                          surf_usm_v(l)%wshf_eb_av(m) =                        &
1998                                          surf_usm_v(l)%wshf_eb_av(m) /        &
1999                                          REAL( average_count_3d, kind=wp )
2000                       ENDDO
2001                    ENDDO
2002                   
2003                CASE ( 'usm_wghf' )
2004!--                 array of heat flux from ground (wall, roof, land)
2005                    DO  m = 1, surf_usm_h%ns
2006                       surf_usm_h%wghf_eb_av(m) =                              &
2007                                          surf_usm_h%wghf_eb_av(m) /           &
2008                                          REAL( average_count_3d, kind=wp )
2009                    ENDDO
2010                    DO  l = 0, 3
2011                       DO  m = 1, surf_usm_v(l)%ns
2012                          surf_usm_v(l)%wghf_eb_av(m) =                        &
2013                                          surf_usm_v(l)%wghf_eb_av(m) /        &
2014                                          REAL( average_count_3d, kind=wp )
2015                       ENDDO
2016                    ENDDO
2017                   
2018                CASE ( 'usm_wghf_window' )
2019!--                 array of heat flux from window ground (wall, roof, land)
2020                    DO  m = 1, surf_usm_h%ns
2021                       surf_usm_h%wghf_eb_window_av(m) =                              &
2022                                          surf_usm_h%wghf_eb_window_av(m) /           &
2023                                          REAL( average_count_3d, kind=wp )
2024                    ENDDO
2025                    DO  l = 0, 3
2026                       DO  m = 1, surf_usm_v(l)%ns
2027                          surf_usm_v(l)%wghf_eb_window_av(m) =                        &
2028                                          surf_usm_v(l)%wghf_eb_window_av(m) /        &
2029                                          REAL( average_count_3d, kind=wp )
2030                       ENDDO
2031                    ENDDO
2032
2033                CASE ( 'usm_wghf_green' )
2034!--                 array of heat flux from green ground (wall, roof, land)
2035                    DO  m = 1, surf_usm_h%ns
2036                       surf_usm_h%wghf_eb_green_av(m) =                              &
2037                                          surf_usm_h%wghf_eb_green_av(m) /           &
2038                                          REAL( average_count_3d, kind=wp )
2039                    ENDDO
2040                    DO  l = 0, 3
2041                       DO  m = 1, surf_usm_v(l)%ns
2042                          surf_usm_v(l)%wghf_eb_green_av(m) =                        &
2043                                          surf_usm_v(l)%wghf_eb_green_av(m) /        &
2044                                          REAL( average_count_3d, kind=wp )
2045                       ENDDO
2046                    ENDDO
2047
2048                CASE ( 'usm_iwghf' )
2049!--                 array of heat flux from indoor ground (wall, roof, land)
2050                    DO  m = 1, surf_usm_h%ns
2051                       surf_usm_h%iwghf_eb_av(m) =                              &
2052                                          surf_usm_h%iwghf_eb_av(m) /           &
2053                                          REAL( average_count_3d, kind=wp )
2054                    ENDDO
2055                    DO  l = 0, 3
2056                       DO  m = 1, surf_usm_v(l)%ns
2057                          surf_usm_v(l)%iwghf_eb_av(m) =                        &
2058                                          surf_usm_v(l)%iwghf_eb_av(m) /        &
2059                                          REAL( average_count_3d, kind=wp )
2060                       ENDDO
2061                    ENDDO
2062                   
2063                CASE ( 'usm_iwghf_window' )
2064!--                 array of heat flux from indoor window ground (wall, roof, land)
2065                    DO  m = 1, surf_usm_h%ns
2066                       surf_usm_h%iwghf_eb_window_av(m) =                              &
2067                                          surf_usm_h%iwghf_eb_window_av(m) /           &
2068                                          REAL( average_count_3d, kind=wp )
2069                    ENDDO
2070                    DO  l = 0, 3
2071                       DO  m = 1, surf_usm_v(l)%ns
2072                          surf_usm_v(l)%iwghf_eb_window_av(m) =                        &
2073                                          surf_usm_v(l)%iwghf_eb_window_av(m) /        &
2074                                          REAL( average_count_3d, kind=wp )
2075                       ENDDO
2076                    ENDDO
2077                   
2078                CASE ( 'usm_t_surf' )
2079!--                 surface temperature for surfaces
2080                    DO  m = 1, surf_usm_h%ns
2081                       surf_usm_h%t_surf_av(m) =                               & 
2082                                          surf_usm_h%t_surf_av(m) /            &
2083                                          REAL( average_count_3d, kind=wp )
2084                    ENDDO
2085                    DO  l = 0, 3
2086                       DO  m = 1, surf_usm_v(l)%ns
2087                          surf_usm_v(l)%t_surf_av(m) =                         &
2088                                          surf_usm_v(l)%t_surf_av(m) /         &
2089                                          REAL( average_count_3d, kind=wp )
2090                       ENDDO
2091                    ENDDO
2092                   
2093                CASE ( 'usm_t_surf_window' )
2094!--                 surface temperature for window surfaces
2095                    DO  m = 1, surf_usm_h%ns
2096                       surf_usm_h%t_surf_window_av(m) =                               & 
2097                                          surf_usm_h%t_surf_window_av(m) /            &
2098                                          REAL( average_count_3d, kind=wp )
2099                    ENDDO
2100                    DO  l = 0, 3
2101                       DO  m = 1, surf_usm_v(l)%ns
2102                          surf_usm_v(l)%t_surf_window_av(m) =                         &
2103                                          surf_usm_v(l)%t_surf_window_av(m) /         &
2104                                          REAL( average_count_3d, kind=wp )
2105                       ENDDO
2106                    ENDDO
2107                   
2108                CASE ( 'usm_t_surf_green' )
2109!--                 surface temperature for green surfaces
2110                    DO  m = 1, surf_usm_h%ns
2111                       surf_usm_h%t_surf_green_av(m) =                               & 
2112                                          surf_usm_h%t_surf_green_av(m) /            &
2113                                          REAL( average_count_3d, kind=wp )
2114                    ENDDO
2115                    DO  l = 0, 3
2116                       DO  m = 1, surf_usm_v(l)%ns
2117                          surf_usm_v(l)%t_surf_green_av(m) =                         &
2118                                          surf_usm_v(l)%t_surf_green_av(m) /         &
2119                                          REAL( average_count_3d, kind=wp )
2120                       ENDDO
2121                    ENDDO
2122                   
2123                CASE ( 'usm_t_surf_10cm' )
2124!--                 near surface temperature for whole surfaces
2125                    DO  m = 1, surf_usm_h%ns
2126                       surf_usm_h%t_surf_10cm_av(m) =                               & 
2127                                          surf_usm_h%t_surf_10cm_av(m) /            &
2128                                          REAL( average_count_3d, kind=wp )
2129                    ENDDO
2130                    DO  l = 0, 3
2131                       DO  m = 1, surf_usm_v(l)%ns
2132                          surf_usm_v(l)%t_surf_10cm_av(m) =                         &
2133                                          surf_usm_v(l)%t_surf_10cm_av(m) /         &
2134                                          REAL( average_count_3d, kind=wp )
2135                       ENDDO
2136                    ENDDO
2137                   
2138                CASE ( 'usm_t_wall' )
2139!--                 wall temperature for  iwl layer of walls and land
2140                    DO  m = 1, surf_usm_h%ns
2141                       surf_usm_h%t_wall_av(iwl,m) =                           &
2142                                          surf_usm_h%t_wall_av(iwl,m) /        &
2143                                          REAL( average_count_3d, kind=wp )
2144                    ENDDO
2145                    DO  l = 0, 3
2146                       DO  m = 1, surf_usm_v(l)%ns
2147                          surf_usm_v(l)%t_wall_av(iwl,m) =                     &
2148                                          surf_usm_v(l)%t_wall_av(iwl,m) /     &
2149                                          REAL( average_count_3d, kind=wp )
2150                       ENDDO
2151                    ENDDO
2152
2153                CASE ( 'usm_t_window' )
2154!--                 window temperature for  iwl layer of walls and land
2155                    DO  m = 1, surf_usm_h%ns
2156                       surf_usm_h%t_window_av(iwl,m) =                           &
2157                                          surf_usm_h%t_window_av(iwl,m) /        &
2158                                          REAL( average_count_3d, kind=wp )
2159                    ENDDO
2160                    DO  l = 0, 3
2161                       DO  m = 1, surf_usm_v(l)%ns
2162                          surf_usm_v(l)%t_window_av(iwl,m) =                     &
2163                                          surf_usm_v(l)%t_window_av(iwl,m) /     &
2164                                          REAL( average_count_3d, kind=wp )
2165                       ENDDO
2166                    ENDDO
2167
2168                CASE ( 'usm_t_green' )
2169!--                 green temperature for  iwl layer of walls and land
2170                    DO  m = 1, surf_usm_h%ns
2171                       surf_usm_h%t_green_av(iwl,m) =                           &
2172                                          surf_usm_h%t_green_av(iwl,m) /        &
2173                                          REAL( average_count_3d, kind=wp )
2174                    ENDDO
2175                    DO  l = 0, 3
2176                       DO  m = 1, surf_usm_v(l)%ns
2177                          surf_usm_v(l)%t_green_av(iwl,m) =                     &
2178                                          surf_usm_v(l)%t_green_av(iwl,m) /     &
2179                                          REAL( average_count_3d, kind=wp )
2180                       ENDDO
2181                    ENDDO
2182
2183
2184           END SELECT
2185
2186        ENDIF
2187
2188    END SUBROUTINE usm_average_3d_data
2189
2190
2191
2192!------------------------------------------------------------------------------!
2193! Description:
2194! ------------
2195!> Set internal Neumann boundary condition at outer soil grid points
2196!> for temperature and humidity.
2197!------------------------------------------------------------------------------!
2198 SUBROUTINE usm_boundary_condition
2199 
2200    IMPLICIT NONE
2201
2202    INTEGER(iwp) :: i      !< grid index x-direction
2203    INTEGER(iwp) :: ioff   !< offset index x-direction indicating location of soil grid point
2204    INTEGER(iwp) :: j      !< grid index y-direction
2205    INTEGER(iwp) :: joff   !< offset index x-direction indicating location of soil grid point
2206    INTEGER(iwp) :: k      !< grid index z-direction
2207    INTEGER(iwp) :: koff   !< offset index x-direction indicating location of soil grid point
2208    INTEGER(iwp) :: l      !< running index surface-orientation
2209    INTEGER(iwp) :: m      !< running index surface elements
2210
2211    koff = surf_usm_h%koff
2212    DO  m = 1, surf_usm_h%ns
2213       i = surf_usm_h%i(m)
2214       j = surf_usm_h%j(m)
2215       k = surf_usm_h%k(m)
2216       pt(k+koff,j,i) = pt(k,j,i)
2217    ENDDO
2218
2219    DO  l = 0, 3
2220       ioff = surf_usm_v(l)%ioff
2221       joff = surf_usm_v(l)%joff
2222       DO  m = 1, surf_usm_v(l)%ns
2223          i = surf_usm_v(l)%i(m)
2224          j = surf_usm_v(l)%j(m)
2225          k = surf_usm_v(l)%k(m)
2226          pt(k,j+joff,i+ioff) = pt(k,j,i)
2227       ENDDO
2228    ENDDO
2229
2230 END SUBROUTINE usm_boundary_condition
2231
2232
2233!------------------------------------------------------------------------------!
2234!
2235! Description:
2236! ------------
2237!> Subroutine checks variables and assigns units.
2238!> It is called out from subroutine check_parameters.
2239!------------------------------------------------------------------------------!
2240    SUBROUTINE usm_check_data_output( variable, unit )
2241       
2242        IMPLICIT NONE
2243 
2244        CHARACTER (len=*),INTENT(IN)    ::  variable !:
2245        CHARACTER (len=*),INTENT(OUT)   ::  unit     !:
2246       
2247        CHARACTER (len=varnamelength)   :: var
2248
2249        var = TRIM(variable)
2250        IF ( var(1:12) == 'usm_rad_net_'  .OR.  var(1:13) == 'usm_rad_insw_'  .OR.        &
2251             var(1:13) == 'usm_rad_inlw_'  .OR.  var(1:16) == 'usm_rad_inswdir_'  .OR.    &
2252             var(1:16) == 'usm_rad_inswdif_'  .OR.  var(1:16) == 'usm_rad_inswref_'  .OR. &
2253             var(1:16) == 'usm_rad_inlwdif_'  .OR.  var(1:16) == 'usm_rad_inlwref_'  .OR. &
2254             var(1:14) == 'usm_rad_outsw_'  .OR.  var(1:14) == 'usm_rad_outlw_'  .OR.     &
2255             var(1:14) == 'usm_rad_ressw_'  .OR.  var(1:14) == 'usm_rad_reslw_'  .OR.     &
2256             var(1:11) == 'usm_rad_hf_'  .OR.                                             &
2257             var(1:9)  == 'usm_wshf_'  .OR.  var(1:9) == 'usm_wghf_' .OR.                 &
2258             var(1:16) == 'usm_wghf_window_' .OR. var(1:15) == 'usm_wghf_green_' .OR.     &
2259             var(1:10) == 'usm_iwghf_' .OR. var(1:17) == 'usm_iwghf_window_' )  THEN
2260            unit = 'W/m2'
[2920]2261        ELSE IF ( var(1:10) == 'usm_t_surf'   .OR.  var(1:10) == 'usm_t_wall' .OR.         &
[2737]2262                  var(1:12) == 'usm_t_window' .OR. var(1:17) == 'usm_t_surf_window' .OR.  &
2263                  var(1:16) == 'usm_t_surf_green'  .OR.                                   &
2264                  var(1:11) == 'usm_t_green' .OR.                                         &
2265                  var(1:15) == 'usm_t_surf_10cm')  THEN
2266            unit = 'K'
2267        ELSE IF ( var(1:9) == 'usm_surfz'  .OR.  var(1:7) == 'usm_svf'  .OR.              & 
2268                  var(1:7) == 'usm_dif'  .OR.  var(1:11) == 'usm_surfcat'  .OR.           &
[2920]2269                  var(1:11) == 'usm_surfalb'  .OR.  var(1:12) == 'usm_surfemis'  .OR.     &
2270                  var(1:9) == 'usm_skyvf' .OR. var(1:9) == 'usm_skyvft' )  THEN
[2737]2271            unit = '1'
2272        ELSE
2273            unit = 'illegal'
2274        ENDIF
2275
2276    END SUBROUTINE usm_check_data_output
2277
2278
2279!------------------------------------------------------------------------------!
2280! Description:
2281! ------------
2282!> Check parameters routine for urban surface model
2283!------------------------------------------------------------------------------!
2284    SUBROUTINE usm_check_parameters
2285   
2286       USE control_parameters,                                                 &
2287           ONLY:  bc_pt_b, bc_q_b, constant_flux_layer, large_scale_forcing,   &
2288                  lsf_surf, topography
2289
2290!
2291!--    Dirichlet boundary conditions are required as the surface fluxes are
2292!--    calculated from the temperature/humidity gradients in the urban surface
2293!--    model
2294       IF ( bc_pt_b == 'neumann'   .OR.   bc_q_b == 'neumann' )  THEN
2295          message_string = 'urban surface model requires setting of '//        &
2296                           'bc_pt_b = "dirichlet" and '//                      &
2297                           'bc_q_b  = "dirichlet"'
2298          CALL message( 'check_parameters', 'PA0590', 1, 2, 0, 6, 0 )
2299       ENDIF
2300
2301       IF ( .NOT.  constant_flux_layer )  THEN
2302          message_string = 'urban surface model requires '//                   &
2303                           'constant_flux_layer = .T.'
2304          CALL message( 'check_parameters', 'PA0591', 1, 2, 0, 6, 0 )
2305       ENDIF
2306!       
2307!--    Surface forcing has to be disabled for LSF in case of enabled
2308!--    urban surface module
2309       IF ( large_scale_forcing )  THEN
2310          lsf_surf = .FALSE.
2311       ENDIF
2312!
2313!--    Topography
2314       IF ( topography == 'flat' )  THEN
2315          message_string = 'topography /= "flat" is required '//               &
2316                           'when using the urban surface model'
2317          CALL message( 'check_parameters', 'PA0592', 1, 2, 0, 6, 0 )
2318       ENDIF
[2920]2319!
2320!--    naheatlayers
2321       IF ( naheatlayers > nzt )  THEN
2322          message_string = 'number of anthropogenic heat layers '//            &
2323                           '"naheatlayers" can not be larger than'//           &
2324                           ' number of domain layers "nzt"'
2325          CALL message( 'check_parameters', 'PA0593', 1, 2, 0, 6, 0 )
2326       ENDIF
[2737]2327
2328    END SUBROUTINE usm_check_parameters
2329
2330
2331!------------------------------------------------------------------------------!
2332!
2333! Description:
2334! ------------
2335!> Output of the 3D-arrays in netCDF and/or AVS format
2336!> for variables of urban_surface model.
2337!> It resorts the urban surface module output quantities from surf style
2338!> indexing into temporary 3D array with indices (i,j,k).
2339!> It is called from subroutine data_output_3d.
2340!------------------------------------------------------------------------------!
2341    SUBROUTINE usm_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
2342       
2343        IMPLICIT NONE
2344
2345        INTEGER(iwp), INTENT(IN)       ::  av        !<
2346        CHARACTER (len=*), INTENT(IN)  ::  variable  !<
2347        INTEGER(iwp), INTENT(IN)       ::  nzb_do    !< lower limit of the data output (usually 0)
2348        INTEGER(iwp), INTENT(IN)       ::  nzt_do    !< vertical upper limit of the data output (usually nz_do3d)
2349        LOGICAL, INTENT(OUT)           ::  found     !<
2350        REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf   !< sp - it has to correspond to module data_output_3d
2351        REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr)     ::  temp_pf    !< temp array for urban surface output procedure
2352       
2353        CHARACTER (len=varnamelength)                          :: var, surfid
2354        INTEGER(iwp), PARAMETER                                :: nd = 5
2355        CHARACTER(len=6), DIMENSION(0:nd-1), PARAMETER         :: dirname = (/ '_roof ', '_south', '_north', '_west ', '_east ' /)
[2920]2356        INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER             :: dirint =  (/    iup_u, isouth_u, inorth_u,  iwest_u,  ieast_u /)
2357        INTEGER(iwp), DIMENSION(0:nd-1), PARAMETER             :: diridx =  (/       -1,        1,        0,        3,        2 /)
2358                                                                     !< index for surf_*_v: 0:3 = (North, South, East, West)
[2737]2359        INTEGER(iwp), DIMENSION(0:nd-1)                        :: dirstart
2360        INTEGER(iwp), DIMENSION(0:nd-1)                        :: dirend
[2920]2361        INTEGER(iwp)                                           :: ids,idsint,idsidx,isurf,isvf,isurfs,isurflt
[2737]2362        INTEGER(iwp)                                           :: is,js,ks,i,j,k,iwl,istat, l, m
2363        INTEGER(iwp)                                           ::  k_topo    !< topography top index
2364
2365        dirstart = (/ startland, startwall, startwall, startwall, startwall /)
2366        dirend = (/ endland, endwall, endwall, endwall, endwall /)
2367
2368        found = .TRUE.
2369        temp_pf = -1._wp
2370       
2371        ids = -1
2372        var = TRIM(variable)
2373        DO i = 0, nd-1
2374            k = len(TRIM(var))
2375            j = len(TRIM(dirname(i)))
2376            IF ( var(k-j+1:k) == dirname(i) )  THEN
2377                ids = i
[2920]2378                idsint = dirint(ids)
2379                idsidx = diridx(ids)
[2737]2380                var = var(:k-j)
2381                EXIT
2382            ENDIF
2383        ENDDO
2384        IF ( ids == -1 )  THEN
2385            var = TRIM(variable)
2386        ENDIF
2387        IF ( var(1:11) == 'usm_t_wall_'  .AND.  len(TRIM(var)) >= 12 )  THEN
2388!--         wall layers
2389            READ(var(12:12), '(I1)', iostat=istat ) iwl
2390            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
2391                var = var(1:10)
2392            ENDIF
2393        ENDIF
2394        IF ( var(1:13) == 'usm_t_window_'  .AND.  len(TRIM(var)) >= 14 )  THEN
2395!--         window layers
2396            READ(var(14:14), '(I1)', iostat=istat ) iwl
2397            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
2398                var = var(1:12)
2399            ENDIF
2400        ENDIF
2401        IF ( var(1:12) == 'usm_t_green_'  .AND.  len(TRIM(var)) >= 13 )  THEN
2402!--         green layers
2403            READ(var(13:13), '(I1)', iostat=istat ) iwl
2404            IF ( istat == 0  .AND.  iwl >= nzb_wall  .AND.  iwl <= nzt_wall )  THEN
2405                var = var(1:11)
2406            ENDIF
2407        ENDIF
2408        IF ( (var(1:8) == 'usm_svf_'  .OR.  var(1:8) == 'usm_dif_')  .AND.  len(TRIM(var)) >= 13 )  THEN
2409!--         svf values to particular surface
2410            surfid = var(9:)
2411            i = index(surfid,'_')
2412            j = index(surfid(i+1:),'_')
2413            READ(surfid(1:i-1),*, iostat=istat ) is
2414            IF ( istat == 0 )  THEN
2415                READ(surfid(i+1:i+j-1),*, iostat=istat ) js
2416            ENDIF
2417            IF ( istat == 0 )  THEN
2418                READ(surfid(i+j+1:),*, iostat=istat ) ks
2419            ENDIF
2420            IF ( istat == 0 )  THEN
2421                var = var(1:7)
2422            ENDIF
2423        ENDIF
2424       
2425        SELECT CASE ( TRIM(var) )
2426
2427          CASE ( 'usm_surfz' )
2428!--           array of lw radiation falling to local surface after i-th reflection
[2920]2429              IF ( idsint == iup_u )  THEN
2430                 DO  m = 1, surf_usm_h%ns
2431                    i = surf_usm_h%i(m)
2432                    j = surf_usm_h%j(m)
2433                    k = surf_usm_h%k(m)
2434                    temp_pf(0,j,i) = MAX( temp_pf(0,j,i), REAL( k, kind=wp) )
2435                 ENDDO
2436              ELSE
2437                 l = idsidx
[2737]2438                 DO  m = 1, surf_usm_v(l)%ns
2439                    i = surf_usm_v(l)%i(m)
2440                    j = surf_usm_v(l)%j(m)
2441                    k = surf_usm_v(l)%k(m)
2442                    temp_pf(0,j,i) = MAX( temp_pf(0,j,i), REAL( k, kind=wp) + 1.0_wp )
2443                 ENDDO
[2920]2444              ENDIF
[2737]2445
2446          CASE ( 'usm_surfcat' )
2447!--           surface category
[2920]2448              IF ( idsint == iup_u )  THEN
2449                 DO  m = 1, surf_usm_h%ns
2450                    i = surf_usm_h%i(m)
2451                    j = surf_usm_h%j(m)
2452                    k = surf_usm_h%k(m)
2453                    temp_pf(k,j,i) = surf_usm_h%surface_types(m)
2454                 ENDDO
2455              ELSE
2456                 l = idsidx
[2737]2457                 DO  m = 1, surf_usm_v(l)%ns
2458                    i = surf_usm_v(l)%i(m)
2459                    j = surf_usm_v(l)%j(m)
2460                    k = surf_usm_v(l)%k(m)
2461                    temp_pf(k,j,i) = surf_usm_v(l)%surface_types(m)
2462                 ENDDO
[2920]2463              ENDIF
[2737]2464             
2465          CASE ( 'usm_surfalb' )
2466!--           surface albedo, weighted average
[2920]2467              IF ( idsint == iup_u )  THEN
2468                 DO  m = 1, surf_usm_h%ns
2469                    i = surf_usm_h%i(m)
2470                    j = surf_usm_h%j(m)
2471                    k = surf_usm_h%k(m)
[2963]2472                    temp_pf(k,j,i) = surf_usm_h%frac(ind_veg_wall,m)     *     &
2473                                     surf_usm_h%albedo(ind_veg_wall,m)  +      &
2474                                     surf_usm_h%frac(ind_pav_green,m)    *     &
2475                                     surf_usm_h%albedo(ind_pav_green,m) +      &
2476                                     surf_usm_h%frac(ind_wat_win,m)      *     &
2477                                     surf_usm_h%albedo(ind_wat_win,m)
[2920]2478                 ENDDO
2479              ELSE
2480                 l = idsidx
[2737]2481                 DO  m = 1, surf_usm_v(l)%ns
2482                    i = surf_usm_v(l)%i(m)
2483                    j = surf_usm_v(l)%j(m)
2484                    k = surf_usm_v(l)%k(m)
[2963]2485                    temp_pf(k,j,i) = surf_usm_v(l)%frac(ind_veg_wall,m)     *  &
2486                                     surf_usm_v(l)%albedo(ind_veg_wall,m)  +   &
2487                                     surf_usm_v(l)%frac(ind_pav_green,m)    *  &
2488                                     surf_usm_v(l)%albedo(ind_pav_green,m) +   &
2489                                     surf_usm_v(l)%frac(ind_wat_win,m)      *  &
2490                                     surf_usm_v(l)%albedo(ind_wat_win,m)
[2737]2491                 ENDDO
[2920]2492              ENDIF
[2737]2493             
2494          CASE ( 'usm_surfemis' )
2495!--           surface emissivity, weighted average
[2920]2496              IF ( idsint == iup_u )  THEN
2497                 DO  m = 1, surf_usm_h%ns
2498                    i = surf_usm_h%i(m)
2499                    j = surf_usm_h%j(m)
2500                    k = surf_usm_h%k(m)
[2963]2501                    temp_pf(k,j,i) =  surf_usm_h%frac(ind_veg_wall,m)      *   &
2502                                      surf_usm_h%emissivity(ind_veg_wall,m)  + &
2503                                      surf_usm_h%frac(ind_pav_green,m)     *   &
2504                                      surf_usm_h%emissivity(ind_pav_green,m) + &
2505                                      surf_usm_h%frac(ind_wat_win,m)       *   &
2506                                      surf_usm_h%emissivity(ind_wat_win,m)
[2920]2507                 ENDDO
2508              ELSE
2509                 l = idsidx
[2737]2510                 DO  m = 1, surf_usm_v(l)%ns
2511                    i = surf_usm_v(l)%i(m)
2512                    j = surf_usm_v(l)%j(m)
2513                    k = surf_usm_v(l)%k(m)
[2963]2514                    temp_pf(k,j,i) = surf_usm_v(l)%frac(ind_veg_wall,m)       *&
2515                                     surf_usm_v(l)%emissivity(ind_veg_wall,m) +&
2516                                     surf_usm_v(l)%frac(ind_pav_green,m)      *&
2517                                     surf_usm_v(l)%emissivity(ind_pav_green,m)+&
2518                                     surf_usm_v(l)%frac(ind_wat_win,m)        *&
2519                                     surf_usm_v(l)%emissivity(ind_wat_win,m)
[2737]2520                 ENDDO
[2920]2521              ENDIF
[2737]2522
2523          CASE ( 'usm_surfwintrans' )
2524!--           transmissivity window tiles
[2920]2525              IF ( idsint == iup_u )  THEN
2526                 DO  m = 1, surf_usm_h%ns
2527                    i = surf_usm_h%i(m)
2528                    j = surf_usm_h%j(m)
2529                    k = surf_usm_h%k(m)
2530                    temp_pf(k,j,i) = surf_usm_h%transmissivity(m)
2531                 ENDDO
2532              ELSE
2533                 l = idsidx
[2737]2534                 DO  m = 1, surf_usm_v(l)%ns
2535                    i = surf_usm_v(l)%i(m)
2536                    j = surf_usm_v(l)%j(m)
2537                    k = surf_usm_v(l)%k(m)
2538                    temp_pf(k,j,i) = surf_usm_v(l)%transmissivity(m)
2539                 ENDDO
[2920]2540              ENDIF
[2737]2541
[2920]2542          CASE ( 'usm_skyvf' )
2543!--           sky view factor
2544              DO isurf = dirstart(ids), dirend(ids)
2545                 IF ( surfl(id,isurf) == idsint )  THEN
2546                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = skyvf(isurf)
2547                 ENDIF
[2737]2548              ENDDO
[2920]2549             
2550          CASE ( 'usm_skyvft' )
2551!--           sky view factor
2552              DO isurf = dirstart(ids), dirend(ids)
2553                 IF ( surfl(id,isurf) == ids )  THEN
2554                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = skyvft(isurf)
2555                 ENDIF
2556              ENDDO
[2737]2557
2558!
2559!-- Not adjusted so far             
2560          CASE ( 'usm_svf', 'usm_dif' )
2561!--           shape view factors or iradiance factors to selected surface
2562              IF ( TRIM(var)=='usm_svf' )  THEN
2563                  k = 1
2564              ELSE
2565                  k = 2
2566              ENDIF
2567              DO isvf = 1, nsvfl
2568                  isurflt = svfsurf(1, isvf)
2569                  isurfs = svfsurf(2, isvf)
2570                             
2571                  IF ( surf(ix,isurfs) == is  .AND.  surf(iy,isurfs) == js  .AND.       &
[2920]2572                       surf(iz,isurfs) == ks  .AND.  surf(id,isurfs) == idsint )  THEN
[2737]2573  !--                 correct source surface
2574                      temp_pf(surfl(iz,isurflt),surfl(iy,isurflt),surfl(ix,isurflt)) = svf(k,isvf)
2575                  ENDIF
2576              ENDDO
2577
2578          CASE ( 'usm_rad_net' )
2579!--           array of complete radiation balance
2580              IF ( av == 0 )  THEN
[2920]2581                 IF ( idsint == iup_u )  THEN
2582                    DO  m = 1, surf_usm_h%ns
2583                       i = surf_usm_h%i(m)
2584                       j = surf_usm_h%j(m)
2585                       k = surf_usm_h%k(m)
2586                       temp_pf(k,j,i) = surf_usm_h%rad_net_l(m)
2587                    ENDDO
2588                 ELSE
2589                    l = idsidx
[2737]2590                    DO  m = 1, surf_usm_v(l)%ns
2591                       i = surf_usm_v(l)%i(m)
2592                       j = surf_usm_v(l)%j(m)
2593                       k = surf_usm_v(l)%k(m)
2594                       temp_pf(k,j,i) = surf_usm_v(l)%rad_net_l(m)
2595                    ENDDO
[2920]2596                 ENDIF
[2737]2597              ELSE
[2920]2598                 IF ( idsint == iup_u )  THEN
2599                    DO  m = 1, surf_usm_h%ns
2600                       i = surf_usm_h%i(m)
2601                       j = surf_usm_h%j(m)
2602                       k = surf_usm_h%k(m)
2603                       temp_pf(k,j,i) = surf_usm_h%rad_net_av(m)
2604                    ENDDO
2605                 ELSE
2606                    l = idsidx
[2737]2607                    DO  m = 1, surf_usm_v(l)%ns
2608                       i = surf_usm_v(l)%i(m)
2609                       j = surf_usm_v(l)%j(m)
2610                       k = surf_usm_v(l)%k(m)
2611                       temp_pf(k,j,i) = surf_usm_v(l)%rad_net_av(m)
2612                    ENDDO
[2920]2613                 ENDIF
[2737]2614              ENDIF
2615
2616          CASE ( 'usm_rad_insw' )
2617!--           array of sw radiation falling to surface after i-th reflection
2618              DO isurf = dirstart(ids), dirend(ids)
[2920]2619                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2620                   IF ( av == 0 )  THEN
2621                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinsw(isurf)
2622                   ELSE
2623                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinsw_av(isurf)
2624                   ENDIF
2625                 ENDIF
2626              ENDDO
2627
2628          CASE ( 'usm_rad_inlw' )
2629!--           array of lw radiation falling to surface after i-th reflection
2630              DO isurf = dirstart(ids), dirend(ids)
[2920]2631                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2632                   IF ( av == 0 )  THEN
2633                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw(isurf)
2634                   ELSE
2635                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw_av(isurf)
2636                   ENDIF
2637                 ENDIF
2638              ENDDO
2639
2640          CASE ( 'usm_rad_inswdir' )
2641!--           array of direct sw radiation falling to surface from sun
2642              DO isurf = dirstart(ids), dirend(ids)
[2920]2643                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2644                   IF ( av == 0 )  THEN
2645                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdir(isurf)
2646                   ELSE
2647                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdir_av(isurf)
2648                   ENDIF
2649                 ENDIF
2650              ENDDO
2651
2652          CASE ( 'usm_rad_inswdif' )
2653!--           array of difusion sw radiation falling to surface from sky and borders of the domain
2654              DO isurf = dirstart(ids), dirend(ids)
[2920]2655                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2656                   IF ( av == 0 )  THEN
2657                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdif(isurf)
2658                   ELSE
2659                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswdif_av(isurf)
2660                   ENDIF
2661                 ENDIF
2662              ENDDO
2663
2664          CASE ( 'usm_rad_inswref' )
2665!--           array of sw radiation falling to surface from reflections
2666              DO isurf = dirstart(ids), dirend(ids)
[2920]2667                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2668                   IF ( av == 0 )  THEN
2669                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = &
2670                       surfinsw(isurf) - surfinswdir(isurf) - surfinswdif(isurf)
2671                   ELSE
2672                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinswref_av(isurf)
2673                   ENDIF
2674                 ENDIF
2675              ENDDO
2676
[2920]2677          CASE ( 'usm_rad_inlwdif' )
2678!--           array of difusion lw radiation falling to surface from sky and borders of the domain
2679              DO isurf = dirstart(ids), dirend(ids)
2680                 IF ( surfl(id,isurf) == idsint )  THEN
2681                   IF ( av == 0 )  THEN
2682                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwdif(isurf)
2683                   ELSE
2684                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwdif_av(isurf)
2685                   ENDIF
2686                 ENDIF
2687              ENDDO
2688
[2737]2689          CASE ( 'usm_rad_inlwref' )
2690!--           array of lw radiation falling to surface from reflections
2691              DO isurf = dirstart(ids), dirend(ids)
[2920]2692                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2693                   IF ( av == 0 )  THEN
2694                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlw(isurf) - surfinlwdif(isurf)
2695                   ELSE
2696                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinlwref_av(isurf)
2697                   ENDIF
2698                 ENDIF
2699              ENDDO
2700
2701          CASE ( 'usm_rad_outsw' )
2702!--           array of sw radiation emitted from surface after i-th reflection
2703              DO isurf = dirstart(ids), dirend(ids)
[2920]2704                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2705                   IF ( av == 0 )  THEN
2706                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutsw(isurf)
2707                   ELSE
2708                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutsw_av(isurf)
2709                   ENDIF
2710                 ENDIF
2711              ENDDO
2712
2713          CASE ( 'usm_rad_outlw' )
2714!--           array of lw radiation emitted from surface after i-th reflection
2715              DO isurf = dirstart(ids), dirend(ids)
[2920]2716                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2717                   IF ( av == 0 )  THEN
2718                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutlw(isurf)
2719                   ELSE
2720                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfoutlw_av(isurf)
2721                   ENDIF
2722                 ENDIF
2723              ENDDO
2724
2725          CASE ( 'usm_rad_ressw' )
2726!--           average of array of residua of sw radiation absorbed in surface after last reflection
2727              DO isurf = dirstart(ids), dirend(ids)
[2920]2728                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2729                   IF ( av == 0 )  THEN
2730                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfins(isurf)
2731                   ELSE
2732                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfins_av(isurf)
2733                   ENDIF
2734                 ENDIF
2735              ENDDO
2736
2737          CASE ( 'usm_rad_reslw' )
2738!--           average of array of residua of lw radiation absorbed in surface after last reflection
2739              DO isurf = dirstart(ids), dirend(ids)
[2920]2740                 IF ( surfl(id,isurf) == idsint )  THEN
[2737]2741                   IF ( av == 0 )  THEN
2742                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinl(isurf)
2743                   ELSE
2744                     temp_pf(surfl(iz,isurf),surfl(iy,isurf),surfl(ix,isurf)) = surfinl_av(isurf)
2745                   ENDIF
2746                 ENDIF
2747              ENDDO
2748 
2749          CASE ( 'usm_rad_hf' )
2750!--           array of heat flux from radiation for surfaces after all reflections
2751              IF ( av == 0 )  THEN
[2920]2752                 IF ( idsint == iup_u )  THEN
2753                    DO  m = 1, surf_usm_h%ns
2754                       i = surf_usm_h%i(m)
2755                       j = surf_usm_h%j(m)
2756                       k = surf_usm_h%k(m)
2757                       temp_pf(k,j,i) = surf_usm_h%surfhf(m)
2758                    ENDDO
2759                 ELSE
2760                    l = idsidx
[2737]2761                    DO  m = 1, surf_usm_v(l)%ns
2762                       i = surf_usm_v(l)%i(m)
2763                       j = surf_usm_v(l)%j(m)
2764                       k = surf_usm_v(l)%k(m)
2765                       temp_pf(k,j,i) = surf_usm_v(l)%surfhf(m)
2766                    ENDDO
[2920]2767                 ENDIF
[2737]2768              ELSE
[2920]2769                 IF ( idsint == iup_u )  THEN
2770                    DO  m = 1, surf_usm_h%ns
2771                       i = surf_usm_h%i(m)
2772                       j = surf_usm_h%j(m)
2773                       k = surf_usm_h%k(m)
2774                       temp_pf(k,j,i) = surf_usm_h%surfhf_av(m)
2775                    ENDDO
2776                 ELSE
2777                    l = idsidx
[2737]2778                    DO  m = 1, surf_usm_v(l)%ns
2779                       i = surf_usm_v(l)%i(m)
2780                       j = surf_usm_v(l)%j(m)
2781                       k = surf_usm_v(l)%k(m)
2782                       temp_pf(k,j,i) = surf_usm_v(l)%surfhf_av(m)
2783                    ENDDO
[2920]2784                 ENDIF
[2737]2785              ENDIF
2786 
2787          CASE ( 'usm_wshf' )
2788!--           array of sensible heat flux from surfaces
2789              IF ( av == 0 )  THEN
[2920]2790                 IF ( idsint == iup_u )  THEN
2791                    DO  m = 1, surf_usm_h%ns
2792                       i = surf_usm_h%i(m)
2793                       j = surf_usm_h%j(m)
2794                       k = surf_usm_h%k(m)
2795                       temp_pf(k,j,i) = surf_usm_h%wshf_eb(m)
2796                    ENDDO
2797                 ELSE
2798                    l = idsidx
[2737]2799                    DO  m = 1, surf_usm_v(l)%ns
2800                       i = surf_usm_v(l)%i(m)
2801                       j = surf_usm_v(l)%j(m)
2802                       k = surf_usm_v(l)%k(m)
2803                       temp_pf(k,j,i) = surf_usm_v(l)%wshf_eb(m)
2804                    ENDDO
[2920]2805                 ENDIF
[2737]2806              ELSE
[2920]2807                 IF ( idsint == iup_u )  THEN
2808                    DO  m = 1, surf_usm_h%ns
2809                       i = surf_usm_h%i(m)
2810                       j = surf_usm_h%j(m)
2811                       k = surf_usm_h%k(m)
2812                       temp_pf(k,j,i) = surf_usm_h%wshf_eb_av(m)
2813                    ENDDO
2814                 ELSE
2815                    l = idsidx
[2737]2816                    DO  m = 1, surf_usm_v(l)%ns
2817                       i = surf_usm_v(l)%i(m)
2818                       j = surf_usm_v(l)%j(m)
2819                       k = surf_usm_v(l)%k(m)
2820                       temp_pf(k,j,i) = surf_usm_v(l)%wshf_eb_av(m)
2821                    ENDDO
[2920]2822                 ENDIF
[2737]2823              ENDIF
2824
2825
2826          CASE ( 'usm_wghf' )
2827!--           array of heat flux from ground (land, wall, roof)
2828              IF ( av == 0 )  THEN
[2920]2829                 IF ( idsint == iup_u )  THEN
2830                    DO  m = 1, surf_usm_h%ns
2831                       i = surf_usm_h%i(m)
2832                       j = surf_usm_h%j(m)
2833                       k = surf_usm_h%k(m)
2834                       temp_pf(k,j,i) = surf_usm_h%wghf_eb(m)
2835                    ENDDO
2836                 ELSE
2837                    l = idsidx
[2737]2838                    DO  m = 1, surf_usm_v(l)%ns
2839                       i = surf_usm_v(l)%i(m)
2840                       j = surf_usm_v(l)%j(m)
2841                       k = surf_usm_v(l)%k(m)
2842                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb(m)
2843                    ENDDO
[2920]2844                 ENDIF
[2737]2845              ELSE
[2920]2846                 IF ( idsint == iup_u )  THEN
2847                    DO  m = 1, surf_usm_h%ns
2848                       i = surf_usm_h%i(m)
2849                       j = surf_usm_h%j(m)
2850                       k = surf_usm_h%k(m)
2851                       temp_pf(k,j,i) = surf_usm_h%wghf_eb_av(m)
2852                    ENDDO
2853                 ELSE
2854                    l = idsidx
[2737]2855                    DO  m = 1, surf_usm_v(l)%ns
2856                       i = surf_usm_v(l)%i(m)
2857                       j = surf_usm_v(l)%j(m)
2858                       k = surf_usm_v(l)%k(m)
2859                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_av(m)
2860                    ENDDO
[2920]2861                 ENDIF
[2737]2862              ENDIF
2863
2864          CASE ( 'usm_wghf_window' )
2865!--           array of heat flux from window ground (land, wall, roof)
2866
2867              IF ( av == 0 )  THEN
[2920]2868                 IF ( idsint == iup_u )  THEN
2869                    DO  m = 1, surf_usm_h%ns
2870                       i = surf_usm_h%i(m)
2871                       j = surf_usm_h%j(m)
2872                       k = surf_usm_h%k(m)
2873                       temp_pf(k,j,i) = surf_usm_h%wghf_eb_window(m)
2874                    ENDDO
2875                 ELSE
2876                    l = idsidx
[2737]2877                    DO  m = 1, surf_usm_v(l)%ns
2878                       i = surf_usm_v(l)%i(m)
2879                       j = surf_usm_v(l)%j(m)
2880                       k = surf_usm_v(l)%k(m)
2881                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_window(m)
2882                    ENDDO
[2920]2883                 ENDIF
[2737]2884              ELSE
[2920]2885                 IF ( idsint == iup_u )  THEN
2886                    DO  m = 1, surf_usm_h%ns
2887                       i = surf_usm_h%i(m)
2888                       j = surf_usm_h%j(m)
2889                       k = surf_usm_h%k(m)
2890                       temp_pf(k,j,i) = surf_usm_h%wghf_eb_window_av(m)
2891                    ENDDO
2892                 ELSE
2893                    l = idsidx
[2737]2894                    DO  m = 1, surf_usm_v(l)%ns
2895                       i = surf_usm_v(l)%i(m)
2896                       j = surf_usm_v(l)%j(m)
2897                       k = surf_usm_v(l)%k(m)
2898                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_window_av(m)
2899                    ENDDO
[2920]2900                 ENDIF
[2737]2901              ENDIF
2902
2903          CASE ( 'usm_wghf_green' )
2904!--           array of heat flux from green ground (land, wall, roof)
2905
2906              IF ( av == 0 )  THEN
[2920]2907                 IF ( idsint == iup_u )  THEN
2908                    DO  m = 1, surf_usm_h%ns
2909                       i = surf_usm_h%i(m)
2910                       j = surf_usm_h%j(m)
2911                       k = surf_usm_h%k(m)
2912                       temp_pf(k,j,i) = surf_usm_h%wghf_eb_green(m)
2913                    ENDDO
2914                 ELSE
2915                    l = idsidx
[2737]2916                    DO  m = 1, surf_usm_v(l)%ns
2917                       i = surf_usm_v(l)%i(m)
2918                       j = surf_usm_v(l)%j(m)
2919                       k = surf_usm_v(l)%k(m)
2920                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_green(m)
2921                    ENDDO
[2920]2922                 ENDIF
[2737]2923              ELSE
[2920]2924                 IF ( idsint == iup_u )  THEN
2925                    DO  m = 1, surf_usm_h%ns
2926                       i = surf_usm_h%i(m)
2927                       j = surf_usm_h%j(m)
2928                       k = surf_usm_h%k(m)
2929                       temp_pf(k,j,i) = surf_usm_h%wghf_eb_green_av(m)
2930                    ENDDO
2931                 ELSE
2932                    l = idsidx
[2737]2933                    DO  m = 1, surf_usm_v(l)%ns
2934                       i = surf_usm_v(l)%i(m)
2935                       j = surf_usm_v(l)%j(m)
2936                       k = surf_usm_v(l)%k(m)
2937                       temp_pf(k,j,i) = surf_usm_v(l)%wghf_eb_green_av(m)
2938                    ENDDO
[2920]2939                 ENDIF
[2737]2940              ENDIF
2941
2942          CASE ( 'usm_iwghf' )
2943!--           array of heat flux from indoor ground (land, wall, roof)
2944              IF ( av == 0 )  THEN
[2920]2945                 IF ( idsint == iup_u )  THEN
2946                    DO  m = 1, surf_usm_h%ns
2947                       i = surf_usm_h%i(m)
2948                       j = surf_usm_h%j(m)
2949                       k = surf_usm_h%k(m)
2950                       temp_pf(k,j,i) = surf_usm_h%iwghf_eb(m)
2951                    ENDDO
2952                 ELSE
2953                    l = idsidx
[2737]2954                    DO  m = 1, surf_usm_v(l)%ns
2955                       i = surf_usm_v(l)%i(m)
2956                       j = surf_usm_v(l)%j(m)
2957                       k = surf_usm_v(l)%k(m)
2958                       temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb(m)
2959                    ENDDO
[2920]2960                 ENDIF
[2737]2961              ELSE
[2920]2962                 IF ( idsint == iup_u )  THEN
2963                    DO  m = 1, surf_usm_h%ns
2964                       i = surf_usm_h%i(m)
2965                       j = surf_usm_h%j(m)
2966                       k = surf_usm_h%k(m)
2967                       temp_pf(k,j,i) = surf_usm_h%iwghf_eb_av(m)
2968                    ENDDO
2969                 ELSE
2970                    l = idsidx
[2737]2971                    DO  m = 1, surf_usm_v(l)%ns
2972                       i = surf_usm_v(l)%i(m)
2973                       j = surf_usm_v(l)%j(m)
2974                       k = surf_usm_v(l)%k(m)
2975                       temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_av(m)
2976                    ENDDO
[2920]2977                 ENDIF
[2737]2978              ENDIF
2979
2980          CASE ( 'usm_iwghf_window' )
2981!--           array of heat flux from indoor window ground (land, wall, roof)
2982
2983              IF ( av == 0 )  THEN
[2920]2984                 IF ( idsint == iup_u )  THEN
2985                    DO  m = 1, surf_usm_h%ns
2986                       i = surf_usm_h%i(m)
2987                       j = surf_usm_h%j(m)
2988                       k = surf_usm_h%k(m)
2989                       temp_pf(k,j,i) = surf_usm_h%iwghf_eb_window(m)
2990                    ENDDO
2991                 ELSE
2992                    l = idsidx
[2737]2993                    DO  m = 1, surf_usm_v(l)%ns
2994                       i = surf_usm_v(l)%i(m)
2995                       j = surf_usm_v(l)%j(m)
2996                       k = surf_usm_v(l)%k(m)
2997                       temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_window(m)
2998                    ENDDO
[2920]2999                 ENDIF
[2737]3000              ELSE
[2920]3001                 IF ( idsint == iup_u )  THEN
3002                    DO  m = 1, surf_usm_h%ns
3003                       i = surf_usm_h%i(m)
3004                       j = surf_usm_h%j(m)
3005                       k = surf_usm_h%k(m)
3006                       temp_pf(k,j,i) = surf_usm_h%iwghf_eb_window_av(m)
3007                    ENDDO
3008                 ELSE
3009                    l = idsidx
[2737]3010                    DO  m = 1, surf_usm_v(l)%ns
3011                       i = surf_usm_v(l)%i(m)
3012                       j = surf_usm_v(l)%j(m)
3013                       k = surf_usm_v(l)%k(m)
3014                       temp_pf(k,j,i) = surf_usm_v(l)%iwghf_eb_window_av(m)
3015                    ENDDO
[2920]3016                 ENDIF
[2737]3017              ENDIF
3018             
3019          CASE ( 'usm_t_surf' )
3020!--           surface temperature for surfaces
3021              IF ( av == 0 )  THEN
[2920]3022                 IF ( idsint == iup_u )  THEN
3023                    DO  m = 1, surf_usm_h%ns
3024                       i = surf_usm_h%i(m)
3025                       j = surf_usm_h%j(m)
3026                       k = surf_usm_h%k(m)
3027                       temp_pf(k,j,i