Ignore:
Timestamp:
Sep 9, 2020 8:27:58 PM (4 years ago)
Author:
pavelkrc
Message:

Radiative transfer model RTM version 4.1

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/land_surface_model_mod.f90

    r4669 r4671  
    2525! -----------------
    2626! $Id$
     27! Implementation of downward facing USM and LSM surfaces for RTM 4.1
     28! Author: J. Resler (Institute of Computer Science, Prague)
     29!
     30! 4669 2020-09-09 13:43:47Z pavelkrc
    2731! Fix calculation of force_radiation_call
    2832!
     
    449453                                   zs_layer = 9999999.9_wp             !< soil layer depths (edge)
    450454
    451     TYPE(surf_type_lsm), POINTER ::  t_soil_h,    & !< Soil temperature (K), horizontal surface elements
    452                                      t_soil_h_p,  & !< Prog. soil temperature (K), horizontal surface elements
    453                                      m_soil_h,    & !< Soil moisture (m3/m3), horizontal surface elements
    454                                      m_soil_h_p     !< Prog. soil moisture (m3/m3), horizontal surface elements
    455 
    456     TYPE(surf_type_lsm), TARGET  ::  t_soil_h_1,  & !<
    457                                      t_soil_h_2,  & !<
    458                                      m_soil_h_1,  & !<
    459                                      m_soil_h_2     !<
    460 
    461     TYPE(surf_type_lsm), DIMENSION(:), POINTER :: &
    462                                      t_soil_v,    & !< Soil temperature (K), vertical surface elements
    463                                      t_soil_v_p,  & !< Prog. soil temperature (K), vertical surface elements
    464                                      m_soil_v,    & !< Soil moisture (m3/m3), vertical surface elements
    465                                      m_soil_v_p     !< Prog. soil moisture (m3/m3), vertical surface elements
    466 
    467     TYPE(surf_type_lsm), DIMENSION(0:3), TARGET ::&
    468                                      t_soil_v_1,  & !<
    469                                      t_soil_v_2,  & !<
    470                                      m_soil_v_1,  & !<
    471                                      m_soil_v_2     !<
    472 
    473     TYPE(surf_type_lsm), POINTER  ::  t_surface_h,    & !< surface temperature (K), horizontal surface elements
     455    TYPE(surf_type_lsm), DIMENSION(:), POINTER ::     &
     456                                     t_soil_h,        & !< Soil temperature (K), horizontal surface elements
     457                                     t_soil_h_p,      & !< Prog. soil temperature (K), horizontal surface elements
     458                                     m_soil_h,        & !< Soil moisture (m3/m3), horizontal surface elements
     459                                     m_soil_h_p         !< Prog. soil moisture (m3/m3), horizontal surface elements
     460
     461    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::    &
     462                                     t_soil_h_1,      & !<
     463                                     t_soil_h_2,      & !<
     464                                     m_soil_h_1,      & !<
     465                                     m_soil_h_2         !<
     466
     467    TYPE(surf_type_lsm), DIMENSION(:), POINTER ::     &
     468                                     t_soil_v,        & !< Soil temperature (K), vertical surface elements
     469                                     t_soil_v_p,      & !< Prog. soil temperature (K), vertical surface elements
     470                                     m_soil_v,        & !< Soil moisture (m3/m3), vertical surface elements
     471                                     m_soil_v_p         !< Prog. soil moisture (m3/m3), vertical surface elements
     472
     473    TYPE(surf_type_lsm), DIMENSION(0:3), TARGET ::    &
     474                                     t_soil_v_1,      & !<
     475                                     t_soil_v_2,      & !<
     476                                     m_soil_v_1,      & !<
     477                                     m_soil_v_2        !<
     478
     479    TYPE(surf_type_lsm), DIMENSION(:), POINTER ::&
     480                                      t_surface_h,    & !< surface temperature (K), horizontal surface elements
    474481                                      t_surface_h_p,  & !< progn. surface temperature (K), horizontal surface elements
    475482                                      m_liq_h,        & !< liquid water reservoir (m), horizontal surface elements
    476483                                      m_liq_h_p         !< progn. liquid water reservoir (m), horizontal surface elements
    477484
    478     TYPE(surf_type_lsm), TARGET   ::  t_surface_h_1,  & !<
     485    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::    &
     486                                      t_surface_h_1,  & !<
    479487                                      t_surface_h_2,  & !<
    480488                                      m_liq_h_1,      & !<
     
    498506                                                        m_soil_av    !< Average of m_soil
    499507
    500     TYPE(surf_type_lsm), TARGET ::  tm_liq_h_m      !< liquid water reservoir tendency (m), horizontal surface elements
    501     TYPE(surf_type_lsm), TARGET ::  tt_surface_h_m  !< surface temperature tendency (K), horizontal surface elements
    502     TYPE(surf_type_lsm), TARGET ::  tt_soil_h_m     !< t_soil storage array, horizontal surface elements
    503     TYPE(surf_type_lsm), TARGET ::  tm_soil_h_m     !< m_soil storage array, horizontal surface elements
     508    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::  tm_liq_h_m      !< liquid water reservoir tendency (m), horizontal surface elements
     509    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::  tt_surface_h_m  !< surface temperature tendency (K), horizontal surface elements
     510    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::  tt_soil_h_m     !< t_soil storage array, horizontal surface elements
     511    TYPE(surf_type_lsm), DIMENSION(0:1), TARGET ::  tm_soil_h_m     !< m_soil storage array, horizontal surface elements
    504512
    505513    TYPE(surf_type_lsm), DIMENSION(0:3), TARGET ::  tm_liq_v_m      !< liquid water reservoir tendency (m), vertical surface elements
     
    871879    INTEGER(iwp) :: m      !< running index surface elements
    872880
    873     koff = surf_lsm_h%koff
    874     DO  m = 1, surf_lsm_h%ns
    875        i = surf_lsm_h%i(m)
    876        j = surf_lsm_h%j(m)
    877        k = surf_lsm_h%k(m)
    878        pt(k+koff,j,i) = pt(k,j,i)
     881    DO l = 0, 1
     882       koff = surf_lsm_h(l)%koff
     883       DO  m = 1, surf_lsm_h(l)%ns
     884          i = surf_lsm_h(l)%i(m)
     885          j = surf_lsm_h(l)%j(m)
     886          k = surf_lsm_h(l)%k(m)
     887          pt(k+koff,j,i) = pt(k,j,i)
     888       ENDDO
    879889    ENDDO
    880890
     
    892902!-- In case of humidity, set boundary conditions also for q and vpt.
    893903    IF ( humidity )  THEN
    894        koff = surf_lsm_h%koff
    895        DO  m = 1, surf_lsm_h%ns
    896           i = surf_lsm_h%i(m)
    897           j = surf_lsm_h%j(m)
    898           k = surf_lsm_h%k(m)
    899           q(k+koff,j,i)   = q(k,j,i)
    900           vpt(k+koff,j,i) = vpt(k,j,i)
     904       DO  l = 0, 1
     905          koff = surf_lsm_h(l)%koff
     906          DO  m = 1, surf_lsm_h(l)%ns
     907             i = surf_lsm_h(l)%i(m)
     908             j = surf_lsm_h(l)%j(m)
     909             k = surf_lsm_h(l)%k(m)
     910             q(k+koff,j,i)   = q(k,j,i)
     911             vpt(k+koff,j,i) = vpt(k,j,i)
     912          ENDDO
    901913       ENDDO
    902 
    903914       DO  l = 0, 3
    904915          ioff = surf_lsm_v(l)%ioff
     
    15971608! Description:
    15981609! ------------
     1610!> Calls energy balance solver and soil energy model the surfaces in all directions.
     1611!------------------------------------------------------------------------------!
     1612 SUBROUTINE lsm_energy_balance ( during_spinup )
     1613
     1614    USE control_parameters,                                                                        &
     1615        ONLY:  calc_soil_moisture_during_spinup
     1616
     1617    LOGICAL      ::  during_spinup  !< flag indicating soil/wall spinup phase
     1618    LOGICAL      ::  calc_sm
     1619    INTEGER(iwp) ::  l              !< loop index for surface types
     1620
     1621    calc_sm = .NOT. during_spinup .OR. calc_soil_moisture_during_spinup
     1622!
     1623!-- Call for horizontal surfaces
     1624    DO l = 0, 1
     1625       CALL lsm_surface_energy_balance( .TRUE., l )
     1626       CALL lsm_soil_model( .TRUE., l, calc_sm )
     1627    ENDDO
     1628!
     1629!--   Call for vertical surfaces
     1630    DO l = 0, 3
     1631       CALL lsm_surface_energy_balance( .FALSE., l )
     1632       CALL lsm_soil_model( .FALSE., l, calc_sm )
     1633    ENDDO
     1634
     1635END SUBROUTINE lsm_energy_balance
     1636
     1637!------------------------------------------------------------------------------!
     1638! Description:
     1639! ------------
    15991640!> Solver for the energy balance at the surface.
    16001641!------------------------------------------------------------------------------!
    1601  SUBROUTINE lsm_energy_balance( horizontal, l )
     1642 SUBROUTINE lsm_surface_energy_balance( horizontal, l )
    16021643
    16031644    USE pegrid
     
    16171658
    16181659    LOGICAL      ::  horizontal !< Flag indicating horizontal or vertical surfaces
     1660    LOGICAL      ::  upward     !< Flag indicating upward horizontal surfaces
    16191661
    16201662    REAL(wp) :: c_surface_tmp,  & !< temporary variable for storing the volumetric heat capacity of the surface
     
    16541696
    16551697    IF ( debug_output_timestep )  THEN
    1656        WRITE( debug_string, * ) 'lsm_energy_balance', horizontal, l
     1698       WRITE( debug_string, * ) 'lsm_surface_energy_balance', horizontal, l
    16571699       CALL debug_message( debug_string, 'start' )
    16581700    ENDIF
    16591701
     1702    upward = .FALSE.
    16601703    IF ( horizontal )  THEN
    1661        surf              => surf_lsm_h
    1662 
    1663        surf_t_surface    => t_surface_h
    1664        surf_t_surface_p  => t_surface_h_p
    1665        surf_tt_surface_m => tt_surface_h_m
    1666        surf_m_liq        => m_liq_h
    1667        surf_m_liq_p      => m_liq_h_p
    1668        surf_tm_liq_m     => tm_liq_h_m
    1669        surf_m_soil       => m_soil_h
    1670        surf_t_soil       => t_soil_h
     1704       surf              => surf_lsm_h(l)
     1705       surf_t_surface    => t_surface_h(l)
     1706       surf_t_surface_p  => t_surface_h_p(l)
     1707       surf_tt_surface_m => tt_surface_h_m(l)
     1708       surf_m_liq        => m_liq_h(l)
     1709       surf_m_liq_p      => m_liq_h_p(l)
     1710       surf_tm_liq_m     => tm_liq_h_m(l)
     1711       surf_m_soil       => m_soil_h(l)
     1712       surf_t_soil       => t_soil_h(l)
     1713       IF ( l == 0 ) upward = .TRUE.
    16711714    ELSE
    16721715       surf              => surf_lsm_v(l)
    1673 
    16741716       surf_t_surface    => t_surface_v(l)
    16751717       surf_t_surface_p  => t_surface_v_p(l)
     
    17771819!        ENDIF
    17781820!
    1779 !--     Calculation of r_a for vertical surfaces
     1821!--     Calculation of r_a for vertical and horizontal downward surfaces
    17801822!--
    17811823!--     heat transfer coefficient for forced convection along vertical walls
     
    18061848!--    holes, the resistance can become negative. For this reason, limit r_a
    18071849!--    to positive values.
    1808        IF ( horizontal .OR.  .NOT. aero_resist_kray )  THEN
     1850       IF ( upward .OR.  .NOT. aero_resist_kray )  THEN
    18091851          surf%r_a(m) = ABS( ( surf%pt1(m) - surf%pt_surface(m) ) /                               &
    18101852                             ( surf%ts(m) * surf%us(m) + 1.0E-20_wp ) )
     
    20552097!
    20562098! update the 3d field of rad_lw_out array to have consistent output
    2057        IF ( horizontal ) THEN
     2099       IF ( upward ) THEN
    20582100          IF ( radiation_scheme == 'rrtmg' ) THEN
    20592101             rad_lw_out(k+k_off,j+j_off,i+i_off) = surf%rad_lw_out(m)
     
    22172259!
    22182260!-- Calculate new roughness lengths (for water surfaces only)
    2219     IF ( horizontal .AND.  .NOT. constant_roughness )  CALL calc_z0_water_surface
     2261    IF ( upward .AND.  .NOT. constant_roughness )  CALL calc_z0_water_surface
    22202262
    22212263    IF ( debug_output_timestep )  THEN
    2222        WRITE( debug_string, * ) 'lsm_energy_balance', horizontal, l
     2264       WRITE( debug_string, * ) 'lsm_surface_energy_balance', horizontal, l
    22232265       CALL debug_message( debug_string, 'end' )
    22242266    ENDIF
     
    22822324    END SUBROUTINE calc_q_surface
    22832325
    2284  END SUBROUTINE lsm_energy_balance
     2326 END SUBROUTINE lsm_surface_energy_balance
    22852327
    22862328
     
    24432485!--    Set initial values for prognostic quantities
    24442486!--    Horizontal surfaces
    2445        tt_surface_h_m%var_1d = 0.0_wp
    2446        tt_soil_h_m%var_2d    = 0.0_wp
    2447        tm_soil_h_m%var_2d    = 0.0_wp
    2448        tm_liq_h_m%var_1d     = 0.0_wp
    2449        surf_lsm_h%c_liq      = 0.0_wp
    2450 
    2451        surf_lsm_h%ghf = 0.0_wp
    2452 
    2453        surf_lsm_h%qsws_liq  = 0.0_wp
    2454        surf_lsm_h%qsws_soil = 0.0_wp
    2455        surf_lsm_h%qsws_veg  = 0.0_wp
    2456 
    2457        surf_lsm_h%r_a        = 50.0_wp
    2458        surf_lsm_h%r_s        = 50.0_wp
    2459        surf_lsm_h%r_canopy   = 0.0_wp
    2460        surf_lsm_h%r_soil     = 0.0_wp
     2487       DO  l = 0, 1
     2488          tt_surface_h_m(l)%var_1d = 0.0_wp
     2489          tt_soil_h_m(l)%var_2d    = 0.0_wp
     2490          tm_soil_h_m(l)%var_2d    = 0.0_wp
     2491          tm_liq_h_m(l)%var_1d     = 0.0_wp
     2492          surf_lsm_h(l)%c_liq      = 0.0_wp
     2493
     2494          surf_lsm_h(l)%ghf = 0.0_wp
     2495
     2496          surf_lsm_h(l)%qsws_liq  = 0.0_wp
     2497          surf_lsm_h(l)%qsws_soil = 0.0_wp
     2498          surf_lsm_h(l)%qsws_veg  = 0.0_wp
     2499
     2500          surf_lsm_h(l)%r_a        = 50.0_wp
     2501          surf_lsm_h(l)%r_s        = 50.0_wp
     2502          surf_lsm_h(l)%r_canopy   = 0.0_wp
     2503          surf_lsm_h(l)%r_soil     = 0.0_wp
     2504       ENDDO
    24612505!
    24622506!--    Do the same for vertical surfaces
     
    24832527!--    Set initial values for prognostic soil quantities
    24842528       IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
    2485           t_soil_h%var_2d = 0.0_wp
    2486           m_soil_h%var_2d = 0.0_wp
    2487           m_liq_h%var_1d  = 0.0_wp
    2488 
     2529          DO  l = 0, 1
     2530             t_soil_h(l)%var_2d = 0.0_wp
     2531             m_soil_h(l)%var_2d = 0.0_wp
     2532             m_liq_h(l)%var_1d  = 0.0_wp
     2533          ENDDO
    24892534          DO  l = 0, 3
    24902535             t_soil_v(l)%var_2d = 0.0_wp
     
    24962541!--    Allocate 3D soil model arrays
    24972542!--    First, for horizontal surfaces
    2498        ALLOCATE ( surf_lsm_h%alpha_vg(nzb_soil:nzt_soil,1:surf_lsm_h%ns)    )
    2499        ALLOCATE ( surf_lsm_h%gamma_w_sat(nzb_soil:nzt_soil,1:surf_lsm_h%ns) )
    2500        ALLOCATE ( surf_lsm_h%lambda_h(nzb_soil:nzt_soil,1:surf_lsm_h%ns)    )
    2501        ALLOCATE ( surf_lsm_h%lambda_h_def(nzb_soil:nzt_soil,1:surf_lsm_h%ns))
    2502        ALLOCATE ( surf_lsm_h%l_vg(nzb_soil:nzt_soil,1:surf_lsm_h%ns)        )
    2503        ALLOCATE ( surf_lsm_h%m_fc(nzb_soil:nzt_soil,1:surf_lsm_h%ns)        )
    2504        ALLOCATE ( surf_lsm_h%m_res(nzb_soil:nzt_soil,1:surf_lsm_h%ns)       )
    2505        ALLOCATE ( surf_lsm_h%m_sat(nzb_soil:nzt_soil,1:surf_lsm_h%ns)       )
    2506        ALLOCATE ( surf_lsm_h%m_wilt(nzb_soil:nzt_soil,1:surf_lsm_h%ns)      )
    2507        ALLOCATE ( surf_lsm_h%n_vg(nzb_soil:nzt_soil,1:surf_lsm_h%ns)        )
    2508        ALLOCATE ( surf_lsm_h%rho_c_total(nzb_soil:nzt_soil,1:surf_lsm_h%ns) )
    2509        ALLOCATE ( surf_lsm_h%rho_c_total_def(nzb_soil:nzt_soil,1:surf_lsm_h%ns) )
    2510        ALLOCATE ( surf_lsm_h%root_fr(nzb_soil:nzt_soil,1:surf_lsm_h%ns)     )
    2511 
    2512        surf_lsm_h%lambda_h     = 0.0_wp
    2513 !
    2514 !--    If required, allocate humidity-related variables for the soil model
    2515        IF ( humidity )  THEN
    2516           ALLOCATE ( surf_lsm_h%lambda_w(nzb_soil:nzt_soil,1:surf_lsm_h%ns) )
    2517           ALLOCATE ( surf_lsm_h%gamma_w(nzb_soil:nzt_soil,1:surf_lsm_h%ns)  )
    2518 
    2519           surf_lsm_h%lambda_w = 0.0_wp
    2520        ENDIF
     2543       DO  l = 0, 1
     2544          ALLOCATE ( surf_lsm_h(l)%alpha_vg(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)    )
     2545          ALLOCATE ( surf_lsm_h(l)%gamma_w_sat(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns) )
     2546          ALLOCATE ( surf_lsm_h(l)%lambda_h(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)    )
     2547          ALLOCATE ( surf_lsm_h(l)%lambda_h_def(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns))
     2548          ALLOCATE ( surf_lsm_h(l)%l_vg(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)        )
     2549          ALLOCATE ( surf_lsm_h(l)%m_fc(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)        )
     2550          ALLOCATE ( surf_lsm_h(l)%m_res(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)       )
     2551          ALLOCATE ( surf_lsm_h(l)%m_sat(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)       )
     2552          ALLOCATE ( surf_lsm_h(l)%m_wilt(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)      )
     2553          ALLOCATE ( surf_lsm_h(l)%n_vg(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)        )
     2554          ALLOCATE ( surf_lsm_h(l)%rho_c_total(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns) )
     2555          ALLOCATE ( surf_lsm_h(l)%rho_c_total_def(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns) )
     2556          ALLOCATE ( surf_lsm_h(l)%root_fr(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)     )
     2557
     2558          surf_lsm_h(l)%lambda_h     = 0.0_wp
     2559!
     2560!--       If required, allocate humidity-related variables for the soil model
     2561          IF ( humidity )  THEN
     2562             ALLOCATE ( surf_lsm_h(l)%lambda_w(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns) )
     2563             ALLOCATE ( surf_lsm_h(l)%gamma_w(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)  )
     2564
     2565             surf_lsm_h(l)%lambda_w = 0.0_wp
     2566          ENDIF
     2567       ENDDO
    25212568!
    25222569!--    For vertical surfaces
     
    25512598!--    fraction.
    25522599!--    Set default values at each surface element.
    2553        ALLOCATE ( surf_lsm_h%albedo_type(1:surf_lsm_h%ns,0:2) )
    2554        ALLOCATE ( surf_lsm_h%emissivity(1:surf_lsm_h%ns,0:2) )
    2555 !
    2556 !--    Initialize albedo type according to its default type, in order to set values
    2557 !--    independent on default albedo_type in radiation model.
    2558        surf_lsm_h%albedo_type(:,ind_veg_wall)  =                               &
    2559                              INT( vegetation_pars(ind_v_at,vegetation_type) )
    2560        surf_lsm_h%albedo_type(:,ind_wat_win)   =                               &
    2561                              INT( water_pars(ind_w_at,water_type)           )
    2562        surf_lsm_h%albedo_type(:,ind_pav_green) =                               &
    2563                              INT( pavement_pars(ind_p_at,pavement_type)     )
    2564        surf_lsm_h%emissivity  = emissivity
     2600       DO  l = 0, 1
     2601          ALLOCATE ( surf_lsm_h(l)%albedo_type(1:surf_lsm_h(l)%ns,0:2) )
     2602          ALLOCATE ( surf_lsm_h(l)%emissivity(1:surf_lsm_h(l)%ns,0:2) )
     2603!
     2604!--       Initialize albedo type according to its default type, in order to set values
     2605!--       independent on default albedo_type in radiation model.
     2606          surf_lsm_h(l)%albedo_type(:,ind_veg_wall)  =                               &
     2607                                INT( vegetation_pars(ind_v_at,vegetation_type) )
     2608          surf_lsm_h(l)%albedo_type(:,ind_wat_win)   =                               &
     2609                                INT( water_pars(ind_w_at,water_type)           )
     2610          surf_lsm_h(l)%albedo_type(:,ind_pav_green) =                               &
     2611                                INT( pavement_pars(ind_p_at,pavement_type)     )
     2612          surf_lsm_h(l)%emissivity  = emissivity
     2613       ENDDO
    25652614       DO  l = 0, 3
    25662615          ALLOCATE ( surf_lsm_v(l)%albedo_type(1:surf_lsm_v(l)%ns,0:2) )
     
    25802629!--    Allocate arrays for relative surface fraction.
    25812630!--    0 - vegetation fraction, 2 - water fraction, 1 - pavement fraction
    2582        ALLOCATE( surf_lsm_h%frac(1:surf_lsm_h%ns,0:2) )
    2583        surf_lsm_h%frac = 0.0_wp
     2631       DO  l = 0, 1
     2632          ALLOCATE( surf_lsm_h(l)%frac(1:surf_lsm_h(l)%ns,0:2) )
     2633          surf_lsm_h(l)%frac = 0.0_wp
     2634       ENDDO
    25842635       DO  l = 0, 3
    25852636          ALLOCATE( surf_lsm_v(l)%frac(1:surf_lsm_v(l)%ns,0:2) )
     
    25952646!--    Allocate arrays for the respective types and their names on the surface
    25962647!--    elements. This will be required to treat deposition of chemical species.
    2597        ALLOCATE( surf_lsm_h%pavement_type(1:surf_lsm_h%ns)   )
    2598        ALLOCATE( surf_lsm_h%vegetation_type(1:surf_lsm_h%ns) )
    2599        ALLOCATE( surf_lsm_h%water_type(1:surf_lsm_h%ns)      )
    2600 
    2601        surf_lsm_h%pavement_type   = 0
    2602        surf_lsm_h%vegetation_type = 0
    2603        surf_lsm_h%water_type      = 0
    2604 
    2605        ALLOCATE( surf_lsm_h%pavement_type_name(1:surf_lsm_h%ns)   )
    2606        ALLOCATE( surf_lsm_h%vegetation_type_name(1:surf_lsm_h%ns) )
    2607        ALLOCATE( surf_lsm_h%water_type_name(1:surf_lsm_h%ns)      )
    2608 
    2609        surf_lsm_h%pavement_type_name   = 'none'
    2610        surf_lsm_h%vegetation_type_name = 'none'
    2611        surf_lsm_h%water_type_name      = 'none'
    2612 
     2648       DO  l = 0, 1
     2649          ALLOCATE( surf_lsm_h(l)%pavement_type(1:surf_lsm_h(l)%ns)   )
     2650          ALLOCATE( surf_lsm_h(l)%vegetation_type(1:surf_lsm_h(l)%ns) )
     2651          ALLOCATE( surf_lsm_h(l)%water_type(1:surf_lsm_h(l)%ns)      )
     2652
     2653          surf_lsm_h(l)%pavement_type   = 0
     2654          surf_lsm_h(l)%vegetation_type = 0
     2655          surf_lsm_h(l)%water_type      = 0
     2656
     2657          ALLOCATE( surf_lsm_h(l)%pavement_type_name(1:surf_lsm_h(l)%ns)   )
     2658          ALLOCATE( surf_lsm_h(l)%vegetation_type_name(1:surf_lsm_h(l)%ns) )
     2659          ALLOCATE( surf_lsm_h(l)%water_type_name(1:surf_lsm_h(l)%ns)      )
     2660
     2661          surf_lsm_h(l)%pavement_type_name   = 'none'
     2662          surf_lsm_h(l)%vegetation_type_name = 'none'
     2663          surf_lsm_h(l)%water_type_name      = 'none'
     2664       ENDDO
    26132665       DO  l = 0, 3
    26142666          ALLOCATE( surf_lsm_v(l)%pavement_type(1:surf_lsm_v(l)%ns)   )
     
    26352687
    26362688          CASE ( 'vegetation' )
    2637 
    2638              surf_lsm_h%vegetation_surface = .TRUE.
    2639              surf_lsm_h%frac(:,ind_veg_wall) = 1.0_wp
     2689             DO  l = 0, 1
     2690                surf_lsm_h(l)%vegetation_surface = .TRUE.
     2691                surf_lsm_h(l)%frac(:,ind_veg_wall) = 1.0_wp
     2692             ENDDO
    26402693             DO  l = 0, 3
    26412694                surf_lsm_v(l)%vegetation_surface = .TRUE.
     
    26442697
    26452698          CASE ( 'water' )
    2646 
    2647              surf_lsm_h%water_surface = .TRUE.
    2648              surf_lsm_h%frac(:,ind_wat_win) = 1.0_wp
    2649 !
    2650 !--          Note, vertical water surface does not really make sense.
     2699!
     2700!--          Note, downward and vertical water surface does not really make sense.
     2701             DO  l = 0, 1
     2702               surf_lsm_h(l)%water_surface = .TRUE.
     2703               surf_lsm_h(l)%frac(:,ind_wat_win) = 1.0_wp
     2704             ENDDO
    26512705             DO  l = 0, 3
    26522706                surf_lsm_v(l)%water_surface   = .TRUE.
     
    26552709
    26562710          CASE ( 'pavement' )
    2657 
    2658              surf_lsm_h%pavement_surface = .TRUE.
    2659                 surf_lsm_h%frac(:,ind_pav_green) = 1.0_wp
     2711             DO  l = 0, 1
     2712                surf_lsm_h(l)%pavement_surface = .TRUE.
     2713                surf_lsm_h(l)%frac(:,ind_pav_green) = 1.0_wp
     2714             ENDDO
    26602715             DO  l = 0, 3
    26612716                surf_lsm_v(l)%pavement_surface   = .TRUE.
     
    26642719
    26652720          CASE ( 'netcdf' )
    2666 
    2667              DO  m = 1, surf_lsm_h%ns
    2668                 i = surf_lsm_h%i(m)
    2669                 j = surf_lsm_h%j(m)
    2670                 IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill )    &
    2671                    surf_lsm_h%vegetation_surface(m) = .TRUE.
    2672                 IF ( pavement_type_f%var(j,i)   /= pavement_type_f%fill )      &
    2673                    surf_lsm_h%pavement_surface(m) = .TRUE.
    2674                 IF ( water_type_f%var(j,i)      /= water_type_f%fill )         &
    2675                    surf_lsm_h%water_surface(m) = .TRUE.
    2676 !
    2677 !--             Check if at least one type is set.
    2678                 IF ( .NOT. surf_lsm_h%vegetation_surface(m)  .AND.             &
    2679                      .NOT. surf_lsm_h%pavement_surface(m)    .AND.             &
    2680                      .NOT. surf_lsm_h%water_surface(m) )  THEN
    2681                    WRITE( message_string, * ) 'Horizontal surface element ' // &
    2682                                        ' at i, j = ',  i, j,                   &
    2683                                        ' is neither a vegetation, ' //         &
    2684                                        'pavement, nor a water surface.'
    2685                    CALL message( 'land_surface_model_mod', 'PA0619',          &
    2686                                   2, 2, myid, 6, 0 )
    2687                 ENDIF
    2688 
     2721             DO  l = 0, 1
     2722                DO  m = 1, surf_lsm_h(l)%ns
     2723                   i = surf_lsm_h(l)%i(m)
     2724                   j = surf_lsm_h(l)%j(m)
     2725                   IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill )    &
     2726                      surf_lsm_h(l)%vegetation_surface(m) = .TRUE.
     2727                   IF ( pavement_type_f%var(j,i)   /= pavement_type_f%fill )      &
     2728                      surf_lsm_h(l)%pavement_surface(m) = .TRUE.
     2729                   IF ( water_type_f%var(j,i)      /= water_type_f%fill )         &
     2730                      surf_lsm_h(l)%water_surface(m) = .TRUE.
     2731   !
     2732   !--             Check if at least one type is set.
     2733                   IF ( .NOT. surf_lsm_h(l)%vegetation_surface(m)  .AND.             &
     2734                        .NOT. surf_lsm_h(l)%pavement_surface(m)    .AND.             &
     2735                        .NOT. surf_lsm_h(l)%water_surface(m) )  THEN
     2736                      WRITE( message_string, * ) 'Horizontal surface element ' // &
     2737                                          ' at i, j = ',  i, j,                   &
     2738                                          ' is neither a vegetation, ' //         &
     2739                                          'pavement, nor a water surface.'
     2740                      CALL message( 'land_surface_model_mod', 'PA0619',          &
     2741                                     2, 2, myid, 6, 0 )
     2742                   ENDIF
     2743                ENDDO
    26892744             ENDDO
    26902745!
     
    27902845!--    location (already checked).
    27912846       IF ( input_pids_static  .AND.  surface_fraction_f%from_file )  THEN
    2792           DO  m = 1, surf_lsm_h%ns
    2793              i = surf_lsm_h%i(m)
    2794              j = surf_lsm_h%j(m)
    2795 !
    2796 !--          0 - vegetation fraction, 1 - pavement fraction, 2 - water fraction
    2797              IF ( surface_fraction_f%frac(ind_veg_wall,j,i) /=                 &
    2798                   surface_fraction_f%fill )  THEN
    2799                 surf_lsm_h%frac(m,ind_veg_wall)  =                             &
    2800                                     surface_fraction_f%frac(ind_veg_wall,j,i)
    2801              ENDIF
    2802              IF ( surface_fraction_f%frac(ind_pav_green,j,i) /=                &
    2803                   surface_fraction_f%fill )  THEN
    2804                 surf_lsm_h%frac(m,ind_pav_green) =                             &
    2805                                     surface_fraction_f%frac(ind_pav_green,j,i)
    2806              ENDIF
    2807              IF ( surface_fraction_f%frac(ind_wat_win,j,i) /=                  &
    2808                   surface_fraction_f%fill )  THEN
    2809                 surf_lsm_h%frac(m,ind_wat_win)   =                             &
    2810                                     surface_fraction_f%frac(ind_wat_win,j,i)
    2811              ENDIF
    2812 !
    2813 !--          Check if sum of relative fractions is zero. This case, give an
    2814 !--          error message.
    2815              IF ( SUM ( surf_lsm_h%frac(m,:) ) == 0.0_wp )  THEN
    2816                 WRITE( message_string, * )                                     &
    2817                                  'surface fractions at grid point (j,i) = (',  &
    2818                                  j, i, ') are all zero.'
    2819                 CALL message( 'land_surface_model_mod', 'PA0688',              &
    2820                                2, 2, myid, 6, 0 )
    2821              ENDIF
    2822 !
    2823 !--          In case the sum of all surfaces is not 1, which may happen
    2824 !--          due to rounding errors or type conversions, normalize the
    2825 !--          fractions to one. Note, at the moment no tile approach is
    2826 !--          implemented, so that relative fractions are either 1 or zero.
    2827              IF ( SUM ( surf_lsm_h%frac(m,:) ) > 1.0_wp  .OR.                  &
    2828                   SUM ( surf_lsm_h%frac(m,:) ) < 1.0_wp  )  THEN
    2829                 surf_lsm_h%frac(m,:) = surf_lsm_h%frac(m,:) /                  &
    2830                                        SUM ( surf_lsm_h%frac(m,:) )
    2831 
    2832              ENDIF
    2833 
     2847          DO  l = 0, 1
     2848             DO  m = 1, surf_lsm_h(l)%ns
     2849                i = surf_lsm_h(l)%i(m)
     2850                j = surf_lsm_h(l)%j(m)
     2851   !
     2852   !--          0 - vegetation fraction, 1 - pavement fraction, 2 - water fraction
     2853                IF ( surface_fraction_f%frac(ind_veg_wall,j,i) /=                 &
     2854                     surface_fraction_f%fill )  THEN
     2855                   surf_lsm_h(l)%frac(m,ind_veg_wall)  =                             &
     2856                                       surface_fraction_f%frac(ind_veg_wall,j,i)
     2857                ENDIF
     2858                IF ( surface_fraction_f%frac(ind_pav_green,j,i) /=                &
     2859                     surface_fraction_f%fill )  THEN
     2860                   surf_lsm_h(l)%frac(m,ind_pav_green) =                             &
     2861                                       surface_fraction_f%frac(ind_pav_green,j,i)
     2862                ENDIF
     2863                IF ( surface_fraction_f%frac(ind_wat_win,j,i) /=                  &
     2864                     surface_fraction_f%fill )  THEN
     2865                   surf_lsm_h(l)%frac(m,ind_wat_win)   =                             &
     2866                                       surface_fraction_f%frac(ind_wat_win,j,i)
     2867                ENDIF
     2868   !
     2869   !--          Check if sum of relative fractions is zero. This case, give an
     2870   !--          error message.
     2871                IF ( SUM ( surf_lsm_h(l)%frac(m,:) ) == 0.0_wp )  THEN
     2872                   WRITE( message_string, * )                                     &
     2873                                    'surface fractions at grid point (j,i) = (',  &
     2874                                    j, i, ') are all zero.'
     2875                   CALL message( 'land_surface_model_mod', 'PA0688',              &
     2876                                  2, 2, myid, 6, 0 )
     2877                ENDIF
     2878   !
     2879   !--          In case the sum of all surfaces is not 1, which may happen
     2880   !--          due to rounding errors or type conversions, normalize the
     2881   !--          fractions to one. Note, at the moment no tile approach is
     2882   !--          implemented, so that relative fractions are either 1 or zero.
     2883                IF ( SUM ( surf_lsm_h(l)%frac(m,:) ) > 1.0_wp  .OR.                  &
     2884                     SUM ( surf_lsm_h(l)%frac(m,:) ) < 1.0_wp  )  THEN
     2885                   surf_lsm_h(l)%frac(m,:) = surf_lsm_h(l)%frac(m,:) /                  &
     2886                                          SUM ( surf_lsm_h(l)%frac(m,:) )
     2887
     2888                ENDIF
     2889             ENDDO
    28342890          ENDDO
    28352891          DO  l = 0, 3
     
    28812937       ELSEIF ( input_pids_static  .AND.  .NOT. surface_fraction_f%from_file ) &
    28822938       THEN
    2883           DO  m = 1, surf_lsm_h%ns
    2884              i = surf_lsm_h%i(m)
    2885              j = surf_lsm_h%j(m)
    2886 
    2887              IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill )       &
    2888                 surf_lsm_h%frac(m,ind_veg_wall)  = 1.0_wp
    2889              IF ( pavement_type_f%var(j,i)   /= pavement_type_f%fill   )       &
    2890                 surf_lsm_h%frac(m,ind_pav_green) = 1.0_wp
    2891              IF ( water_type_f%var(j,i)      /= water_type_f%fill      )       &
    2892                 surf_lsm_h%frac(m,ind_wat_win)   = 1.0_wp
     2939          DO  l = 0, 1
     2940             DO  m = 1, surf_lsm_h(l)%ns
     2941                i = surf_lsm_h(l)%i(m)
     2942                j = surf_lsm_h(l)%j(m)
     2943
     2944                IF ( vegetation_type_f%var(j,i) /= vegetation_type_f%fill )       &
     2945                   surf_lsm_h(l)%frac(m,ind_veg_wall)  = 1.0_wp
     2946                IF ( pavement_type_f%var(j,i)   /= pavement_type_f%fill   )       &
     2947                   surf_lsm_h(l)%frac(m,ind_pav_green) = 1.0_wp
     2948                IF ( water_type_f%var(j,i)      /= water_type_f%fill      )       &
     2949                   surf_lsm_h(l)%frac(m,ind_wat_win)   = 1.0_wp
     2950             ENDDO
    28932951          ENDDO
    28942952          DO  l = 0, 3
     
    29503008!--    Map values to the respective 2D/3D arrays
    29513009!--    Horizontal surfaces
    2952        surf_lsm_h%alpha_vg      = alpha_vangenuchten
    2953        surf_lsm_h%l_vg          = l_vangenuchten
    2954        surf_lsm_h%n_vg          = n_vangenuchten
    2955        surf_lsm_h%gamma_w_sat   = hydraulic_conductivity
    2956        surf_lsm_h%m_sat         = saturation_moisture
    2957        surf_lsm_h%m_fc          = field_capacity
    2958        surf_lsm_h%m_wilt        = wilting_point
    2959        surf_lsm_h%m_res         = residual_moisture
    2960        surf_lsm_h%r_soil_min    = min_soil_resistance
     3010       DO  l = 0, 1
     3011          surf_lsm_h(l)%alpha_vg      = alpha_vangenuchten
     3012          surf_lsm_h(l)%l_vg          = l_vangenuchten
     3013          surf_lsm_h(l)%n_vg          = n_vangenuchten
     3014          surf_lsm_h(l)%gamma_w_sat   = hydraulic_conductivity
     3015          surf_lsm_h(l)%m_sat         = saturation_moisture
     3016          surf_lsm_h(l)%m_fc          = field_capacity
     3017          surf_lsm_h(l)%m_wilt        = wilting_point
     3018          surf_lsm_h(l)%m_res         = residual_moisture
     3019          surf_lsm_h(l)%r_soil_min    = min_soil_resistance
     3020       ENDDO
    29613021!
    29623022!--    Vertical surfaces
     
    29843044!
    29853045!--          Horizontal surfaces
    2986              DO  m = 1, surf_lsm_h%ns
    2987                 i = surf_lsm_h%i(m)
    2988                 j = surf_lsm_h%j(m)
    2989 
    2990                 st = soil_type_f%var_2d(j,i)
    2991                 IF ( st /= soil_type_f%fill )  THEN
    2992                    surf_lsm_h%alpha_vg(:,m)    = soil_pars(0,st)
    2993                    surf_lsm_h%l_vg(:,m)        = soil_pars(1,st)
    2994                    surf_lsm_h%n_vg(:,m)        = soil_pars(2,st)
    2995                    surf_lsm_h%gamma_w_sat(:,m) = soil_pars(3,st)
    2996                    surf_lsm_h%m_sat(:,m)       = soil_pars(4,st)
    2997                    surf_lsm_h%m_fc(:,m)        = soil_pars(5,st)
    2998                    surf_lsm_h%m_wilt(:,m)      = soil_pars(6,st)
    2999                    surf_lsm_h%m_res(:,m)       = soil_pars(7,st)
    3000                 ENDIF
     3046             DO  l = 0, 1
     3047                DO  m = 1, surf_lsm_h(l)%ns
     3048                   i = surf_lsm_h(l)%i(m)
     3049                   j = surf_lsm_h(l)%j(m)
     3050
     3051                   st = soil_type_f%var_2d(j,i)
     3052                   IF ( st /= soil_type_f%fill )  THEN
     3053                      surf_lsm_h(l)%alpha_vg(:,m)    = soil_pars(0,st)
     3054                      surf_lsm_h(l)%l_vg(:,m)        = soil_pars(1,st)
     3055                      surf_lsm_h(l)%n_vg(:,m)        = soil_pars(2,st)
     3056                      surf_lsm_h(l)%gamma_w_sat(:,m) = soil_pars(3,st)
     3057                      surf_lsm_h(l)%m_sat(:,m)       = soil_pars(4,st)
     3058                      surf_lsm_h(l)%m_fc(:,m)        = soil_pars(5,st)
     3059                      surf_lsm_h(l)%m_wilt(:,m)      = soil_pars(6,st)
     3060                      surf_lsm_h(l)%m_res(:,m)       = soil_pars(7,st)
     3061                   ENDIF
     3062                ENDDO
    30013063             ENDDO
    30023064!
     
    30283090!
    30293091!--          Horizontal surfaces
    3030              DO  m = 1, surf_lsm_h%ns
    3031                 i = surf_lsm_h%i(m)
    3032                 j = surf_lsm_h%j(m)
    3033 
    3034                 DO  k = nzb_soil, nzt_soil
    3035                    st = soil_type_f%var_3d(k,j,i)
    3036                    IF ( st /= soil_type_f%fill )  THEN
    3037                       surf_lsm_h%alpha_vg(k,m)    = soil_pars(0,st)
    3038                       surf_lsm_h%l_vg(k,m)        = soil_pars(1,st)
    3039                       surf_lsm_h%n_vg(k,m)        = soil_pars(2,st)
    3040                       surf_lsm_h%gamma_w_sat(k,m) = soil_pars(3,st)
    3041                       surf_lsm_h%m_sat(k,m)       = soil_pars(4,st)
    3042                       surf_lsm_h%m_fc(k,m)        = soil_pars(5,st)
    3043                       surf_lsm_h%m_wilt(k,m)      = soil_pars(6,st)
    3044                       surf_lsm_h%m_res(k,m)       = soil_pars(7,st)
    3045                    ENDIF
     3092             DO  l = 0, 1
     3093                DO  m = 1, surf_lsm_h(l)%ns
     3094                   i = surf_lsm_h(l)%i(m)
     3095                   j = surf_lsm_h(l)%j(m)
     3096
     3097                   DO  k = nzb_soil, nzt_soil
     3098                      st = soil_type_f%var_3d(k,j,i)
     3099                      IF ( st /= soil_type_f%fill )  THEN
     3100                         surf_lsm_h(l)%alpha_vg(k,m)    = soil_pars(0,st)
     3101                         surf_lsm_h(l)%l_vg(k,m)        = soil_pars(1,st)
     3102                         surf_lsm_h(l)%n_vg(k,m)        = soil_pars(2,st)
     3103                         surf_lsm_h(l)%gamma_w_sat(k,m) = soil_pars(3,st)
     3104                         surf_lsm_h(l)%m_sat(k,m)       = soil_pars(4,st)
     3105                         surf_lsm_h(l)%m_fc(k,m)        = soil_pars(5,st)
     3106                         surf_lsm_h(l)%m_wilt(k,m)      = soil_pars(6,st)
     3107                         surf_lsm_h(l)%m_res(k,m)       = soil_pars(7,st)
     3108                      ENDIF
     3109                   ENDDO
    30463110                ENDDO
    30473111             ENDDO
     
    30833147!
    30843148!--          Horizontal surfaces
    3085              DO  m = 1, surf_lsm_h%ns
    3086                 i = surf_lsm_h%i(m)
    3087                 j = surf_lsm_h%j(m)
    3088 
    3089                 IF ( soil_pars_f%pars_xy(0,j,i) /= soil_pars_f%fill )              &
    3090                    surf_lsm_h%alpha_vg(:,m)    = soil_pars_f%pars_xy(0,j,i)
    3091                 IF ( soil_pars_f%pars_xy(1,j,i) /= soil_pars_f%fill )              &
    3092                    surf_lsm_h%l_vg(:,m)        = soil_pars_f%pars_xy(1,j,i)
    3093                 IF ( soil_pars_f%pars_xy(2,j,i) /= soil_pars_f%fill )              &
    3094                    surf_lsm_h%n_vg(:,m)        = soil_pars_f%pars_xy(2,j,i)
    3095                 IF ( soil_pars_f%pars_xy(3,j,i) /= soil_pars_f%fill )              &
    3096                    surf_lsm_h%gamma_w_sat(:,m) = soil_pars_f%pars_xy(3,j,i)
    3097                 IF ( soil_pars_f%pars_xy(4,j,i) /= soil_pars_f%fill )              &
    3098                    surf_lsm_h%m_sat(:,m)       = soil_pars_f%pars_xy(4,j,i)
    3099                 IF ( soil_pars_f%pars_xy(5,j,i) /= soil_pars_f%fill )              &
    3100                    surf_lsm_h%m_fc(:,m)        = soil_pars_f%pars_xy(5,j,i)
    3101                 IF ( soil_pars_f%pars_xy(6,j,i) /= soil_pars_f%fill )              &
    3102                    surf_lsm_h%m_wilt(:,m)      = soil_pars_f%pars_xy(6,j,i)
    3103                 IF ( soil_pars_f%pars_xy(7,j,i) /= soil_pars_f%fill )              &
    3104                    surf_lsm_h%m_res(:,m)       = soil_pars_f%pars_xy(7,j,i)
    3105 
     3149             DO  l = 0, 1
     3150                DO  m = 1, surf_lsm_h(l)%ns
     3151                   i = surf_lsm_h(l)%i(m)
     3152                   j = surf_lsm_h(l)%j(m)
     3153
     3154                   IF ( soil_pars_f%pars_xy(0,j,i) /= soil_pars_f%fill )              &
     3155                      surf_lsm_h(l)%alpha_vg(:,m)    = soil_pars_f%pars_xy(0,j,i)
     3156                   IF ( soil_pars_f%pars_xy(1,j,i) /= soil_pars_f%fill )              &
     3157                      surf_lsm_h(l)%l_vg(:,m)        = soil_pars_f%pars_xy(1,j,i)
     3158                   IF ( soil_pars_f%pars_xy(2,j,i) /= soil_pars_f%fill )              &
     3159                      surf_lsm_h(l)%n_vg(:,m)        = soil_pars_f%pars_xy(2,j,i)
     3160                   IF ( soil_pars_f%pars_xy(3,j,i) /= soil_pars_f%fill )              &
     3161                      surf_lsm_h(l)%gamma_w_sat(:,m) = soil_pars_f%pars_xy(3,j,i)
     3162                   IF ( soil_pars_f%pars_xy(4,j,i) /= soil_pars_f%fill )              &
     3163                      surf_lsm_h(l)%m_sat(:,m)       = soil_pars_f%pars_xy(4,j,i)
     3164                   IF ( soil_pars_f%pars_xy(5,j,i) /= soil_pars_f%fill )              &
     3165                      surf_lsm_h(l)%m_fc(:,m)        = soil_pars_f%pars_xy(5,j,i)
     3166                   IF ( soil_pars_f%pars_xy(6,j,i) /= soil_pars_f%fill )              &
     3167                      surf_lsm_h(l)%m_wilt(:,m)      = soil_pars_f%pars_xy(6,j,i)
     3168                   IF ( soil_pars_f%pars_xy(7,j,i) /= soil_pars_f%fill )              &
     3169                      surf_lsm_h(l)%m_res(:,m)       = soil_pars_f%pars_xy(7,j,i)
     3170                ENDDO
    31063171             ENDDO
    31073172!
     
    31393204!
    31403205!--          Horizontal surfaces
    3141              DO  m = 1, surf_lsm_h%ns
    3142                 i = surf_lsm_h%i(m)
    3143                 j = surf_lsm_h%j(m)
    3144 
    3145                 DO  k = nzb_soil, nzt_soil
    3146                    IF ( soil_pars_f%pars_xyz(0,k,j,i) /= soil_pars_f%fill )        &
    3147                       surf_lsm_h%alpha_vg(k,m)    = soil_pars_f%pars_xyz(0,k,j,i)
    3148                    IF ( soil_pars_f%pars_xyz(1,k,j,i) /= soil_pars_f%fill )        &
    3149                       surf_lsm_h%l_vg(k,m)        = soil_pars_f%pars_xyz(1,k,j,i)
    3150                    IF ( soil_pars_f%pars_xyz(2,k,j,i) /= soil_pars_f%fill )        &
    3151                       surf_lsm_h%n_vg(k,m)        = soil_pars_f%pars_xyz(2,k,j,i)
    3152                    IF ( soil_pars_f%pars_xyz(3,k,j,i) /= soil_pars_f%fill )        &
    3153                       surf_lsm_h%gamma_w_sat(k,m) = soil_pars_f%pars_xyz(3,k,j,i)
    3154                    IF ( soil_pars_f%pars_xyz(4,k,j,i) /= soil_pars_f%fill )        &
    3155                       surf_lsm_h%m_sat(k,m)       = soil_pars_f%pars_xyz(4,k,j,i)
    3156                    IF ( soil_pars_f%pars_xyz(5,k,j,i) /= soil_pars_f%fill )        &
    3157                       surf_lsm_h%m_fc(k,m)        = soil_pars_f%pars_xyz(5,k,j,i)
    3158                    IF ( soil_pars_f%pars_xyz(6,k,j,i) /= soil_pars_f%fill )        &
    3159                       surf_lsm_h%m_wilt(k,m)      = soil_pars_f%pars_xyz(6,k,j,i)
    3160                    IF ( soil_pars_f%pars_xyz(7,k,j,i) /= soil_pars_f%fill )        &
    3161                       surf_lsm_h%m_res(k,m)       = soil_pars_f%pars_xyz(7,k,j,i)
     3206             DO  l = 0, 1
     3207                DO  m = 1, surf_lsm_h(l)%ns
     3208                   i = surf_lsm_h(l)%i(m)
     3209                   j = surf_lsm_h(l)%j(m)
     3210
     3211                   DO  k = nzb_soil, nzt_soil
     3212                      IF ( soil_pars_f%pars_xyz(0,k,j,i) /= soil_pars_f%fill )        &
     3213                         surf_lsm_h(l)%alpha_vg(k,m)    = soil_pars_f%pars_xyz(0,k,j,i)
     3214                      IF ( soil_pars_f%pars_xyz(1,k,j,i) /= soil_pars_f%fill )        &
     3215                         surf_lsm_h(l)%l_vg(k,m)        = soil_pars_f%pars_xyz(1,k,j,i)
     3216                      IF ( soil_pars_f%pars_xyz(2,k,j,i) /= soil_pars_f%fill )        &
     3217                         surf_lsm_h(l)%n_vg(k,m)        = soil_pars_f%pars_xyz(2,k,j,i)
     3218                      IF ( soil_pars_f%pars_xyz(3,k,j,i) /= soil_pars_f%fill )        &
     3219                         surf_lsm_h(l)%gamma_w_sat(k,m) = soil_pars_f%pars_xyz(3,k,j,i)
     3220                      IF ( soil_pars_f%pars_xyz(4,k,j,i) /= soil_pars_f%fill )        &
     3221                         surf_lsm_h(l)%m_sat(k,m)       = soil_pars_f%pars_xyz(4,k,j,i)
     3222                      IF ( soil_pars_f%pars_xyz(5,k,j,i) /= soil_pars_f%fill )        &
     3223                         surf_lsm_h(l)%m_fc(k,m)        = soil_pars_f%pars_xyz(5,k,j,i)
     3224                      IF ( soil_pars_f%pars_xyz(6,k,j,i) /= soil_pars_f%fill )        &
     3225                         surf_lsm_h(l)%m_wilt(k,m)      = soil_pars_f%pars_xyz(6,k,j,i)
     3226                      IF ( soil_pars_f%pars_xyz(7,k,j,i) /= soil_pars_f%fill )        &
     3227                         surf_lsm_h(l)%m_res(k,m)       = soil_pars_f%pars_xyz(7,k,j,i)
     3228                   ENDDO
    31623229                ENDDO
    3163 
    31643230             ENDDO
    31653231!
     
    32573323!
    32583324!--    Map values onto horizontal elemements
    3259        DO  m = 1, surf_lsm_h%ns
    3260           IF ( surf_lsm_h%vegetation_surface(m) )  THEN
    3261              surf_lsm_h%r_canopy_min(m)     = min_canopy_resistance
    3262              surf_lsm_h%lai(m)              = leaf_area_index
    3263              surf_lsm_h%c_veg(m)            = vegetation_coverage
    3264              surf_lsm_h%g_d(m)              = canopy_resistance_coefficient
    3265              surf_lsm_h%z0(m)               = z0_vegetation
    3266              surf_lsm_h%z0h(m)              = z0h_vegetation
    3267              surf_lsm_h%z0q(m)              = z0q_vegetation
    3268              surf_lsm_h%lambda_surface_s(m) = lambda_surface_stable
    3269              surf_lsm_h%lambda_surface_u(m) = lambda_surface_unstable
    3270              surf_lsm_h%f_sw_in(m)          = f_shortwave_incoming
    3271              surf_lsm_h%c_surface(m)        = c_surface
    3272              surf_lsm_h%albedo_type(m,ind_veg_wall) = albedo_type
    3273              surf_lsm_h%emissivity(m,ind_veg_wall)  = emissivity
    3274 
    3275              surf_lsm_h%vegetation_type(m)      = vegetation_type
    3276              surf_lsm_h%vegetation_type_name(m) = vegetation_type_name(vegetation_type)
    3277           ELSE
    3278              surf_lsm_h%lai(m)   = 0.0_wp
    3279              surf_lsm_h%c_veg(m) = 0.0_wp
    3280              surf_lsm_h%g_d(m)   = 0.0_wp
    3281           ENDIF
    3282 
     3325       DO  l = 0, 1
     3326          DO  m = 1, surf_lsm_h(l)%ns
     3327             IF ( surf_lsm_h(l)%vegetation_surface(m) )  THEN
     3328                surf_lsm_h(l)%r_canopy_min(m)     = min_canopy_resistance
     3329                surf_lsm_h(l)%lai(m)              = leaf_area_index
     3330                surf_lsm_h(l)%c_veg(m)            = vegetation_coverage
     3331                surf_lsm_h(l)%g_d(m)              = canopy_resistance_coefficient
     3332                surf_lsm_h(l)%z0(m)               = z0_vegetation
     3333                surf_lsm_h(l)%z0h(m)              = z0h_vegetation
     3334                surf_lsm_h(l)%z0q(m)              = z0q_vegetation
     3335                surf_lsm_h(l)%lambda_surface_s(m) = lambda_surface_stable
     3336                surf_lsm_h(l)%lambda_surface_u(m) = lambda_surface_unstable
     3337                surf_lsm_h(l)%f_sw_in(m)          = f_shortwave_incoming
     3338                surf_lsm_h(l)%c_surface(m)        = c_surface
     3339                surf_lsm_h(l)%albedo_type(m,ind_veg_wall) = albedo_type
     3340                surf_lsm_h(l)%emissivity(m,ind_veg_wall)  = emissivity
     3341
     3342                surf_lsm_h(l)%vegetation_type(m)      = vegetation_type
     3343                surf_lsm_h(l)%vegetation_type_name(m) = vegetation_type_name(vegetation_type)
     3344             ELSE
     3345                surf_lsm_h(l)%lai(m)   = 0.0_wp
     3346                surf_lsm_h(l)%c_veg(m) = 0.0_wp
     3347                surf_lsm_h(l)%g_d(m)   = 0.0_wp
     3348             ENDIF
     3349          ENDDO
    32833350       ENDDO
    32843351!
     
    33203387!
    33213388!--       Horizontal surfaces
    3322           DO  m = 1, surf_lsm_h%ns
    3323              i = surf_lsm_h%i(m)
    3324              j = surf_lsm_h%j(m)
    3325 
    3326              st = vegetation_type_f%var(j,i)
    3327              IF ( st /= vegetation_type_f%fill  .AND.  st /= 0 )  THEN
    3328                 surf_lsm_h%r_canopy_min(m)     = vegetation_pars(ind_v_rc_min,st)
    3329                 surf_lsm_h%lai(m)              = vegetation_pars(ind_v_rc_lai,st)
    3330                 surf_lsm_h%c_veg(m)            = vegetation_pars(ind_v_c_veg,st)
    3331                 surf_lsm_h%g_d(m)              = vegetation_pars(ind_v_gd,st)
    3332                 surf_lsm_h%z0(m)               = vegetation_pars(ind_v_z0,st)
    3333                 surf_lsm_h%z0h(m)              = vegetation_pars(ind_v_z0qh,st)
    3334                 surf_lsm_h%z0q(m)              = vegetation_pars(ind_v_z0qh,st)
    3335                 surf_lsm_h%lambda_surface_s(m) = vegetation_pars(ind_v_lambda_s,st)
    3336                 surf_lsm_h%lambda_surface_u(m) = vegetation_pars(ind_v_lambda_u,st)
    3337                 surf_lsm_h%f_sw_in(m)          = vegetation_pars(ind_v_f_sw_in,st)
    3338                 surf_lsm_h%c_surface(m)        = vegetation_pars(ind_v_c_surf,st)
    3339                 surf_lsm_h%albedo_type(m,ind_veg_wall) = INT( vegetation_pars(ind_v_at,st) )
    3340                 surf_lsm_h%emissivity(m,ind_veg_wall)  = vegetation_pars(ind_v_emis,st)
    3341 
    3342                 surf_lsm_h%vegetation_type(m)      = st
    3343                 surf_lsm_h%vegetation_type_name(m) = vegetation_type_name(st)
    3344              ENDIF
     3389          DO  l = 0, 1
     3390             DO  m = 1, surf_lsm_h(l)%ns
     3391                i = surf_lsm_h(l)%i(m)
     3392                j = surf_lsm_h(l)%j(m)
     3393
     3394                st = vegetation_type_f%var(j,i)
     3395                IF ( st /= vegetation_type_f%fill  .AND.  st /= 0 )  THEN
     3396                   surf_lsm_h(l)%r_canopy_min(m)     = vegetation_pars(ind_v_rc_min,st)
     3397                   surf_lsm_h(l)%lai(m)              = vegetation_pars(ind_v_rc_lai,st)
     3398                   surf_lsm_h(l)%c_veg(m)            = vegetation_pars(ind_v_c_veg,st)
     3399                   surf_lsm_h(l)%g_d(m)              = vegetation_pars(ind_v_gd,st)
     3400                   surf_lsm_h(l)%z0(m)               = vegetation_pars(ind_v_z0,st)
     3401                   surf_lsm_h(l)%z0h(m)              = vegetation_pars(ind_v_z0qh,st)
     3402                   surf_lsm_h(l)%z0q(m)              = vegetation_pars(ind_v_z0qh,st)
     3403                   surf_lsm_h(l)%lambda_surface_s(m) = vegetation_pars(ind_v_lambda_s,st)
     3404                   surf_lsm_h(l)%lambda_surface_u(m) = vegetation_pars(ind_v_lambda_u,st)
     3405                   surf_lsm_h(l)%f_sw_in(m)          = vegetation_pars(ind_v_f_sw_in,st)
     3406                   surf_lsm_h(l)%c_surface(m)        = vegetation_pars(ind_v_c_surf,st)
     3407                   surf_lsm_h(l)%albedo_type(m,ind_veg_wall) = INT( vegetation_pars(ind_v_at,st) )
     3408                   surf_lsm_h(l)%emissivity(m,ind_veg_wall)  = vegetation_pars(ind_v_emis,st)
     3409
     3410                   surf_lsm_h(l)%vegetation_type(m)      = st
     3411                   surf_lsm_h(l)%vegetation_type_name(m) = vegetation_type_name(st)
     3412                ENDIF
     3413             ENDDO
    33453414          ENDDO
    33463415!
     
    33813450!
    33823451!--       Horizontal surfaces
    3383           DO  m = 1, surf_lsm_h%ns
    3384 
    3385              i = surf_lsm_h%i(m)
    3386              j = surf_lsm_h%j(m)
    3387 !
    3388 !--          If surface element is not a vegetation surface and any value in
    3389 !--          vegetation_pars is given, neglect this information and give an
    3390 !--          informative message that this value will not be used.
    3391              IF ( .NOT. surf_lsm_h%vegetation_surface(m)  .AND.                &
    3392                    ANY( vegetation_pars_f%pars_xy(:,j,i) /=                    &
    3393                    vegetation_pars_f%fill ) )  THEN
    3394                 WRITE( message_string, * )                                     &
    3395                                  'surface element at grid point (j,i) = (',    &
    3396                                  j, i, ') is not a vegetation surface, ',      &
    3397                                  'so that information given in ',              &
    3398                                  'vegetation_pars at this point is neglected.'
    3399                 CALL message( 'land_surface_model_mod', 'PA0436', 0, 0, myid, 6, 0 )
    3400              ELSE
    3401 
    3402                 IF ( vegetation_pars_f%pars_xy(ind_v_rc_min,j,i) /=            &
    3403                      vegetation_pars_f%fill )                                  &
    3404                    surf_lsm_h%r_canopy_min(m)  =                               &
    3405                                    vegetation_pars_f%pars_xy(ind_v_rc_min,j,i)
    3406                 IF ( vegetation_pars_f%pars_xy(ind_v_rc_lai,j,i) /=            &
    3407                      vegetation_pars_f%fill )                                  &
    3408                    surf_lsm_h%lai(m)           =                               &
    3409                                    vegetation_pars_f%pars_xy(ind_v_rc_lai,j,i)
    3410                 IF ( vegetation_pars_f%pars_xy(ind_v_c_veg,j,i) /=             &
    3411                      vegetation_pars_f%fill )                                  &
    3412                    surf_lsm_h%c_veg(m)         =                               &
    3413                                    vegetation_pars_f%pars_xy(ind_v_c_veg,j,i)
    3414                 IF ( vegetation_pars_f%pars_xy(ind_v_gd,j,i) /=                &
    3415                      vegetation_pars_f%fill )                                  &
    3416                    surf_lsm_h%g_d(m)           =                               &
    3417                                    vegetation_pars_f%pars_xy(ind_v_gd,j,i)
    3418                 IF ( vegetation_pars_f%pars_xy(ind_v_z0,j,i) /=                &
    3419                      vegetation_pars_f%fill )                                  &
    3420                    surf_lsm_h%z0(m)            =                               &
    3421                                    vegetation_pars_f%pars_xy(ind_v_z0,j,i)
    3422                 IF ( vegetation_pars_f%pars_xy(ind_v_z0qh,j,i) /=              &
    3423                      vegetation_pars_f%fill )  THEN
    3424                    surf_lsm_h%z0h(m)           =                               &
    3425                                    vegetation_pars_f%pars_xy(ind_v_z0qh,j,i)
    3426                    surf_lsm_h%z0q(m)           =                               &
    3427                                    vegetation_pars_f%pars_xy(ind_v_z0qh,j,i)
     3452          DO  l = 0, 1
     3453             DO  m = 1, surf_lsm_h(l)%ns
     3454                i = surf_lsm_h(l)%i(m)
     3455                j = surf_lsm_h(l)%j(m)
     3456!
     3457!--             If surface element is not a vegetation surface and any value in
     3458!--             vegetation_pars is given, neglect this information and give an
     3459!--             informative message that this value will not be used.
     3460                IF ( .NOT. surf_lsm_h(l)%vegetation_surface(m)  .AND.                &
     3461                      ANY( vegetation_pars_f%pars_xy(:,j,i) /=                    &
     3462                      vegetation_pars_f%fill ) )  THEN
     3463                   WRITE( message_string, * )                                     &
     3464                                    'surface element at grid point (j,i) = (',    &
     3465                                    j, i, ') is not a vegetation surface, ',      &
     3466                                    'so that information given in ',              &
     3467                                    'vegetation_pars at this point is neglected.'
     3468                   CALL message( 'land_surface_model_mod', 'PA0436', 0, 0, myid, 6, 0 )
     3469                ELSE
     3470
     3471                   IF ( vegetation_pars_f%pars_xy(ind_v_rc_min,j,i) /=            &
     3472                        vegetation_pars_f%fill )                                  &
     3473                      surf_lsm_h(l)%r_canopy_min(m)  =                               &
     3474                                      vegetation_pars_f%pars_xy(ind_v_rc_min,j,i)
     3475                   IF ( vegetation_pars_f%pars_xy(ind_v_rc_lai,j,i) /=            &
     3476                        vegetation_pars_f%fill )                                  &
     3477                      surf_lsm_h(l)%lai(m)           =                               &
     3478                                      vegetation_pars_f%pars_xy(ind_v_rc_lai,j,i)
     3479                   IF ( vegetation_pars_f%pars_xy(ind_v_c_veg,j,i) /=             &
     3480                        vegetation_pars_f%fill )                                  &
     3481                      surf_lsm_h(l)%c_veg(m)         =                               &
     3482                                      vegetation_pars_f%pars_xy(ind_v_c_veg,j,i)
     3483                   IF ( vegetation_pars_f%pars_xy(ind_v_gd,j,i) /=                &
     3484                        vegetation_pars_f%fill )                                  &
     3485                      surf_lsm_h(l)%g_d(m)           =                               &
     3486                                      vegetation_pars_f%pars_xy(ind_v_gd,j,i)
     3487                   IF ( vegetation_pars_f%pars_xy(ind_v_z0,j,i) /=                &
     3488                        vegetation_pars_f%fill )                                  &
     3489                      surf_lsm_h(l)%z0(m)            =                               &
     3490                                      vegetation_pars_f%pars_xy(ind_v_z0,j,i)
     3491                   IF ( vegetation_pars_f%pars_xy(ind_v_z0qh,j,i) /=              &
     3492                        vegetation_pars_f%fill )  THEN
     3493                      surf_lsm_h(l)%z0h(m)           =                               &
     3494                                      vegetation_pars_f%pars_xy(ind_v_z0qh,j,i)
     3495                      surf_lsm_h(l)%z0q(m)           =                               &
     3496                                      vegetation_pars_f%pars_xy(ind_v_z0qh,j,i)
     3497                   ENDIF
     3498                   IF ( vegetation_pars_f%pars_xy(ind_v_lambda_s,j,i) /=          &
     3499                        vegetation_pars_f%fill )                                  &
     3500                      surf_lsm_h(l)%lambda_surface_s(m) =                            &
     3501                                      vegetation_pars_f%pars_xy(ind_v_lambda_s,j,i)
     3502                   IF ( vegetation_pars_f%pars_xy(ind_v_lambda_u,j,i) /=          &
     3503                        vegetation_pars_f%fill )                                  &
     3504                      surf_lsm_h(l)%lambda_surface_u(m) =                            &
     3505                                      vegetation_pars_f%pars_xy(ind_v_lambda_u,j,i)
     3506                   IF ( vegetation_pars_f%pars_xy(ind_v_f_sw_in,j,i) /=           &
     3507                        vegetation_pars_f%fill )                                  &
     3508                      surf_lsm_h(l)%f_sw_in(m)          =                            &
     3509                                      vegetation_pars_f%pars_xy(ind_v_f_sw_in,j,i)
     3510                   IF ( vegetation_pars_f%pars_xy(ind_v_c_surf,j,i) /=            &
     3511                        vegetation_pars_f%fill )                                  &
     3512                      surf_lsm_h(l)%c_surface(m)        =                            &
     3513                                      vegetation_pars_f%pars_xy(ind_v_c_surf,j,i)
     3514                   IF ( vegetation_pars_f%pars_xy(ind_v_at,j,i) /=                &
     3515                        vegetation_pars_f%fill )                                  &
     3516                      surf_lsm_h(l)%albedo_type(m,ind_veg_wall) =                    &
     3517                                      INT( vegetation_pars_f%pars_xy(ind_v_at,j,i) )
     3518                   IF ( vegetation_pars_f%pars_xy(ind_v_emis,j,i) /=              &
     3519                        vegetation_pars_f%fill )                                  &
     3520                      surf_lsm_h(l)%emissivity(m,ind_veg_wall)  =                    &
     3521                                      vegetation_pars_f%pars_xy(ind_v_emis,j,i)
    34283522                ENDIF
    3429                 IF ( vegetation_pars_f%pars_xy(ind_v_lambda_s,j,i) /=          &
    3430                      vegetation_pars_f%fill )                                  &
    3431                    surf_lsm_h%lambda_surface_s(m) =                            &
    3432                                    vegetation_pars_f%pars_xy(ind_v_lambda_s,j,i)
    3433                 IF ( vegetation_pars_f%pars_xy(ind_v_lambda_u,j,i) /=          &
    3434                      vegetation_pars_f%fill )                                  &
    3435                    surf_lsm_h%lambda_surface_u(m) =                            &
    3436                                    vegetation_pars_f%pars_xy(ind_v_lambda_u,j,i)
    3437                 IF ( vegetation_pars_f%pars_xy(ind_v_f_sw_in,j,i) /=           &
    3438                      vegetation_pars_f%fill )                                  &
    3439                    surf_lsm_h%f_sw_in(m)          =                            &
    3440                                    vegetation_pars_f%pars_xy(ind_v_f_sw_in,j,i)
    3441                 IF ( vegetation_pars_f%pars_xy(ind_v_c_surf,j,i) /=            &
    3442                      vegetation_pars_f%fill )                                  &
    3443                    surf_lsm_h%c_surface(m)        =                            &
    3444                                    vegetation_pars_f%pars_xy(ind_v_c_surf,j,i)
    3445                 IF ( vegetation_pars_f%pars_xy(ind_v_at,j,i) /=                &
    3446                      vegetation_pars_f%fill )                                  &
    3447                    surf_lsm_h%albedo_type(m,ind_veg_wall) =                    &
    3448                                    INT( vegetation_pars_f%pars_xy(ind_v_at,j,i) )
    3449                 IF ( vegetation_pars_f%pars_xy(ind_v_emis,j,i) /=              &
    3450                      vegetation_pars_f%fill )                                  &
    3451                    surf_lsm_h%emissivity(m,ind_veg_wall)  =                    &
    3452                                    vegetation_pars_f%pars_xy(ind_v_emis,j,i)
    3453              ENDIF
     3523             ENDDO
    34543524          ENDDO
    34553525!
     
    35653635!
    35663636!--    Map values onto horizontal elemements
    3567        DO  m = 1, surf_lsm_h%ns
    3568           IF ( surf_lsm_h%water_surface(m) )  THEN
    3569              IF ( TRIM( initializing_actions ) /= 'read_restart_data' )        &
    3570                 t_soil_h%var_2d(:,m)        = water_temperature
    3571              surf_lsm_h%z0(m)               = z0_water
    3572              surf_lsm_h%z0h(m)              = z0h_water
    3573              surf_lsm_h%z0q(m)              = z0q_water
    3574              surf_lsm_h%lambda_surface_s(m) = 1.0E10_wp
    3575              surf_lsm_h%lambda_surface_u(m) = 1.0E10_wp
    3576              surf_lsm_h%c_surface(m)        = 0.0_wp
    3577              surf_lsm_h%albedo_type(m,ind_wat_win) = albedo_type
    3578              surf_lsm_h%emissivity(m,ind_wat_win)  = emissivity
    3579 
    3580              surf_lsm_h%water_type(m)      = water_type
    3581              surf_lsm_h%water_type_name(m) = water_type_name(water_type)
    3582           ENDIF
     3637       DO  l = 0, 1
     3638          DO  m = 1, surf_lsm_h(l)%ns
     3639             IF ( surf_lsm_h(l)%water_surface(m) )  THEN
     3640                IF ( TRIM( initializing_actions ) /= 'read_restart_data' )        &
     3641                   t_soil_h(l)%var_2d(:,m)        = water_temperature
     3642                surf_lsm_h(l)%z0(m)               = z0_water
     3643                surf_lsm_h(l)%z0h(m)              = z0h_water
     3644                surf_lsm_h(l)%z0q(m)              = z0q_water
     3645                surf_lsm_h(l)%lambda_surface_s(m) = 1.0E10_wp
     3646                surf_lsm_h(l)%lambda_surface_u(m) = 1.0E10_wp
     3647                surf_lsm_h(l)%c_surface(m)        = 0.0_wp
     3648                surf_lsm_h(l)%albedo_type(m,ind_wat_win) = albedo_type
     3649                surf_lsm_h(l)%emissivity(m,ind_wat_win)  = emissivity
     3650
     3651                surf_lsm_h(l)%water_type(m)      = water_type
     3652                surf_lsm_h(l)%water_type_name(m) = water_type_name(water_type)
     3653             ENDIF
     3654          ENDDO
    35833655       ENDDO
    35843656!
     
    35893661             IF ( surf_lsm_v(l)%water_surface(m) )  THEN
    35903662                IF ( TRIM( initializing_actions ) /= 'read_restart_data' )     &
    3591                    t_soil_v(l)%var_2d(:,m)           = water_temperature
     3663                   t_soil_v(l)%var_2d(:,m)        = water_temperature
    35923664                surf_lsm_v(l)%z0(m)               = z0_water
    35933665                surf_lsm_v(l)%z0h(m)              = z0h_water
     
    36153687!
    36163688!--       Horizontal surfaces
    3617           DO  m = 1, surf_lsm_h%ns
    3618              i = surf_lsm_h%i(m)
    3619              j = surf_lsm_h%j(m)
    3620 
    3621              st = water_type_f%var(j,i)
    3622              IF ( st /= water_type_f%fill  .AND.  st /= 0 )  THEN
    3623                 IF ( TRIM( initializing_actions ) /= 'read_restart_data' )     &
    3624                    t_soil_h%var_2d(:,m) = water_pars(ind_w_temp,st)
    3625                 surf_lsm_h%z0(m)     = water_pars(ind_w_z0,st)
    3626                 surf_lsm_h%z0h(m)    = water_pars(ind_w_z0h,st)
    3627                 surf_lsm_h%z0q(m)    = water_pars(ind_w_z0h,st)
    3628                 surf_lsm_h%lambda_surface_s(m) = water_pars(ind_w_lambda_s,st)
    3629                 surf_lsm_h%lambda_surface_u(m) = water_pars(ind_w_lambda_u,st)
    3630                 surf_lsm_h%c_surface(m)        = 0.0_wp
    3631                 surf_lsm_h%albedo_type(m,ind_wat_win) = INT( water_pars(ind_w_at,st) )
    3632                 surf_lsm_h%emissivity(m,ind_wat_win)  = water_pars(ind_w_emis,st)
    3633 
    3634                 surf_lsm_h%water_type(m)      = st
    3635                 surf_lsm_h%water_type_name(m) = water_type_name(st)
    3636              ENDIF
     3689          DO  l = 0, 1
     3690             DO  m = 1, surf_lsm_h(l)%ns
     3691                i = surf_lsm_h(l)%i(m)
     3692                j = surf_lsm_h(l)%j(m)
     3693
     3694                st = water_type_f%var(j,i)
     3695                IF ( st /= water_type_f%fill  .AND.  st /= 0 )  THEN
     3696                   IF ( TRIM( initializing_actions ) /= 'read_restart_data' )     &
     3697                      t_soil_h(l)%var_2d(:,m) = water_pars(ind_w_temp,st)
     3698                   surf_lsm_h(l)%z0(m)     = water_pars(ind_w_z0,st)
     3699                   surf_lsm_h(l)%z0h(m)    = water_pars(ind_w_z0h,st)
     3700                   surf_lsm_h(l)%z0q(m)    = water_pars(ind_w_z0h,st)
     3701                   surf_lsm_h(l)%lambda_surface_s(m) = water_pars(ind_w_lambda_s,st)
     3702                   surf_lsm_h(l)%lambda_surface_u(m) = water_pars(ind_w_lambda_u,st)
     3703                   surf_lsm_h(l)%c_surface(m)        = 0.0_wp
     3704                   surf_lsm_h(l)%albedo_type(m,ind_wat_win) = INT( water_pars(ind_w_at,st) )
     3705                   surf_lsm_h(l)%emissivity(m,ind_wat_win)  = water_pars(ind_w_emis,st)
     3706
     3707                   surf_lsm_h(l)%water_type(m)      = st
     3708                   surf_lsm_h(l)%water_type_name(m) = water_type_name(st)
     3709                ENDIF
     3710             ENDDO
    36373711          ENDDO
    36383712!
     
    36753749!
    36763750!--       Horizontal surfaces
    3677           DO  m = 1, surf_lsm_h%ns
    3678              i = surf_lsm_h%i(m)
    3679              j = surf_lsm_h%j(m)
    3680 !
    3681 !--          If surface element is not a water surface and any value in
    3682 !--          water_pars is given, neglect this information and give an
    3683 !--          informative message that this value will not be used.
    3684              IF ( .NOT. surf_lsm_h%water_surface(m)  .AND.                     &
    3685                    ANY( water_pars_f%pars_xy(:,j,i) /= water_pars_f%fill ) )  THEN
    3686                 WRITE( message_string, * )                                     &
    3687                               'surface element at grid point (j,i) = (',       &
    3688                               j, i, ') is not a water surface, ',              &
    3689                               'so that information given in ',                 &
    3690                               'water_pars at this point is neglected.'
    3691                 CALL message( 'land_surface_model_mod', 'PA0645', 0, 0, myid, 6, 0 )
    3692              ELSE
    3693                 IF ( water_pars_f%pars_xy(ind_w_temp,j,i) /=                   &
    3694                      water_pars_f%fill  .AND.                                  &
    3695                      TRIM( initializing_actions ) /= 'read_restart_data' )     &
    3696                       t_soil_h%var_2d(:,m) = water_pars_f%pars_xy(ind_w_temp,j,i)
    3697 
    3698                 IF ( water_pars_f%pars_xy(ind_w_z0,j,i) /= water_pars_f%fill ) &
    3699                    surf_lsm_h%z0(m)     = water_pars_f%pars_xy(ind_w_z0,j,i)
    3700 
    3701                 IF ( water_pars_f%pars_xy(ind_w_z0h,j,i) /= water_pars_f%fill )&
    3702                 THEN
    3703                    surf_lsm_h%z0h(m)    = water_pars_f%pars_xy(ind_w_z0h,j,i)
    3704                    surf_lsm_h%z0q(m)    = water_pars_f%pars_xy(ind_w_z0h,j,i)
     3751          DO  l = 0, 1
     3752             DO  m = 1, surf_lsm_h(l)%ns
     3753                i = surf_lsm_h(l)%i(m)
     3754                j = surf_lsm_h(l)%j(m)
     3755   !
     3756   !--          If surface element is not a water surface and any value in
     3757   !--          water_pars is given, neglect this information and give an
     3758   !--          informative message that this value will not be used.
     3759                IF ( .NOT. surf_lsm_h(l)%water_surface(m)  .AND.                     &
     3760                      ANY( water_pars_f%pars_xy(:,j,i) /= water_pars_f%fill ) )  THEN
     3761                   WRITE( message_string, * )                                     &
     3762                                 'surface element at grid point (j,i) = (',       &
     3763                                 j, i, ') is not a water surface, ',              &
     3764                                 'so that information given in ',                 &
     3765                                 'water_pars at this point is neglected.'
     3766                   CALL message( 'land_surface_model_mod', 'PA0645', 0, 0, myid, 6, 0 )
     3767                ELSE
     3768                   IF ( water_pars_f%pars_xy(ind_w_temp,j,i) /=                   &
     3769                        water_pars_f%fill  .AND.                                  &
     3770                        TRIM( initializing_actions ) /= 'read_restart_data' )     &
     3771                         t_soil_h(l)%var_2d(:,m) = water_pars_f%pars_xy(ind_w_temp,j,i)
     3772
     3773                   IF ( water_pars_f%pars_xy(ind_w_z0,j,i) /= water_pars_f%fill ) &
     3774                      surf_lsm_h(l)%z0(m)     = water_pars_f%pars_xy(ind_w_z0,j,i)
     3775
     3776                   IF ( water_pars_f%pars_xy(ind_w_z0h,j,i) /= water_pars_f%fill )&
     3777                   THEN
     3778                      surf_lsm_h(l)%z0h(m)    = water_pars_f%pars_xy(ind_w_z0h,j,i)
     3779                      surf_lsm_h(l)%z0q(m)    = water_pars_f%pars_xy(ind_w_z0h,j,i)
     3780                   ENDIF
     3781                   IF ( water_pars_f%pars_xy(ind_w_lambda_s,j,i) /=               &
     3782                        water_pars_f%fill )                                       &
     3783                      surf_lsm_h(l)%lambda_surface_s(m) =                            &
     3784                                           water_pars_f%pars_xy(ind_w_lambda_s,j,i)
     3785
     3786                   IF ( water_pars_f%pars_xy(ind_w_lambda_u,j,i) /=               &
     3787                         water_pars_f%fill )                                      &
     3788                      surf_lsm_h(l)%lambda_surface_u(m) =                            &
     3789                                           water_pars_f%pars_xy(ind_w_lambda_u,j,i)
     3790
     3791                   IF ( water_pars_f%pars_xy(ind_w_at,j,i) /=                     &
     3792                        water_pars_f%fill )                                       &
     3793                      surf_lsm_h(l)%albedo_type(m,ind_wat_win) =                     &
     3794                                          INT( water_pars_f%pars_xy(ind_w_at,j,i) )
     3795
     3796                   IF ( water_pars_f%pars_xy(ind_w_emis,j,i) /=                   &
     3797                        water_pars_f%fill )                                       &
     3798                      surf_lsm_h(l)%emissivity(m,ind_wat_win) =                      &
     3799                                             water_pars_f%pars_xy(ind_w_emis,j,i)
    37053800                ENDIF
    3706                 IF ( water_pars_f%pars_xy(ind_w_lambda_s,j,i) /=               &
    3707                      water_pars_f%fill )                                       &
    3708                    surf_lsm_h%lambda_surface_s(m) =                            &
    3709                                         water_pars_f%pars_xy(ind_w_lambda_s,j,i)
    3710 
    3711                 IF ( water_pars_f%pars_xy(ind_w_lambda_u,j,i) /=               &
    3712                       water_pars_f%fill )                                      &
    3713                    surf_lsm_h%lambda_surface_u(m) =                            &
    3714                                         water_pars_f%pars_xy(ind_w_lambda_u,j,i)
    3715 
    3716                 IF ( water_pars_f%pars_xy(ind_w_at,j,i) /=                     &
    3717                      water_pars_f%fill )                                       &
    3718                    surf_lsm_h%albedo_type(m,ind_wat_win) =                     &
    3719                                        INT( water_pars_f%pars_xy(ind_w_at,j,i) )
    3720 
    3721                 IF ( water_pars_f%pars_xy(ind_w_emis,j,i) /=                   &
    3722                      water_pars_f%fill )                                       &
    3723                    surf_lsm_h%emissivity(m,ind_wat_win) =                      &
    3724                                           water_pars_f%pars_xy(ind_w_emis,j,i)
    3725              ENDIF
     3801             ENDDO
    37263802          ENDDO
    37273803!
     
    38443920!--    Level 1 initialization of pavement type surfaces. Horizontally
    38453921!--    homogeneous characteristics are assumed
    3846        surf_lsm_h%nzt_pavement = pavement_depth_level
    3847        DO  m = 1, surf_lsm_h%ns
    3848           IF ( surf_lsm_h%pavement_surface(m) )  THEN
    3849              surf_lsm_h%nzt_pavement(m)        = nzt_pavement
    3850              surf_lsm_h%z0(m)                  = z0_pavement
    3851              surf_lsm_h%z0h(m)                 = z0h_pavement
    3852              surf_lsm_h%z0q(m)                 = z0q_pavement
    3853              surf_lsm_h%lambda_surface_s(m)    = pavement_heat_conduct         &
    3854                                                   * ddz_soil(nzb_soil)         &
    3855                                                   * 2.0_wp
    3856              surf_lsm_h%lambda_surface_u(m)    = pavement_heat_conduct         &
    3857                                                   * ddz_soil(nzb_soil)         &
    3858                                                   * 2.0_wp
    3859              surf_lsm_h%c_surface(m)           = pavement_heat_capacity        &
    3860                                                         * dz_soil(nzb_soil)    &
    3861                                                         * 0.25_wp
    3862 
    3863              surf_lsm_h%albedo_type(m,ind_pav_green) = albedo_type
    3864              surf_lsm_h%emissivity(m,ind_pav_green)  = emissivity
    3865 
    3866              surf_lsm_h%pavement_type(m)      = pavement_type
    3867              surf_lsm_h%pavement_type_name(m) = pavement_type_name(pavement_type)
    3868 
    3869              IF ( pavement_type /= 0 )  THEN
    3870                 DO  k = nzb_soil, surf_lsm_h%nzt_pavement(m)
    3871                    surf_lsm_h%lambda_h_def(k,m)    =                           &
    3872                                      pavement_subsurface_pars_1(k,pavement_type)
    3873                    surf_lsm_h%rho_c_total_def(k,m) =                           &
    3874                                      pavement_subsurface_pars_2(k,pavement_type)
    3875                 ENDDO
    3876              ELSE
    3877                 surf_lsm_h%lambda_h_def(:,m)     = pavement_heat_conduct
    3878                 surf_lsm_h%rho_c_total_def(:,m)  = pavement_heat_capacity
     3922       DO  l = 0, 1
     3923          surf_lsm_h(l)%nzt_pavement = pavement_depth_level
     3924          DO  m = 1, surf_lsm_h(l)%ns
     3925             IF ( surf_lsm_h(l)%pavement_surface(m) )  THEN
     3926                surf_lsm_h(l)%nzt_pavement(m)        = nzt_pavement
     3927                surf_lsm_h(l)%z0(m)                  = z0_pavement
     3928                surf_lsm_h(l)%z0h(m)                 = z0h_pavement
     3929                surf_lsm_h(l)%z0q(m)                 = z0q_pavement
     3930                surf_lsm_h(l)%lambda_surface_s(m)    = pavement_heat_conduct         &
     3931                                                     * ddz_soil(nzb_soil)         &
     3932                                                     * 2.0_wp
     3933                surf_lsm_h(l)%lambda_surface_u(m)    = pavement_heat_conduct         &
     3934                                                     * ddz_soil(nzb_soil)         &
     3935                                                     * 2.0_wp
     3936                surf_lsm_h(l)%c_surface(m)           = pavement_heat_capacity        &
     3937                                                           * dz_soil(nzb_soil)    &
     3938                                                           * 0.25_wp
     3939
     3940                surf_lsm_h(l)%albedo_type(m,ind_pav_green) = albedo_type
     3941                surf_lsm_h(l)%emissivity(m,ind_pav_green)  = emissivity
     3942
     3943                surf_lsm_h(l)%pavement_type(m)      = pavement_type
     3944                surf_lsm_h(l)%pavement_type_name(m) = pavement_type_name(pavement_type)
     3945
     3946                IF ( pavement_type /= 0 )  THEN
     3947                   DO  k = nzb_soil, surf_lsm_h(l)%nzt_pavement(m)
     3948                      surf_lsm_h(l)%lambda_h_def(k,m)    =                           &
     3949                                        pavement_subsurface_pars_1(k,pavement_type)
     3950                      surf_lsm_h(l)%rho_c_total_def(k,m) =                           &
     3951                                        pavement_subsurface_pars_2(k,pavement_type)
     3952                   ENDDO
     3953                ELSE
     3954                   surf_lsm_h(l)%lambda_h_def(:,m)     = pavement_heat_conduct
     3955                   surf_lsm_h(l)%rho_c_total_def(:,m)  = pavement_heat_capacity
     3956                ENDIF
    38793957             ENDIF
    3880           ENDIF
     3958          ENDDO
    38813959       ENDDO
    38823960
     
    39264004!
    39274005!--       Horizontal surfaces
    3928           DO  m = 1, surf_lsm_h%ns
    3929              i = surf_lsm_h%i(m)
    3930              j = surf_lsm_h%j(m)
    3931 
    3932              st = pavement_type_f%var(j,i)
    3933              IF ( st /= pavement_type_f%fill  .AND.  st /= 0 )  THEN
    3934 !
    3935 !--             Determine deepmost index of pavement layer
    3936                 DO  k = nzb_soil, nzt_soil
    3937                    IF ( pavement_subsurface_pars_1(k,st) == 9999999.9_wp       &
    3938                    .OR. pavement_subsurface_pars_2(k,st) == 9999999.9_wp)      &
    3939                    THEN
    3940                       surf_lsm_h%nzt_pavement(m) = k-1
    3941                       EXIT
    3942                    ENDIF
    3943                 ENDDO
    3944 
    3945                 surf_lsm_h%z0(m)                = pavement_pars(ind_p_z0,st)
    3946                 surf_lsm_h%z0h(m)               = pavement_pars(ind_p_z0h,st)
    3947                 surf_lsm_h%z0q(m)               = pavement_pars(ind_p_z0h,st)
    3948 
    3949                 surf_lsm_h%lambda_surface_s(m)  =                              &
    3950                                               pavement_subsurface_pars_1(0,st) &
    3951                                                   * ddz_soil(nzb_soil)         &
    3952                                                   * 2.0_wp
    3953                 surf_lsm_h%lambda_surface_u(m)  =                              &
    3954                                               pavement_subsurface_pars_1(0,st) &
    3955                                                   * ddz_soil(nzb_soil)         &
    3956                                                   * 2.0_wp
    3957                 surf_lsm_h%c_surface(m)         =                              &
    3958                                                pavement_subsurface_pars_2(0,st)&
    3959                                                         * dz_soil(nzb_soil)    &
    3960                                                         * 0.25_wp
    3961                 surf_lsm_h%albedo_type(m,ind_pav_green) = INT( pavement_pars(ind_p_at,st) )
    3962                 surf_lsm_h%emissivity(m,ind_pav_green)  = pavement_pars(ind_p_emis,st)
    3963 
    3964                 surf_lsm_h%pavement_type(m)      = st
    3965                 surf_lsm_h%pavement_type_name(m) = pavement_type_name(st)
    3966 
    3967                 DO  k = nzb_soil, surf_lsm_h%nzt_pavement(m)
    3968                    surf_lsm_h%lambda_h_def(k,m)    =                           &
    3969                                      pavement_subsurface_pars_1(k,pavement_type)
    3970                    surf_lsm_h%rho_c_total_def(k,m) =                           &
    3971                                      pavement_subsurface_pars_2(k,pavement_type)
    3972                 ENDDO
    3973              ENDIF
     4006          DO  l = 0, 1
     4007             DO  m = 1, surf_lsm_h(l)%ns
     4008                i = surf_lsm_h(l)%i(m)
     4009                j = surf_lsm_h(l)%j(m)
     4010
     4011                st = pavement_type_f%var(j,i)
     4012                IF ( st /= pavement_type_f%fill  .AND.  st /= 0 )  THEN
     4013   !
     4014   !--             Determine deepmost index of pavement layer
     4015                   DO  k = nzb_soil, nzt_soil
     4016                      IF ( pavement_subsurface_pars_1(k,st) == 9999999.9_wp       &
     4017                      .OR. pavement_subsurface_pars_2(k,st) == 9999999.9_wp)      &
     4018                      THEN
     4019                         surf_lsm_h(l)%nzt_pavement(m) = k-1
     4020                         EXIT
     4021                      ENDIF
     4022                   ENDDO
     4023
     4024                   surf_lsm_h(l)%z0(m)                = pavement_pars(ind_p_z0,st)
     4025                   surf_lsm_h(l)%z0h(m)               = pavement_pars(ind_p_z0h,st)
     4026                   surf_lsm_h(l)%z0q(m)               = pavement_pars(ind_p_z0h,st)
     4027
     4028                   surf_lsm_h(l)%lambda_surface_s(m)  =                              &
     4029                                                 pavement_subsurface_pars_1(0,st) &
     4030                                                     * ddz_soil(nzb_soil)         &
     4031                                                     * 2.0_wp
     4032                   surf_lsm_h(l)%lambda_surface_u(m)  =                              &
     4033                                                 pavement_subsurface_pars_1(0,st) &
     4034                                                     * ddz_soil(nzb_soil)         &
     4035                                                     * 2.0_wp
     4036                   surf_lsm_h(l)%c_surface(m)         =                              &
     4037                                                  pavement_subsurface_pars_2(0,st)&
     4038                                                           * dz_soil(nzb_soil)    &
     4039                                                           * 0.25_wp
     4040                   surf_lsm_h(l)%albedo_type(m,ind_pav_green) = INT( pavement_pars(ind_p_at,st) )
     4041                   surf_lsm_h(l)%emissivity(m,ind_pav_green)  = pavement_pars(ind_p_emis,st)
     4042
     4043                   surf_lsm_h(l)%pavement_type(m)      = st
     4044                   surf_lsm_h(l)%pavement_type_name(m) = pavement_type_name(st)
     4045
     4046                   DO  k = nzb_soil, surf_lsm_h(l)%nzt_pavement(m)
     4047                      surf_lsm_h(l)%lambda_h_def(k,m)    =                           &
     4048                                        pavement_subsurface_pars_1(k,pavement_type)
     4049                      surf_lsm_h(l)%rho_c_total_def(k,m) =                           &
     4050                                        pavement_subsurface_pars_2(k,pavement_type)
     4051                   ENDDO
     4052                ENDIF
     4053             ENDDO
    39744054          ENDDO
    39754055!
     
    40364116!
    40374117!--       Horizontal surfaces
    4038           DO  m = 1, surf_lsm_h%ns
    4039              i = surf_lsm_h%i(m)
    4040              j = surf_lsm_h%j(m)
    4041 !
    4042 !--          If surface element is not a pavement surface and any value in
    4043 !--          pavement_pars is given, neglect this information and give an
    4044 !--          informative message that this value will not be used.
    4045              IF ( .NOT. surf_lsm_h%pavement_surface(m)  .AND.                  &
    4046                    ANY( pavement_pars_f%pars_xy(:,j,i) /=                      &
    4047                    pavement_pars_f%fill ) )  THEN
    4048                 WRITE( message_string, * )                                     &
    4049                               'surface element at grid point (j,i) = (',       &
    4050                               j, i, ') is not a pavement surface, ',           &
    4051                               'so that information given in ',                 &
    4052                               'pavement_pars at this point is neglected.'
    4053                 CALL message( 'land_surface_model_mod', 'PA0647', 0, 0, myid, 6, 0 )
    4054              ELSE
    4055                 IF ( pavement_pars_f%pars_xy(ind_p_z0,j,i) /=                  &
    4056                      pavement_pars_f%fill )                                    &
    4057                    surf_lsm_h%z0(m)  = pavement_pars_f%pars_xy(ind_p_z0,j,i)
    4058                 IF ( pavement_pars_f%pars_xy(ind_p_z0h,j,i) /=                 &
    4059                      pavement_pars_f%fill )  THEN
    4060                    surf_lsm_h%z0h(m) = pavement_pars_f%pars_xy(ind_p_z0h,j,i)
    4061                    surf_lsm_h%z0q(m) = pavement_pars_f%pars_xy(ind_p_z0h,j,i)
     4118          DO  l = 0, 1
     4119             DO  m = 1, surf_lsm_h(l)%ns
     4120                i = surf_lsm_h(l)%i(m)
     4121                j = surf_lsm_h(l)%j(m)
     4122   !
     4123   !--          If surface element is not a pavement surface and any value in
     4124   !--          pavement_pars is given, neglect this information and give an
     4125   !--          informative message that this value will not be used.
     4126                IF ( .NOT. surf_lsm_h(l)%pavement_surface(m)  .AND.                  &
     4127                      ANY( pavement_pars_f%pars_xy(:,j,i) /=                      &
     4128                      pavement_pars_f%fill ) )  THEN
     4129                   WRITE( message_string, * )                                     &
     4130                                 'surface element at grid point (j,i) = (',       &
     4131                                 j, i, ') is not a pavement surface, ',           &
     4132                                 'so that information given in ',                 &
     4133                                 'pavement_pars at this point is neglected.'
     4134                   CALL message( 'land_surface_model_mod', 'PA0647', 0, 0, myid, 6, 0 )
     4135                ELSE
     4136                   IF ( pavement_pars_f%pars_xy(ind_p_z0,j,i) /=                  &
     4137                        pavement_pars_f%fill )                                    &
     4138                      surf_lsm_h(l)%z0(m)  = pavement_pars_f%pars_xy(ind_p_z0,j,i)
     4139                   IF ( pavement_pars_f%pars_xy(ind_p_z0h,j,i) /=                 &
     4140                        pavement_pars_f%fill )  THEN
     4141                      surf_lsm_h(l)%z0h(m) = pavement_pars_f%pars_xy(ind_p_z0h,j,i)
     4142                      surf_lsm_h(l)%z0q(m) = pavement_pars_f%pars_xy(ind_p_z0h,j,i)
     4143                   ENDIF
     4144                   IF ( pavement_pars_f%pars_xy(ind_p_at,j,i) /=                  &
     4145                        pavement_pars_f%fill )                                    &
     4146                      surf_lsm_h(l)%albedo_type(m,ind_pav_green) =                   &
     4147                                      INT( pavement_pars_f%pars_xy(ind_p_at,j,i) )
     4148                   IF ( pavement_pars_f%pars_xy(ind_p_emis,j,i) /=                &
     4149                        pavement_pars_f%fill )                                    &
     4150                      surf_lsm_h(l)%emissivity(m,ind_pav_green)  =                   &
     4151                                      pavement_pars_f%pars_xy(ind_p_emis,j,i)
    40624152                ENDIF
    4063                 IF ( pavement_pars_f%pars_xy(ind_p_at,j,i) /=                  &
    4064                      pavement_pars_f%fill )                                    &
    4065                    surf_lsm_h%albedo_type(m,ind_pav_green) =                   &
    4066                                    INT( pavement_pars_f%pars_xy(ind_p_at,j,i) )
    4067                 IF ( pavement_pars_f%pars_xy(ind_p_emis,j,i) /=                &
    4068                      pavement_pars_f%fill )                                    &
    4069                    surf_lsm_h%emissivity(m,ind_pav_green)  =                   &
    4070                                    pavement_pars_f%pars_xy(ind_p_emis,j,i)
    4071              ENDIF
    4072 
     4153             ENDDO
    40734154          ENDDO
    40744155!
     
    41244205!--       Set pavement depth to nzt_soil. Please note, this is just a
    41254206!--       workaround at the moment.
    4126           DO  m = 1, surf_lsm_h%ns
    4127              IF ( surf_lsm_h%pavement_surface(m) )  THEN
    4128 
    4129                 i = surf_lsm_h%i(m)
    4130                 j = surf_lsm_h%j(m)
    4131 
    4132                 surf_lsm_h%nzt_pavement(m) = nzt_soil
    4133 
    4134                 IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i) &
    4135                      /= pavement_subsurface_pars_f%fill )  THEN
    4136                    surf_lsm_h%lambda_surface_s(m)  =                           &
    4137                       pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i)&
    4138                                                   * ddz_soil(nzb_soil)         &
    4139                                                   * 2.0_wp
    4140                    surf_lsm_h%lambda_surface_u(m)  =                           &
    4141                       pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i)&
    4142                                                   * ddz_soil(nzb_soil)         &
    4143                                                   * 2.0_wp
     4207          DO  l = 0, 1
     4208             DO  m = 1, surf_lsm_h(l)%ns
     4209                IF ( surf_lsm_h(l)%pavement_surface(m) )  THEN
     4210
     4211                   i = surf_lsm_h(l)%i(m)
     4212                   j = surf_lsm_h(l)%j(m)
     4213
     4214                   surf_lsm_h(l)%nzt_pavement(m) = nzt_soil
     4215
     4216                   IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i) &
     4217                        /= pavement_subsurface_pars_f%fill )  THEN
     4218                      surf_lsm_h(l)%lambda_surface_s(m)  =                           &
     4219                         pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i)&
     4220                                                     * ddz_soil(nzb_soil)         &
     4221                                                     * 2.0_wp
     4222                      surf_lsm_h(l)%lambda_surface_u(m)  =                           &
     4223                         pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,0,j,i)&
     4224                                                     * ddz_soil(nzb_soil)         &
     4225                                                     * 2.0_wp
     4226                   ENDIF
     4227                   IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,0,j,i) /= &
     4228                        pavement_subsurface_pars_f%fill )  THEN
     4229                      surf_lsm_h(l)%c_surface(m)     =                               &
     4230                         pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,0,j,i)   &
     4231                                                     * dz_soil(nzb_soil)          &
     4232                                                     * 0.25_wp
     4233                   ENDIF
     4234
     4235                   DO  k = nzb_soil, nzt_soil
     4236                      IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,k,j,i) /= &
     4237                           pavement_subsurface_pars_f%fill )  THEN
     4238                         surf_lsm_h(l)%lambda_h_def(k,m) =                              &
     4239                             pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,k,j,i)
     4240                      ENDIF
     4241                      IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,k,j,i) /= &
     4242                           pavement_subsurface_pars_f%fill )  THEN
     4243                         surf_lsm_h(l)%rho_c_total_def(k,m) =                           &
     4244                             pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,k,j,i)
     4245                      ENDIF
     4246                   ENDDO
     4247
    41444248                ENDIF
    4145                 IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,0,j,i) /= &
    4146                      pavement_subsurface_pars_f%fill )  THEN
    4147                    surf_lsm_h%c_surface(m)     =                               &
    4148                       pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,0,j,i)   &
    4149                                                   * dz_soil(nzb_soil)          &
    4150                                                   * 0.25_wp
    4151                 ENDIF
    4152 
    4153                 DO  k = nzb_soil, nzt_soil
    4154                    IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,k,j,i) /= &
    4155                         pavement_subsurface_pars_f%fill )  THEN
    4156                       surf_lsm_h%lambda_h_def(k,m) =                              &
    4157                           pavement_subsurface_pars_f%pars_xyz(ind_p_lambda_h,k,j,i)
    4158                    ENDIF
    4159                    IF ( pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,k,j,i) /= &
    4160                         pavement_subsurface_pars_f%fill )  THEN
    4161                       surf_lsm_h%rho_c_total_def(k,m) =                           &
    4162                           pavement_subsurface_pars_f%pars_xyz(ind_p_rho_c,k,j,i)
    4163                    ENDIF
    4164                 ENDDO
    4165 
    4166              ENDIF
     4249             ENDDO
    41674250          ENDDO
    41684251          DO  l = 0, 3
     
    42174300!--    the albedo type has been already given by the pars, albedo_type overwrites these values.
    42184301       IF ( albedo_type_f%from_file )  THEN
    4219           DO  m = 1, surf_lsm_h%ns
    4220              i = surf_lsm_h%i(m)
    4221              j = surf_lsm_h%j(m)
    4222              IF ( albedo_type_f%var(j,i) /= albedo_type_f%fill )                                   &
    4223                 surf_lsm_h%albedo_type(m,:) = albedo_type_f%var(j,i)
     4302          DO  l = 0, 1
     4303             DO  m = 1, surf_lsm_h(l)%ns
     4304                i = surf_lsm_h(l)%i(m)
     4305                j = surf_lsm_h(l)%j(m)
     4306                IF ( albedo_type_f%var(j,i) /= albedo_type_f%fill )                                   &
     4307                   surf_lsm_h(l)%albedo_type(m,:) = albedo_type_f%var(j,i)
     4308             ENDDO
    42244309          ENDDO
    42254310          DO  l = 0, 3
     
    42444329!--       data provided from input file contains fill values at some locations.
    42454330!--       Level 1, initialization via profiles given in parameter file
    4246           DO  m = 1, surf_lsm_h%ns
    4247              IF ( surf_lsm_h%vegetation_surface(m)  .OR.                       &
    4248                   surf_lsm_h%pavement_surface(m) )  THEN
    4249                 DO  k = nzb_soil, nzt_soil
    4250                    t_soil_h%var_2d(k,m) = soil_temperature(k)
    4251                    m_soil_h%var_2d(k,m) = soil_moisture(k)
    4252                 ENDDO
    4253                 t_soil_h%var_2d(nzt_soil+1,m) = deep_soil_temperature
    4254              ENDIF
     4331          DO  l = 0, 1
     4332             DO  m = 1, surf_lsm_h(l)%ns
     4333                IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                       &
     4334                     surf_lsm_h(l)%pavement_surface(m) )  THEN
     4335                   DO  k = nzb_soil, nzt_soil
     4336                      t_soil_h(l)%var_2d(k,m) = soil_temperature(k)
     4337                      m_soil_h(l)%var_2d(k,m) = soil_moisture(k)
     4338                   ENDDO
     4339                   t_soil_h(l)%var_2d(nzt_soil+1,m) = deep_soil_temperature
     4340                ENDIF
     4341             ENDDO
    42554342          ENDDO
    42564343          DO  l = 0, 3
     
    45434630
    45444631             IF ( init_3d%lod_msoil == 1 )  THEN
    4545                 DO  m = 1, surf_lsm_h%ns
    4546                    IF ( surf_lsm_h%vegetation_surface(m)  .OR.                 &
    4547                         surf_lsm_h%pavement_surface(m) )  THEN
    4548 
    4549                       CALL interpolate_soil_profile(                           &
    4550                                        m_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    4551                                        init_3d%msoil_1d(:),                    &
    4552                                        zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
    4553                                        nzb_soil, nzt_soil,                     &
    4554                                        nzb_soil, init_3d%nzs-1 )
    4555                    ENDIF
     4632                DO  l = 0, 1
     4633                   DO  m = 1, surf_lsm_h(l)%ns
     4634                      IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                 &
     4635                           surf_lsm_h(l)%pavement_surface(m) )  THEN
     4636
     4637                         CALL interpolate_soil_profile(                           &
     4638                                          m_soil_h(l)%var_2d(nzb_soil:nzt_soil,m),&
     4639                                          init_3d%msoil_1d(:),                    &
     4640                                          zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
     4641                                          nzb_soil, nzt_soil,                     &
     4642                                          nzb_soil, init_3d%nzs-1 )
     4643                      ENDIF
     4644                   ENDDO
    45564645                ENDDO
    45574646                DO  l = 0, 3
     
    45694658                ENDDO
    45704659             ELSE
    4571 
    4572                 DO  m = 1, surf_lsm_h%ns
    4573                    IF ( surf_lsm_h%vegetation_surface(m)  .OR.                 &
    4574                         surf_lsm_h%pavement_surface(m) )  THEN
    4575                       i = surf_lsm_h%i(m)
    4576                       j = surf_lsm_h%j(m)
    4577 
    4578                       IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil )     &
    4579                          CALL interpolate_soil_profile(                        &
    4580                                        m_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    4581                                        init_3d%msoil_3d(:,j,i),                &
    4582                                        zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
    4583                                        nzb_soil, nzt_soil,                     &
    4584                                        nzb_soil, init_3d%nzs-1 )
    4585                    ENDIF
     4660                DO  l = 0, 1
     4661                   DO  m = 1, surf_lsm_h(l)%ns
     4662                      IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                 &
     4663                           surf_lsm_h(l)%pavement_surface(m) )  THEN
     4664                         i = surf_lsm_h(l)%i(m)
     4665                         j = surf_lsm_h(l)%j(m)
     4666
     4667                         IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil )     &
     4668                            CALL interpolate_soil_profile(                        &
     4669                                          m_soil_h(l)%var_2d(nzb_soil:nzt_soil,m),&
     4670                                          init_3d%msoil_3d(:,j,i),                &
     4671                                          zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
     4672                                          nzb_soil, nzt_soil,                     &
     4673                                          nzb_soil, init_3d%nzs-1 )
     4674                      ENDIF
     4675                   ENDDO
    45864676                ENDDO
    45874677                DO  l = 0, 3
     
    46144704
    46154705             IF ( init_3d%lod_tsoil == 1 )  THEN ! change to 1 if provided correctly by INIFOR
    4616                 DO  m = 1, surf_lsm_h%ns
    4617                    IF ( surf_lsm_h%vegetation_surface(m)  .OR.                 &
    4618                         surf_lsm_h%pavement_surface(m) )  THEN
    4619                       CALL interpolate_soil_profile(                           &
    4620                                        t_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    4621                                        init_3d%tsoil_1d(:),                    &
    4622                                        zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
    4623                                        nzb_soil, nzt_soil,                     &
    4624                                        nzb_soil, init_3d%nzs-1 )
    4625 !
    4626 !--                   Set boundary condition, i.e. deep soil temperature
    4627                       t_soil_h%var_2d(nzt_soil+1,m) = t_soil_h%var_2d(nzt_soil,m)
    4628                    ENDIF
     4706                DO  l = 0, 1
     4707                   DO  m = 1, surf_lsm_h(l)%ns
     4708                      IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                 &
     4709                           surf_lsm_h(l)%pavement_surface(m) )  THEN
     4710                         CALL interpolate_soil_profile(                           &
     4711                                          t_soil_h(l)%var_2d(nzb_soil:nzt_soil,m),   &
     4712                                          init_3d%tsoil_1d(:),                    &
     4713                                          zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
     4714                                          nzb_soil, nzt_soil,                     &
     4715                                          nzb_soil, init_3d%nzs-1 )
     4716   !
     4717   !--                   Set boundary condition, i.e. deep soil temperature
     4718                         t_soil_h(l)%var_2d(nzt_soil+1,m) = t_soil_h(l)%var_2d(nzt_soil,m)
     4719                      ENDIF
     4720                   ENDDO
    46294721                ENDDO
    46304722                DO  l = 0, 3
     
    46464738                ENDDO
    46474739             ELSE
    4648 
    4649                 DO  m = 1, surf_lsm_h%ns
    4650                    IF ( surf_lsm_h%vegetation_surface(m)  .OR.                 &
    4651                         surf_lsm_h%pavement_surface(m) )  THEN
    4652                       i = surf_lsm_h%i(m)
    4653                       j = surf_lsm_h%j(m)
    4654 
    4655                       IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil )     &
    4656                          CALL interpolate_soil_profile(                        &
    4657                                        t_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    4658                                        init_3d%tsoil_3d(:,j,i),                &
    4659                                        zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
    4660                                        nzb_soil, nzt_soil,                     &
    4661                                        nzb_soil, init_3d%nzs-1 )
    4662 !
    4663 !--                   Set boundary condition, i.e. deep soil temperature
    4664                       t_soil_h%var_2d(nzt_soil+1,m) = t_soil_h%var_2d(nzt_soil,m)
    4665                    ENDIF
     4740                DO  l = 0, 1
     4741                   DO  m = 1, surf_lsm_h(l)%ns
     4742                      IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                 &
     4743                           surf_lsm_h(l)%pavement_surface(m) )  THEN
     4744                         i = surf_lsm_h(l)%i(m)
     4745                         j = surf_lsm_h(l)%j(m)
     4746
     4747                         IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil )     &
     4748                            CALL interpolate_soil_profile(                        &
     4749                                          t_soil_h(l)%var_2d(nzb_soil:nzt_soil,m),   &
     4750                                          init_3d%tsoil_3d(:,j,i),                &
     4751                                          zs(nzb_soil:nzt_soil), init_3d%z_soil,  &
     4752                                          nzb_soil, nzt_soil,                     &
     4753                                          nzb_soil, init_3d%nzs-1 )
     4754   !
     4755   !--                   Set boundary condition, i.e. deep soil temperature
     4756                         t_soil_h(l)%var_2d(nzt_soil+1,m) = t_soil_h(l)%var_2d(nzt_soil,m)
     4757                      ENDIF
     4758                   ENDDO
    46664759                ENDDO
    46674760                DO  l = 0, 3
     
    46994792!--       point errors. Hence, limit the soil moisture to its saturation value
    47004793!--       and give a warning.
    4701           DO  m = 1, surf_lsm_h%ns
    4702              IF ( surf_lsm_h%vegetation_surface(m)  .OR.                       &
    4703                   surf_lsm_h%pavement_surface(m) )  THEN
    4704                 DO  k = nzb_soil, nzt_soil
    4705                    IF ( m_soil_h%var_2d(k,m) > surf_lsm_h%m_sat(k,m) )  THEN
    4706                       m_soil_h%var_2d(k,m) = surf_lsm_h%m_sat(k,m)
    4707                       WRITE( message_string, * ) 'soil moisture is higher '//  &
    4708                             'than its saturation value at (k,j,i) ', k,        &
    4709                             surf_lsm_h%i(m), surf_lsm_h%j(m), ' and is ' //    &
    4710                             'thus limited to this value to maintain stability.'
    4711                       CALL message( 'lsm_init', 'PA0458', 0, 1, myid, 6, 0 )
    4712                    ENDIF
    4713                 ENDDO
    4714              ENDIF
     4794          DO l = 0, 1
     4795             DO  m = 1, surf_lsm_h(l)%ns
     4796                IF ( surf_lsm_h(l)%vegetation_surface(m)  .OR.                       &
     4797                     surf_lsm_h(l)%pavement_surface(m) )  THEN
     4798                   DO  k = nzb_soil, nzt_soil
     4799                      IF ( m_soil_h(l)%var_2d(k,m) > surf_lsm_h(l)%m_sat(k,m) )  THEN
     4800                         m_soil_h(l)%var_2d(k,m) = surf_lsm_h(l)%m_sat(k,m)
     4801                         WRITE( message_string, * ) 'soil moisture is higher '//     &
     4802                               'than its saturation value at (k,j,i) ', k,           &
     4803                               surf_lsm_h(l)%i(m), surf_lsm_h(l)%j(m), ' and is ' // &
     4804                               'thus limited to this value to maintain stability.'
     4805                         CALL message( 'lsm_init', 'PA0458', 0, 1, myid, 6, 0 )
     4806                      ENDIF
     4807                   ENDDO
     4808                ENDIF
     4809             ENDDO
    47154810          ENDDO
    47164811          DO  l = 0, 3
     
    47374832!
    47384833!--       Further initialization
    4739           DO  m = 1, surf_lsm_h%ns
    4740 
    4741              i   = surf_lsm_h%i(m)
    4742              j   = surf_lsm_h%j(m)
    4743              k   = surf_lsm_h%k(m)
    4744 !
    4745 !--          Initialize surface temperature with soil temperature in the uppermost
    4746 !--          uppermost layer
    4747              t_surface_h%var_1d(m)    = t_soil_h%var_2d(nzb_soil,m)
    4748              surf_lsm_h%pt_surface(m) = t_soil_h%var_2d(nzb_soil,m) / exner(nzb)
    4749 
    4750              IF ( bulk_cloud_model  .OR. cloud_droplets ) THEN
    4751                 surf_lsm_h%pt1(m) = pt(k,j,i) + lv_d_cp * d_exner(k) * ql(k,j,i)
    4752              ELSE
    4753                 surf_lsm_h%pt1(m) = pt(k,j,i)
    4754              ENDIF
    4755 !
    4756 !--          Assure that r_a cannot be zero at model start
    4757              IF ( surf_lsm_h%pt1(m) == surf_lsm_h%pt_surface(m) )              &
    4758                 surf_lsm_h%pt1(m) = surf_lsm_h%pt1(m) + 1.0E-20_wp
    4759 
    4760              surf_lsm_h%us(m)   = 0.1_wp
    4761              surf_lsm_h%ts(m)   = ( surf_lsm_h%pt1(m) - surf_lsm_h%pt_surface(m) )&
    4762                                   / surf_lsm_h%r_a(m)
    4763              surf_lsm_h%shf(m)  = - surf_lsm_h%us(m) * surf_lsm_h%ts(m)        &
    4764                                   * rho_surface
     4834          DO l = 0, 1
     4835             DO  m = 1, surf_lsm_h(l)%ns
     4836
     4837                i   = surf_lsm_h(l)%i(m)
     4838                j   = surf_lsm_h(l)%j(m)
     4839                k   = surf_lsm_h(l)%k(m)
     4840   !
     4841   !--          Initialize surface temperature with soil temperature in the uppermost
     4842   !--          uppermost layer
     4843                t_surface_h(l)%var_1d(m)    = t_soil_h(l)%var_2d(nzb_soil,m)
     4844                surf_lsm_h(l)%pt_surface(m) = t_soil_h(l)%var_2d(nzb_soil,m) / exner(nzb)
     4845
     4846                IF ( bulk_cloud_model  .OR. cloud_droplets ) THEN
     4847                   surf_lsm_h(l)%pt1(m) = pt(k,j,i) + lv_d_cp * d_exner(k) * ql(k,j,i)
     4848                ELSE
     4849                   surf_lsm_h(l)%pt1(m) = pt(k,j,i)
     4850                ENDIF
     4851   !
     4852   !--          Assure that r_a cannot be zero at model start
     4853                IF ( surf_lsm_h(l)%pt1(m) == surf_lsm_h(l)%pt_surface(m) )              &
     4854                   surf_lsm_h(l)%pt1(m) = surf_lsm_h(l)%pt1(m) + 1.0E-20_wp
     4855
     4856                surf_lsm_h(l)%us(m)   = 0.1_wp
     4857                surf_lsm_h(l)%ts(m)   = ( surf_lsm_h(l)%pt1(m) - surf_lsm_h(l)%pt_surface(m) )&
     4858                                     / surf_lsm_h(l)%r_a(m)
     4859                surf_lsm_h(l)%shf(m)  = - surf_lsm_h(l)%us(m) * surf_lsm_h(l)%ts(m)        &
     4860                                     * rho_surface
     4861             ENDDO
    47654862          ENDDO
    47664863!
     
    48024899!--    Level 1 initialization of root distribution - provided by the user via
    48034900!--    via namelist.
    4804        DO  m = 1, surf_lsm_h%ns
    4805           DO  k = nzb_soil, nzt_soil
    4806              surf_lsm_h%root_fr(k,m) = root_fraction(k)
     4901       DO l = 0, 1
     4902          DO  m = 1, surf_lsm_h(l)%ns
     4903             DO  k = nzb_soil, nzt_soil
     4904                surf_lsm_h(l)%root_fr(k,m) = root_fraction(k)
     4905             ENDDO
    48074906          ENDDO
    48084907       ENDDO
     
    48744973!
    48754974!--       Map calculated root fractions
    4876           DO  m = 1, surf_lsm_h%ns
    4877              DO  k = nzb_soil, nzt_soil
    4878                 IF ( surf_lsm_h%pavement_surface(m)  .AND.                     &
    4879                      k <= surf_lsm_h%nzt_pavement(m) )  THEN
    4880                    surf_lsm_h%root_fr(k,m) = 0.0_wp
    4881                 ELSE
    4882                    surf_lsm_h%root_fr(k,m) = root_fraction(k)
     4975          DO  l = 0, 1
     4976                DO  m = 1, surf_lsm_h(l)%ns
     4977                   DO  k = nzb_soil, nzt_soil
     4978                      IF ( surf_lsm_h(l)%pavement_surface(m)  .AND.                     &
     4979                           k <= surf_lsm_h(l)%nzt_pavement(m) )  THEN
     4980                         surf_lsm_h(l)%root_fr(k,m) = 0.0_wp
     4981                      ELSE
     4982                         surf_lsm_h(l)%root_fr(k,m) = root_fraction(k)
     4983                      ENDIF
     4984
     4985                   ENDDO
     4986!
     4987!--             Normalize so that the sum = 1. Only relevant when the root
     4988!--             distribution was set to zero due to pavement at some layers.
     4989                IF ( SUM( surf_lsm_h(l)%root_fr(:,m) ) > 0.0_wp )  THEN
     4990                   DO k = nzb_soil, nzt_soil
     4991                      surf_lsm_h(l)%root_fr(k,m) = surf_lsm_h(l)%root_fr(k,m)           &
     4992                      / SUM( surf_lsm_h(l)%root_fr(:,m) )
     4993                   ENDDO
    48834994                ENDIF
    4884 
    4885              ENDDO
    4886 !
    4887 !--          Normalize so that the sum = 1. Only relevant when the root
    4888 !--          distribution was set to zero due to pavement at some layers.
    4889              IF ( SUM( surf_lsm_h%root_fr(:,m) ) > 0.0_wp )  THEN
    4890                 DO k = nzb_soil, nzt_soil
    4891                    surf_lsm_h%root_fr(k,m) = surf_lsm_h%root_fr(k,m)           &
    4892                    / SUM( surf_lsm_h%root_fr(:,m) )
    4893                 ENDDO
    4894              ENDIF
     4995             ENDDO
    48954996          ENDDO
    48964997          DO  l = 0, 3
     
    49205021!--    Take value from file
    49215022       IF ( root_area_density_lsm_f%from_file )  THEN
    4922           DO  m = 1, surf_lsm_h%ns
    4923              IF ( surf_lsm_h%vegetation_surface(m) )  THEN
    4924                 i = surf_lsm_h%i(m)
    4925                 j = surf_lsm_h%j(m)
    4926                 DO  k = nzb_soil, nzt_soil
    4927                    surf_lsm_h%root_fr(k,m) = root_area_density_lsm_f%var(k,j,i)
    4928                 ENDDO
    4929 
    4930              ENDIF
     5023          DO  l = 0, 1
     5024             DO  m = 1, surf_lsm_h(l)%ns
     5025                IF ( surf_lsm_h(l)%vegetation_surface(m) )  THEN
     5026                   i = surf_lsm_h(l)%i(m)
     5027                   j = surf_lsm_h(l)%j(m)
     5028                   DO  k = nzb_soil, nzt_soil
     5029                      surf_lsm_h(l)%root_fr(k,m) = root_area_density_lsm_f%var(k,j,i)
     5030                   ENDDO
     5031
     5032                ENDIF
     5033             ENDDO
    49315034          ENDDO
    4932 
    49335035          DO  l = 0, 3
    49345036             DO  m = 1, surf_lsm_v(l)%ns
     
    49835085!--    Finally, make some consistency checks.
    49845086!--    Ceck for illegal combination of LAI and vegetation coverage.
    4985        IF ( ANY( .NOT. surf_lsm_h%pavement_surface  .AND.                      &
    4986                  surf_lsm_h%lai == 0.0_wp  .AND.  surf_lsm_h%c_veg == 1.0_wp ) &
    4987           )  THEN
    4988           message_string = 'For non-pavement surfaces the combination ' //     &
    4989                            ' lai = 0.0 and c_veg = 1.0 is not allowed.'
    4990           CALL message( 'lsm_init', 'PA0671', 2, 2, myid, 6, 0 )
    4991        ENDIF
    4992 
     5087       DO  l = 0, 1
     5088          IF ( ANY( .NOT. surf_lsm_h(l)%pavement_surface  .AND.                      &
     5089                    surf_lsm_h(l)%lai == 0.0_wp  .AND.  surf_lsm_h(l)%c_veg == 1.0_wp ) &
     5090             )  THEN
     5091             message_string = 'For non-pavement surfaces the combination ' //     &
     5092                              ' lai = 0.0 and c_veg = 1.0 is not allowed.'
     5093             CALL message( 'lsm_init', 'PA0671', 2, 2, myid, 6, 0 )
     5094          ENDIF
     5095       ENDDO
    49935096       DO  l = 0, 3
    49945097          IF ( ANY( .NOT. surf_lsm_v(l)%pavement_surface  .AND.                &
     
    50055108!--    necessary. This case, give an informative message. Note, to avoid
    50065109!--    that the job-protocoll is messed-up, this message is only given once.
     5110!--    ???what about the vertical surfaces???
    50075111       flag_exceed_z0  = .FALSE.
    50085112       flag_exceed_z0h = .FALSE.
    5009        DO  m = 1, surf_lsm_h%ns
    5010           IF ( surf_lsm_h%z0(m) > 0.5_wp * surf_lsm_h%z_mo(m) )  THEN
    5011              surf_lsm_h%z0(m) = 0.5_wp * surf_lsm_h%z_mo(m)
    5012              flag_exceed_z0   = .TRUE.
    5013           ENDIF
    5014           IF ( surf_lsm_h%z0h(m) > 0.5_wp * surf_lsm_h%z_mo(m) )  THEN
    5015              surf_lsm_h%z0h(m) = 0.5_wp * surf_lsm_h%z_mo(m)
    5016              surf_lsm_h%z0q(m) = 0.5_wp * surf_lsm_h%z_mo(m)
    5017              flag_exceed_z0h   = .TRUE.
    5018           ENDIF
     5113       DO l = 0, 1
     5114          DO  m = 1, surf_lsm_h(l)%ns
     5115             IF ( surf_lsm_h(l)%z0(m) > 0.5_wp * surf_lsm_h(l)%z_mo(m) )  THEN
     5116                surf_lsm_h(l)%z0(m) = 0.5_wp * surf_lsm_h(l)%z_mo(m)
     5117                flag_exceed_z0   = .TRUE.
     5118             ENDIF
     5119             IF ( surf_lsm_h(l)%z0h(m) > 0.5_wp * surf_lsm_h(l)%z_mo(m) )  THEN
     5120                surf_lsm_h(l)%z0h(m) = 0.5_wp * surf_lsm_h(l)%z_mo(m)
     5121                surf_lsm_h(l)%z0q(m) = 0.5_wp * surf_lsm_h(l)%z_mo(m)
     5122                flag_exceed_z0h   = .TRUE.
     5123             ENDIF
     5124          ENDDO
    50195125       ENDDO
    50205126#if defined( __parallel )
     
    51025208!
    51035209!--    Horizontal surfaces
    5104        ALLOCATE ( m_liq_h_1%var_1d(1:surf_lsm_h%ns)                      )
    5105        ALLOCATE ( m_liq_h_2%var_1d(1:surf_lsm_h%ns)                      )
    5106        ALLOCATE ( t_surface_h_1%var_1d(1:surf_lsm_h%ns)                  )
    5107        ALLOCATE ( t_surface_h_2%var_1d(1:surf_lsm_h%ns)                  )
    5108        ALLOCATE ( m_soil_h_1%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h%ns)   )
    5109        ALLOCATE ( m_soil_h_2%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h%ns)   )
    5110        ALLOCATE ( t_soil_h_1%var_2d(nzb_soil:nzt_soil+1,1:surf_lsm_h%ns) )
    5111        ALLOCATE ( t_soil_h_2%var_2d(nzb_soil:nzt_soil+1,1:surf_lsm_h%ns) )
     5210       DO l = 0, 1
     5211          ALLOCATE ( m_liq_h_1(l)%var_1d(1:surf_lsm_h(l)%ns)                      )
     5212          ALLOCATE ( m_liq_h_2(l)%var_1d(1:surf_lsm_h(l)%ns)                      )
     5213          ALLOCATE ( t_surface_h_1(l)%var_1d(1:surf_lsm_h(l)%ns)                  )
     5214          ALLOCATE ( t_surface_h_2(l)%var_1d(1:surf_lsm_h(l)%ns)                  )
     5215          ALLOCATE ( m_soil_h_1(l)%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)   )
     5216          ALLOCATE ( m_soil_h_2(l)%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)   )
     5217          ALLOCATE ( t_soil_h_1(l)%var_2d(nzb_soil:nzt_soil+1,1:surf_lsm_h(l)%ns) )
     5218          ALLOCATE ( t_soil_h_2(l)%var_2d(nzb_soil:nzt_soil+1,1:surf_lsm_h(l)%ns) )
     5219       ENDDO
    51125220!
    51135221!--    Vertical surfaces
     
    51265234!--    Allocate array for heat flux in W/m2, required for radiation?
    51275235!--    Consider to remove this array
    5128        ALLOCATE( surf_lsm_h%surfhf(1:surf_lsm_h%ns) )
     5236       DO l = 0, 1
     5237          ALLOCATE( surf_lsm_h(l)%surfhf(1:surf_lsm_h(l)%ns) )
     5238       ENDDO
    51295239       DO  l = 0, 3
    51305240          ALLOCATE( surf_lsm_v(l)%surfhf(1:surf_lsm_v(l)%ns) )
     
    51355245!--    Allocate intermediate timestep arrays
    51365246!--    Horizontal surfaces
    5137        ALLOCATE ( tm_liq_h_m%var_1d(1:surf_lsm_h%ns)                     )
    5138        ALLOCATE ( tt_surface_h_m%var_1d(1:surf_lsm_h%ns)                 )
    5139        ALLOCATE ( tm_soil_h_m%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h%ns)  )
    5140        ALLOCATE ( tt_soil_h_m%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h%ns)  )
     5247       DO l = 0, 1
     5248          ALLOCATE ( tm_liq_h_m(l)%var_1d(1:surf_lsm_h(l)%ns)                     )
     5249          ALLOCATE ( tt_surface_h_m(l)%var_1d(1:surf_lsm_h(l)%ns)                 )
     5250          ALLOCATE ( tm_soil_h_m(l)%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)  )
     5251          ALLOCATE ( tt_soil_h_m(l)%var_2d(nzb_soil:nzt_soil,1:surf_lsm_h(l)%ns)  )
     5252       ENDDO
    51415253!
    51425254!--    Horizontal surfaces
     
    51515263!--    Allocate 2D vegetation model arrays
    51525264!--    Horizontal surfaces
    5153        ALLOCATE ( surf_lsm_h%building_surface(1:surf_lsm_h%ns)    )
    5154        ALLOCATE ( surf_lsm_h%c_liq(1:surf_lsm_h%ns)               )
    5155        ALLOCATE ( surf_lsm_h%c_surface(1:surf_lsm_h%ns)           )
    5156        ALLOCATE ( surf_lsm_h%c_veg(1:surf_lsm_h%ns)               )
    5157        ALLOCATE ( surf_lsm_h%f_sw_in(1:surf_lsm_h%ns)             )
    5158        ALLOCATE ( surf_lsm_h%ghf(1:surf_lsm_h%ns)                 )
    5159        ALLOCATE ( surf_lsm_h%g_d(1:surf_lsm_h%ns)                 )
    5160        ALLOCATE ( surf_lsm_h%lai(1:surf_lsm_h%ns)                 )
    5161        ALLOCATE ( surf_lsm_h%lambda_surface_u(1:surf_lsm_h%ns)    )
    5162        ALLOCATE ( surf_lsm_h%lambda_surface_s(1:surf_lsm_h%ns)    )
    5163        ALLOCATE ( surf_lsm_h%nzt_pavement(1:surf_lsm_h%ns)        )
    5164        ALLOCATE ( surf_lsm_h%pavement_surface(1:surf_lsm_h%ns)    )
    5165        ALLOCATE ( surf_lsm_h%qsws_soil(1:surf_lsm_h%ns)           )
    5166        ALLOCATE ( surf_lsm_h%qsws_liq(1:surf_lsm_h%ns)            )
    5167        ALLOCATE ( surf_lsm_h%qsws_veg(1:surf_lsm_h%ns)            )
    5168        ALLOCATE ( surf_lsm_h%rad_net_l(1:surf_lsm_h%ns)           )
    5169        ALLOCATE ( surf_lsm_h%r_a(1:surf_lsm_h%ns)                 )
    5170        ALLOCATE ( surf_lsm_h%r_canopy(1:surf_lsm_h%ns)            )
    5171        ALLOCATE ( surf_lsm_h%r_soil(1:surf_lsm_h%ns)              )
    5172        ALLOCATE ( surf_lsm_h%r_soil_min(1:surf_lsm_h%ns)          )
    5173        ALLOCATE ( surf_lsm_h%r_s(1:surf_lsm_h%ns)                 )
    5174        ALLOCATE ( surf_lsm_h%r_canopy_min(1:surf_lsm_h%ns)        )
    5175        ALLOCATE ( surf_lsm_h%vegetation_surface(1:surf_lsm_h%ns)  )
    5176        ALLOCATE ( surf_lsm_h%water_surface(1:surf_lsm_h%ns)       )
    5177 
    5178        surf_lsm_h%water_surface        = .FALSE.
    5179        surf_lsm_h%pavement_surface     = .FALSE.
    5180        surf_lsm_h%vegetation_surface   = .FALSE.
    5181 
    5182 !
    5183 !--    Set default values
    5184        surf_lsm_h%r_canopy_min = 0.0_wp
    5185 
     5265       DO l = 0, 1
     5266          ALLOCATE ( surf_lsm_h(l)%building_surface(1:surf_lsm_h(l)%ns)    )
     5267          ALLOCATE ( surf_lsm_h(l)%c_liq(1:surf_lsm_h(l)%ns)               )
     5268          ALLOCATE ( surf_lsm_h(l)%c_surface(1:surf_lsm_h(l)%ns)           )
     5269          ALLOCATE ( surf_lsm_h(l)%c_veg(1:surf_lsm_h(l)%ns)               )
     5270          ALLOCATE ( surf_lsm_h(l)%f_sw_in(1:surf_lsm_h(l)%ns)             )
     5271          ALLOCATE ( surf_lsm_h(l)%ghf(1:surf_lsm_h(l)%ns)                 )
     5272          ALLOCATE ( surf_lsm_h(l)%g_d(1:surf_lsm_h(l)%ns)                 )
     5273          ALLOCATE ( surf_lsm_h(l)%lai(1:surf_lsm_h(l)%ns)                 )
     5274          ALLOCATE ( surf_lsm_h(l)%lambda_surface_u(1:surf_lsm_h(l)%ns)    )
     5275          ALLOCATE ( surf_lsm_h(l)%lambda_surface_s(1:surf_lsm_h(l)%ns)    )
     5276          ALLOCATE ( surf_lsm_h(l)%nzt_pavement(1:surf_lsm_h(l)%ns)        )
     5277          ALLOCATE ( surf_lsm_h(l)%pavement_surface(1:surf_lsm_h(l)%ns)    )
     5278          ALLOCATE ( surf_lsm_h(l)%qsws_soil(1:surf_lsm_h(l)%ns)           )
     5279          ALLOCATE ( surf_lsm_h(l)%qsws_liq(1:surf_lsm_h(l)%ns)            )
     5280          ALLOCATE ( surf_lsm_h(l)%qsws_veg(1:surf_lsm_h(l)%ns)            )
     5281          ALLOCATE ( surf_lsm_h(l)%rad_net_l(1:surf_lsm_h(l)%ns)           )
     5282          ALLOCATE ( surf_lsm_h(l)%r_a(1:surf_lsm_h(l)%ns)                 )
     5283          ALLOCATE ( surf_lsm_h(l)%r_canopy(1:surf_lsm_h(l)%ns)            )
     5284          ALLOCATE ( surf_lsm_h(l)%r_soil(1:surf_lsm_h(l)%ns)              )
     5285          ALLOCATE ( surf_lsm_h(l)%r_soil_min(1:surf_lsm_h(l)%ns)          )
     5286          ALLOCATE ( surf_lsm_h(l)%r_s(1:surf_lsm_h(l)%ns)                 )
     5287          ALLOCATE ( surf_lsm_h(l)%r_canopy_min(1:surf_lsm_h(l)%ns)        )
     5288          ALLOCATE ( surf_lsm_h(l)%vegetation_surface(1:surf_lsm_h(l)%ns)  )
     5289          ALLOCATE ( surf_lsm_h(l)%water_surface(1:surf_lsm_h(l)%ns)       )
     5290
     5291          surf_lsm_h(l)%water_surface        = .FALSE.
     5292          surf_lsm_h(l)%pavement_surface     = .FALSE.
     5293          surf_lsm_h(l)%vegetation_surface   = .FALSE.
     5294
     5295!
     5296!--       Set default values
     5297          surf_lsm_h(l)%r_canopy_min = 0.0_wp
     5298       ENDDO
    51865299!
    51875300!--    Vertical surfaces
     
    54105523
    54115524       IF ( horizontal )  THEN
    5412           surf           => surf_lsm_h
    5413 
    5414           surf_m_soil    => m_soil_h
    5415           surf_m_soil_p  => m_soil_h_p
    5416           surf_t_soil    => t_soil_h
    5417           surf_t_soil_p  => t_soil_h_p
    5418           surf_tm_soil_m => tm_soil_h_m
    5419           surf_tt_soil_m => tt_soil_h_m
     5525          surf           => surf_lsm_h(l)
     5526
     5527          surf_m_soil    => m_soil_h(l)
     5528          surf_m_soil_p  => m_soil_h_p(l)
     5529          surf_t_soil    => t_soil_h(l)
     5530          surf_t_soil_p  => t_soil_h_p(l)
     5531          surf_tm_soil_m => tm_soil_h_m(l)
     5532          surf_tt_soil_m => tt_soil_h_m(l)
    54205533       ELSE
    54215534          surf           => surf_lsm_v(l)
     
    58805993          CASE ( 'c_liq*' )
    58815994             IF ( ALLOCATED( c_liq_av ) ) THEN
    5882                 DO  m = 1, surf_lsm_h%ns
    5883                    i   = surf_lsm_h%i(m)
    5884                    j   = surf_lsm_h%j(m)
    5885                    c_liq_av(j,i) = c_liq_av(j,i) + surf_lsm_h%c_liq(m)
     5995                DO  m = 1, surf_lsm_h(0)%ns
     5996                   i   = surf_lsm_h(0)%i(m)
     5997                   j   = surf_lsm_h(0)%j(m)
     5998                   c_liq_av(j,i) = c_liq_av(j,i) + surf_lsm_h(0)%c_liq(m)
    58865999                ENDDO
    58876000             ENDIF
     
    58896002          CASE ( 'c_soil*' )
    58906003             IF ( ALLOCATED( c_soil_av ) ) THEN
    5891                 DO  m = 1, surf_lsm_h%ns
    5892                    i   = surf_lsm_h%i(m)
    5893                    j   = surf_lsm_h%j(m)
    5894                    c_soil_av(j,i) = c_soil_av(j,i) + (1.0 - surf_lsm_h%c_veg(m))
     6004                DO  m = 1, surf_lsm_h(0)%ns
     6005                   i   = surf_lsm_h(0)%i(m)
     6006                   j   = surf_lsm_h(0)%j(m)
     6007                   c_soil_av(j,i) = c_soil_av(j,i) + (1.0 - surf_lsm_h(0)%c_veg(m))
    58956008                ENDDO
    58966009             ENDIF
     
    58986011          CASE ( 'c_veg*' )
    58996012             IF ( ALLOCATED( c_veg_av ) ) THEN
    5900                 DO  m = 1, surf_lsm_h%ns
    5901                    i   = surf_lsm_h%i(m)
    5902                    j   = surf_lsm_h%j(m)
    5903                    c_veg_av(j,i) = c_veg_av(j,i) + surf_lsm_h%c_veg(m)
     6013                DO  m = 1, surf_lsm_h(0)%ns
     6014                   i   = surf_lsm_h(0)%i(m)
     6015                   j   = surf_lsm_h(0)%j(m)
     6016                   c_veg_av(j,i) = c_veg_av(j,i) + surf_lsm_h(0)%c_veg(m)
    59046017                ENDDO
    59056018             ENDIF
     
    59076020          CASE ( 'lai*' )
    59086021             IF ( ALLOCATED( lai_av ) ) THEN
    5909                 DO  m = 1, surf_lsm_h%ns
    5910                    i   = surf_lsm_h%i(m)
    5911                    j   = surf_lsm_h%j(m)
    5912                    lai_av(j,i) = lai_av(j,i) + surf_lsm_h%lai(m)
     6022                DO  m = 1, surf_lsm_h(0)%ns
     6023                   i   = surf_lsm_h(0)%i(m)
     6024                   j   = surf_lsm_h(0)%j(m)
     6025                   lai_av(j,i) = lai_av(j,i) + surf_lsm_h(0)%lai(m)
    59136026                ENDDO
    59146027             ENDIF
     
    59166029          CASE ( 'm_liq*' )
    59176030             IF ( ALLOCATED( m_liq_av ) ) THEN
    5918                 DO  m = 1, surf_lsm_h%ns
    5919                    i   = surf_lsm_h%i(m)
    5920                    j   = surf_lsm_h%j(m)
    5921                    m_liq_av(j,i) = m_liq_av(j,i) + m_liq_h%var_1d(m)
     6031                DO  m = 1, surf_lsm_h(0)%ns
     6032                   i   = surf_lsm_h(0)%i(m)
     6033                   j   = surf_lsm_h(0)%j(m)
     6034                   m_liq_av(j,i) = m_liq_av(j,i) + m_liq_h(0)%var_1d(m)
    59226035                ENDDO
    59236036             ENDIF
     
    59256038          CASE ( 'm_soil' )
    59266039             IF ( ALLOCATED( m_soil_av ) ) THEN
    5927                 DO  m = 1, surf_lsm_h%ns
    5928                    i   = surf_lsm_h%i(m)
    5929                    j   = surf_lsm_h%j(m)
     6040                DO  m = 1, surf_lsm_h(0)%ns
     6041                   i   = surf_lsm_h(0)%i(m)
     6042                   j   = surf_lsm_h(0)%j(m)
    59306043                   DO  k = nzb_soil, nzt_soil
    5931                       m_soil_av(k,j,i) = m_soil_av(k,j,i) + m_soil_h%var_2d(k,m)
     6044                      m_soil_av(k,j,i) = m_soil_av(k,j,i) + m_soil_h(0)%var_2d(k,m)
    59326045                   ENDDO
    59336046                ENDDO
     
    59366049          CASE ( 'qsws_liq*' )
    59376050             IF ( ALLOCATED( qsws_liq_av ) ) THEN
    5938                 DO  m = 1, surf_lsm_h%ns
    5939                    i   = surf_lsm_h%i(m)
    5940                    j   = surf_lsm_h%j(m)
     6051                DO  m = 1, surf_lsm_h(0)%ns
     6052                   i   = surf_lsm_h(0)%i(m)
     6053                   j   = surf_lsm_h(0)%j(m)
    59416054                   qsws_liq_av(j,i) = qsws_liq_av(j,i) +                       &
    5942                                          surf_lsm_h%qsws_liq(m)
     6055                                         surf_lsm_h(0)%qsws_liq(m)
    59436056                ENDDO
    59446057             ENDIF
     
    59466059          CASE ( 'qsws_soil*' )
    59476060             IF ( ALLOCATED( qsws_soil_av ) ) THEN
    5948                 DO  m = 1, surf_lsm_h%ns
    5949                    i   = surf_lsm_h%i(m)
    5950                    j   = surf_lsm_h%j(m)
     6061                DO  m = 1, surf_lsm_h(0)%ns
     6062                   i   = surf_lsm_h(0)%i(m)
     6063                   j   = surf_lsm_h(0)%j(m)
    59516064                   qsws_soil_av(j,i) = qsws_soil_av(j,i) +                     &
    5952                                           surf_lsm_h%qsws_soil(m)
     6065                                          surf_lsm_h(0)%qsws_soil(m)
    59536066                ENDDO
    59546067             ENDIF
     
    59566069          CASE ( 'qsws_veg*' )
    59576070             IF ( ALLOCATED(qsws_veg_av ) ) THEN
    5958                 DO  m = 1, surf_lsm_h%ns
    5959                    i   = surf_lsm_h%i(m)
    5960                    j   = surf_lsm_h%j(m)
     6071                DO  m = 1, surf_lsm_h(0)%ns
     6072                   i   = surf_lsm_h(0)%i(m)
     6073                   j   = surf_lsm_h(0)%j(m)
    59616074                   qsws_veg_av(j,i) = qsws_veg_av(j,i) +                       &
    5962                                          surf_lsm_h%qsws_veg(m)
     6075                                         surf_lsm_h(0)%qsws_veg(m)
    59636076                ENDDO
    59646077             ENDIF
     
    59666079          CASE ( 'r_s*' )
    59676080             IF ( ALLOCATED( r_s_av) ) THEN
    5968                 DO  m = 1, surf_lsm_h%ns
    5969                    i   = surf_lsm_h%i(m)
    5970                    j   = surf_lsm_h%j(m)
    5971                    r_s_av(j,i) = r_s_av(j,i) + surf_lsm_h%r_s(m)
     6081                DO  m = 1, surf_lsm_h(0)%ns
     6082                   i   = surf_lsm_h(0)%i(m)
     6083                   j   = surf_lsm_h(0)%j(m)
     6084                   r_s_av(j,i) = r_s_av(j,i) + surf_lsm_h(0)%r_s(m)
    59726085                ENDDO
    59736086             ENDIF
     
    59756088          CASE ( 't_soil' )
    59766089             IF ( ALLOCATED( t_soil_av ) ) THEN
    5977                 DO  m = 1, surf_lsm_h%ns
    5978                    i   = surf_lsm_h%i(m)
    5979                    j   = surf_lsm_h%j(m)
     6090                DO  m = 1, surf_lsm_h(0)%ns
     6091                   i   = surf_lsm_h(0)%i(m)
     6092                   j   = surf_lsm_h(0)%j(m)
    59806093                   DO  k = nzb_soil, nzt_soil
    5981                       t_soil_av(k,j,i) = t_soil_av(k,j,i) + t_soil_h%var_2d(k,m)
     6094                      t_soil_av(k,j,i) = t_soil_av(k,j,i) + t_soil_h(0)%var_2d(k,m)
    59826095                   ENDDO
    59836096                ENDDO
     
    61966309       CASE ( 'c_liq*_xy' )        ! 2d-array
    61976310          IF ( av == 0 )  THEN
    6198              DO  m = 1, surf_lsm_h%ns
    6199                 i                   = surf_lsm_h%i(m)
    6200                 j                   = surf_lsm_h%j(m)
    6201                 local_pf(i,j,nzb+1) = surf_lsm_h%c_liq(m) * surf_lsm_h%c_veg(m)
     6311             DO  m = 1, surf_lsm_h(0)%ns
     6312                i                   = surf_lsm_h(0)%i(m)
     6313                j                   = surf_lsm_h(0)%j(m)
     6314                local_pf(i,j,nzb+1) = surf_lsm_h(0)%c_liq(m) * surf_lsm_h(0)%c_veg(m)
    62026315             ENDDO
    62036316          ELSE
     
    62186331       CASE ( 'c_soil*_xy' )        ! 2d-array
    62196332          IF ( av == 0 )  THEN
    6220              DO  m = 1, surf_lsm_h%ns
    6221                 i                   = surf_lsm_h%i(m)
    6222                 j                   = surf_lsm_h%j(m)
    6223                 local_pf(i,j,nzb+1) = 1.0_wp - surf_lsm_h%c_veg(m)
     6333             DO  m = 1, surf_lsm_h(0)%ns
     6334                i                   = surf_lsm_h(0)%i(m)
     6335                j                   = surf_lsm_h(0)%j(m)
     6336                local_pf(i,j,nzb+1) = 1.0_wp - surf_lsm_h(0)%c_veg(m)
    62246337             ENDDO
    62256338          ELSE
     
    62406353       CASE ( 'c_veg*_xy' )        ! 2d-array
    62416354          IF ( av == 0 )  THEN
    6242              DO  m = 1, surf_lsm_h%ns
    6243                 i                   = surf_lsm_h%i(m)
    6244                 j                   = surf_lsm_h%j(m)
    6245                 local_pf(i,j,nzb+1) = surf_lsm_h%c_veg(m)
     6355             DO  m = 1, surf_lsm_h(0)%ns
     6356                i                   = surf_lsm_h(0)%i(m)
     6357                j                   = surf_lsm_h(0)%j(m)
     6358                local_pf(i,j,nzb+1) = surf_lsm_h(0)%c_veg(m)
    62466359             ENDDO
    62476360          ELSE
     
    62626375       CASE ( 'lai*_xy' )        ! 2d-array
    62636376          IF ( av == 0 )  THEN
    6264              DO  m = 1, surf_lsm_h%ns
    6265                 i                   = surf_lsm_h%i(m)
    6266                 j                   = surf_lsm_h%j(m)
    6267                 local_pf(i,j,nzb+1) = surf_lsm_h%lai(m)
     6377             DO  m = 1, surf_lsm_h(0)%ns
     6378                i                   = surf_lsm_h(0)%i(m)
     6379                j                   = surf_lsm_h(0)%j(m)
     6380                local_pf(i,j,nzb+1) = surf_lsm_h(0)%lai(m)
    62686381             ENDDO
    62696382          ELSE
     
    62846397       CASE ( 'm_liq*_xy' )        ! 2d-array
    62856398          IF ( av == 0 )  THEN
    6286              DO  m = 1, surf_lsm_h%ns
    6287                 i                   = surf_lsm_h%i(m)
    6288                 j                   = surf_lsm_h%j(m)
    6289                 local_pf(i,j,nzb+1) = m_liq_h%var_1d(m)
     6399             DO  m = 1, surf_lsm_h(0)%ns
     6400                i                   = surf_lsm_h(0)%i(m)
     6401                j                   = surf_lsm_h(0)%j(m)
     6402                local_pf(i,j,nzb+1) = m_liq_h(0)%var_1d(m)
    62906403             ENDDO
    62916404          ELSE
     
    63066419       CASE ( 'm_soil_xy', 'm_soil_xz', 'm_soil_yz' )
    63076420          IF ( av == 0 )  THEN
    6308              DO  m = 1, surf_lsm_h%ns
    6309                 i   = surf_lsm_h%i(m)
    6310                 j   = surf_lsm_h%j(m)
     6421             DO  m = 1, surf_lsm_h(0)%ns
     6422                i   = surf_lsm_h(0)%i(m)
     6423                j   = surf_lsm_h(0)%j(m)
    63116424                DO k = nzb_soil, nzt_soil
    6312                    local_pf(i,j,k) = m_soil_h%var_2d(k,m)
     6425                   local_pf(i,j,k) = m_soil_h(0)%var_2d(k,m)
    63136426                ENDDO
    63146427             ENDDO
     
    63346447       CASE ( 'qsws_liq*_xy' )        ! 2d-array
    63356448          IF ( av == 0 ) THEN
    6336              DO  m = 1, surf_lsm_h%ns
    6337                 i                   = surf_lsm_h%i(m)
    6338                 j                   = surf_lsm_h%j(m)
    6339                 local_pf(i,j,nzb+1) = surf_lsm_h%qsws_liq(m)
     6449             DO  m = 1, surf_lsm_h(0)%ns
     6450                i                   = surf_lsm_h(0)%i(m)
     6451                j                   = surf_lsm_h(0)%j(m)
     6452                local_pf(i,j,nzb+1) = surf_lsm_h(0)%qsws_liq(m)
    63406453             ENDDO
    63416454          ELSE
     
    63566469       CASE ( 'qsws_soil*_xy' )        ! 2d-array
    63576470          IF ( av == 0 ) THEN
    6358              DO  m = 1, surf_lsm_h%ns
    6359                 i                   = surf_lsm_h%i(m)
    6360                 j                   = surf_lsm_h%j(m)
    6361                 local_pf(i,j,nzb+1) =  surf_lsm_h%qsws_soil(m)
     6471             DO  m = 1, surf_lsm_h(0)%ns
     6472                i                   = surf_lsm_h(0)%i(m)
     6473                j                   = surf_lsm_h(0)%j(m)
     6474                local_pf(i,j,nzb+1) =  surf_lsm_h(0)%qsws_soil(m)
    63626475             ENDDO
    63636476          ELSE
     
    63786491       CASE ( 'qsws_veg*_xy' )        ! 2d-array
    63796492          IF ( av == 0 ) THEN
    6380              DO  m = 1, surf_lsm_h%ns
    6381                 i                   = surf_lsm_h%i(m)
    6382                 j                   = surf_lsm_h%j(m)
    6383                 local_pf(i,j,nzb+1) =  surf_lsm_h%qsws_veg(m)
     6493             DO  m = 1, surf_lsm_h(0)%ns
     6494                i                   = surf_lsm_h(0)%i(m)
     6495                j                   = surf_lsm_h(0)%j(m)
     6496                local_pf(i,j,nzb+1) =  surf_lsm_h(0)%qsws_veg(m)
    63846497             ENDDO
    63856498          ELSE
     
    64016514       CASE ( 'r_s*_xy' )        ! 2d-array
    64026515          IF ( av == 0 )  THEN
    6403              DO  m = 1, surf_lsm_h%ns
    6404                 i                   = surf_lsm_h%i(m)
    6405                 j                   = surf_lsm_h%j(m)
    6406                 local_pf(i,j,nzb+1) = surf_lsm_h%r_s(m)
     6516             DO  m = 1, surf_lsm_h(0)%ns
     6517                i                   = surf_lsm_h(0)%i(m)
     6518                j                   = surf_lsm_h(0)%j(m)
     6519                local_pf(i,j,nzb+1) = surf_lsm_h(0)%r_s(m)
    64076520             ENDDO
    64086521          ELSE
     
    64236536       CASE ( 't_soil_xy', 't_soil_xz', 't_soil_yz' )
    64246537          IF ( av == 0 )  THEN
    6425              DO  m = 1, surf_lsm_h%ns
    6426                 i   = surf_lsm_h%i(m)
    6427                 j   = surf_lsm_h%j(m)
     6538             DO  m = 1, surf_lsm_h(0)%ns
     6539                i   = surf_lsm_h(0)%i(m)
     6540                j   = surf_lsm_h(0)%j(m)
    64286541                DO k = nzb_soil, nzt_soil
    6429                    local_pf(i,j,k) = t_soil_h%var_2d(k,m)
     6542                   local_pf(i,j,k) = t_soil_h(0)%var_2d(k,m)
    64306543                ENDDO
    64316544             ENDDO
     
    64986611
    64996612         IF ( av == 0 )  THEN
    6500             DO  m = 1, surf_lsm_h%ns
    6501                 i   = surf_lsm_h%i(m)
    6502                 j   = surf_lsm_h%j(m)
     6613            DO  m = 1, surf_lsm_h(0)%ns
     6614                i   = surf_lsm_h(0)%i(m)
     6615                j   = surf_lsm_h(0)%j(m)
    65036616                DO  k = nzb_soil, nzt_soil
    6504                    local_pf(i,j,k) = m_soil_h%var_2d(k,m)
     6617                   local_pf(i,j,k) = m_soil_h(0)%var_2d(k,m)
    65056618                ENDDO
    65066619            ENDDO
     
    65226635
    65236636         IF ( av == 0 )  THEN
    6524             DO  m = 1, surf_lsm_h%ns
    6525                i   = surf_lsm_h%i(m)
    6526                j   = surf_lsm_h%j(m)
     6637            DO  m = 1, surf_lsm_h(0)%ns
     6638               i   = surf_lsm_h(0)%i(m)
     6639               j   = surf_lsm_h(0)%j(m)
    65276640               DO  k = nzb_soil, nzt_soil
    6528                   local_pf(i,j,k) = t_soil_h%var_2d(k,m)
     6641                  local_pf(i,j,k) = t_soil_h(0)%var_2d(k,m)
    65296642               ENDDO
    65306643            ENDDO
     
    65776690
    65786691       CALL wrd_write_string( 'ns_h_on_file_lsm' )
    6579        WRITE ( 14 )  surf_lsm_h%ns
     6692       WRITE ( 14 )  surf_lsm_h(0:1)%ns
    65806693
    65816694       CALL wrd_write_string( 'ns_v_on_file_lsm' )
     
    66336746       ENDIF
    66346747
    6635        CALL wrd_write_string( 'lsm_start_index_h' )
    6636        WRITE ( 14 )  surf_lsm_h%start_index
    6637 
    6638        CALL wrd_write_string( 'lsm_end_index_h' )
    6639        WRITE ( 14 )  surf_lsm_h%end_index
    6640 
    6641        CALL wrd_write_string( 't_soil_h' )
    6642        WRITE ( 14 )  t_soil_h%var_2d
     6748       DO  l = 0, 1
     6749          CALL wrd_write_string( 'lsm_start_index_h' )
     6750          WRITE ( 14 )  surf_lsm_h(l)%start_index
     6751
     6752          CALL wrd_write_string( 'lsm_end_index_h' )
     6753          WRITE ( 14 )  surf_lsm_h(l)%end_index
     6754
     6755          WRITE( dum, '(I1)')  l
     6756          CALL wrd_write_string( 't_soil_h(' // dum // ')' )
     6757          WRITE ( 14 )  t_soil_h(l)%var_2d
     6758       ENDDO
    66436759
    66446760       DO  l = 0, 3
     
    66566772       ENDDO
    66576773
    6658        CALL wrd_write_string( 'lsm_start_index_h' )
    6659        WRITE ( 14 )  surf_lsm_h%start_index
    6660 
    6661        CALL wrd_write_string( 'lsm_end_index_h' )
    6662        WRITE ( 14 )  surf_lsm_h%end_index
    6663 
    6664        CALL wrd_write_string( 'm_soil_h' )
    6665        WRITE ( 14 )  m_soil_h%var_2d
     6774       DO  l = 0, 1
     6775          CALL wrd_write_string( 'lsm_start_index_h' )
     6776          WRITE ( 14 )  surf_lsm_h(l)%start_index
     6777
     6778          CALL wrd_write_string( 'lsm_end_index_h' )
     6779          WRITE ( 14 )  surf_lsm_h(l)%end_index
     6780
     6781          WRITE( dum, '(I1)')  l
     6782          CALL wrd_write_string( 'm_soil_h(' // dum // ')' )
     6783          WRITE ( 14 )  m_soil_h(l)%var_2d
     6784       ENDDO
    66666785
    66676786       DO  l = 0, 3
     
    66796798       ENDDO
    66806799
    6681        CALL wrd_write_string( 'lsm_start_index_h' )
    6682        WRITE ( 14 )  surf_lsm_h%start_index
    6683 
    6684        CALL wrd_write_string( 'lsm_end_index_h' )
    6685        WRITE ( 14 )  surf_lsm_h%end_index
    6686 
    6687        CALL wrd_write_string( 'm_liq_h' )
    6688        WRITE ( 14 )  m_liq_h%var_1d
     6800       DO  l = 0, 1
     6801          CALL wrd_write_string( 'lsm_start_index_h' )
     6802          WRITE ( 14 )  surf_lsm_h(l)%start_index
     6803
     6804          CALL wrd_write_string( 'lsm_end_index_h' )
     6805          WRITE ( 14 )  surf_lsm_h(l)%end_index
     6806
     6807          WRITE( dum, '(I1)')  l
     6808          CALL wrd_write_string( 'm_liq_h(' // dum // ')' )
     6809          WRITE ( 14 )  m_liq_h(l)%var_1d
     6810       ENDDO
    66896811
    66906812       DO  l = 0, 3
     
    67026824       ENDDO
    67036825
    6704        CALL wrd_write_string( 'lsm_start_index_h' )
    6705        WRITE ( 14 )  surf_lsm_h%start_index
    6706 
    6707        CALL wrd_write_string( 'lsm_end_index_h' )
    6708        WRITE ( 14 )  surf_lsm_h%end_index
    6709 
    6710        CALL wrd_write_string( 't_surface_h' )
    6711        WRITE ( 14 )  t_surface_h%var_1d
     6826       DO  l = 0, 1
     6827          CALL wrd_write_string( 'lsm_start_index_h' )
     6828          WRITE ( 14 )  surf_lsm_h(l)%start_index
     6829
     6830          CALL wrd_write_string( 'lsm_end_index_h' )
     6831          WRITE ( 14 )  surf_lsm_h(l)%end_index
     6832
     6833          WRITE( dum, '(I1)')  l
     6834          CALL wrd_write_string( 't_surface_h(' // dum // ')' )
     6835          WRITE ( 14 )  t_surface_h(l)%var_1d
     6836       ENDDO
    67126837
    67136838       DO  l = 0, 3
     
    67386863       IF ( ALLOCATED( t_soil_av ) )  CALL wrd_mpi_io( 't_soil_av', t_soil_av, nzb_soil, nzt_soil )
    67396864
    6740        CALL rd_mpi_io_surface_filetypes( surf_lsm_h%start_index, surf_lsm_h%end_index,             &
    6741                                          surface_data_to_write, global_start_index )
    6742 
    6743        CALL wrd_mpi_io( 'lsm_start_index_h',  surf_lsm_h%start_index )
    6744        CALL wrd_mpi_io( 'lsm_end_index_h',  surf_lsm_h%end_index )
    6745        CALL wrd_mpi_io( 'lsm_global_start_index_h', global_start_index )
    6746        CALL wrd_mpi_io_surface( 't_soil_h', t_soil_h%var_2d )
    6747        CALL wrd_mpi_io_surface( 'm_soil_h',  m_soil_h%var_2d )
    6748        CALL wrd_mpi_io_surface( 'm_liq_h', m_liq_h%var_1d )
    6749        CALL wrd_mpi_io_surface( 't_surface_h', t_surface_h%var_1d )
     6865       DO  l = 0, 1
     6866
     6867          WRITE( dum, '(I1)' )  l
     6868
     6869          CALL rd_mpi_io_surface_filetypes( surf_lsm_h(l)%start_index, surf_lsm_h(l)%end_index,             &
     6870                                            surface_data_to_write, global_start_index )
     6871
     6872          CALL wrd_mpi_io( 'lsm_start_index_h_' // dum,  surf_lsm_h(l)%start_index )
     6873          CALL wrd_mpi_io( 'lsm_end_index_h_' // dum,  surf_lsm_h(l)%end_index )
     6874          CALL wrd_mpi_io( 'lsm_global_start_index_h_' // dum, global_start_index )
     6875
     6876          IF ( .NOT. surface_data_to_write )  CYCLE
     6877
     6878          CALL wrd_mpi_io_surface( 't_soil_h(' // dum // ')', t_soil_h(l)%var_2d )
     6879          CALL wrd_mpi_io_surface( 'm_soil_h(' // dum // ')',  m_soil_h(l)%var_2d )
     6880          CALL wrd_mpi_io_surface( 'm_liq_h(' // dum // ')', m_liq_h(l)%var_1d )
     6881          CALL wrd_mpi_io_surface( 't_surface_h(' // dum // ')', t_surface_h(l)%var_1d )
     6882       ENDDO
    67506883
    67516884       DO  l = 0, 3
     
    67966929    INTEGER(iwp) ::  k                 !<
    67976930    INTEGER(iwp) ::  l                 !< running index surface orientation
    6798     INTEGER(iwp) ::  ns_h_on_file_lsm  !< number of horizontal surface elements (natural type) on file
    67996931    INTEGER(iwp) ::  nxlc              !<
    68006932    INTEGER(iwp) ::  nxlf              !<
     
    68106942    INTEGER(iwp) ::  nys_on_file       !< index of south boundary on former local domain
    68116943
     6944    INTEGER(iwp) ::  ns_h_on_file_lsm(0:1) !< number of horizontal surface elements (natural type) on file
    68126945    INTEGER(iwp) ::  ns_v_on_file_lsm(0:3) !< number of vertical surface elements (natural type) on file
    68136946
     
    68216954    REAL(wp), DIMENSION(nzb_soil:nzt_soil,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d   !<
    68226955
    6823     TYPE(surf_type_lsm), SAVE :: tmp_walltype_h_1d   !< temporary 1D array containing the respective surface variable stored on file, horizontal surfaces
    6824     TYPE(surf_type_lsm), SAVE :: tmp_walltype_h_2d   !< temporary 2D array containing the respective surface variable stored on file, horizontal surfaces
    6825     TYPE(surf_type_lsm), SAVE :: tmp_walltype_h_2d2  !< temporary 2D array containing the respective surface variable stored on file, horizontal surfaces
     6956    TYPE(surf_type_lsm), DIMENSION(0:1), SAVE :: tmp_walltype_h_1d   !< temporary 1D array containing the respective surface variable stored on file, horizontal surfaces
     6957    TYPE(surf_type_lsm), DIMENSION(0:1), SAVE :: tmp_walltype_h_2d   !< temporary 2D array containing the respective surface variable stored on file, horizontal surfaces
     6958    TYPE(surf_type_lsm), DIMENSION(0:1), SAVE :: tmp_walltype_h_2d2  !< temporary 2D array containing the respective surface variable stored on file, horizontal surfaces
    68266959
    68276960    TYPE(surf_type_lsm), DIMENSION(0:3), SAVE :: tmp_walltype_v_1d   !< temporary 1D array containing the respective surface variable stored on file, vertical surfaces
     
    68396972             READ ( 13 ) ns_h_on_file_lsm
    68406973
    6841              IF ( ALLOCATED( tmp_walltype_h_1d%var_1d ) )                      &
    6842                 DEALLOCATE( tmp_walltype_h_1d%var_1d )
    6843              IF ( ALLOCATED( tmp_walltype_h_2d%var_2d ) )                      &
    6844                 DEALLOCATE( tmp_walltype_h_2d%var_2d )
    6845              IF ( ALLOCATED( tmp_walltype_h_2d2%var_2d ) )                     &
    6846                 DEALLOCATE( tmp_walltype_h_2d2%var_2d )
     6974             DO l = 0, 1
     6975                IF ( ALLOCATED( tmp_walltype_h_1d(l)%var_1d ) )                   &
     6976                   DEALLOCATE( tmp_walltype_h_1d(l)%var_1d )
     6977                IF ( ALLOCATED( tmp_walltype_h_2d(l)%var_2d ) )                   &
     6978                   DEALLOCATE( tmp_walltype_h_2d(l)%var_2d )
     6979                IF ( ALLOCATED( tmp_walltype_h_2d2(l)%var_2d ) )                  &
     6980                   DEALLOCATE( tmp_walltype_h_2d2(l)%var_2d )
     6981             ENDDO
    68476982
    68486983!
    68496984!--          Allocate temporary arrays to store surface data
    6850              ALLOCATE( tmp_walltype_h_1d%var_1d(1:ns_h_on_file_lsm) )
    6851              ALLOCATE( tmp_walltype_h_2d%var_2d(nzb_soil:nzt_soil+1,           &
    6852                                                 1:ns_h_on_file_lsm) )
    6853              ALLOCATE( tmp_walltype_h_2d2%var_2d(nzb_soil:nzt_soil,            &
    6854                        1:ns_h_on_file_lsm)  )
     6985             DO l = 0, 1
     6986                ALLOCATE( tmp_walltype_h_1d(l)                                    &
     6987                             %var_1d(1:ns_h_on_file_lsm(l)) )
     6988                ALLOCATE( tmp_walltype_h_2d(l)                                    &
     6989                             %var_2d(nzb_soil:nzt_soil+1, 1:ns_h_on_file_lsm(l)) )
     6990                ALLOCATE( tmp_walltype_h_2d2(l)                                   &
     6991                             %var_2d(nzb_soil:nzt_soil,1:ns_h_on_file_lsm(l))  )
     6992             ENDDO
    68556993
    68566994          ENDIF
     
    69917129            ENDIF
    69927130
    6993        CASE ( 't_soil_h' )
     7131       CASE ( 't_soil_h(0)' )
    69947132
    69957133          IF ( k == 1 )  THEN
    6996              IF ( .NOT.  ALLOCATED( t_soil_h%var_2d ) )                        &
    6997                 ALLOCATE( t_soil_h%var_2d(nzb_soil:nzt_soil+1,                 &
    6998                                           1:surf_lsm_h%ns) )
    6999              READ ( 13 )  tmp_walltype_h_2d%var_2d
     7134             IF ( .NOT.  ALLOCATED( t_soil_h(0)%var_2d ) )                     &
     7135                ALLOCATE( t_soil_h(0)%var_2d(nzb_soil:nzt_soil+1,              &
     7136                                          1:surf_lsm_h(0)%ns) )
     7137             READ ( 13 )  tmp_walltype_h_2d(0)%var_2d
    70007138          ENDIF
    70017139          CALL surface_restore_elements(                                       &
    7002                                      t_soil_h%var_2d,                          &
    7003                                      tmp_walltype_h_2d%var_2d,                 &
    7004                                      surf_lsm_h%start_index,                   &
     7140                                     t_soil_h(0)%var_2d,                       &
     7141                                     tmp_walltype_h_2d(0)%var_2d,              &
     7142                                     surf_lsm_h(0)%start_index,                &
     7143                                     start_index_on_file,                      &
     7144                                     end_index_on_file,                        &
     7145                                     nxlc, nysc,                               &
     7146                                     nxlf, nxrf, nysf, nynf,                   &
     7147                                     nys_on_file, nyn_on_file,                 &
     7148                                     nxl_on_file,nxr_on_file )
     7149
     7150       CASE ( 't_soil_h(1)' )
     7151
     7152          IF ( k == 1 )  THEN
     7153             IF ( .NOT.  ALLOCATED( t_soil_h(1)%var_2d ) )                     &
     7154                ALLOCATE( t_soil_h(1)%var_2d(nzb_soil:nzt_soil+1,              &
     7155                                          1:surf_lsm_h(1)%ns) )
     7156             write(9,*) surf_lsm_h(1)%ns
     7157             READ ( 13 )  tmp_walltype_h_2d(1)%var_2d
     7158          ENDIF
     7159          CALL surface_restore_elements(                                       &
     7160                                     t_soil_h(1)%var_2d,                       &
     7161                                     tmp_walltype_h_2d(1)%var_2d,              &
     7162                                     surf_lsm_h(1)%start_index,                &
    70057163                                     start_index_on_file,                      &
    70067164                                     end_index_on_file,                        &
     
    70867244                                  nxl_on_file,nxr_on_file )
    70877245
    7088        CASE ( 'm_soil_h' )
     7246       CASE ( 'm_soil_h(0)' )
    70897247
    70907248          IF ( k == 1 )  THEN
    7091              IF ( .NOT.  ALLOCATED( m_soil_h%var_2d ) )                        &
    7092                 ALLOCATE( m_soil_h%var_2d(nzb_soil:nzt_soil+1,                 &
    7093                                           1:surf_lsm_h%ns) )
    7094              READ ( 13 )  tmp_walltype_h_2d2%var_2d
     7249             IF ( .NOT.  ALLOCATED( m_soil_h(0)%var_2d ) )                     &
     7250                ALLOCATE( m_soil_h(0)%var_2d(nzb_soil:nzt_soil+1,              &
     7251                                          1:surf_lsm_h(0)%ns) )
     7252             READ ( 13 )  tmp_walltype_h_2d2(0)%var_2d
    70957253          ENDIF
    70967254          CALL surface_restore_elements(                                       &
    7097                                     m_soil_h%var_2d,                           &
    7098                                     tmp_walltype_h_2d2%var_2d,                 &
    7099                                     surf_lsm_h%start_index,                    &
     7255                                    m_soil_h(0)%var_2d,                        &
     7256                                    tmp_walltype_h_2d2(0)%var_2d,              &
     7257                                    surf_lsm_h(0)%start_index,                 &
     7258                                    start_index_on_file,                       &
     7259                                    end_index_on_file,                         &
     7260                                    nxlc, nysc,                                &
     7261                                    nxlf, nxrf, nysf, nynf,                    &
     7262                                    nys_on_file, nyn_on_file,                  &
     7263                                    nxl_on_file,nxr_on_file )
     7264
     7265       CASE ( 'm_soil_h(1)' )
     7266
     7267          IF ( k == 1 )  THEN
     7268             IF ( .NOT.  ALLOCATED( m_soil_h(1)%var_2d ) )                     &
     7269                ALLOCATE( m_soil_h(1)%var_2d(nzb_soil:nzt_soil+1,              &
     7270                                          1:surf_lsm_h(1)%ns) )
     7271             READ ( 13 )  tmp_walltype_h_2d2(1)%var_2d
     7272          ENDIF
     7273          CALL surface_restore_elements(                                       &
     7274                                    m_soil_h(1)%var_2d,                        &
     7275                                    tmp_walltype_h_2d2(1)%var_2d,              &
     7276                                    surf_lsm_h(1)%start_index,                 &
    71007277                                    start_index_on_file,                       &
    71017278                                    end_index_on_file,                         &
     
    71847361
    71857362
    7186        CASE ( 'm_liq_h' )
     7363       CASE ( 'm_liq_h(0)' )
    71877364
    71887365          IF ( k == 1 )  THEN
    7189              IF ( .NOT.  ALLOCATED( m_liq_h%var_1d ) )                         &
    7190                 ALLOCATE( m_liq_h%var_1d(1:surf_lsm_h%ns) )
    7191              READ ( 13 )  tmp_walltype_h_1d%var_1d
     7366             IF ( .NOT.  ALLOCATED( m_liq_h(0)%var_1d ) )                      &
     7367                ALLOCATE( m_liq_h(0)%var_1d(1:surf_lsm_h(0)%ns) )
     7368             READ ( 13 )  tmp_walltype_h_1d(0)%var_1d
    71927369          ENDIF
    71937370          CALL surface_restore_elements(                                       &
    7194                                      m_liq_h%var_1d,                           &
    7195                                      tmp_walltype_h_1d%var_1d,                 &
    7196                                      surf_lsm_h%start_index,                   &
     7371                                     m_liq_h(0)%var_1d,                        &
     7372                                     tmp_walltype_h_1d(0)%var_1d,              &
     7373                                     surf_lsm_h(0)%start_index,                &
     7374                                     start_index_on_file,                      &
     7375                                     end_index_on_file,                        &
     7376                                     nxlc, nysc,                               &
     7377                                     nxlf, nxrf, nysf, nynf,                   &
     7378                                     nys_on_file, nyn_on_file,                 &
     7379                                     nxl_on_file,nxr_on_file )
     7380
     7381       CASE ( 'm_liq_h(1)' )
     7382
     7383          IF ( k == 1 )  THEN
     7384             IF ( .NOT.  ALLOCATED( m_liq_h(1)%var_1d ) )                      &
     7385                ALLOCATE( m_liq_h(1)%var_1d(1:surf_lsm_h(1)%ns) )
     7386             READ ( 13 )  tmp_walltype_h_1d(1)%var_1d
     7387          ENDIF
     7388          CALL surface_restore_elements(                                       &
     7389                                     m_liq_h(1)%var_1d,                        &
     7390                                     tmp_walltype_h_1d(1)%var_1d,              &
     7391                                     surf_lsm_h(1)%start_index,                &
    71977392                                     start_index_on_file,                      &
    71987393                                     end_index_on_file,                        &
     
    72787473
    72797474
    7280        CASE ( 't_surface_h' )
     7475       CASE ( 't_surface_h(0)' )
    72817476
    72827477          IF ( k == 1 )  THEN
    7283              IF ( .NOT.  ALLOCATED( t_surface_h%var_1d ) )                     &
    7284                 ALLOCATE( t_surface_h%var_1d(1:surf_lsm_h%ns) )
    7285              READ ( 13 )  tmp_walltype_h_1d%var_1d
     7478             IF ( .NOT.  ALLOCATED( t_surface_h(0)%var_1d ) )                  &
     7479                ALLOCATE( t_surface_h(0)%var_1d(1:surf_lsm_h(0)%ns) )
     7480             READ ( 13 )  tmp_walltype_h_1d(0)%var_1d
    72867481          ENDIF
    72877482          CALL surface_restore_elements(                                       &
    7288                                      t_surface_h%var_1d,                       &
    7289                                      tmp_walltype_h_1d%var_1d,                 &
    7290                                      surf_lsm_h%start_index,                   &
     7483                                     t_surface_h(0)%var_1d,                    &
     7484                                     tmp_walltype_h_1d(0)%var_1d,              &
     7485                                     surf_lsm_h(0)%start_index,                &
     7486                                     start_index_on_file,                      &
     7487                                     end_index_on_file,                        &
     7488                                     nxlc, nysc,                               &
     7489                                     nxlf, nxrf, nysf, nynf,                   &
     7490                                     nys_on_file, nyn_on_file,                 &
     7491                                     nxl_on_file,nxr_on_file )
     7492
     7493       CASE ( 't_surface_h(1)' )
     7494
     7495          IF ( k == 1 )  THEN
     7496             IF ( .NOT.  ALLOCATED( t_surface_h(1)%var_1d ) )                  &
     7497                ALLOCATE( t_surface_h(1)%var_1d(1:surf_lsm_h(1)%ns) )
     7498             READ ( 13 )  tmp_walltype_h_1d(1)%var_1d
     7499          ENDIF
     7500          CALL surface_restore_elements(                                       &
     7501                                     t_surface_h(1)%var_1d,                    &
     7502                                     tmp_walltype_h_1d(1)%var_1d,              &
     7503                                     surf_lsm_h(1)%start_index,                &
    72917504                                     start_index_on_file,                      &
    72927505                                     end_index_on_file,                        &
     
    74667679    ENDIF
    74677680
    7468     CALL rrd_mpi_io( 'lsm_start_index_h',  surf_lsm_h%start_index )
    7469     CALL rrd_mpi_io( 'lsm_end_index_h',  surf_lsm_h%end_index )
    7470     CALL rrd_mpi_io( 'lsm_global_start_index_h', global_start )
    7471 
    7472     CALL rd_mpi_io_surface_filetypes ( surf_lsm_h%start_index, surf_lsm_h%end_index, ldum,         &
    7473                                        global_start )
    7474 
    7475     CALL rrd_mpi_io_surface( 't_soil_h', t_soil_h%var_2d )
    7476     CALL rrd_mpi_io_surface( 'm_soil_h',  m_soil_h%var_2d )
    7477     CALL rrd_mpi_io_surface( 'm_liq_h', m_liq_h%var_1d )
    7478     CALL rrd_mpi_io_surface( 't_surface_h', t_surface_h%var_1d )
    7479 
     7681    DO l = 0, 1
     7682
     7683       WRITE( dum, '(I1)')  l
     7684
     7685       CALL rrd_mpi_io( 'lsm_start_index_h_' // dum,  surf_lsm_h(l)%start_index )
     7686       CALL rrd_mpi_io( 'lsm_end_index_h_' // dum,  surf_lsm_h(l)%end_index )
     7687       CALL rrd_mpi_io( 'lsm_global_start_index_h_' // dum, global_start )
     7688
     7689       CALL rd_mpi_io_surface_filetypes ( surf_lsm_h(l)%start_index, surf_lsm_h(l)%end_index, ldum,         &
     7690                                          global_start )
     7691
     7692       CALL rrd_mpi_io_surface( 't_soil_h(' // dum // ')', t_soil_h(l)%var_2d )
     7693       CALL rrd_mpi_io_surface( 'm_soil_h(' // dum // ')',  m_soil_h(l)%var_2d )
     7694       CALL rrd_mpi_io_surface( 'm_liq_h(' // dum // ')', m_liq_h(l)%var_1d )
     7695       CALL rrd_mpi_io_surface( 't_surface_h(' // dum // ')', t_surface_h(l)%var_1d )
     7696    ENDDO
    74807697
    74817698    DO  l = 0, 3
     
    75317748!       REAL(wp) :: re_0 !< near-surface roughness Reynolds number
    75327749
    7533        DO  m = 1, surf_lsm_h%ns
    7534 
    7535           i   = surf_lsm_h%i(m)
    7536           j   = surf_lsm_h%j(m)
    7537 
    7538           IF ( surf_lsm_h%water_surface(m) )  THEN
     7750       DO  m = 1, surf_lsm_h(0)%ns
     7751!--       only upward facin horizontal surfaces are considered for water surface processing
     7752          i   = surf_lsm_h(0)%i(m)
     7753          j   = surf_lsm_h(0)%j(m)
     7754
     7755          IF ( surf_lsm_h(0)%water_surface(m) )  THEN
    75397756
    75407757!
     
    75637780!--          Charnock relation is used. Add a security factor of 1E-8 to avoid
    75647781!--          divisions by zero.
    7565              surf_lsm_h%z0(m)  = ( 0.11_wp * molecular_viscosity /             &
    7566                                  ( surf_lsm_h%us(m) + 1E-8_wp ) )              &
    7567                                + ( alpha_ch * surf_lsm_h%us(m)**2 / g )
    7568 
    7569              surf_lsm_h%z0h(m) = 0.40_wp * molecular_viscosity /               &
    7570                                  ( surf_lsm_h%us(m) + 1E-8_wp )
    7571              surf_lsm_h%z0q(m) = 0.62_wp * molecular_viscosity /               &
    7572                                  ( surf_lsm_h%us(m) + 1E-8_wp )
    7573 
    7574 
    7575              IF ( surf_lsm_h%z0(m) > 0.1_wp * surf_lsm_h%z_mo(m) )  THEN
    7576                 surf_lsm_h%z0(m) = 0.1_wp * surf_lsm_h%z_mo(m)
     7782             surf_lsm_h(0)%z0(m)  = ( 0.11_wp * molecular_viscosity /             &
     7783                                 ( surf_lsm_h(0)%us(m) + 1E-8_wp ) )              &
     7784                               + ( alpha_ch * surf_lsm_h(0)%us(m)**2 / g )
     7785
     7786             surf_lsm_h(0)%z0h(m) = 0.40_wp * molecular_viscosity /               &
     7787                                 ( surf_lsm_h(0)%us(m) + 1E-8_wp )
     7788             surf_lsm_h(0)%z0q(m) = 0.62_wp * molecular_viscosity /               &
     7789                                 ( surf_lsm_h(0)%us(m) + 1E-8_wp )
     7790
     7791
     7792             IF ( surf_lsm_h(0)%z0(m) > 0.1_wp * surf_lsm_h(0)%z_mo(m) )  THEN
     7793                surf_lsm_h(0)%z0(m) = 0.1_wp * surf_lsm_h(0)%z_mo(m)
    75777794                flag_exceed_z0   = .TRUE.
    75787795             ENDIF
    75797796
    7580              IF ( surf_lsm_h%z0h(m) >= 0.1_wp * surf_lsm_h%z_mo(m) )  THEN
    7581                 surf_lsm_h%z0h(m) = 0.1_wp * surf_lsm_h%z_mo(m)
     7797             IF ( surf_lsm_h(0)%z0h(m) >= 0.1_wp * surf_lsm_h(0)%z_mo(m) )  THEN
     7798                surf_lsm_h(0)%z0h(m) = 0.1_wp * surf_lsm_h(0)%z_mo(m)
    75827799                flag_exceed_z0h   = .TRUE.
    75837800             ENDIF
    75847801
    7585              IF ( surf_lsm_h%z0q(m) >= 0.1_wp * surf_lsm_h%z_mo(m) )  THEN
    7586                 surf_lsm_h%z0q(m) = 0.1_wp * surf_lsm_h%z_mo(m)
     7802             IF ( surf_lsm_h(0)%z0q(m) >= 0.1_wp * surf_lsm_h(0)%z_mo(m) )  THEN
     7803                surf_lsm_h(0)%z0q(m) = 0.1_wp * surf_lsm_h(0)%z_mo(m)
    75877804                flag_exceed_z0h   = .TRUE.
    75887805             ENDIF
Note: See TracChangeset for help on using the changeset viewer.