Changeset 4452 for palm


Ignore:
Timestamp:
Mar 10, 2020 8:15:32 PM (11 months ago)
Author:
suehring
Message:

Bugfix in calc_albedo

File:
1 edited

Legend:

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

    r4442 r4452  
    2828! -----------------
    2929! $Id$
    30 ! - Change order of dimension in surface arrays %frac, %emissivity and %albedo
     30! Bugfix in calc_albedo
     31!
     32! 4442 2020-03-04 19:21:13Z suehring
     33! - Change order of dimension in surface arrays %frac, %emissivity and %albedo
    3134!   to allow for better vectorization in the radiation interactions.
    3235! - Minor formatting issues
    33 ! 
     36!
    3437! 4441 2020-03-04 19:20:35Z suehring
    3538! bugfixes: cpp-directives for serial mode moved, small changes to get serial mode compiled
    36 ! 
     39!
    3740! 4400 2020-02-10 20:32:41Z suehring
    3841! Initialize radiation arrays with zero
    39 ! 
     42!
    4043! 4392 2020-01-31 16:14:57Z pavelkrc
    4144! - Add debug tracing of large radiative fluxes (option trace_fluxes_above)
    4245! - Print exact counts of SVF and CSF if debut_output is enabled
    4346! - Update descriptions of RTM 3.0 and related comments
    44 ! 
     47!
    4548! 4360 2020-01-07 11:25:50Z suehring
    4649! Renamed pc_heating_rate, pc_transpiration_rate, pc_transpiration_rate to
    4750! pcm_heating_rate, pcm_latent_rate, pcm_transpiration_rate
    48 ! 
     51!
    4952! 4340 2019-12-16 08:17:03Z Giersch
    50 ! Albedo indices for building_surface_pars are now declared as parameters to 
     53! Albedo indices for building_surface_pars are now declared as parameters to
    5154! prevent an error if the gfortran compiler with -Werror=unused-value is used
    52 ! 
     55!
    5356! 4291 2019-11-11 12:36:54Z moh.hefny
    54 ! Enabled RTM in case of biometeorology even if there is no vertical 
     57! Enabled RTM in case of biometeorology even if there is no vertical
    5558! surfaces or 3D vegetation in the domain
    56 ! 
     59!
    5760! 4286 2019-10-30 16:01:14Z resler
    5861! - Fix wrong treating of time_rad during interpolation in radiation_model_mod
     
    6265! 4271 2019-10-23 10:46:41Z maronga
    6366! Bugfix: missing parentheses in calculation of snow albedo
    64 ! 
     67!
    6568! 4245 2019-09-30 08:40:37Z pavelkrc
    6669! Initialize explicit per-surface albedos from building_surface_pars
    67 ! 
     70!
    6871! 4238 2019-09-25 16:06:01Z suehring
    6972! Modify check in order to avoid equality comparisons of floating points
    70 ! 
     73!
    7174! 4227 2019-09-10 18:04:34Z gronemeier
    7275! implement new palm_date_time_mod
    73 ! 
     76!
    7477! 4226 2019-09-10 17:03:24Z suehring
    7578! - Netcdf input routine for dimension length renamed
    76 ! - Define time variable for external radiation input relative to time_utc_init 
    77 ! 
     79! - Define time variable for external radiation input relative to time_utc_init
     80!
    7881! 4210 2019-09-02 13:07:09Z suehring
    7982! - Revise steering of splitting diffuse and direct radiation
     
    8184! - Optimize mapping of radiation components onto 2D arrays, avoid unnecessary
    8285!   operations
    83 ! 
     86!
    8487! 4208 2019-09-02 09:01:07Z suehring
    8588! Bugfix in accessing albedo_pars in the clear-sky branch
    86 ! (merge from branch resler) 
    87 ! 
     89! (merge from branch resler)
     90!
    8891! 4198 2019-08-29 15:17:48Z gronemeier
    8992! Prohibit execution of radiation model if rotation_angle is not zero
    90 ! 
     93!
    9194! 4197 2019-08-29 14:33:32Z suehring
    9295! Revise steering of surface albedo initialization when albedo_pars is provided
    93 ! 
     96!
    9497! 4190 2019-08-27 15:42:37Z suehring
    95 ! Implement external radiation forcing also for level-of-detail = 2 
     98! Implement external radiation forcing also for level-of-detail = 2
    9699! (horizontally 2D radiation)
    97 ! 
     100!
    98101! 4188 2019-08-26 14:15:47Z suehring
    99102! Minor adjustment in error message
    100 ! 
     103!
    101104! 4187 2019-08-26 12:43:15Z suehring
    102 ! - Take external radiation from root domain dynamic input if not provided for 
     105! - Take external radiation from root domain dynamic input if not provided for
    103106!   each nested domain
    104107! - Combine MPI_ALLREDUCE calls to reduce mpi overhead
    105 ! 
     108!
    106109! 4182 2019-08-22 15:20:23Z scharf
    107110! Corrected "Former revisions" section
    108 ! 
     111!
    109112! 4179 2019-08-21 11:16:12Z suehring
    110113! Remove debug prints
    111 ! 
     114!
    112115! 4178 2019-08-21 11:13:06Z suehring
    113116! External radiation forcing implemented.
    114 ! 
     117!
    115118! 4168 2019-08-16 13:50:17Z suehring
    116119! Replace function get_topography_top_index by topo_top_ind
    117 ! 
     120!
    118121! 4157 2019-08-14 09:19:12Z suehring
    119122! Give informative message on raytracing distance only by core zero
    120 ! 
     123!
    121124! 4148 2019-08-08 11:26:00Z suehring
    122125! Comments added
    123 ! 
     126!
    124127! 4134 2019-08-02 18:39:57Z suehring
    125128! Bugfix in formatted write statement
    126 ! 
     129!
    127130! 4127 2019-07-30 14:47:10Z suehring
    128131! Remove unused pch_index (merge from branch resler)
    129 ! 
     132!
    130133! 4089 2019-07-11 14:30:27Z suehring
    131134! - Correct level 2 initialization of spectral albedos in rrtmg branch, long- and
    132 !   shortwave albedos were mixed-up. 
     135!   shortwave albedos were mixed-up.
    133136! - Change order of albedo_pars so that it is now consistent with the defined
    134 !   order of albedo_pars in PIDS 
    135 ! 
     137!   order of albedo_pars in PIDS
     138!
    136139! 4069 2019-07-01 14:05:51Z Giersch
    137 ! Masked output running index mid has been introduced as a local variable to 
     140! Masked output running index mid has been introduced as a local variable to
    138141! avoid runtime error (Loop variable has been modified) in time_integration
    139 ! 
     142!
    140143! 4067 2019-07-01 13:29:25Z suehring
    141144! Bugfix, pass dummy string to MPI_INFO_SET (J. Resler)
    142 ! 
     145!
    143146! 4039 2019-06-18 10:32:41Z suehring
    144147! Bugfix for masked data output
    145 ! 
     148!
    146149! 4008 2019-05-30 09:50:11Z moh.hefny
    147150! Bugfix in check variable when a variable's string is less than 3
     
    150153!
    151154! 3992 2019-05-22 16:49:38Z suehring
    152 ! Bugfix in rrtmg radiation branch in a nested run when the lowest prognistic 
    153 ! grid points in a child domain are all inside topography 
    154 ! 
     155! Bugfix in rrtmg radiation branch in a nested run when the lowest prognistic
     156! grid points in a child domain are all inside topography
     157!
    155158! 3987 2019-05-22 09:52:13Z kanani
    156159! Introduce alternative switch for debug output during timestepping
    157 ! 
     160!
    158161! 3943 2019-05-02 09:50:41Z maronga
    159162! Missing blank characteer added.
    160 ! 
     163!
    161164! 3900 2019-04-16 15:17:43Z suehring
    162165! Fixed initialization problem
    163 ! 
     166!
    164167! 3885 2019-04-11 11:29:34Z kanani
    165 ! Changes related to global restructuring of location messages and introduction 
     168! Changes related to global restructuring of location messages and introduction
    166169! of additional debug messages
    167 ! 
     170!
    168171! 3881 2019-04-10 09:31:22Z suehring
    169172! Output of albedo and emissivity moved from USM, bugfixes in initialization
    170173! of albedo
    171 ! 
     174!
    172175! 3861 2019-04-04 06:27:41Z maronga
    173176! Bugfix: factor of 4.0 required instead of 3.0 in calculation of rad_lw_out_change_0
    174 ! 
     177!
    175178! 3859 2019-04-03 20:30:31Z maronga
    176179! Added some descriptions
    177 ! 
     180!
    178181! 3847 2019-04-01 14:51:44Z suehring
    179182! Implement check for dt_radiation (must be > 0)
    180 ! 
     183!
    181184! 3846 2019-04-01 13:55:30Z suehring
    182185! unused variable removed
    183 ! 
     186!
    184187! 3814 2019-03-26 08:40:31Z pavelkrc
    185188! Change zenith(0:0) and others to scalar.
    186189! Code review.
    187190! Rename exported nzu, nzp and related variables due to name conflict
    188 ! 
     191!
    189192! 3771 2019-02-28 12:19:33Z raasch
    190193! rrtmg preprocessor for directives moved/added, save attribute added to temporary
    191194! pointers to avoid compiler warnings about outlived pointer targets,
    192195! statement added to avoid compiler warning about unused variable
    193 ! 
     196!
    194197! 3769 2019-02-28 10:16:49Z moh.hefny
    195198! removed unused variables and subroutine radiation_radflux_gridbox
     
    197200! 3767 2019-02-27 08:18:02Z raasch
    198201! unused variable for file index removed from rrd-subroutines parameter list
    199 ! 
     202!
    200203! 3760 2019-02-21 18:47:35Z moh.hefny
    201204! Bugfix: initialized simulated_time before calculating solar position
     
    205208! (resler, pavelkrc)
    206209! Bugfixes: add further required MRT factors to read/write_svf,
    207 ! fix for aggregating view factors to eliminate local noise in reflected 
    208 ! irradiance at mutually close surfaces (corners, presence of trees) in the 
     210! fix for aggregating view factors to eliminate local noise in reflected
     211! irradiance at mutually close surfaces (corners, presence of trees) in the
    209212! angular discretization scheme.
    210213!
     
    213216!
    214217! 3705 2019-01-29 19:56:39Z suehring
    215 ! Make variables that are sampled in virtual measurement module public 
    216 ! 
     218! Make variables that are sampled in virtual measurement module public
     219!
    217220! 3704 2019-01-29 19:51:41Z suehring
    218221! Some interface calls moved to module_interface + cleanup
    219 ! 
     222!
    220223! 3667 2019-01-10 14:26:24Z schwenkel
    221224! Modified check for rrtmg input files
    222 ! 
     225!
    223226! 3655 2019-01-07 16:51:22Z knoop
    224227! nopointer option removed
    225 ! 
     228!
    226229! 1496 2014-12-02 17:25:50Z maronga
    227230! Initial revision
    228 ! 
     231!
    229232!
    230233! Description:
     
    247250!------------------------------------------------------------------------------!
    248251 MODULE radiation_model_mod
    249  
     252
    250253    USE arrays_3d,                                                             &
    251254        ONLY:  dzw, hyp, nc, pt, p, q, ql, u, v, w, zu, zw, exner, d_exner
     
    276279
    277280    USE grid_variables,                                                        &
    278          ONLY:  ddx, ddy, dx, dy 
     281         ONLY:  ddx, ddy, dx, dy
    279282
    280283    USE indices,                                                               &
     
    366369!-- Predefined Land surface classes (albedo_type) after Briegleb (1992)
    367370    CHARACTER(37), DIMENSION(0:33), PARAMETER :: albedo_type_name = (/      &
    368                                    'user defined                         ', & !  0 
     371                                   'user defined                         ', & !  0
    369372                                   'ocean                                ', & !  1
    370373                                   'mixed farming, tall grassland        ', & !  2
    371                                    'tall/medium grassland                ', & !  3 
     374                                   'tall/medium grassland                ', & !  3
    372375                                   'evergreen shrubland                  ', & !  4
    373376                                   'short grassland/meadow/shrubland     ', & !  5
     
    377380                                   'tropical evergreen broadleaved forest', & !  9
    378381                                   'medium/tall grassland/woodland       ', & ! 10
    379                                    'desert, sandy                        ', & ! 11 
    380                                    'desert, rocky                        ', & ! 12 
     382                                   'desert, sandy                        ', & ! 11
     383                                   'desert, rocky                        ', & ! 12
    381384                                   'tundra                               ', & ! 13
    382                                    'land ice                             ', & ! 14 
    383                                    'sea ice                              ', & ! 15 
     385                                   'land ice                             ', & ! 14
     386                                   'sea ice                              ', & ! 15
    384387                                   'snow                                 ', & ! 16
    385388                                   'bare soil                            ', & ! 17
     
    415418
    416419    INTEGER(iwp) :: albedo_type  = 9999999_iwp, &     !< Albedo surface type
    417                     dots_rad     = 0_iwp              !< starting index for timeseries output 
     420                    dots_rad     = 0_iwp              !< starting index for timeseries output
    418421
    419422    LOGICAL ::  unscheduled_radiation_calls = .TRUE., & !< flag parameter indicating whether additional calls of the radiation code are allowed
     
    467470    REAL(wp), PARAMETER :: emissivity_atm_clsky = 0.8_wp       !< emissivity of the clear-sky atmosphere
    468471!
    469 !-- Land surface albedos for solar zenith angle of 60degree after Briegleb (1992)     
     472!-- Land surface albedos for solar zenith angle of 60degree after Briegleb (1992)
    470473!-- (broadband, longwave, shortwave ):   bb,      lw,      sw,
    471     REAL(wp), DIMENSION(0:2,1:33), PARAMETER :: albedo_pars = RESHAPE( (/& 
     474    REAL(wp), DIMENSION(0:2,1:33), PARAMETER :: albedo_pars = RESHAPE( (/&
    472475                                   0.06_wp, 0.06_wp, 0.06_wp,            & !  1
    473476                                   0.19_wp, 0.28_wp, 0.09_wp,            & !  2
     
    541544
    542545!
    543 !-- The following variables should be only changed with care, as this will 
     546!-- The following variables should be only changed with care, as this will
    544547!-- require further setting of some variables, which is currently not
    545548!-- implemented (aerosols, ice phase).
     
    570573                                             rrtm_cliqwp,    & !< in-cloud liquid water path (g/m2)
    571574                                             rrtm_co2vmr,    & !< CO2 volume mixing ratio (g/mol)
    572                                              rrtm_emis,      & !< surface emissivity (0-1) 
     575                                             rrtm_emis,      & !< surface emissivity (0-1)
    573576                                             rrtm_h2ovmr,    & !< H2O volume mixing ratio
    574577                                             rrtm_n2ovmr,    & !< N2O volume mixing ratio
     
    582585                                             rrtm_tlev,      & !< actual temperature (K, zw-grid)
    583586                                             rrtm_lwdflx,    & !< RRTM output of incoming longwave radiation flux (W/m2)
    584                                              rrtm_lwdflxc,   & !< RRTM output of outgoing clear sky longwave radiation flux (W/m2) 
     587                                             rrtm_lwdflxc,   & !< RRTM output of outgoing clear sky longwave radiation flux (W/m2)
    585588                                             rrtm_lwuflx,    & !< RRTM output of outgoing longwave radiation flux (W/m2)
    586589                                             rrtm_lwuflxc,   & !< RRTM output of incoming clear sky longwave radiation flux (W/m2)
     
    590593                                             rrtm_lwhrc,     & !< RRTM output of incoming longwave clear sky radiation heating rate (K/d)
    591594                                             rrtm_swdflx,    & !< RRTM output of incoming shortwave radiation flux (W/m2)
    592                                              rrtm_swdflxc,   & !< RRTM output of outgoing clear sky shortwave radiation flux (W/m2) 
     595                                             rrtm_swdflxc,   & !< RRTM output of outgoing clear sky shortwave radiation flux (W/m2)
    593596                                             rrtm_swuflx,    & !< RRTM output of outgoing shortwave radiation flux (W/m2)
    594597                                             rrtm_swuflxc,   & !< RRTM output of incoming clear sky shortwave radiation flux (W/m2)
    595598                                             rrtm_swhr,      & !< RRTM output of shortwave radiation heating rate (K/d)
    596                                              rrtm_swhrc,     & !< RRTM output of incoming shortwave clear sky radiation heating rate (K/d) 
    597                                              rrtm_dirdflux,  & !< RRTM output of incoming direct shortwave (W/m2) 
    598                                              rrtm_difdflux     !< RRTM output of incoming diffuse shortwave (W/m2) 
     599                                             rrtm_swhrc,     & !< RRTM output of incoming shortwave clear sky radiation heating rate (K/d)
     600                                             rrtm_dirdflux,  & !< RRTM output of incoming direct shortwave (W/m2)
     601                                             rrtm_difdflux     !< RRTM output of incoming diffuse shortwave (W/m2)
    599602
    600603    REAL(wp), DIMENSION(1) ::                rrtm_aldif,     & !< surface albedo for longwave diffuse radiation
     
    734737!-- type for calculation of svf
    735738    TYPE t_svf
    736         INTEGER(iwp)                               :: isurflt           !< 
    737         INTEGER(iwp)                               :: isurfs            !< 
    738         REAL(wp)                                   :: rsvf              !< 
    739         REAL(wp)                                   :: rtransp           !< 
     739        INTEGER(iwp)                               :: isurflt           !<
     740        INTEGER(iwp)                               :: isurfs            !<
     741        REAL(wp)                                   :: rsvf              !<
     742        REAL(wp)                                   :: rtransp           !<
    740743    END TYPE
    741744
    742745!-- type for calculation of csf
    743746    TYPE t_csf
    744         INTEGER(iwp)                               :: ip                !< 
    745         INTEGER(iwp)                               :: itx               !< 
    746         INTEGER(iwp)                               :: ity               !< 
    747         INTEGER(iwp)                               :: itz               !< 
     747        INTEGER(iwp)                               :: ip                !<
     748        INTEGER(iwp)                               :: itx               !<
     749        INTEGER(iwp)                               :: ity               !<
     750        INTEGER(iwp)                               :: itz               !<
    748751        INTEGER(iwp)                               :: isurfs            !< Idx of source face / -1 for sky
    749752        REAL(wp)                                   :: rcvf              !< Canopy view factor for faces /
     
    788791    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinswdif      !< array of diffuse sw radiation from sky and model boundary falling to local surface
    789792    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfinlwdif      !< array of diffuse lw radiation from sky and model boundary falling to local surface
    790    
     793
    791794                                                                        !< Outward radiation is only valid for nonvirtual surfaces
    792795    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  surfoutsl        !< array of reflected sw radiation for local surface in i-th reflection
     
    835838    INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE      ::  boxes            !< coordinates of gridboxes being crossed by ray
    836839    REAL(wp), DIMENSION(:), ALLOCATABLE            ::  crlens           !< array of crossing lengths of ray for particular grid boxes
    837     INTEGER(iwp), DIMENSION(:), ALLOCATABLE        ::  lad_ip           !< array of numbers of process where lad is stored 
     840    INTEGER(iwp), DIMENSION(:), ALLOCATABLE        ::  lad_ip           !< array of numbers of process where lad is stored
    838841#if defined( __parallel )
    839842    INTEGER(kind=MPI_ADDRESS_KIND), &
     
    876879    REAL(wp), DIMENSION(:), ALLOCATABLE            :: emiss_surf         !< emissivity of the wall surface
    877880!
    878 !-- External radiation. Depending on the given level of detail either a 1D or 
     881!-- External radiation. Depending on the given level of detail either a 1D or
    879882!-- a 3D array will be allocated.
    880     TYPE( real_1d_3d ) ::  rad_lw_in_f     !< external incoming longwave radiation, from observation or model 
    881     TYPE( real_1d_3d ) ::  rad_sw_in_f     !< external incoming shortwave radiation, from observation or model 
     883    TYPE( real_1d_3d ) ::  rad_lw_in_f     !< external incoming longwave radiation, from observation or model
     884    TYPE( real_1d_3d ) ::  rad_sw_in_f     !< external incoming shortwave radiation, from observation or model
    882885    TYPE( real_1d_3d ) ::  rad_sw_in_dif_f !< external incoming shortwave radiation, diffuse part, from observation or model
    883     TYPE( real_1d_3d ) ::  time_rad_f      !< time dimension for external radiation, from observation or model 
     886    TYPE( real_1d_3d ) ::  time_rad_f      !< time dimension for external radiation, from observation or model
    884887
    885888    INTERFACE radiation_check_data_output
     
    894897       MODULE PROCEDURE radiation_check_data_output_pr
    895898    END INTERFACE radiation_check_data_output_pr
    896  
     899
    897900    INTERFACE radiation_check_parameters
    898901       MODULE PROCEDURE radiation_check_parameters
    899902    END INTERFACE radiation_check_parameters
    900  
     903
    901904    INTERFACE radiation_clearsky
    902905       MODULE PROCEDURE radiation_clearsky
    903906    END INTERFACE radiation_clearsky
    904  
     907
    905908    INTERFACE radiation_constant
    906909       MODULE PROCEDURE radiation_constant
    907910    END INTERFACE radiation_constant
    908  
     911
    909912    INTERFACE radiation_control
    910913       MODULE PROCEDURE radiation_control
     
    933936    INTERFACE radiation_header
    934937       MODULE PROCEDURE radiation_header
    935     END INTERFACE radiation_header 
    936  
     938    END INTERFACE radiation_header
     939
    937940    INTERFACE radiation_init
    938941       MODULE PROCEDURE radiation_init
     
    942945       MODULE PROCEDURE radiation_parin
    943946    END INTERFACE radiation_parin
    944    
     947
    945948    INTERFACE radiation_rrtmg
    946949       MODULE PROCEDURE radiation_rrtmg
     
    969972       MODULE PROCEDURE radiation_interaction_init
    970973    END INTERFACE radiation_interaction_init
    971  
     974
    972975    INTERFACE radiation_presimulate_solar_pos
    973976       MODULE PROCEDURE radiation_presimulate_solar_pos
     
    10051008           radiation_read_svf, radiation_presimulate_solar_pos
    10061009
    1007    
     1010
    10081011!
    10091012!-- Public variables and constants / NEEDS SORTING
     
    10421045!------------------------------------------------------------------------------!
    10431046    SUBROUTINE radiation_control
    1044  
    1045  
     1047
     1048
    10461049       IMPLICIT NONE
    10471050
     
    10541057          CASE ( 'constant' )
    10551058             CALL radiation_constant
    1056          
    1057           CASE ( 'clear-sky' ) 
     1059
     1060          CASE ( 'clear-sky' )
    10581061             CALL radiation_clearsky
    1059        
     1062
    10601063          CASE ( 'rrtmg' )
    10611064             CALL radiation_rrtmg
    1062              
     1065
    10631066          CASE ( 'external' )
    10641067!
     
    10841087!------------------------------------------------------------------------------!
    10851088    SUBROUTINE radiation_check_data_output( variable, unit, i, ilen, k )
    1086  
    1087  
     1089
     1090
    10881091       USE control_parameters,                                                 &
    10891092           ONLY: data_output, message_string
     
    12621265! ------------
    12631266!> Check data output of profiles for radiation model
    1264 !------------------------------------------------------------------------------! 
     1267!------------------------------------------------------------------------------!
    12651268    SUBROUTINE radiation_check_data_output_pr( variable, var_count, unit,      &
    12661269               dopr_unit )
    1267  
     1270
    12681271       USE arrays_3d,                                                          &
    12691272           ONLY: zu
     
    12791282
    12801283       IMPLICIT NONE
    1281    
    1282        CHARACTER (LEN=*) ::  unit      !< 
    1283        CHARACTER (LEN=*) ::  variable  !< 
     1284
     1285       CHARACTER (LEN=*) ::  unit      !<
     1286       CHARACTER (LEN=*) ::  variable  !<
    12841287       CHARACTER (LEN=*) ::  dopr_unit !< local value of dopr_unit
    1285  
    1286        INTEGER(iwp) ::  var_count     !< 
     1288
     1289       INTEGER(iwp) ::  var_count     !<
    12871290
    12881291       SELECT CASE ( TRIM( variable ) )
    1289        
     1292
    12901293         CASE ( 'rad_net' )
    12911294             IF ( (  .NOT.  radiation )  .OR.  radiation_scheme == 'constant' )&
     
    13151318                dopr_unit  = 'W/m2'
    13161319                hom(:,2,100,:)  = SPREAD( zw, 2, statistic_regions+1 )
    1317                 unit = dopr_unit 
     1320                unit = dopr_unit
    13181321             ENDIF
    13191322
     
    13301333                dopr_unit  = 'W/m2'
    13311334                hom(:,2,101,:)  = SPREAD( zw, 2, statistic_regions+1 )
    1332                 unit = dopr_unit   
     1335                unit = dopr_unit
    13331336             ENDIF
    13341337
     
    14311434
    14321435    END SUBROUTINE radiation_check_data_output_pr
    1433  
    1434  
     1436
     1437
    14351438!------------------------------------------------------------------------------!
    14361439! Description:
     
    14441447
    14451448       USE netcdf_data_input_mod,                                              &
    1446            ONLY:  input_pids_static                 
    1447    
     1449           ONLY:  input_pids_static
     1450
    14481451       IMPLICIT NONE
    1449        
    1450 !
    1451 !--    In case no urban-surface or land-surface model is applied, usage of 
    1452 !--    a radiation model make no sense.         
     1452
     1453!
     1454!--    In case no urban-surface or land-surface model is applied, usage of
     1455!--    a radiation model make no sense.
    14531456       IF ( .NOT. land_surface  .AND.  .NOT. urban_surface )  THEN
    14541457          message_string = 'Usage of radiation module is only allowed if ' //  &
     
    14661469       ELSEIF ( radiation_scheme == 'rrtmg' )  THEN
    14671470#if ! defined ( __rrtmg )
    1468           message_string = 'radiation_scheme = "rrtmg" requires ' //           & 
     1471          message_string = 'radiation_scheme = "rrtmg" requires ' //           &
    14691472                           'compilation of PALM with pre-processor ' //        &
    14701473                           'directive -D__rrtmg'
     
    14721475#endif
    14731476#if defined ( __rrtmg ) && ! defined( __netcdf )
    1474           message_string = 'radiation_scheme = "rrtmg" requires ' //           & 
     1477          message_string = 'radiation_scheme = "rrtmg" requires ' //           &
    14751478                           'the use of NetCDF (preprocessor directive ' //     &
    14761479                           '-D__netcdf'
     
    14801483       ENDIF
    14811484!
    1482 !--    Checks performed only if data is given via namelist only. 
     1485!--    Checks performed only if data is given via namelist only.
    14831486       IF ( .NOT. input_pids_static )  THEN
    14841487          IF ( albedo_type == 0  .AND.  albedo == 9999999.9_wp  .AND.          &
    14851488               radiation_scheme == 'clear-sky')  THEN
    1486              message_string = 'radiation_scheme = "clear-sky" in combination'//& 
     1489             message_string = 'radiation_scheme = "clear-sky" in combination'//&
    14871490                              'with albedo_type = 0 requires setting of'//     &
    14881491                              'albedo /= 9999999.9'
     
    14921495          IF ( albedo_type == 0  .AND.  radiation_scheme == 'rrtmg'  .AND.     &
    14931496             ( albedo_lw_dif == 9999999.9_wp .OR. albedo_lw_dir == 9999999.9_wp&
    1494           .OR. albedo_sw_dif == 9999999.9_wp .OR. albedo_sw_dir == 9999999.9_wp& 
     1497          .OR. albedo_sw_dif == 9999999.9_wp .OR. albedo_sw_dir == 9999999.9_wp&
    14951498             ) ) THEN
    1496              message_string = 'radiation_scheme = "rrtmg" in combination' //   & 
     1499             message_string = 'radiation_scheme = "rrtmg" in combination' //   &
    14971500                              'with albedo_type = 0 requires setting of ' //   &
    14981501                              'albedo_lw_dif /= 9999999.9' //                  &
     
    15061509!--    Parallel rad_angular_discretization without raytrace_mpi_rma is not implemented
    15071510!--    Serial mode does not allow mpi_rma
    1508 #if defined( __parallel )     
     1511#if defined( __parallel )
    15091512       IF ( rad_angular_discretization  .AND.  .NOT. raytrace_mpi_rma )  THEN
    15101513          message_string = 'rad_angular_discretization can only be used ' //  &
     
    15221525       IF ( cloud_droplets  .AND.   radiation_scheme == 'rrtmg'  .AND.         &
    15231526            average_radiation ) THEN
    1524           message_string = 'average_radiation = .T. with radiation_scheme'//   & 
     1527          message_string = 'average_radiation = .T. with radiation_scheme'//   &
    15251528                           '= "rrtmg" in combination cloud_droplets = .T.'//   &
    15261529                           'is not implementd'
     
    15391542!--    Check for dt_radiation
    15401543       IF ( dt_radiation <= 0.0 )  THEN
    1541           message_string = 'dt_radiation must be > 0.0' 
    1542           CALL message( 'check_parameters', 'PA0591', 1, 2, 0, 6, 0 ) 
     1544          message_string = 'dt_radiation must be > 0.0'
     1545          CALL message( 'check_parameters', 'PA0591', 1, 2, 0, 6, 0 )
    15431546       ENDIF
    15441547!
     
    15491552                           'model.&Using rotation_angle /= 0.0 is not allowed in combination ' //  &
    15501553                           'with the radiation model at the moment!'
    1551           CALL message( 'check_parameters', 'PA0675', 1, 2, 0, 6, 0 ) 
     1554          CALL message( 'check_parameters', 'PA0675', 1, 2, 0, 6, 0 )
    15521555       ENDIF
    1553  
    1554     END SUBROUTINE radiation_check_parameters 
    1555  
    1556  
     1556
     1557    END SUBROUTINE radiation_check_parameters
     1558
     1559
    15571560!------------------------------------------------------------------------------!
    15581561! Description:
     
    15611564!------------------------------------------------------------------------------!
    15621565    SUBROUTINE radiation_init
    1563    
     1566
    15641567       IMPLICIT NONE
    15651568
    1566        INTEGER(iwp) ::  i         !< running index x-direction 
     1569       INTEGER(iwp) ::  i         !< running index x-direction
    15671570       INTEGER(iwp) ::  is        !< running index for input surface elements
    1568        INTEGER(iwp) ::  ioff      !< offset in x between surface element reference grid point in atmosphere and actual surface 
    1569        INTEGER(iwp) ::  j         !< running index y-direction 
    1570        INTEGER(iwp) ::  joff      !< offset in y between surface element reference grid point in atmosphere and actual surface 
     1571       INTEGER(iwp) ::  ioff      !< offset in x between surface element reference grid point in atmosphere and actual surface
     1572       INTEGER(iwp) ::  j         !< running index y-direction
     1573       INTEGER(iwp) ::  joff      !< offset in y between surface element reference grid point in atmosphere and actual surface
    15711574       INTEGER(iwp) ::  k         !< running index z-direction
    1572        INTEGER(iwp) ::  l         !< running index for orientation of vertical surfaces 
     1575       INTEGER(iwp) ::  l         !< running index for orientation of vertical surfaces
    15731576       INTEGER(iwp) ::  m         !< running index for surface elements
    15741577       INTEGER(iwp) ::  ntime = 0 !< number of available external radiation timesteps
     
    15841587!      or if biometeorology output is required for flat surfaces.
    15851588!--    The namelist parameter radiation_interactions_on can override this behavior.
    1586 !--    (This check cannot be performed in check_parameters, because vertical_surfaces_exist is first set in 
     1589!--    (This check cannot be performed in check_parameters, because vertical_surfaces_exist is first set in
    15871590!--    init_surface_arrays.)
    15881591       IF ( radiation_interactions_on )  THEN
     
    15911594             average_radiation         = .TRUE.
    15921595          ELSE
    1593              radiation_interactions_on = .FALSE.   !< reset namelist parameter: no interactions 
     1596             radiation_interactions_on = .FALSE.   !< reset namelist parameter: no interactions
    15941597                                                   !< calculations necessary in case of flat surface
    15951598          ENDIF
    15961599       ELSEIF ( vertical_surfaces_exist  .OR.  plant_canopy  .OR.  biometeorology )  THEN
    15971600          message_string = 'radiation_interactions_on is set to .FALSE. although '     // &
    1598                            'vertical surfaces and/or trees or biometeorology exist '   // & 
     1601                           'vertical surfaces and/or trees or biometeorology exist '   // &
    15991602                           'is ON. The model will run without RTM (no shadows, no '    // &
    16001603                           'radiation reflections)'
     
    16071610
    16081611!
    1609 !--    If required, initialize radiation interactions between surfaces 
     1612!--    If required, initialize radiation interactions between surfaces
    16101613!--    via sky-view factors. This must be done before radiation is initialized.
    16111614       IF ( radiation_interactions )  CALL radiation_interaction_init
     
    16151618                  surf_lsm_h%ns > 0  )   THEN
    16161619          ALLOCATE( surf_lsm_h%rad_net(1:surf_lsm_h%ns) )
    1617           surf_lsm_h%rad_net = 0.0_wp 
     1620          surf_lsm_h%rad_net = 0.0_wp
    16181621       ENDIF
    16191622       IF ( .NOT. ALLOCATED ( surf_usm_h%rad_net )  .AND.                      &
    16201623                  surf_usm_h%ns > 0  )  THEN
    16211624          ALLOCATE( surf_usm_h%rad_net(1:surf_usm_h%ns) )
    1622           surf_usm_h%rad_net = 0.0_wp 
     1625          surf_usm_h%rad_net = 0.0_wp
    16231626       ENDIF
    16241627       DO  l = 0, 3
     
    16261629                     surf_lsm_v(l)%ns > 0  )  THEN
    16271630             ALLOCATE( surf_lsm_v(l)%rad_net(1:surf_lsm_v(l)%ns) )
    1628              surf_lsm_v(l)%rad_net = 0.0_wp 
     1631             surf_lsm_v(l)%rad_net = 0.0_wp
    16291632          ENDIF
    16301633          IF ( .NOT. ALLOCATED ( surf_usm_v(l)%rad_net )  .AND.                &
    16311634                     surf_usm_v(l)%ns > 0  )  THEN
    16321635             ALLOCATE( surf_usm_v(l)%rad_net(1:surf_usm_v(l)%ns) )
    1633              surf_usm_v(l)%rad_net = 0.0_wp 
     1636             surf_usm_v(l)%rad_net = 0.0_wp
    16341637          ENDIF
    16351638       ENDDO
     
    16411644                  surf_lsm_h%ns > 0  )   THEN
    16421645          ALLOCATE( surf_lsm_h%rad_lw_out_change_0(1:surf_lsm_h%ns) )
    1643           surf_lsm_h%rad_lw_out_change_0 = 0.0_wp 
     1646          surf_lsm_h%rad_lw_out_change_0 = 0.0_wp
    16441647       ENDIF
    16451648       IF ( .NOT. ALLOCATED ( surf_usm_h%rad_lw_out_change_0 )  .AND.          &
    16461649                  surf_usm_h%ns > 0  )  THEN
    16471650          ALLOCATE( surf_usm_h%rad_lw_out_change_0(1:surf_usm_h%ns) )
    1648           surf_usm_h%rad_lw_out_change_0 = 0.0_wp 
     1651          surf_usm_h%rad_lw_out_change_0 = 0.0_wp
    16491652       ENDIF
    16501653       DO  l = 0, 3
     
    16521655                     surf_lsm_v(l)%ns > 0  )  THEN
    16531656             ALLOCATE( surf_lsm_v(l)%rad_lw_out_change_0(1:surf_lsm_v(l)%ns) )
    1654              surf_lsm_v(l)%rad_lw_out_change_0 = 0.0_wp 
     1657             surf_lsm_v(l)%rad_lw_out_change_0 = 0.0_wp
    16551658          ENDIF
    16561659          IF ( .NOT. ALLOCATED ( surf_usm_v(l)%rad_lw_out_change_0 )  .AND.    &
    16571660                     surf_usm_v(l)%ns > 0  )  THEN
    16581661             ALLOCATE( surf_usm_v(l)%rad_lw_out_change_0(1:surf_usm_v(l)%ns) )
    1659              surf_usm_v(l)%rad_lw_out_change_0 = 0.0_wp 
     1662             surf_usm_v(l)%rad_lw_out_change_0 = 0.0_wp
    16601663          ENDIF
    16611664       ENDDO
     
    16761679          ALLOCATE( surf_lsm_h%rad_lw_ref(1:surf_lsm_h%ns) )
    16771680          ALLOCATE( surf_lsm_h%rad_lw_res(1:surf_lsm_h%ns) )
    1678           surf_lsm_h%rad_sw_in  = 0.0_wp 
    1679           surf_lsm_h%rad_sw_out = 0.0_wp 
    1680           surf_lsm_h%rad_sw_dir = 0.0_wp 
    1681           surf_lsm_h%rad_sw_dif = 0.0_wp 
    1682           surf_lsm_h%rad_sw_ref = 0.0_wp 
    1683           surf_lsm_h%rad_sw_res = 0.0_wp 
    1684           surf_lsm_h%rad_lw_in  = 0.0_wp 
    1685           surf_lsm_h%rad_lw_out = 0.0_wp 
    1686           surf_lsm_h%rad_lw_dif = 0.0_wp 
    1687           surf_lsm_h%rad_lw_ref = 0.0_wp 
    1688           surf_lsm_h%rad_lw_res = 0.0_wp 
     1681          surf_lsm_h%rad_sw_in  = 0.0_wp
     1682          surf_lsm_h%rad_sw_out = 0.0_wp
     1683          surf_lsm_h%rad_sw_dir = 0.0_wp
     1684          surf_lsm_h%rad_sw_dif = 0.0_wp
     1685          surf_lsm_h%rad_sw_ref = 0.0_wp
     1686          surf_lsm_h%rad_sw_res = 0.0_wp
     1687          surf_lsm_h%rad_lw_in  = 0.0_wp
     1688          surf_lsm_h%rad_lw_out = 0.0_wp
     1689          surf_lsm_h%rad_lw_dif = 0.0_wp
     1690          surf_lsm_h%rad_lw_ref = 0.0_wp
     1691          surf_lsm_h%rad_lw_res = 0.0_wp
    16891692       ENDIF
    16901693       IF ( .NOT. ALLOCATED ( surf_usm_h%rad_sw_in )  .AND.                    &
     
    17011704          ALLOCATE( surf_usm_h%rad_lw_ref(1:surf_usm_h%ns) )
    17021705          ALLOCATE( surf_usm_h%rad_lw_res(1:surf_usm_h%ns) )
    1703           surf_usm_h%rad_sw_in  = 0.0_wp 
    1704           surf_usm_h%rad_sw_out = 0.0_wp 
    1705           surf_usm_h%rad_sw_dir = 0.0_wp 
    1706           surf_usm_h%rad_sw_dif = 0.0_wp 
    1707           surf_usm_h%rad_sw_ref = 0.0_wp 
    1708           surf_usm_h%rad_sw_res = 0.0_wp 
    1709           surf_usm_h%rad_lw_in  = 0.0_wp 
    1710           surf_usm_h%rad_lw_out = 0.0_wp 
    1711           surf_usm_h%rad_lw_dif = 0.0_wp 
    1712           surf_usm_h%rad_lw_ref = 0.0_wp 
    1713           surf_usm_h%rad_lw_res = 0.0_wp 
     1706          surf_usm_h%rad_sw_in  = 0.0_wp
     1707          surf_usm_h%rad_sw_out = 0.0_wp
     1708          surf_usm_h%rad_sw_dir = 0.0_wp
     1709          surf_usm_h%rad_sw_dif = 0.0_wp
     1710          surf_usm_h%rad_sw_ref = 0.0_wp
     1711          surf_usm_h%rad_sw_res = 0.0_wp
     1712          surf_usm_h%rad_lw_in  = 0.0_wp
     1713          surf_usm_h%rad_lw_out = 0.0_wp
     1714          surf_usm_h%rad_lw_dif = 0.0_wp
     1715          surf_usm_h%rad_lw_ref = 0.0_wp
     1716          surf_usm_h%rad_lw_res = 0.0_wp
    17141717       ENDIF
    17151718       DO  l = 0, 3
     
    17291732             ALLOCATE( surf_lsm_v(l)%rad_lw_res(1:surf_lsm_v(l)%ns) )
    17301733
    1731              surf_lsm_v(l)%rad_sw_in  = 0.0_wp 
     1734             surf_lsm_v(l)%rad_sw_in  = 0.0_wp
    17321735             surf_lsm_v(l)%rad_sw_out = 0.0_wp
    17331736             surf_lsm_v(l)%rad_sw_dir = 0.0_wp
     
    17361739             surf_lsm_v(l)%rad_sw_res = 0.0_wp
    17371740
    1738              surf_lsm_v(l)%rad_lw_in  = 0.0_wp 
    1739              surf_lsm_v(l)%rad_lw_out = 0.0_wp 
    1740              surf_lsm_v(l)%rad_lw_dif = 0.0_wp 
    1741              surf_lsm_v(l)%rad_lw_ref = 0.0_wp 
    1742              surf_lsm_v(l)%rad_lw_res = 0.0_wp 
     1741             surf_lsm_v(l)%rad_lw_in  = 0.0_wp
     1742             surf_lsm_v(l)%rad_lw_out = 0.0_wp
     1743             surf_lsm_v(l)%rad_lw_dif = 0.0_wp
     1744             surf_lsm_v(l)%rad_lw_ref = 0.0_wp
     1745             surf_lsm_v(l)%rad_lw_res = 0.0_wp
    17431746          ENDIF
    17441747          IF ( .NOT. ALLOCATED ( surf_usm_v(l)%rad_sw_in )  .AND.              &
     
    17551758             ALLOCATE( surf_usm_v(l)%rad_lw_ref(1:surf_usm_v(l)%ns) )
    17561759             ALLOCATE( surf_usm_v(l)%rad_lw_res(1:surf_usm_v(l)%ns) )
    1757              surf_usm_v(l)%rad_sw_in  = 0.0_wp 
     1760             surf_usm_v(l)%rad_sw_in  = 0.0_wp
    17581761             surf_usm_v(l)%rad_sw_out = 0.0_wp
    17591762             surf_usm_v(l)%rad_sw_dir = 0.0_wp
     
    17611764             surf_usm_v(l)%rad_sw_ref = 0.0_wp
    17621765             surf_usm_v(l)%rad_sw_res = 0.0_wp
    1763              surf_usm_v(l)%rad_lw_in  = 0.0_wp 
    1764              surf_usm_v(l)%rad_lw_out = 0.0_wp 
    1765              surf_usm_v(l)%rad_lw_dif = 0.0_wp 
    1766              surf_usm_v(l)%rad_lw_ref = 0.0_wp 
    1767              surf_usm_v(l)%rad_lw_res = 0.0_wp 
     1766             surf_usm_v(l)%rad_lw_in  = 0.0_wp
     1767             surf_usm_v(l)%rad_lw_out = 0.0_wp
     1768             surf_usm_v(l)%rad_lw_dif = 0.0_wp
     1769             surf_usm_v(l)%rad_lw_ref = 0.0_wp
     1770             surf_usm_v(l)%rad_lw_res = 0.0_wp
    17681771          ENDIF
    17691772       ENDDO
     
    18551858          ENDDO
    18561859!
    1857 !--       Level 2 initialization of broadband albedo via given albedo_type. 
    1858 !--       Only if albedo_type is non-zero. In case of urban surface and 
     1860!--       Level 2 initialization of broadband albedo via given albedo_type.
     1861!--       Only if albedo_type is non-zero. In case of urban surface and
    18591862!--       input data is read from ASCII file, albedo_type will be zero, so that
    1860 !--       albedo won't be overwritten. 
     1863!--       albedo won't be overwritten.
    18611864          DO  m = 1, surf_lsm_h%ns
    18621865             IF ( surf_lsm_h%albedo_type(m,ind_veg_wall) /= 0 )                &
     
    19311934                   surf_usm_h%albedo(m,ind_wat_win)   = albedo_pars_f%pars_xy(0,j,i)
    19321935                ENDIF
    1933              ENDDO 
    1934 !
    1935 !--          Vertical surfaces           
     1936             ENDDO
     1937!
     1938!--          Vertical surfaces
    19361939             DO  l = 0, 3
    19371940
     
    19611964             ENDDO
    19621965
    1963           ENDIF 
     1966          ENDIF
    19641967!
    19651968!--       Read explicit albedo values from building surface pars. If present,
     
    20512054!--       Allocate albedos for short/longwave radiation, horizontal surfaces
    20522055!--       for wall/green/window (USM) or vegetation/pavement/water surfaces
    2053 !--       (LSM). 
     2056!--       (LSM).
    20542057          ALLOCATE ( surf_lsm_h%aldif(1:surf_lsm_h%ns,0:2)       )
    20552058          ALLOCATE ( surf_lsm_h%aldir(1:surf_lsm_h%ns,0:2)       )
     
    20712074
    20722075!
    2073 !--       Allocate broadband albedo (temporary for the current radiation 
     2076!--       Allocate broadband albedo (temporary for the current radiation
    20742077!--       implementations)
    20752078          IF ( .NOT. ALLOCATED(surf_lsm_h%albedo) )                            &
     
    20792082
    20802083!
    2081 !--       Allocate albedos for short/longwave radiation, vertical surfaces 
     2084!--       Allocate albedos for short/longwave radiation, vertical surfaces
    20822085          DO  l = 0, 3
    20832086
     
    21112114          ENDDO
    21122115!
    2113 !--       Level 1 initialization of spectral albedos via namelist 
    2114 !--       paramters. Please note, this case all surface tiles are initialized 
    2115 !--       the same. 
     2116!--       Level 1 initialization of spectral albedos via namelist
     2117!--       paramters. Please note, this case all surface tiles are initialized
     2118!--       the same.
    21162119          IF ( surf_lsm_h%ns > 0 )  THEN
    21172120             surf_lsm_h%aldif  = albedo_lw_dif
     
    21622165
    21632166!
    2164 !--       Level 2 initialization of spectral albedos via albedo_type. 
    2165 !--       Please note, for natural- and urban-type surfaces, a tile approach 
    2166 !--       is applied so that the resulting albedo is calculated via the weighted 
    2167 !--       average of respective surface fractions. 
     2167!--       Level 2 initialization of spectral albedos via albedo_type.
     2168!--       Please note, for natural- and urban-type surfaces, a tile approach
     2169!--       is applied so that the resulting albedo is calculated via the weighted
     2170!--       average of respective surface fractions.
    21682171          DO  m = 1, surf_lsm_h%ns
    21692172!
     
    21862189          ENDDO
    21872190!
    2188 !--       For urban surface only if albedo has not been already initialized 
     2191!--       For urban surface only if albedo has not been already initialized
    21892192!--       in the urban-surface model via the ASCII file.
    21902193          IF ( .NOT. surf_usm_h%albedo_from_ascii )  THEN
     
    22312234             ENDDO
    22322235!
    2233 !--          For urban surface only if albedo has not been already initialized 
     2236!--          For urban surface only if albedo has not been already initialized
    22342237!--          in the urban-surface model via the ASCII file.
    22352238             IF ( .NOT. surf_usm_v(l)%albedo_from_ascii )  THEN
     
    22692272                   IF ( albedo_pars_f%pars_xy(0,j,i) /= albedo_pars_f%fill )   &
    22702273                      surf_lsm_h%albedo(m,ind_type) =                          &
    2271                                              albedo_pars_f%pars_xy(0,j,i)     
     2274                                             albedo_pars_f%pars_xy(0,j,i)
    22722275                   IF ( albedo_pars_f%pars_xy(1,j,i) /= albedo_pars_f%fill )   &
    22732276                      surf_lsm_h%aldir(m,ind_type) =                           &
    2274                                              albedo_pars_f%pars_xy(1,j,i)     
     2277                                             albedo_pars_f%pars_xy(1,j,i)
    22752278                   IF ( albedo_pars_f%pars_xy(1,j,i) /= albedo_pars_f%fill )   &
    22762279                      surf_lsm_h%aldif(m,ind_type) =                           &
    2277                                              albedo_pars_f%pars_xy(1,j,i)     
     2280                                             albedo_pars_f%pars_xy(1,j,i)
    22782281                   IF ( albedo_pars_f%pars_xy(2,j,i) /= albedo_pars_f%fill )   &
    22792282                      surf_lsm_h%asdir(m,ind_type) =                           &
    2280                                              albedo_pars_f%pars_xy(2,j,i)     
     2283                                             albedo_pars_f%pars_xy(2,j,i)
    22812284                   IF ( albedo_pars_f%pars_xy(2,j,i) /= albedo_pars_f%fill )   &
    22822285                      surf_lsm_h%asdif(m,ind_type) =                           &
     
    22852288             ENDDO
    22862289!
    2287 !--          For urban surface only if albedo has not been already initialized 
     2290!--          For urban surface only if albedo has not been already initialized
    22882291!--          in the urban-surface model via the ASCII file.
    22892292             IF ( .NOT. surf_usm_h%albedo_from_ascii )  THEN
     
    23782381                ENDDO
    23792382!
    2380 !--             For urban surface only if albedo has not been already initialized 
     2383!--             For urban surface only if albedo has not been already initialized
    23812384!--             in the urban-surface model via the ASCII file.
    23822385                IF ( .NOT. surf_usm_v(l)%albedo_from_ascii )  THEN
     
    24112414                                         albedo_pars_f%pars_xy(2,j+joff,i+ioff)
    24122415                      ENDIF
    2413 !                     
     2416!
    24142417!--                   Spectral albedos especially for building green surfaces
    24152418                      IF ( albedo_pars_f%pars_xy(3,j+joff,i+ioff) /=           &
     
    24272430                                         albedo_pars_f%pars_xy(4,j+joff,i+ioff)
    24282431                      ENDIF
    2429 !                     
     2432!
    24302433!--                   Spectral albedos especially for building window surfaces
    24312434                      IF ( albedo_pars_f%pars_xy(5,j+joff,i+ioff) /=           &
     
    26402643
    26412644!
    2642 !--       Calculate initial values of current (cosine of) the zenith angle and 
     2645!--       Calculate initial values of current (cosine of) the zenith angle and
    26432646!--       whether the sun is up
    26442647          CALL get_date_time( time_since_reference_point, &
     
    27852788
    27862789!
    2787 !--       Allocate 1-element array for surface temperature 
     2790!--       Allocate 1-element array for surface temperature
    27882791!--       (RRTMG anticipates an array as passed argument).
    27892792          ALLOCATE ( rrtm_tsfc(1) )
    27902793!
    2791 !--       Allocate surface emissivity. 
     2794!--       Allocate surface emissivity.
    27922795!--       Values will be given directly before calling rrtm_lw.
    27932796          ALLOCATE ( rrtm_emis(0:0,1:nbndlw+1) )
     
    28012804                            '&Please provide <jobname>_lsw file in the INPUT directory.'
    28022805             CALL message( 'radiation_init', 'PA0583', 1, 2, 0, 6, 0 )
    2803           ENDIF         
     2806          ENDIF
    28042807          INQUIRE( FILE='rrtmg_sw.nc', EXIST=sw_exists )
    28052808          IF ( .NOT. sw_exists )  THEN
     
    28082811                            '&Please provide <jobname>_rsw file in the INPUT directory.'
    28092812             CALL message( 'radiation_init', 'PA0584', 1, 2, 0, 6, 0 )
    2810           ENDIF         
    2811          
     2813          ENDIF
     2814
    28122815          IF ( lw_radiation )  CALL rrtmg_lw_ini ( c_p )
    28132816          IF ( sw_radiation )  CALL rrtmg_sw_ini ( c_p )
    2814          
     2817
    28152818!
    28162819!--       Set input files for RRTMG
    2817           INQUIRE(FILE="RAD_SND_DATA", EXIST=snd_exists) 
     2820          INQUIRE(FILE="RAD_SND_DATA", EXIST=snd_exists)
    28182821          IF ( .NOT. snd_exists )  THEN
    28192822             rrtm_input_file = "rrtmg_lw.nc"
     
    28232826!--       Read vertical layers for RRTMG from sounding data
    28242827!--       The routine provides nzt_rad, hyp_snd(1:nzt_rad),
    2825 !--       t_snd(nzt+2:nzt_rad), rrtm_play(1:nzt_rad), rrtm_plev(1_nzt_rad+1), 
     2828!--       t_snd(nzt+2:nzt_rad), rrtm_play(1:nzt_rad), rrtm_plev(1_nzt_rad+1),
    28262829!--       rrtm_tlay(nzt+2:nzt_rad), rrtm_tlev(nzt+2:nzt_rad+1)
    28272830          CALL read_sounding_data
     
    28342837       ENDIF
    28352838!
    2836 !--    Initializaion actions exclusively required for external 
     2839!--    Initializaion actions exclusively required for external
    28372840!--    radiation forcing
    28382841       IF ( radiation_scheme == 'external' )  THEN
    28392842!
    2840 !--       Open the radiation input file. Note, for child domain, a dynamic 
    2841 !--       input file is often not provided. In order to do not need to 
     2843!--       Open the radiation input file. Note, for child domain, a dynamic
     2844!--       input file is often not provided. In order to do not need to
    28422845!--       duplicate the dynamic input file just for the radiation input, take
    2843 !--       it from the dynamic file for the parent if not available for the 
    2844 !--       child domain(s). In this case this is possible because radiation 
     2846!--       it from the dynamic file for the parent if not available for the
     2847!--       child domain(s). In this case this is possible because radiation
    28452848!--       input should be the same for each model.
    28462849          INQUIRE( FILE = TRIM( input_file_dynamic ),                          &
    28472850                   EXIST = radiation_input_root_domain  )
    2848                    
     2851
    28492852          IF ( .NOT. input_pids_dynamic  .AND.                                 &
    28502853               .NOT. radiation_input_root_domain )  THEN
     
    28592862!
    28602863!--       Open dynamic input file for child domain if available, else, open
    2861 !--       dynamic input file for the root domain. 
     2864!--       dynamic input file for the root domain.
    28622865          IF ( input_pids_dynamic )  THEN
    28632866             CALL open_read_file( TRIM( input_file_dynamic ) //                &
     
    28682871                                  pids_id )
    28692872          ENDIF
    2870                                
     2873
    28712874          CALL inquire_num_variables( pids_id, num_var_pids )
    2872 !         
     2875!
    28732876!--       Allocate memory to store variable names and read them
    28742877          ALLOCATE( vars_pids(1:num_var_pids) )
    28752878          CALL inquire_variable_names( pids_id, vars_pids )
    2876 !         
     2879!
    28772880!--       Input time dimension.
    28782881          IF ( check_existence( vars_pids, 'time_rad' ) )  THEN
    28792882             CALL get_dimension_length( pids_id, ntime, 'time_rad' )
    2880          
     2883
    28812884             ALLOCATE( time_rad_f%var1d(0:ntime-1) )
    2882 !                                                                                 
    2883 !--          Read variable                   
     2885!
     2886!--          Read variable
    28842887             CALL get_variable( pids_id, 'time_rad', time_rad_f%var1d )
    2885                                
     2888
    28862889             time_rad_f%from_file = .TRUE.
    2887           ENDIF           
    2888 !         
     2890          ENDIF
     2891!
    28892892!--       Input shortwave downwelling.
    28902893          IF ( check_existence( vars_pids, 'rad_sw_in' ) )  THEN
    2891 !         
     2894!
    28922895!--          Get _FillValue attribute
    28932896             CALL get_attribute( pids_id, char_fill, rad_sw_in_f%fill,         &
    28942897                                 .FALSE., 'rad_sw_in' )
    2895 !         
     2898!
    28962899!--          Get level-of-detail
    28972900             CALL get_attribute( pids_id, char_lod, rad_sw_in_f%lod,           &
     
    29052908!
    29062909!--          Level-of-detail 2 - radiation depends on time_rad, y, x
    2907              ELSEIF ( rad_sw_in_f%lod == 2 )  THEN 
     2910             ELSEIF ( rad_sw_in_f%lod == 2 )  THEN
    29082911                ALLOCATE( rad_sw_in_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) )
    2909                
     2912
    29102913                CALL get_variable( pids_id, 'rad_sw_in', rad_sw_in_f%var3d,    &
    29112914                                   nxl, nxr, nys, nyn, 0, ntime-1 )
    2912                                    
     2915
    29132916                rad_sw_in_f%from_file = .TRUE.
    2914              ELSE 
     2917             ELSE
    29152918                message_string = '"rad_sw_in" has no valid lod attribute'
    29162919                CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 )
    29172920             ENDIF
    29182921          ENDIF
    2919 !         
     2922!
    29202923!--       Input longwave downwelling.
    29212924          IF ( check_existence( vars_pids, 'rad_lw_in' ) )  THEN
    2922 !         
     2925!
    29232926!--          Get _FillValue attribute
    29242927             CALL get_attribute( pids_id, char_fill, rad_lw_in_f%fill,         &
    29252928                                 .FALSE., 'rad_lw_in' )
    2926 !         
     2929!
    29272930!--          Get level-of-detail
    29282931             CALL get_attribute( pids_id, char_lod, rad_lw_in_f%lod,           &
     
    29362939!
    29372940!--          Level-of-detail 2 - radiation depends on time_rad, y, x
    2938              ELSEIF ( rad_lw_in_f%lod == 2 )  THEN 
     2941             ELSEIF ( rad_lw_in_f%lod == 2 )  THEN
    29392942                ALLOCATE( rad_lw_in_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) )
    2940                
     2943
    29412944                CALL get_variable( pids_id, 'rad_lw_in', rad_lw_in_f%var3d,    &
    29422945                                   nxl, nxr, nys, nyn, 0, ntime-1 )
    2943                                    
     2946
    29442947                rad_lw_in_f%from_file = .TRUE.
    2945              ELSE 
     2948             ELSE
    29462949                message_string = '"rad_lw_in" has no valid lod attribute'
    29472950                CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 )
    29482951             ENDIF
    29492952          ENDIF
    2950 !         
     2953!
    29512954!--       Input shortwave downwelling, diffuse part.
    29522955          IF ( check_existence( vars_pids, 'rad_sw_in_dif' ) )  THEN
    2953 !         
     2956!
    29542957!--          Read _FillValue attribute
    29552958             CALL get_attribute( pids_id, char_fill, rad_sw_in_dif_f%fill,     &
    29562959                                 .FALSE., 'rad_sw_in_dif' )
    2957 !         
     2960!
    29582961!--          Get level-of-detail
    29592962             CALL get_attribute( pids_id, char_lod, rad_sw_in_dif_f%lod,       &
     
    29682971!
    29692972!--          Level-of-detail 2 - radiation depends on time_rad, y, x
    2970              ELSEIF ( rad_sw_in_dif_f%lod == 2 )  THEN 
     2973             ELSEIF ( rad_sw_in_dif_f%lod == 2 )  THEN
    29712974                ALLOCATE( rad_sw_in_dif_f%var3d(0:ntime-1,nys:nyn,nxl:nxr) )
    2972                
     2975
    29732976                CALL get_variable( pids_id, 'rad_sw_in_dif',                   &
    29742977                                   rad_sw_in_dif_f%var3d,                      &
    29752978                                   nxl, nxr, nys, nyn, 0, ntime-1 )
    2976                                    
     2979
    29772980                rad_sw_in_dif_f%from_file = .TRUE.
    2978              ELSE 
     2981             ELSE
    29792982                message_string = '"rad_sw_in_dif" has no valid lod attribute'
    29802983                CALL message( 'radiation_init', 'PA0646', 1, 2, 0, 6, 0 )
    29812984             ENDIF
    29822985          ENDIF
    2983 !         
     2986!
    29842987!--       Finally, close the input file and deallocate temporary arrays
    29852988          DEALLOCATE( vars_pids )
    2986          
     2989
    29872990          CALL close_input_file( pids_id )
    29882991#endif
     
    29952998             CALL message( 'radiation_init', 'PA0195', 1, 2, 0, 6, 0 )
    29962999          ENDIF
    2997          
     3000
    29983001          IF ( .NOT. time_rad_f%from_file )  THEN
    29993002             message_string = 'In case of external radiation forcing ' //      &
     
    30013004             CALL message( 'radiation_init', 'PA0196', 1, 2, 0, 6, 0 )
    30023005          ENDIF
    3003          
     3006
    30043007          CALL get_date_time( 0.0_wp, second_of_day=second_of_day )
    30053008
     
    30183021             ENDIF
    30193022          ENDIF
    3020          
     3023
    30213024          IF ( ALLOCATED( rad_lw_in_f%var1d ) )  THEN
    30223025             IF ( ANY( rad_lw_in_f%var1d == rad_lw_in_f%fill ) )  THEN
     
    30263029             ENDIF
    30273030          ENDIF
    3028          
     3031
    30293032          IF ( ALLOCATED( rad_sw_in_dif_f%var1d ) )  THEN
    30303033             IF ( ANY( rad_sw_in_dif_f%var1d == rad_sw_in_dif_f%fill ) )  THEN
     
    30343037             ENDIF
    30353038          ENDIF
    3036          
     3039
    30373040          IF ( ALLOCATED( rad_sw_in_f%var3d ) )  THEN
    30383041             IF ( ANY( rad_sw_in_f%var3d == rad_sw_in_f%fill ) )  THEN
     
    30423045             ENDIF
    30433046          ENDIF
    3044          
     3047
    30453048          IF ( ALLOCATED( rad_lw_in_f%var3d ) )  THEN
    30463049             IF ( ANY( rad_lw_in_f%var3d == rad_lw_in_f%fill ) )  THEN
     
    30503053             ENDIF
    30513054          ENDIF
    3052          
     3055
    30533056          IF ( ALLOCATED( rad_sw_in_dif_f%var3d ) )  THEN
    30543057             IF ( ANY( rad_sw_in_dif_f%var3d == rad_sw_in_dif_f%fill ) )  THEN
     
    30683071          ENDIF
    30693072!
    3070 !--       All radiation input should have the same level of detail. The sum 
     3073!--       All radiation input should have the same level of detail. The sum
    30713074!--       of lods divided by the number of available radiation arrays must be
    3072 !--       1 (if all are lod = 1) or 2 (if all are lod = 2). 
     3075!--       1 (if all are lod = 1) or 2 (if all are lod = 2).
    30733076          IF ( REAL( MERGE( rad_lw_in_f%lod, 0, rad_lw_in_f%from_file ) +       &
    30743077                     MERGE( rad_sw_in_f%lod, 0, rad_sw_in_f%from_file ) +       &
     
    31093112          CASE ( 'constant' )
    31103113             CALL radiation_constant
    3111              
     3114
    31123115          CASE ( 'external' )
    31133116!
     
    31283131
    31293132!
    3130 !--    If required, read or calculate and write out the SVF 
     3133!--    If required, read or calculate and write out the SVF
    31313134       IF ( radiation_interactions .AND. read_svf)  THEN
    31323135!
     
    31473150
    31483151!
    3149 !--    Adjust radiative fluxes. In case of urban and land surfaces, also 
     3152!--    Adjust radiative fluxes. In case of urban and land surfaces, also
    31503153!--    call an initial interaction.
    31513154       IF ( radiation_interactions )  THEN
     
    31713174       INTEGER(iwp) ::  l   !< running index for surface orientation
    31723175       INTEGER(iwp) ::  t   !< index of current timestep
    3173        INTEGER(iwp) ::  tm  !< index of previous timestep 
    3174        
     3176       INTEGER(iwp) ::  tm  !< index of previous timestep
     3177
    31753178       LOGICAL      ::  horizontal !< flag indicating treatment of horinzontal surfaces
    3176        
    3177        REAL(wp) ::  fac_dt               !< interpolation factor 
     3179
     3180       REAL(wp) ::  fac_dt               !< interpolation factor
    31783181       REAL(wp) ::  second_of_day_init   !< second of the day at model start
    31793182
    3180        TYPE(surf_type), POINTER ::  surf !< pointer on respective surface type, used to generalize routine   
     3183       TYPE(surf_type), POINTER ::  surf !< pointer on respective surface type, used to generalize routine
    31813184
    31823185!
     
    31873190       CALL calc_zenith( day_of_year, second_of_day )
    31883191!
    3189 !--    Interpolate external radiation on current timestep 
     3192!--    Interpolate external radiation on current timestep
    31903193       IF ( time_since_reference_point  <= 0.0_wp )  THEN
    31913194          t      = 0
     
    31983201             t = t + 1
    31993202          ENDDO
    3200          
     3203
    32013204          tm = MAX( t-1, 0 )
    3202          
     3205
    32033206          fac_dt = ( time_since_reference_point                                &
    32043207                   - time_rad_f%var1d(tm) + dt_3d )                            &
     
    32073210       ENDIF
    32083211!
    3209 !--    Call clear-sky calculation for each surface orientation. 
     3212!--    Call clear-sky calculation for each surface orientation.
    32103213!--    First, horizontal surfaces
    32113214       horizontal = .TRUE.
     
    32233226          CALL radiation_external_surf
    32243227       ENDDO
    3225        
     3228
    32263229       CONTAINS
    32273230
     
    32293232
    32303233             USE control_parameters
    3231          
     3234
    32323235             IMPLICIT NONE
    32333236
    3234              INTEGER(iwp) ::  i    !< grid index along x-dimension   
    3235              INTEGER(iwp) ::  j    !< grid index along y-dimension 
    3236              INTEGER(iwp) ::  k    !< grid index along z-dimension   
     3237             INTEGER(iwp) ::  i    !< grid index along x-dimension
     3238             INTEGER(iwp) ::  j    !< grid index along y-dimension
     3239             INTEGER(iwp) ::  k    !< grid index along z-dimension
    32373240             INTEGER(iwp) ::  m    !< running index for surface elements
    3238              
    3239              REAL(wp) ::  lw_in     !< downwelling longwave radiation, interpolated value     
     3241
     3242             REAL(wp) ::  lw_in     !< downwelling longwave radiation, interpolated value
    32403243             REAL(wp) ::  sw_in     !< downwelling shortwave radiation, interpolated value
    3241              REAL(wp) ::  sw_in_dif !< downwelling diffuse shortwave radiation, interpolated value   
     3244             REAL(wp) ::  sw_in_dif !< downwelling diffuse shortwave radiation, interpolated value
    32423245
    32433246             IF ( surf%ns < 1 )  RETURN
     
    32453248!--          level-of-detail = 1. Note, here it must be distinguished between
    32463249!--          averaged radiation and non-averaged radiation for the upwelling
    3247 !--          fluxes. 
     3250!--          fluxes.
    32483251             IF ( rad_sw_in_f%lod == 1 )  THEN
    3249              
     3252
    32503253                sw_in = ( 1.0_wp - fac_dt ) * rad_sw_in_f%var1d(tm)            &
    32513254                                   + fac_dt * rad_sw_in_f%var1d(t)
    3252                                          
     3255
    32533256                lw_in = ( 1.0_wp - fac_dt ) * rad_lw_in_f%var1d(tm)            &
    32543257                                   + fac_dt * rad_lw_in_f%var1d(t)
    32553258!
    3256 !--             Limit shortwave incoming radiation to positive values, in order 
     3259!--             Limit shortwave incoming radiation to positive values, in order
    32573260!--             to overcome possible observation errors.
    32583261                sw_in = MAX( 0.0_wp, sw_in )
    32593262                sw_in = MERGE( sw_in, 0.0_wp, sun_up )
    3260                          
    3261                 surf%rad_sw_in = sw_in                                         
     3263
     3264                surf%rad_sw_in = sw_in
    32623265                surf%rad_lw_in = lw_in
    3263              
     3266
    32643267                IF ( average_radiation )  THEN
    32653268                   surf%rad_sw_out = albedo_urb * surf%rad_sw_in
    3266                    
     3269
    32673270                   surf%rad_lw_out = emissivity_urb * sigma_sb * t_rad_urb**4  &
    32683271                                  + ( 1.0_wp - emissivity_urb ) * surf%rad_lw_in
    3269                    
     3272
    32703273                   surf%rad_net = surf%rad_sw_in - surf%rad_sw_out             &
    32713274                                + surf%rad_lw_in - surf%rad_lw_out
    3272                    
     3275
    32733276                   surf%rad_lw_out_change_0 = 4.0_wp * emissivity_urb          &
    32743277                                                     * sigma_sb                &
     
    32843287                                             surf%albedo(m,ind_wat_win) )      &
    32853288                                           * surf%rad_sw_in(m)
    3286                    
     3289
    32873290                      surf%rad_lw_out(m) = ( surf%frac(m,ind_veg_wall)  *      &
    32883291                                             surf%emissivity(m,ind_veg_wall)   &
     
    32943297                                           * sigma_sb                          &
    32953298                                           * ( surf%pt_surface(m) * exner(k) )**4
    3296                    
     3299
    32973300                      surf%rad_lw_out_change_0(m) =                            &
    32983301                                         ( surf%frac(m,ind_veg_wall)  *        &
     
    33053308                                         * ( surf%pt_surface(m) * exner(k) )**3
    33063309                   ENDDO
    3307                
     3310
    33083311                ENDIF
    33093312!
    3310 !--             If diffuse shortwave radiation is available, store it on 
     3313!--             If diffuse shortwave radiation is available, store it on
    33113314!--             the respective files.
    3312                 IF ( rad_sw_in_dif_f%from_file )  THEN 
     3315                IF ( rad_sw_in_dif_f%from_file )  THEN
    33133316                   sw_in_dif= ( 1.0_wp - fac_dt ) * rad_sw_in_dif_f%var1d(tm)  &
    33143317                                         + fac_dt * rad_sw_in_dif_f%var1d(t)
    3315                                          
     3318
    33163319                   IF ( ALLOCATED( rad_sw_in_diff ) )  rad_sw_in_diff = sw_in_dif
    33173320                   IF ( ALLOCATED( rad_sw_in_dir  ) )  rad_sw_in_dir  = sw_in  &
    33183321                                                                    - sw_in_dif
    3319 !             
    3320 !--                Diffuse longwave radiation equals the total downwelling 
     3322!
     3323!--                Diffuse longwave radiation equals the total downwelling
    33213324!--                longwave radiation
    33223325                   IF ( ALLOCATED( rad_lw_in_diff ) )  rad_lw_in_diff = lw_in
     
    33303333                   j = surf%j(m)
    33313334                   k = surf%k(m)
    3332                    
     3335
    33333336                   surf%rad_sw_in(m) = ( 1.0_wp - fac_dt )                     &
    33343337                                            * rad_sw_in_f%var3d(tm,j,i)        &
    3335                                    + fac_dt * rad_sw_in_f%var3d(t,j,i)                                   
    3336 !
    3337 !--                Limit shortwave incoming radiation to positive values, in 
     3338                                   + fac_dt * rad_sw_in_f%var3d(t,j,i)
     3339!
     3340!--                Limit shortwave incoming radiation to positive values, in
    33383341!--                order to overcome possible observation errors.
    33393342                   surf%rad_sw_in(m) = MAX( 0.0_wp, surf%rad_sw_in(m) )
    33403343                   surf%rad_sw_in(m) = MERGE( surf%rad_sw_in(m), 0.0_wp, sun_up )
    3341                                          
     3344
    33423345                   surf%rad_lw_in(m) = ( 1.0_wp - fac_dt )                     &
    33433346                                            * rad_lw_in_f%var3d(tm,j,i)        &
    3344                                    + fac_dt * rad_lw_in_f%var3d(t,j,i) 
    3345 !
    3346 !--                Weighted average according to surface fraction. 
     3347                                   + fac_dt * rad_lw_in_f%var3d(t,j,i)
     3348!
     3349!--                Weighted average according to surface fraction.
    33473350                   surf%rad_sw_out(m) = ( surf%frac(m,ind_veg_wall)  *         &
    33483351                                          surf%albedo(m,ind_veg_wall)          &
     
    33763379                                   + surf%rad_lw_in(m) - surf%rad_lw_out(m)
    33773380!
    3378 !--                If diffuse shortwave radiation is available, store it on 
    3379 !--                the respective files. 
     3381!--                If diffuse shortwave radiation is available, store it on
     3382!--                the respective files.
    33803383                   IF ( rad_sw_in_dif_f%from_file )  THEN
    33813384                      IF ( ALLOCATED( rad_sw_in_diff ) )                       &
     
    33843387                                     + fac_dt * rad_sw_in_dif_f%var3d(t,j,i)
    33853388!
    3386 !--                   dir = sw_in - sw_in_dif. 
     3389!--                   dir = sw_in - sw_in_dif.
    33873390                      IF ( ALLOCATED( rad_sw_in_dir  ) )                       &
    33883391                         rad_sw_in_dir(j,i)  = surf%rad_sw_in(m) -             &
    33893392                                               rad_sw_in_diff(j,i)
    3390 !                 
    3391 !--                   Diffuse longwave radiation equals the total downwelling 
     3393!
     3394!--                   Diffuse longwave radiation equals the total downwelling
    33923395!--                   longwave radiation
    33933396                      IF ( ALLOCATED( rad_lw_in_diff ) )                       &
     
    33993402             ENDIF
    34003403!
    3401 !--          Store radiation also on 2D arrays, which are still used for 
     3404!--          Store radiation also on 2D arrays, which are still used for
    34023405!--          direct-diffuse splitting. Note, this is only required
    34033406!--          for horizontal surfaces, which covers all x,y position.
     
    34063409                   i = surf%i(m)
    34073410                   j = surf%j(m)
    3408                    
     3411
    34093412                   rad_sw_in(0,j,i)  = surf%rad_sw_in(m)
    34103413                   rad_lw_in(0,j,i)  = surf%rad_lw_in(m)
     
    34133416                ENDDO
    34143417             ENDIF
    3415  
     3418
    34163419          END SUBROUTINE radiation_external_surf
    34173420
    3418     END SUBROUTINE radiation_external   
    3419    
     3421    END SUBROUTINE radiation_external
     3422
    34203423!------------------------------------------------------------------------------!
    34213424! Description:
     
    34313434       LOGICAL      ::  horizontal !< flag indicating treatment of horinzontal surfaces
    34323435
    3433        REAL(wp)     ::  pt1       !< potential temperature at first grid level or mean value at urban layer top 
    3434        REAL(wp)     ::  pt1_l     !< potential temperature at first grid level or mean value at urban layer top at local subdomain 
    3435        REAL(wp)     ::  ql1       !< liquid water mixing ratio at first grid level or mean value at urban layer top 
    3436        REAL(wp)     ::  ql1_l     !< liquid water mixing ratio at first grid level or mean value at urban layer top at local subdomain 
    3437 
    3438        TYPE(surf_type), POINTER ::  surf !< pointer on respective surface type, used to generalize routine   
     3436       REAL(wp)     ::  pt1       !< potential temperature at first grid level or mean value at urban layer top
     3437       REAL(wp)     ::  pt1_l     !< potential temperature at first grid level or mean value at urban layer top at local subdomain
     3438       REAL(wp)     ::  ql1       !< liquid water mixing ratio at first grid level or mean value at urban layer top
     3439       REAL(wp)     ::  ql1_l     !< liquid water mixing ratio at first grid level or mean value at urban layer top at local subdomain
     3440
     3441       TYPE(surf_type), POINTER ::  surf !< pointer on respective surface type, used to generalize routine
    34393442
    34403443!
     
    34523455!--    Calculate value of the Exner function at model surface
    34533456!
    3454 !--    In case averaged radiation is used, calculate mean temperature and 
     3457!--    In case averaged radiation is used, calculate mean temperature and
    34553458!--    liquid water mixing ratio at the urban-layer top.
    34563459       IF ( average_radiation ) THEN
     
    34613464          IF ( bulk_cloud_model  .OR.  cloud_droplets  )  ql1_l = SUM( ql(nz_urban_t,nys:nyn,nxl:nxr) )
    34623465
    3463 #if defined( __parallel )     
     3466#if defined( __parallel )
    34643467          IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
    34653468          CALL MPI_ALLREDUCE( pt1_l, pt1, 1, MPI_REAL, MPI_SUM, comm2d, ierr )
     
    34773480          ENDIF
    34783481#else
    3479           pt1 = pt1_l 
     3482          pt1 = pt1_l
    34803483          IF ( bulk_cloud_model  .OR.  cloud_droplets )  ql1 = ql1_l
    34813484#endif
     
    34873490       ENDIF
    34883491!
    3489 !--    Call clear-sky calculation for each surface orientation. 
     3492!--    Call clear-sky calculation for each surface orientation.
    34903493!--    First, horizontal surfaces
    34913494       horizontal = .TRUE.
     
    35183521
    35193522!
    3520 !--          Calculate radiation fluxes and net radiation (rad_net) assuming 
    3521 !--          homogeneous urban radiation conditions. 
    3522              IF ( average_radiation ) THEN       
     3523!--          Calculate radiation fluxes and net radiation (rad_net) assuming
     3524!--          homogeneous urban radiation conditions.
     3525             IF ( average_radiation ) THEN
    35233526
    35243527                k = nz_urban_t
     
    35263529                surf%rad_sw_in  = solar_constant * sky_trans * cos_zenith
    35273530                surf%rad_sw_out = albedo_urb * surf%rad_sw_in
    3528                
     3531
    35293532                surf%rad_lw_in  = emissivity_atm_clsky * sigma_sb * (pt1 * exner(k+1))**4
    35303533
     
    35393542
    35403543!
    3541 !--          Calculate radiation fluxes and net radiation (rad_net) for each surface 
     3544!--          Calculate radiation fluxes and net radiation (rad_net) for each surface
    35423545!--          element.
    35433546             ELSE
     
    35513554
    35523555!
    3553 !--                Weighted average according to surface fraction. 
    3554 !--                ATTENTION: when radiation interactions are switched on the 
    3555 !--                calculated fluxes below are not actually used as they are 
     3556!--                Weighted average according to surface fraction.
     3557!--                ATTENTION: when radiation interactions are switched on the
     3558!--                calculated fluxes below are not actually used as they are
    35563559!--                overwritten in radiation_interaction.
    35573560                   surf%rad_sw_out(m) = ( surf%frac(m,ind_veg_wall)  *         &
     
    36113614                ENDDO
    36123615             ENDIF
    3613  
     3616
    36143617          END SUBROUTINE radiation_clearsky_surf
    36153618
     
    36283631
    36293632       INTEGER(iwp) ::  l         !< running index for surface orientation
    3630        
     3633
    36313634       LOGICAL      ::  horizontal !< flag indicating treatment of horinzontal surfaces
    36323635
    3633        REAL(wp)     ::  pt1       !< potential temperature at first grid level or mean value at urban layer top 
    3634        REAL(wp)     ::  pt1_l     !< potential temperature at first grid level or mean value at urban layer top at local subdomain 
    3635        REAL(wp)     ::  ql1       !< liquid water mixing ratio at first grid level or mean value at urban layer top 
    3636        REAL(wp)     ::  ql1_l     !< liquid water mixing ratio at first grid level or mean value at urban layer top at local subdomain 
     3636       REAL(wp)     ::  pt1       !< potential temperature at first grid level or mean value at urban layer top
     3637       REAL(wp)     ::  pt1_l     !< potential temperature at first grid level or mean value at urban layer top at local subdomain
     3638       REAL(wp)     ::  ql1       !< liquid water mixing ratio at first grid level or mean value at urban layer top
     3639       REAL(wp)     ::  ql1_l     !< liquid water mixing ratio at first grid level or mean value at urban layer top at local subdomain
    36373640
    36383641       TYPE(surf_type), POINTER ::  surf !< pointer on respective surface type, used to generalize routine
    36393642
    36403643!
    3641 !--    In case averaged radiation is used, calculate mean temperature and 
     3644!--    In case averaged radiation is used, calculate mean temperature and
    36423645!--    liquid water mixing ratio at the urban-layer top.
    3643        IF ( average_radiation ) THEN   
     3646       IF ( average_radiation ) THEN
    36443647          pt1   = 0.0_wp
    36453648          IF ( bulk_cloud_model  .OR.  cloud_droplets )  ql1   = 0.0_wp
     
    36483651          IF ( bulk_cloud_model  .OR.  cloud_droplets )  ql1_l = SUM( ql(nz_urban_t,nys:nyn,nxl:nxr) )
    36493652
    3650 #if defined( __parallel )     
     3653#if defined( __parallel )
    36513654          IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
    36523655          CALL MPI_ALLREDUCE( pt1_l, pt1, 1, MPI_REAL, MPI_SUM, comm2d, ierr )
     
    37293732             ELSE
    37303733!
    3731 !--             Determine index offset between surface element and adjacent 
     3734!--             Determine index offset between surface element and adjacent
    37323735!--             atmospheric grid point
    37333736                ioff = surf%ioff
     
    37533756
    37543757!
    3755 !--                Weighted average according to surface fraction. 
     3758!--                Weighted average according to surface fraction.
    37563759                   surf%rad_lw_out(m) = ( surf%frac(m,ind_veg_wall)  *         &
    37573760                                          surf%emissivity(m,ind_veg_wall)      &
     
    38023805
    38033806          END SUBROUTINE radiation_constant_surf
    3804          
     3807
    38053808
    38063809    END SUBROUTINE radiation_constant
     
    38153818
    38163819       IMPLICIT NONE
    3817  
     3820
    38183821       INTEGER(iwp), INTENT(IN) ::  io            !< Unit of the output file
    3819    
    3820 
    3821        
     3822
     3823
     3824
    38223825!
    38233826!--    Write radiation model header
     
    38343837       ELSEIF ( radiation_scheme == "external" )  THEN
    38353838          WRITE( io, 14 )
    3836        ENDIF 
     3839       ENDIF
    38373840
    38383841       IF ( albedo_type_f%from_file  .OR.  vegetation_type_f%from_file  .OR.   &
     
    38403843            building_type_f%from_file )  THEN
    38413844             WRITE( io, 13 )
    3842        ELSE 
     3845       ELSE
    38433846          IF ( albedo_type == 0 )  THEN
    38443847             WRITE( io, 7 ) albedo
     
    38503853          WRITE( io, 9 )
    38513854       ENDIF
    3852        
     3855
    38533856       WRITE( io, 12 ) dt_radiation
    3854  
     3857
    38553858
    38563859 3 FORMAT (//' Radiation model information:'/                                  &
     
    38733876
    38743877    END SUBROUTINE radiation_header
    3875    
     3878
    38763879
    38773880!------------------------------------------------------------------------------!
     
    38853888       IMPLICIT NONE
    38863889
    3887        CHARACTER (LEN=80) ::  line  !< dummy string that contains the current line of the parameter file 
    3888        
     3890       CHARACTER (LEN=80) ::  line  !< dummy string that contains the current line of the parameter file
     3891
    38893892       NAMELIST /radiation_par/   albedo, albedo_lw_dif, albedo_lw_dir,         &
    38903893                                  albedo_sw_dif, albedo_sw_dir, albedo_type,    &
     
    39033906                                  unscheduled_radiation_calls
    39043907
    3905    
     3908
    39063909       NAMELIST /radiation_parameters/ albedo, albedo_lw_dif, albedo_lw_dir,    &
    39073910                                  albedo_sw_dif, albedo_sw_dir, albedo_type,    &
     
    39193922                                  svfnorm_report_thresh, sw_radiation,          &
    39203923                                  unscheduled_radiation_calls
    3921    
     3924
    39223925       line = ' '
    3923        
     3926
    39243927!
    39253928!--    Try to find radiation model namelist
     
    39793982
    39803983 14    CONTINUE
    3981        
     3984
    39823985    END SUBROUTINE radiation_parin
    39833986
     
    40394042       zenith(0) = cos_zenith
    40404043!
    4041 !--    Calculate surface albedo. In case average radiation is applied, 
     4044!--    Calculate surface albedo. In case average radiation is applied,
    40424045!--    this is not required.
    40434046#if defined( __netcdf )
     
    40714074       IF ( average_radiation ) THEN
    40724075!
    4073 !--       Determine minimum topography top index. 
     4076!--       Determine minimum topography top index.
    40744077          k_topo_l = MINVAL( topo_top_ind(nys:nyn,nxl:nxr,0) )
    40754078#if defined( __parallel )
     
    40794082          k_topo = k_topo_l
    40804083#endif
    4081        
     4084
    40824085          rrtm_asdir(1)  = albedo_urb
    40834086          rrtm_asdif(1)  = albedo_urb
     
    40874090          rrtm_emis = emissivity_urb
    40884091!
    4089 !--       Calculate mean pt profile. 
     4092!--       Calculate mean pt profile.
    40904093          CALL calc_mean_profile( pt, 4 )
    40914094          pt_av = hom(:, 1, 4, 0)
    4092          
     4095
    40934096          IF ( humidity )  THEN
    40944097             CALL calc_mean_profile( q, 41 )
     
    41044107             ! average ql is now in hom(:, 1, 54, 0)
    41054108             ql_av = hom(:, 1, 54, 0)
    4106              
     4109
    41074110             DO k = nzb+1, nzt+1
    41084111                rrtm_tlay(0,k) = pt_av(k) * ( (hyp(k) ) / 100000._wp       &
     
    41264129
    41274130!
    4128 !--       Avoid temperature/humidity jumps at the top of the PALM domain by 
     4131!--       Avoid temperature/humidity jumps at the top of the PALM domain by
    41294132!--       linear interpolation from nzt+2 to nzt+7. Jumps are induced by
    41304133!--       discrepancies between the values in the  domain and those above that
     
    41634166                rrtm_cliqwp(0,k) =  ql_av(k) * 1000._wp *                   &
    41644167                                    (rrtm_plev(0,k) - rrtm_plev(0,k+1))     &
    4165                                     * 100._wp / g 
     4168                                    * 100._wp / g
    41664169
    41674170                IF ( rrtm_cliqwp(0,k) > 0._wp )  THEN
     
    41894192!--       Set surface temperature
    41904193          rrtm_tsfc = t_rad_urb
    4191          
    4192           IF ( lw_radiation )  THEN 
     4194
     4195          IF ( lw_radiation )  THEN
    41934196!
    41944197!--          Due to technical reasons, copy optical depth to dummy arguments
    41954198!--          which are allocated on the exact size as the rrtmg_lw is called.
    41964199!--          As one dimesion is allocated with zero size, compiler complains
    4197 !--          that rank of the array does not match that of the 
    4198 !--          assumed-shaped arguments in the RRTMG library. In order to 
    4199 !--          avoid this, write to dummy arguments and give pass the entire 
    4200 !--          dummy array. Seems to be the only existing work-around. 
     4200!--          that rank of the array does not match that of the
     4201!--          assumed-shaped arguments in the RRTMG library. In order to
     4202!--          avoid this, write to dummy arguments and give pass the entire
     4203!--          dummy array. Seems to be the only existing work-around.
    42014204             ALLOCATE( rrtm_lw_taucld_dum(1:nbndlw+1,0:0,k_topo+1:nzt_rad+1) )
    42024205             ALLOCATE( rrtm_lw_tauaer_dum(0:0,k_topo+1:nzt_rad+1,1:nbndlw+1) )
     
    42064209             rrtm_lw_tauaer_dum =                                              &
    42074210                             rrtm_lw_tauaer(0:0,k_topo+1:nzt_rad+1,1:nbndlw+1)
    4208          
    4209              CALL rrtmg_lw( 1,                                                 &                                       
     4211
     4212             CALL rrtmg_lw( 1,                                                 &
    42104213                            nzt_rad-k_topo,                                    &
    42114214                            rrtm_icld,                                         &
     
    42344237                            rrtm_cicewp(:,k_topo+1:),                 &
    42354238                            rrtm_cliqwp(:,k_topo+1:),                 &
    4236                             rrtm_reice(:,k_topo+1:),                  & 
     4239                            rrtm_reice(:,k_topo+1:),                  &
    42374240                            rrtm_reliq(:,k_topo+1:),                  &
    42384241                            rrtm_lw_tauaer_dum,                                &
     
    42454248                            rrtm_lwuflx_dt(:,k_topo:),                &
    42464249                            rrtm_lwuflxc_dt(:,k_topo:) )
    4247                            
     4250
    42484251             DEALLOCATE ( rrtm_lw_taucld_dum )
    42494252             DEALLOCATE ( rrtm_lw_tauaer_dum )
     
    42564259             rad_lw_in_diff(:,:) = rad_lw_in(k_topo,:,:)
    42574260!
    4258 !--          Save heating rates (convert from K/d to K/h). 
     4261!--          Save heating rates (convert from K/d to K/h).
    42594262!--          Further, even though an aggregated radiation is computed, map
    4260 !--          signle-column profiles on top of any topography, in order to 
     4263!--          signle-column profiles on top of any topography, in order to
    42614264!--          obtain correct near surface radiation heating/cooling rates.
    42624265             DO  i = nxl, nxr
     
    42744277          IF ( sw_radiation .AND. sun_up )  THEN
    42754278!
    4276 !--          Due to technical reasons, copy optical depths and other 
    4277 !--          to dummy arguments which are allocated on the exact size as the 
     4279!--          Due to technical reasons, copy optical depths and other
     4280!--          to dummy arguments which are allocated on the exact size as the
    42784281!--          rrtmg_sw is called.
    42794282!--          As one dimesion is allocated with zero size, compiler complains
    4280 !--          that rank of the array does not match that of the 
    4281 !--          assumed-shaped arguments in the RRTMG library. In order to 
    4282 !--          avoid this, write to dummy arguments and give pass the entire 
    4283 !--          dummy array. Seems to be the only existing work-around. 
     4283!--          that rank of the array does not match that of the
     4284!--          assumed-shaped arguments in the RRTMG library. In order to
     4285!--          avoid this, write to dummy arguments and give pass the entire
     4286!--          dummy array. Seems to be the only existing work-around.
    42844287             ALLOCATE( rrtm_sw_taucld_dum(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1) )
    42854288             ALLOCATE( rrtm_sw_ssacld_dum(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1) )
     
    42904293             ALLOCATE( rrtm_sw_asmaer_dum(0:0,k_topo+1:nzt_rad+1,1:nbndsw+1) )
    42914294             ALLOCATE( rrtm_sw_ecaer_dum(0:0,k_topo+1:nzt_rad+1,1:naerec+1)  )
    4292      
     4295
    42934296             rrtm_sw_taucld_dum = rrtm_sw_taucld(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1)
    42944297             rrtm_sw_ssacld_dum = rrtm_sw_ssacld(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1)
     
    43094312                            rrtm_tlev(:,k_topo+1:nzt_rad+2),                   &
    43104313                            rrtm_tsfc,                                         &
    4311                             rrtm_h2ovmr(:,k_topo+1:nzt_rad+1),                 &                               
    4312                             rrtm_o3vmr(:,k_topo+1:nzt_rad+1),                  &       
     4314                            rrtm_h2ovmr(:,k_topo+1:nzt_rad+1),                 &
     4315                            rrtm_o3vmr(:,k_topo+1:nzt_rad+1),                  &
    43134316                            rrtm_co2vmr(:,k_topo+1:nzt_rad+1),                 &
    43144317                            rrtm_ch4vmr(:,k_topo+1:nzt_rad+1),                 &
    43154318                            rrtm_n2ovmr(:,k_topo+1:nzt_rad+1),                 &
    43164319                            rrtm_o2vmr(:,k_topo+1:nzt_rad+1),                  &
    4317                             rrtm_asdir,                                        & 
     4320                            rrtm_asdir,                                        &
    43184321                            rrtm_asdif,                                        &
    43194322                            rrtm_aldir,                                        &
     
    43394342                            rrtm_sw_asmaer_dum,                                &
    43404343                            rrtm_sw_ecaer_dum,                                 &
    4341                             rrtm_swuflx(:,k_topo:nzt_rad+1),                   & 
    4342                             rrtm_swdflx(:,k_topo:nzt_rad+1),                   & 
    4343                             rrtm_swhr(:,k_topo+1:nzt_rad+1),                   & 
    4344                             rrtm_swuflxc(:,k_topo:nzt_rad+1),                  & 
     4344                            rrtm_swuflx(:,k_topo:nzt_rad+1),                   &
     4345                            rrtm_swdflx(:,k_topo:nzt_rad+1),                   &
     4346                            rrtm_swhr(:,k_topo+1:nzt_rad+1),                   &
     4347                            rrtm_swuflxc(:,k_topo:nzt_rad+1),                  &
    43454348                            rrtm_swdflxc(:,k_topo:nzt_rad+1),                  &
    43464349                            rrtm_swhrc(:,k_topo+1:nzt_rad+1),                  &
    43474350                            rrtm_dirdflux(:,k_topo:nzt_rad+1),                 &
    43484351                            rrtm_difdflux(:,k_topo:nzt_rad+1) )
    4349                            
     4352
    43504353             DEALLOCATE( rrtm_sw_taucld_dum )
    43514354             DEALLOCATE( rrtm_sw_ssacld_dum )
     
    43564359             DEALLOCATE( rrtm_sw_asmaer_dum )
    43574360             DEALLOCATE( rrtm_sw_ecaer_dum )
    4358  
     4361
    43594362!
    43604363!--          Save radiation fluxes for the entire depth of the model domain
     
    43824385          ENDIF
    43834386!
    4384 !--    RRTMG is called for each (j,i) grid point separately, starting at the 
     4387!--    RRTMG is called for each (j,i) grid point separately, starting at the
    43854388!--    highest topography level. Here no RTM is used since average_radiation is false
    43864389       ELSE
     
    44104413                      rrtm_tlay(0,k) = pt(k,j,i) * exner(k)                     &
    44114414                                        + lv_d_cp * ql(k,j,i)
    4412                       rrtm_h2ovmr(0,k) = mol_mass_air_d_wv * q(k,j,i) 
     4415                      rrtm_h2ovmr(0,k) = mol_mass_air_d_wv * q(k,j,i)
    44134416                   ENDDO
    44144417                ELSE
     
    44204423                      DO k = nzb+1, nzt+1
    44214424                         rrtm_h2ovmr(0,k) =  mol_mass_air_d_wv * q(k,j,i)
    4422                       ENDDO   
     4425                      ENDDO
    44234426                   ELSE
    44244427                      rrtm_h2ovmr(0,nzb+1:nzt+1) = 0.0_wp
     
    44274430
    44284431!
    4429 !--             Avoid temperature/humidity jumps at the top of the LES domain by 
     4432!--             Avoid temperature/humidity jumps at the top of the LES domain by
    44304433!--             linear interpolation from nzt+2 to nzt+7
    44314434                DO k = nzt+2, nzt+7
     
    44634466                      rrtm_cliqwp(0,k) =  ql(k,j,i) * 1000.0_wp *              &
    44644467                                          (rrtm_plev(0,k) - rrtm_plev(0,k+1))  &
    4465                                           * 100.0_wp / g 
     4468                                          * 100.0_wp / g
    44664469
    44674470                      IF ( rrtm_cliqwp(0,k) > 0.0_wp )  THEN
     
    44734476                         IF ( bulk_cloud_model )  THEN
    44744477!
    4475 !--                         Calculete effective droplet radius. In case of using 
    4476 !--                         cloud_scheme = 'morrison' and a non reasonable number 
    4477 !--                         of cloud droplets the inital aerosol number 
     4478!--                         Calculete effective droplet radius. In case of using
     4479!--                         cloud_scheme = 'morrison' and a non reasonable number
     4480!--                         of cloud droplets the inital aerosol number
    44784481!--                         concentration is considered.
    44794482                            IF ( microphysics_morrison )  THEN
     
    44824485                               ELSE
    44834486                                  nc_rad = na_init
    4484                                ENDIF 
     4487                               ENDIF
    44854488                            ELSE
    44864489                               nc_rad = nc_const
    4487                             ENDIF 
     4490                            ENDIF
    44884491
    44894492                            rrtm_reliq(0,k) = 1.0E6_wp * ( 3.0_wp * ql(k,j,i)     &
     
    45254528
    45264529!
    4527 !--             Write surface emissivity and surface temperature at current 
     4530!--             Write surface emissivity and surface temperature at current
    45284531!--             surface element on RRTMG-shaped array.
    45294532!--             Please note, as RRTMG is a single column model, surface attributes
    4530 !--             are only obtained from horizontally aligned surfaces (for 
    4531 !--             simplicity). Taking surface attributes from horizontal and 
    4532 !--             vertical walls would lead to multiple solutions. 
     4533!--             are only obtained from horizontally aligned surfaces (for
     4534!--             simplicity). Taking surface attributes from horizontal and
     4535!--             vertical walls would lead to multiple solutions.
    45334536!--             Moreover, for natural- and urban-type surfaces, several surface
    4534 !--             classes can exist at a surface element next to each other. 
    4535 !--             To obtain bulk parameters, apply a weighted average for these 
    4536 !--             surfaces. 
     4537!--             classes can exist at a surface element next to each other.
     4538!--             To obtain bulk parameters, apply a weighted average for these
     4539!--             surfaces.
    45374540                DO  m = surf_lsm_h%start_index(j,i), surf_lsm_h%end_index(j,i)
    45384541                   rrtm_emis = surf_lsm_h%frac(m,ind_veg_wall)  *              &
    45394542                               surf_lsm_h%emissivity(m,ind_veg_wall)  +        &
    45404543                               surf_lsm_h%frac(m,ind_pav_green) *              &
    4541                                surf_lsm_h%emissivity(m,ind_pav_green) +        & 
     4544                               surf_lsm_h%emissivity(m,ind_pav_green) +        &
    45424545                               surf_lsm_h%frac(m,ind_wat_win)   *              &
    45434546                               surf_lsm_h%emissivity(m,ind_wat_win)
    45444547                   rrtm_tsfc = surf_lsm_h%pt_surface(m) * exner(nzb)
    4545                 ENDDO             
     4548                ENDDO
    45464549                DO  m = surf_usm_h%start_index(j,i), surf_usm_h%end_index(j,i)
    45474550                   rrtm_emis = surf_usm_h%frac(m,ind_veg_wall)  *              &
    45484551                               surf_usm_h%emissivity(m,ind_veg_wall)  +        &
    45494552                               surf_usm_h%frac(m,ind_pav_green) *              &
    4550                                surf_usm_h%emissivity(m,ind_pav_green) +        & 
     4553                               surf_usm_h%emissivity(m,ind_pav_green) +        &
    45514554                               surf_usm_h%frac(m,ind_wat_win)   *              &
    45524555                               surf_usm_h%emissivity(m,ind_wat_win)
     
    45624565!--                which are allocated on the exact size as the rrtmg_lw is called.
    45634566!--                As one dimesion is allocated with zero size, compiler complains
    4564 !--                that rank of the array does not match that of the 
    4565 !--                assumed-shaped arguments in the RRTMG library. In order to 
    4566 !--                avoid this, write to dummy arguments and give pass the entire 
    4567 !--                dummy array. Seems to be the only existing work-around. 
     4567!--                that rank of the array does not match that of the
     4568!--                assumed-shaped arguments in the RRTMG library. In order to
     4569!--                avoid this, write to dummy arguments and give pass the entire
     4570!--                dummy array. Seems to be the only existing work-around.
    45684571                   ALLOCATE( rrtm_lw_taucld_dum(1:nbndlw+1,0:0,k_topo+1:nzt_rad+1) )
    45694572                   ALLOCATE( rrtm_lw_tauaer_dum(0:0,k_topo+1:nzt_rad+1,1:nbndlw+1) )
     
    45744577                               rrtm_lw_tauaer(0:0,k_topo+1:nzt_rad+1,1:nbndlw+1)
    45754578
    4576                    CALL rrtmg_lw( 1,                                           &                                       
     4579                   CALL rrtmg_lw( 1,                                           &
    45774580                                  nzt_rad-k_topo,                              &
    45784581                                  rrtm_icld,                                   &
     
    46014604                                  rrtm_cicewp(:,k_topo+1:nzt_rad+1),           &
    46024605                                  rrtm_cliqwp(:,k_topo+1:nzt_rad+1),           &
    4603                                   rrtm_reice(:,k_topo+1:nzt_rad+1),            & 
     4606                                  rrtm_reice(:,k_topo+1:nzt_rad+1),            &
    46044607                                  rrtm_reliq(:,k_topo+1:nzt_rad+1),            &
    46054608                                  rrtm_lw_tauaer_dum,                          &
     
    46304633
    46314634!
    4632 !--                Save surface radiative fluxes and change in LW heating rate 
     4635!--                Save surface radiative fluxes and change in LW heating rate
    46334636!--                onto respective surface elements
    46344637!--                Horizontal surfaces
     
    46384641                      surf_lsm_h%rad_lw_out(m)          = rrtm_lwuflx(0,k_topo)
    46394642                      surf_lsm_h%rad_lw_out_change_0(m) = rrtm_lwuflx_dt(0,k_topo)
    4640                    ENDDO             
     4643                   ENDDO
    46414644                   DO  m = surf_usm_h%start_index(j,i),                        &
    46424645                           surf_usm_h%end_index(j,i)
     
    46444647                      surf_usm_h%rad_lw_out(m)          = rrtm_lwuflx(0,k_topo)
    46454648                      surf_usm_h%rad_lw_out_change_0(m) = rrtm_lwuflx_dt(0,k_topo)
    4646                    ENDDO 
    4647 !
    4648 !--                Vertical surfaces. Fluxes are obtain at vertical level of the 
     4649                   ENDDO
     4650!
     4651!--                Vertical surfaces. Fluxes are obtain at vertical level of the
    46494652!--                respective surface element
    46504653                   DO  l = 0, 3
     
    46554658                         surf_lsm_v(l)%rad_lw_out(m)          = rrtm_lwuflx(0,k)
    46564659                         surf_lsm_v(l)%rad_lw_out_change_0(m) = rrtm_lwuflx_dt(0,k)
    4657                       ENDDO             
     4660                      ENDDO
    46584661                      DO  m = surf_usm_v(l)%start_index(j,i),                  &
    46594662                              surf_usm_v(l)%end_index(j,i)
     
    46624665                         surf_usm_v(l)%rad_lw_out(m)          = rrtm_lwuflx(0,k)
    46634666                         surf_usm_v(l)%rad_lw_out_change_0(m) = rrtm_lwuflx_dt(0,k)
    4664                       ENDDO 
     4667                      ENDDO
    46654668                   ENDDO
    46664669
     
    46694672                IF ( sw_radiation .AND. sun_up )  THEN
    46704673!
    4671 !--                Get albedo for direct/diffusive long/shortwave radiation at 
    4672 !--                current (y,x)-location from surface variables. 
    4673 !--                Only obtain it from horizontal surfaces, as RRTMG is a single 
     4674!--                Get albedo for direct/diffusive long/shortwave radiation at
     4675!--                current (y,x)-location from surface variables.
     4676!--                Only obtain it from horizontal surfaces, as RRTMG is a single
    46744677!--                column model
    4675 !--                (Please note, only one loop will entered, controlled by 
     4678!--                (Please note, only one loop will entered, controlled by
    46764679!--                start-end index.)
    46774680                   DO  m = surf_lsm_h%start_index(j,i),                        &
    46784681                           surf_lsm_h%end_index(j,i)
    46794682                      rrtm_asdir(1)  = SUM( surf_lsm_h%frac(m,:) *             &
    4680                                             surf_lsm_h%rrtm_asdir(:,m) )
     4683                                            surf_lsm_h%rrtm_asdir(m,:) )
    46814684                      rrtm_asdif(1)  = SUM( surf_lsm_h%frac(m,:) *             &
    4682                                             surf_lsm_h%rrtm_asdif(:,m) )
     4685                                            surf_lsm_h%rrtm_asdif(m,:) )
    46834686                      rrtm_aldir(1)  = SUM( surf_lsm_h%frac(m,:) *             &
    4684                                             surf_lsm_h%rrtm_aldir(:,m) )
     4687                                            surf_lsm_h%rrtm_aldir(m,:) )
    46854688                      rrtm_aldif(1)  = SUM( surf_lsm_h%frac(m,:) *             &
    4686                                             surf_lsm_h%rrtm_aldif(:,m) )
    4687                    ENDDO             
     4689                                            surf_lsm_h%rrtm_aldif(m,:) )
     4690                   ENDDO
    46884691                   DO  m = surf_usm_h%start_index(j,i),                        &
    46894692                           surf_usm_h%end_index(j,i)
    46904693                      rrtm_asdir(1)  = SUM( surf_usm_h%frac(m,:) *             &
    4691                                             surf_usm_h%rrtm_asdir(:,m) )
     4694                                            surf_usm_h%rrtm_asdir(m,:) )
    46924695                      rrtm_asdif(1)  = SUM( surf_usm_h%frac(m,:) *             &
    4693                                             surf_usm_h%rrtm_asdif(:,m) )
     4696                                            surf_usm_h%rrtm_asdif(m,:) )
    46944697                      rrtm_aldir(1)  = SUM( surf_usm_h%frac(m,:) *             &
    4695                                             surf_usm_h%rrtm_aldir(:,m) )
     4698                                            surf_usm_h%rrtm_aldir(m,:) )
    46964699                      rrtm_aldif(1)  = SUM( surf_usm_h%frac(m,:) *             &
    4697                                             surf_usm_h%rrtm_aldif(:,m) )
     4700                                            surf_usm_h%rrtm_aldif(m,:) )
    46984701                   ENDDO
    46994702!
    4700 !--                Due to technical reasons, copy optical depths and other 
    4701 !--                to dummy arguments which are allocated on the exact size as the 
     4703!--                Due to technical reasons, copy optical depths and other
     4704!--                to dummy arguments which are allocated on the exact size as the
    47024705!--                rrtmg_sw is called.
    47034706!--                As one dimesion is allocated with zero size, compiler complains
    4704 !--                that rank of the array does not match that of the 
    4705 !--                assumed-shaped arguments in the RRTMG library. In order to 
    4706 !--                avoid this, write to dummy arguments and give pass the entire 
    4707 !--                dummy array. Seems to be the only existing work-around. 
     4707!--                that rank of the array does not match that of the
     4708!--                assumed-shaped arguments in the RRTMG library. In order to
     4709!--                avoid this, write to dummy arguments and give pass the entire
     4710!--                dummy array. Seems to be the only existing work-around.
    47084711                   ALLOCATE( rrtm_sw_taucld_dum(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1) )
    47094712                   ALLOCATE( rrtm_sw_ssacld_dum(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1) )
     
    47144717                   ALLOCATE( rrtm_sw_asmaer_dum(0:0,k_topo+1:nzt_rad+1,1:nbndsw+1) )
    47154718                   ALLOCATE( rrtm_sw_ecaer_dum(0:0,k_topo+1:nzt_rad+1,1:naerec+1)  )
    4716      
     4719
    47174720                   rrtm_sw_taucld_dum = rrtm_sw_taucld(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1)
    47184721                   rrtm_sw_ssacld_dum = rrtm_sw_ssacld(1:nbndsw+1,0:0,k_topo+1:nzt_rad+1)
     
    47334736                                  rrtm_tlev(:,k_topo+1:nzt_rad+2),             &
    47344737                                  rrtm_tsfc,                                   &
    4735                                   rrtm_h2ovmr(:,k_topo+1:nzt_rad+1),           &                               
    4736                                   rrtm_o3vmr(:,k_topo+1:nzt_rad+1),            &       
     4738                                  rrtm_h2ovmr(:,k_topo+1:nzt_rad+1),           &
     4739                                  rrtm_o3vmr(:,k_topo+1:nzt_rad+1),            &
    47374740                                  rrtm_co2vmr(:,k_topo+1:nzt_rad+1),           &
    47384741                                  rrtm_ch4vmr(:,k_topo+1:nzt_rad+1),           &
    47394742                                  rrtm_n2ovmr(:,k_topo+1:nzt_rad+1),           &
    47404743                                  rrtm_o2vmr(:,k_topo+1:nzt_rad+1),            &
    4741                                   rrtm_asdir,                                  & 
     4744                                  rrtm_asdir,                                  &
    47424745                                  rrtm_asdif,                                  &
    47434746                                  rrtm_aldir,                                  &
     
    47634766                                  rrtm_sw_asmaer_dum,                          &
    47644767                                  rrtm_sw_ecaer_dum,                           &
    4765                                   rrtm_swuflx(:,k_topo:nzt_rad+1),             & 
    4766                                   rrtm_swdflx(:,k_topo:nzt_rad+1),             & 
    4767                                   rrtm_swhr(:,k_topo+1:nzt_rad+1),             & 
    4768                                   rrtm_swuflxc(:,k_topo:nzt_rad+1),            & 
     4768                                  rrtm_swuflx(:,k_topo:nzt_rad+1),             &
     4769                                  rrtm_swdflx(:,k_topo:nzt_rad+1),             &
     4770                                  rrtm_swhr(:,k_topo+1:nzt_rad+1),             &
     4771                                  rrtm_swuflxc(:,k_topo:nzt_rad+1),            &
    47694772                                  rrtm_swdflxc(:,k_topo:nzt_rad+1),            &
    47704773                                  rrtm_swhrc(:,k_topo+1:nzt_rad+1),            &
     
    48004803                      surf_lsm_h%rad_sw_in(m)     = rrtm_swdflx(0,k_topo)
    48014804                      surf_lsm_h%rad_sw_out(m)    = rrtm_swuflx(0,k_topo)
    4802                    ENDDO             
     4805                   ENDDO
    48034806                   DO  m = surf_usm_h%start_index(j,i),                        &
    48044807                           surf_usm_h%end_index(j,i)
    48054808                      surf_usm_h%rad_sw_in(m)     = rrtm_swdflx(0,k_topo)
    48064809                      surf_usm_h%rad_sw_out(m)    = rrtm_swuflx(0,k_topo)
    4807                    ENDDO 
    4808 !
    4809 !--                Vertical surfaces. Fluxes are obtain at respective vertical 
     4810                   ENDDO
     4811!
     4812!--                Vertical surfaces. Fluxes are obtain at respective vertical
    48104813!--                level of the surface element
    48114814                   DO  l = 0, 3
     
    48154818                         surf_lsm_v(l)%rad_sw_in(m)  = rrtm_swdflx(0,k)
    48164819                         surf_lsm_v(l)%rad_sw_out(m) = rrtm_swuflx(0,k)
    4817                       ENDDO             
     4820                      ENDDO
    48184821                      DO  m = surf_usm_v(l)%start_index(j,i),                  &
    48194822                              surf_usm_v(l)%end_index(j,i)
     
    48214824                         surf_usm_v(l)%rad_sw_in(m)  = rrtm_swdflx(0,k)
    48224825                         surf_usm_v(l)%rad_sw_out(m) = rrtm_swuflx(0,k)
    4823                       ENDDO 
     4826                      ENDDO
    48244827                   ENDDO
    48254828!
     
    48294832                   rad_sw_out = 0.0_wp
    48304833!--             !!!!!!!! ATTENSION !!!!!!!!!!!!!!!
    4831 !--             Surface radiative fluxes should be also set to zero here                 
     4834!--             Surface radiative fluxes should be also set to zero here
    48324835!--                Save surface radiative fluxes onto respective surface elements
    48334836!--                Horizontal surfaces
     
    48364839                      surf_lsm_h%rad_sw_in(m)     = 0.0_wp
    48374840                      surf_lsm_h%rad_sw_out(m)    = 0.0_wp
    4838                    ENDDO             
     4841                   ENDDO
    48394842                   DO  m = surf_usm_h%start_index(j,i),                        &
    48404843                           surf_usm_h%end_index(j,i)
    48414844                      surf_usm_h%rad_sw_in(m)     = 0.0_wp
    48424845                      surf_usm_h%rad_sw_out(m)    = 0.0_wp
    4843                    ENDDO 
    4844 !
    4845 !--                Vertical surfaces. Fluxes are obtain at respective vertical 
     4846                   ENDDO
     4847!
     4848!--                Vertical surfaces. Fluxes are obtain at respective vertical
    48464849!--                level of the surface element
    48474850                   DO  l = 0, 3
     
    48514854                         surf_lsm_v(l)%rad_sw_in(m)  = 0.0_wp
    48524855                         surf_lsm_v(l)%rad_sw_out(m) = 0.0_wp
    4853                       ENDDO             
     4856                      ENDDO
    48544857                      DO  m = surf_usm_v(l)%start_index(j,i),                  &
    48554858                              surf_usm_v(l)%end_index(j,i)
     
    48574860                         surf_usm_v(l)%rad_sw_in(m)  = 0.0_wp
    48584861                         surf_usm_v(l)%rad_sw_out(m) = 0.0_wp
    4859                       ENDDO 
     4862                      ENDDO
    48604863                   ENDDO
    48614864                ENDIF
     
    48684871!--    Finally, calculate surface net radiation for surface elements.
    48694872       IF (  .NOT.  radiation_interactions  ) THEN
    4870 !--       First, for horizontal surfaces   
     4873!--       First, for horizontal surfaces
    48714874          DO  m = 1, surf_lsm_h%ns
    48724875             surf_lsm_h%rad_net(m) = surf_lsm_h%rad_sw_in(m)                   &
     
    48824885          ENDDO
    48834886!
    4884 !--       Vertical surfaces. 
     4887!--       Vertical surfaces.
    48854888!--       Todo: weight with azimuth and zenith angle according to their orientation!
    4886           DO  l = 0, 3     
     4889          DO  l = 0, 3
    48874890             DO  m = 1, surf_lsm_v(l)%ns
    48884891                surf_lsm_v(l)%rad_net(m) = surf_lsm_v(l)%rad_sw_in(m)          &
     
    49074910
    49084911       CALL exchange_horiz( rad_sw_in,  nbgp )
    4909        CALL exchange_horiz( rad_sw_out, nbgp ) 
     4912       CALL exchange_horiz( rad_sw_out, nbgp )
    49104913       CALL exchange_horiz( rad_sw_hr,    nbgp )
    49114914       CALL exchange_horiz( rad_sw_cs_hr, nbgp )
     
    49214924!> Calculate the cosine of the zenith angle (variable is called zenith)
    49224925!------------------------------------------------------------------------------!
    4923     SUBROUTINE calc_zenith( day_of_year, second_of_day ) 
     4926    SUBROUTINE calc_zenith( day_of_year, second_of_day )
    49244927
    49254928       USE palm_date_time_mod,                                                 &
     
    49734976! Description:
    49744977! ------------
    4975 !> Calculates surface albedo components based on Briegleb (1992) and 
     4978!> Calculates surface albedo components based on Briegleb (1992) and
    49764979!> Briegleb et al. (1986)
    49774980!------------------------------------------------------------------------------!
     
    49904993!
    49914994!--           Loop over surface elements
    4992               DO  ind_type = 0, SIZE( surf%albedo_type, 1 ) - 1
    4993            
     4995              DO  ind_type = 0, SIZE( surf%albedo_type, 2 ) - 1
     4996
    49944997!
    49954998!--              Ocean
     
    50815084           surf%rrtm_asdif = albedo_urb
    50825085           surf%rrtm_aldir = albedo_urb
    5083            surf%rrtm_aldif = albedo_urb 
     5086           surf%rrtm_aldif = albedo_urb
    50845087!
    50855088!--     Darkness
     
    51505153          DEALLOCATE ( rrtm_sw_tauaer )
    51515154          DEALLOCATE ( rrtm_sw_ssaaer )
    5152           DEALLOCATE ( rrtm_sw_asmaer ) 
    5153           DEALLOCATE ( rrtm_sw_ecaer )   
    5154  
     5155          DEALLOCATE ( rrtm_sw_asmaer )
     5156          DEALLOCATE ( rrtm_sw_ecaer )
     5157
    51555158          DEALLOCATE ( rrtm_swdflx  )
    51565159          DEALLOCATE ( rrtm_swdflxc )
     
    52365239
    52375240!
    5238 !--    Calculate pressure levels on zu and zw grid. Sounding data is added at 
    5239 !--    top of the LES domain. This routine does not consider horizontal or 
     5241!--    Calculate pressure levels on zu and zw grid. Sounding data is added at
     5242!--    top of the LES domain. This routine does not consider horizontal or
    52405243!--    vertical variability of pressure and temperature
    52415244       ALLOCATE ( rrtm_play(0:0,nzb+1:nzt_rad+1)   )
     
    52635266
    52645267!
    5265 !--    Calculate temperature/humidity levels at top of the LES domain. 
    5266 !--    Currently, the temperature is taken from sounding data (might lead to a 
    5267 !--    temperature jump at interface. To do: Humidity is currently not 
     5268!--    Calculate temperature/humidity levels at top of the LES domain.
     5269!--    Currently, the temperature is taken from sounding data (might lead to a
     5270!--    temperature jump at interface. To do: Humidity is currently not
    52685271!--    calculated above the LES domain.
    52695272       ALLOCATE ( rrtm_tlay(0:0,nzb+1:nzt_rad+1)   )
     
    52995302       ALLOCATE ( rrtm_sw_tauaer(0:0,nzb+1:nzt_rad+1,1:nbndsw+1) )
    53005303       ALLOCATE ( rrtm_sw_ssaaer(0:0,nzb+1:nzt_rad+1,1:nbndsw+1) )
    5301        ALLOCATE ( rrtm_sw_asmaer(0:0,nzb+1:nzt_rad+1,1:nbndsw+1) ) 
    5302        ALLOCATE ( rrtm_sw_ecaer(0:0,nzb+1:nzt_rad+1,1:naerec+1) )   
     5304       ALLOCATE ( rrtm_sw_asmaer(0:0,nzb+1:nzt_rad+1,1:nbndsw+1) )
     5305       ALLOCATE ( rrtm_sw_ecaer(0:0,nzb+1:nzt_rad+1,1:naerec+1) )
    53035306
    53045307!
     
    53325335       rrtm_swdflx  = 0.0_wp
    53335336       rrtm_swuflx  = 0.0_wp
    5334        rrtm_swhr    = 0.0_wp 
     5337       rrtm_swhr    = 0.0_wp
    53355338       rrtm_swuflxc = 0.0_wp
    53365339       rrtm_swdflxc = 0.0_wp
     
    53485351       rrtm_lwdflx  = 0.0_wp
    53495352       rrtm_lwuflx  = 0.0_wp
    5350        rrtm_lwhr    = 0.0_wp 
     5353       rrtm_lwhr    = 0.0_wp
    53515354       rrtm_lwuflxc = 0.0_wp
    53525355       rrtm_lwdflxc = 0.0_wp
     
    53675370!> Read trace gas data from file and convert into trace gas paths / volume
    53685371!> mixing ratios. If a user-defined input file is provided it needs to follow
    5369 !> the convections used in RRTMG (see respective netCDF files shipped with 
     5372!> the convections used in RRTMG (see respective netCDF files shipped with
    53705373!> RRTMG)
    53715374!------------------------------------------------------------------------------!
     
    54225425          DEALLOCATE ( rrtm_cfc22vmr )
    54235426          DEALLOCATE ( rrtm_ccl4vmr  )
    5424           DEALLOCATE ( rrtm_h2ovmr  )     
     5427          DEALLOCATE ( rrtm_h2ovmr  )
    54255428       ENDIF
    54265429
     
    54465449       nc_stat = NF90_INQ_DIMID( id, "Pressure", id_dim )
    54475450       CALL netcdf_handle_error_rad( 'read_trace_gas_data', 550 )
    5448        nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, len = np) 
     5451       nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, len = np)
    54495452       CALL netcdf_handle_error_rad( 'read_trace_gas_data', 550 )
    54505453
    54515454       nc_stat = NF90_INQ_DIMID( id, "Absorber", id_dim )
    54525455       CALL netcdf_handle_error_rad( 'read_trace_gas_data', 550 )
    5453        nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, len = nabs ) 
     5456       nc_stat = NF90_INQUIRE_DIMENSION( id, id_dim, len = nabs )
    54545457       CALL netcdf_handle_error_rad( 'read_trace_gas_data', 550 )
    5455    
    5456 
    5457 !
    5458 !--    Allocate pressure, and trace gas arrays     
     5458
     5459
     5460!
     5461!--    Allocate pressure, and trace gas arrays
    54595462       ALLOCATE( p_mls(1:np) )
    5460        ALLOCATE( trace_mls(1:num_trace_gases,1:np) ) 
    5461        ALLOCATE( trace_mls_tmp(1:nabs,1:np) ) 
     5463       ALLOCATE( trace_mls(1:num_trace_gases,1:np) )
     5464       ALLOCATE( trace_mls_tmp(1:nabs,1:np) )
    54625465
    54635466
     
    54825485!
    54835486!--       Replace missing values by zero
    5484           WHERE ( trace_mls(n,:) > 2.0_wp ) 
     5487          WHERE ( trace_mls(n,:) > 2.0_wp )
    54855488             trace_mls(n,:) = 0.0_wp
    54865489          END WHERE
     
    54975500       ALLOCATE ( rrtm_plev_tmp(1:nzt_rad+2) )
    54985501
    5499        rrtm_play_tmp(1:nzt_rad)   = rrtm_play(0,1:nzt_rad) 
     5502       rrtm_play_tmp(1:nzt_rad)   = rrtm_play(0,1:nzt_rad)
    55005503       rrtm_plev_tmp(1:nzt_rad+1) = rrtm_plev(0,1:nzt_rad+1)
    55015504       rrtm_play_tmp(nzt_rad+1)   = rrtm_plev(0,nzt_rad+1) * 0.5_wp
    55025505       rrtm_plev_tmp(nzt_rad+2)   = MIN( 1.0E-4_wp, 0.25_wp                    &
    55035506                                         * rrtm_plev(0,nzt_rad+1) )
    5504  
     5507
    55055508!
    55065509!--    Calculate trace gas path (zero at surface) with interpolation to the
     
    55095512
    55105513       trace_mls_path(nzb+1,:) = 0.0_wp
    5511        
     5514
    55125515       DO k = nzb+2, nzt_rad+2
    55135516          DO m = 1, num_trace_gases
     
    55165519!
    55175520!--          When the pressure level is higher than the trace gas pressure
    5518 !--          level, assume that 
    5519              IF ( rrtm_plev_tmp(k-1) > p_mls(1) )  THEN             
    5520                
     5521!--          level, assume that
     5522             IF ( rrtm_plev_tmp(k-1) > p_mls(1) )  THEN
     5523
    55215524                trace_mls_path(k,m) = trace_mls_path(k,m) + trace_mls(m,1)     &
    55225525                                      * ( rrtm_plev_tmp(k-1)                   &
     
    55265529
    55275530!
    5528 !--          Integrate for each sounding level from the contributing p_mls 
     5531!--          Integrate for each sounding level from the contributing p_mls
    55295532!--          levels
    55305533             DO n = 2, np
     
    55535556             ENDDO
    55545557
    5555              IF ( rrtm_plev_tmp(k) < p_mls(np) )  THEN 
     5558             IF ( rrtm_plev_tmp(k) < p_mls(np) )  THEN
    55565559                trace_mls_path(k,m) = trace_mls_path(k,m) + trace_mls(m,np)    &
    55575560                                      * ( MIN( rrtm_plev_tmp(k-1), p_mls(np) ) &
    55585561                                          - rrtm_plev_tmp(k)                   &
    5559                                         ) / g 
    5560              ENDIF 
     5562                                        ) / g
     5563             ENDIF
    55615564          ENDDO
    55625565       ENDDO
     
    56175620
    56185621                rrtm_h2ovmr(0,:) = trace_path_tmp(:)
    5619                
     5622
    56205623             CASE DEFAULT
    56215624
     
    57725775     INTEGER(iwp)                      ::  pc_box_dimshift    !< transform for best accuracy
    57735776     INTEGER(iwp), DIMENSION(0:3)      ::  reorder = (/ 1, 0, 3, 2 /)
    5774                                            
     5777
    57755778     REAL(wp), DIMENSION(3,3)          ::  mrot               !< grid rotation matrix (zyx)
    57765779     REAL(wp), DIMENSION(3,0:nsurf_type)::  vnorm             !< face direction normal vectors (zyx)
     
    57975800     REAL(wp)                          ::  area_surf          !< total area of surfaces in all processor
    57985801     REAL(wp)                          ::  area_hor           !< total horizontal area of domain in all processor
    5799 #if defined( __parallel )     
     5802#if defined( __parallel )
    58005803     REAL(wp), DIMENSION(1:7)          ::  combine_allreduce   !< dummy array used to combine several MPI_ALLREDUCE calls
    58015804     REAL(wp), DIMENSION(1:7)          ::  combine_allreduce_l !< dummy array used to combine several MPI_ALLREDUCE calls
     
    58505853     ENDIF
    58515854!
    5852 !--  Split downwelling shortwave radiation into a diffuse and a direct part. 
     5855!--  Split downwelling shortwave radiation into a diffuse and a direct part.
    58535856!--  Note, if radiation scheme is RRTMG or diffuse radiation is externally
    58545857!--  prescribed, this is not required. Please note, in case of external
    5855 !--  radiation, the clear-sky model is applied during spinup, so that 
    5856 !--  radiation need to be split also in this case. 
     5858!--  radiation, the clear-sky model is applied during spinup, so that
     5859!--  radiation need to be split also in this case.
    58575860     IF ( radiation_scheme == 'constant'   .OR.                                &
    58585861          radiation_scheme == 'clear-sky'  .OR.                                &
     
    61416144     IF ( .NOT.  surface_reflections )  THEN
    61426145!
    6143 !--     Set nrefsteps to 0 to disable reflections       
     6146!--     Set nrefsteps to 0 to disable reflections
    61446147        nrefsteps = 0
    61456148        surfoutsl = albedo_surf * surfins
     
    65496552     combine_allreduce_l(6) = emiss_sum_surfl
    65506553     combine_allreduce_l(7) = area_surfl
    6551      
     6554
    65526555     CALL MPI_ALLREDUCE( combine_allreduce_l,                                  &
    65536556                         combine_allreduce,                                    &
     
    65576560                         comm2d,                                               &
    65586561                         ierr )
    6559      
     6562
    65606563     pinsw          = combine_allreduce(1)
    65616564     pinlw          = combine_allreduce(2)
     
    65806583     IF ( area_surf /= 0.0_wp )  emissivity_urb = emiss_sum_surf / area_surf
    65816584!
    6582 !--  Temporally comment out calculation of effective radiative temperature. 
     6585!--  Temporally comment out calculation of effective radiative temperature.
    65836586!--  See below for more explanation.
    65846587!--  (3) temperature
    65856588!--   first we calculate an effective horizontal area to account for
    65866589!--   the effect of vertical surfaces (which contributes to LW emission)
    6587 !--   We simply use the ratio of the total LW to the incoming LW flux 
     6590!--   We simply use the ratio of the total LW to the incoming LW flux
    65886591      area_hor = pinlw / rad_lw_in_diff(nyn,nxl)
    65896592      t_rad_urb = ( ( pemitlw - pabslw + emissivity_urb * pinlw ) / &
    6590            (emissivity_urb * sigma_sb * area_hor) )**0.25_wp     
     6593           (emissivity_urb * sigma_sb * area_hor) )**0.25_wp
    65916594
    65926595     IF ( debug_output_timestep )  CALL debug_message( 'radiation_interaction', 'end' )
     
    66756678!> It follows Boland, Ridley & Brown (2008)
    66766679!------------------------------------------------------------------------------!
    6677     SUBROUTINE calc_diffusion_radiation 
     6680    SUBROUTINE calc_diffusion_radiation
    66786681
    66796682       USE palm_date_time_mod,                                                 &
     
    67066709                          0.000719_wp * cos(2.0_wp * year_angle) +             &
    67076710                          0.000077_wp * sin(2.0_wp * year_angle))
    6708        
    6709 !--   
     6711
     6712!--
    67106713!--     Under a very low angle, we keep extraterestrial radiation at
    67116714!--     the last small value, therefore the clearness index will be pushed
    6712 !--     towards 0 while keeping full continuity.   
     6715!--     towards 0 while keeping full continuity.
    67136716        IF ( cos_zenith <= lowest_solarUp )  THEN
    67146717            corrected_solarUp = lowest_solarUp
     
    67166719            corrected_solarUp = cos_zenith
    67176720        ENDIF
    6718        
     6721
    67196722        horizontalETR = etr * corrected_solarUp
    6720        
     6723
    67216724        DO i = nxl, nxr
    67226725            DO j = nys, nyn
     
    67286731            ENDDO
    67296732        ENDDO
    6730        
     6733
    67316734    END SUBROUTINE calc_diffusion_radiation
    67326735
     
    68266829
    68276830 END SUBROUTINE radiation_interaction
    6828    
     6831
    68296832!------------------------------------------------------------------------------!
    68306833! Description:
     
    69346937#endif
    69356938!
    6936 !--    Stretching (non-uniform grid spacing) is not considered in the radiation 
    6937 !--    model. Therefore, vertical stretching has to be applied above the area 
     6939!--    Stretching (non-uniform grid spacing) is not considered in the radiation
     6940!--    model. Therefore, vertical stretching has to be applied above the area
    69386941!--    where the parts of the radiation model which assume constant grid spacing
    6939 !--    are active. ABS (...) is required because the default value of 
    6940 !--    dz_stretch_level_start is -9999999.9_wp (negative). 
     6942!--    are active. ABS (...) is required because the default value of
     6943!--    dz_stretch_level_start is -9999999.9_wp (negative).
    69416944       IF ( ABS( dz_stretch_level_start(1) ) <= zw(nz_urban_t) ) THEN
    69426945          WRITE( message_string, * ) 'The lowest level where vertical ',       &
     
    69446947                                     'greater than ', zw(nz_urban_t)
    69456948          CALL message( 'radiation_interaction_init', 'PA0496', 1, 2, 0, 6, 0 )
    6946        ENDIF 
     6949       ENDIF
    69476950!
    69486951!--    global number of urban and plant layers
     
    73657368!------------------------------------------------------------------------------!
    73667369    SUBROUTINE radiation_calc_svf
    7367    
     7370
    73687371        IMPLICIT NONE
    7369        
     7372
    73707373        INTEGER(iwp)                                  :: i, j, k, d, ip, jp
    73717374        INTEGER(iwp)                                  :: isvf, ksvf, icsf, kcsf, npcsfl, isvf_surflt, imrt, imrtf, ipcgb
     
    74117414        INTEGER(kind=MPI_ADDRESS_KIND)                :: size_lad_rma
    74127415#endif
    7413 !   
     7416!
    74147417        INTEGER(iwp), DIMENSION(0:svfnorm_report_num) :: svfnorm_counts
    74157418
     
    74407443        ray_skip_maxdist = 0
    74417444        ray_skip_minval = 0
    7442        
     7445
    74437446!--     initialize temporary terrain and plant canopy height arrays (global 2D array!)
    74447447        ALLOCATE( nzterr(0:(nx+1)*(ny+1)-1) )
     
    74827485            IF ( raytrace_mpi_rma )  THEN
    74837486                ALLOCATE( lad_s_ray(maxboxesg) )
    7484                
     7487
    74857488                ! set conditions for RMA communication
    74867489                CALL MPI_Info_create(minfo, ierr)
     
    75567559                    FLUSH(9)
    75577560                ENDIF
    7558                
     7561
    75597562            ELSE
    75607563                ALLOCATE( sub_lad_g(0:(nx+1)*(ny+1)*nz_plant-1) )
     
    75897592        !--
    75907593        dsitrans(:,:) = 0._wp
    7591        
     7594
    75927595        DO isurflt = 1, nsurfl
    75937596!--         determine face centers
     
    76947697                                 SUM(ztransp(itarg0:itarg0+lowest_free_ray-1) &
    76957698                                             * vffrac(itarg0:itarg0+lowest_free_ray-1))
    7696  
     7699
    76977700!--            Save direct solar transparency
    76987701               j = MODULO(NINT(azmid/                                          &
     
    77607763                        CALL debug_message( debug_string, 'info' )
    77617764                     ENDIF
    7762                      
     7765
    77637766                     nsvfla = k
    77647767                  ENDIF
     
    77857788                     CYCLE
    77867789                  ENDIF
    7787                  
     7790
    77887791                  sd = surf(id, isurfs)
    77897792                  sa = (/ REAL(surf(iz, isurfs), wp) - 0.5_wp * kdir(sd),  &
     
    78017804                     CYCLE
    78027805                  ENDIF
    7803                  
     7806
    78047807                  difvf = dot_product((/ kdir(sd), jdir(sd), idir(sd) /), uv) & ! cosine of source normal and direction
    78057808                      * dot_product((/ kdir(td), jdir(td), idir(td) /), -uv) &  ! cosine of target normal and reverse direction
     
    78457848                        CALL debug_message( debug_string, 'info' )
    78467849                     ENDIF
    7847                      
     7850
    78487851                     nsvfla = k
    78497852                  ENDIF
     
    80828085!--     deallocate temporary global arrays
    80838086        DEALLOCATE(nzterr)
    8084        
     8087
    80858088        IF ( plant_canopy )  THEN
    80868089!--         finalize mpi_rma communication and deallocate temporary arrays
     
    82028205!--         sort and merge csf for the last time, keeping the array size to minimum
    82038206            CALL merge_and_grow_csf(-1)
    8204            
     8207
    82058208!--         aggregate csb among processors
    82068209!--         allocate necessary arrays
     
    82148217            ALLOCATE( ipcsflt(0:numprocs-1) )
    82158218            ALLOCATE( dpcsflt(0:numprocs-1) )
    8216            
    8217 !--         fill out arrays of csf values and 
     8219
     8220!--         fill out arrays of csf values and
    82188221!--         arrays of number of elements and displacements
    82198222!--         for particular precessors
     
    82578260                dcsflt(jp) = d
    82588261            ENDDO
    8259            
     8262
    82608263!--         deallocate temporary acsf array
    82618264!--         DEALLOCATE(acsf) - ifort has a problem with deallocation of allocatable target
     
    82678270                DEALLOCATE(acsf2)
    82688271            ENDIF
    8269                    
     8272
    82708273#if defined( __parallel )
    82718274!--         scatter and gather the number of elements to and from all processor
     
    83098312                FLUSH(9)
    83108313            ENDIF
    8311            
     8314
    83128315#else
    83138316            npcsfl = ncsfl
     
    83708373                ENDDO
    83718374            ENDIF
    8372            
     8375
    83738376!--         deallocation of temporary arrays
    83748377            IF ( npcbl > 0 )  DEALLOCATE( gridpcbl )
     
    83798382               CALL debug_message( debug_string, 'info' )
    83808383            ENDIF
    8381            
     8384
    83828385        ENDIF
    83838386
     
    83888391
    83898392        RETURN  !todo: remove
    8390        
     8393
    83918394!        WRITE( message_string, * )  &
    83928395!            'I/O error when processing shape view factors / ',  &
    83938396!            'plant canopy sink factors / direct irradiance factors.'
    83948397!        CALL message( 'init_urban_surface', 'PA0502', 2, 2, 0, 6, 0 )
    8395        
     8398
    83968399    END SUBROUTINE radiation_calc_svf
    83978400
    8398    
     8401
    83998402!------------------------------------------------------------------------------!
    84008403! Description:
     
    84408443        INTEGER(iwp), DIMENSION(3)             :: dimnext      !< next dimension increments along path
    84418444        INTEGER(iwp), DIMENSION(3)             :: dimdelta     !< dimension direction = +- 1
    8442         INTEGER(iwp)                           :: px, py       !< number of processors in x and y dir before 
     8445        INTEGER(iwp)                           :: px, py       !< number of processors in x and y dir before
    84438446                                                               !< the processor in the question
    84448447        INTEGER(iwp)                           :: ip           !< number of processor where gridbox reside
    84458448        INTEGER(iwp)                           :: ig           !< 1D index of gridbox in global 2D array
    8446        
     8449
    84478450        REAL(wp)                               :: eps = 1E-10_wp !< epsilon for value comparison
    84488451        REAL(wp)                               :: lad_s_target   !< recieved lad_s of particular grid box
     
    84618464            CALL merge_and_grow_csf(k)
    84628465        ENDIF
    8463        
     8466
    84648467        transparency = 1._wp
    84658468        ncsb = 0
     
    85348537            dimnextdist(seldim) = (dimnext(seldim) - .5_wp - src(seldim)) / uvect(seldim)
    85358538        ENDDO
    8536        
     8539
    85378540        IF ( plant_canopy )  THEN
    85388541#if defined( __parallel )
     
    85498552                    ENDIF
    85508553                ENDDO
    8551                
     8554
    85528555!--             wait for all pending local requests complete
    85538556                CALL MPI_Win_flush_local_all(win_lad, ierr)
     
    85578560                ENDIF
    85588561                CALL cpu_log( log_point_s(77), 'rad_rma_lad', 'stop' )
    8559                
     8562
    85608563            ENDIF
    85618564#endif
     
    85868589
    85878590                transparency = transparency * (1._wp - cursink)
    8588                
     8591
    85898592            ENDDO
    85908593        ENDIF
    8591        
     8594
    85928595        visible = .TRUE.
    85938596
    85948597    END SUBROUTINE raytrace
    8595    
    8596  
     8598
     8599
    85978600!------------------------------------------------------------------------------!
    85988601! Description:
     
    86408643      INTEGER(iwp), DIMENSION(2)             ::  dimnext      !< next dimension increments along path
    86418644      INTEGER(iwp), DIMENSION(2)             ::  dimdelta     !< dimension direction = +- 1
    8642       INTEGER(iwp)                           ::  px, py       !< number of processors in x and y dir before 
     8645      INTEGER(iwp)                           ::  px, py       !< number of processors in x and y dir before
    86438646                                                              !< the processor in the question
    86448647      INTEGER(iwp)                           ::  ip           !< number of processor where gridbox reside
    86458648      INTEGER(iwp)                           ::  ig           !< 1D index of gridbox in global 2D array
    86468649      INTEGER(iwp)                           ::  maxboxes     !< max no of CSF created
    8647       INTEGER(iwp)                           ::  nly          !< maximum  plant canopy height 
     8650      INTEGER(iwp)                           ::  nly          !< maximum  plant canopy height
    86488651      INTEGER(iwp)                           ::  ntrack
    8649      
     8652
    86508653      INTEGER(iwp)                           ::  zb0
    86518654      INTEGER(iwp)                           ::  zb1
     
    86618664      INTEGER(MPI_ADDRESS_KIND)              ::  wdisp        !< RMA window displacement
    86628665#endif
    8663      
     8666
    86648667      REAL(wp)                               ::  eps = 1E-10_wp !< epsilon for value comparison
    86658668      REAL(wp)                               ::  zbottom, ztop !< urban surface boundary in real numbers
     
    86698672      REAL(wp)                               ::  dxxyy         !< square of real horizontal distance
    86708673      REAL(wp)                               ::  curtrans      !< transparency of current PC box crossing
    8671      
    8672 
    8673      
     8674
     8675
     8676
    86748677      yxorigin(:) = origin(2:3)
    86758678      transparency(:) = 1._wp !-- Pre-set the all rays to transparent before reducing
     
    88168819      IF ( plant_canopy )  THEN
    88178820!--      Request LAD WHERE applicable
    8818 !--     
     8821!--
    88198822#if defined( __parallel )
    88208823         IF ( raytrace_mpi_rma )  THEN
     
    88988901            ENDIF
    88998902         ENDDO
    8900          
     8903
    89018904         DEALLOCATE( target_surfl )
    8902          
     8905
    89038906      ELSE
    89048907         itarget(:) = -1
     
    89078910      IF ( plant_canopy )  THEN
    89088911!--      Skip the PCB around origin if requested (for MRT, the PCB might not be there)
    8909 !--     
     8912!--
    89108913         IF ( skip_1st_pcb  .AND.  NINT(origin(1)) <= plantt_max )  THEN
    89118914            rt2_track_lad(NINT(origin(1), iwp), 1) = 0._wp
     
    89138916
    89148917!--      Assert that we have space allocated for CSFs
    8915 !--     
     8918!--
    89168919         maxboxes = (ntrack + MAX(CEILING(origin(1)-.5_wp) - nz_urban_b,          &
    89178920                                  nz_urban_t - CEILING(origin(1)-.5_wp))) * nrays
     
    89268929
    89278930!--      Calculate transparencies and store new CSFs
    8928 !--     
     8931!--
    89298932         zbottom = REAL(nz_urban_b, wp) - .5_wp
    89308933         ztop = REAL(plantt_max, wp) + .5_wp
    89318934
    89328935!--      Reverse direction of radiation (face->sky), only when calc_svf
    8933 !--     
     8936!--
    89348937         IF ( calc_svf )  THEN
    89358938            DO  i = 1, ntrack ! for each column
     
    90029005
    90039006!--      Forward direction of radiation (sky->face), always
    9004 !--     
     9007!--
    90059008         DO  i = ntrack, 1, -1 ! for each column backwards
    90069009            dxxyy = ((dy*yxdir(1))**2 + (dx*yxdir(2))**2) * (rt2_track_dist(i)-rt2_track_dist(i-1))**2
     
    91079110
    91089111   END SUBROUTINE raytrace_2d
    9109  
     9112
    91109113
    91119114!------------------------------------------------------------------------------!
     
    91329135      ndsidir = 0
    91339136      sun_direction = .TRUE.
    9134      
     9137
    91359138!
    91369139!--   Process spinup time if configured
     
    91769179         CALL calc_zenith( day_of_year, second_of_day )
    91779180         IF ( cos_zenith > 0 )  THEN
    9178 !--         
     9181!--
    91799182!--         Identify solar direction vector (discretized number) 1)
    91809183            i = MODULO(NINT(ATAN2(sun_dir_lon, sun_dir_lat)               &
     
    92109213        IMPLICIT NONE
    92119214        INTEGER(iwp),   INTENT(in)  :: x, y, z, d, x2, y2, z2, d2
    9212      
     9215
    92139216        surface_facing = .FALSE.
    92149217
     
    92559258
    92569259        surface_facing = .TRUE.
    9257        
     9260
    92589261    END FUNCTION surface_facing
    92599262
     
    92709273
    92719274       IMPLICIT NONE
    9272        
     9275
    92739276       CHARACTER(rad_version_len)   :: rad_version_field
    9274        
     9277
    92759278       INTEGER(iwp)                 :: i
    92769279       INTEGER(iwp)                 :: ndsidir_from_file = 0
     
    92969299                                              'that has written the svf data ',&
    92979300                                              'and the one that will read it ',&
    9298                                               'is not allowed' 
     9301                                              'is not allowed'
    92999302                   CALL message( 'check_open', 'PA0491', 1, 2, 0, 6, 0 )
    93009303                ENDIF
    93019304
    93029305             ENDIF
    9303              
    9304 !
    9305 !--          Open binary file 
     9306
     9307!
     9308!--          Open binary file
    93069309             CALL check_open( 88 )
    93079310
     
    93159318                 CALL message( 'radiation_read_svf', 'PA0482', 1, 2, 0, 6, 0 )
    93169319             ENDIF
    9317              
     9320
    93189321!
    93199322!--          read nsvfl, ncsfl, nsurfl, nmrtf
    93209323             READ ( 88 ) nsvfl, ncsfl, nsurfl_from_file, npcbl_from_file,      &
    93219324                         ndsidir_from_file, nmrtbl_from_file, nmrtf
    9322              
     9325
    93239326             IF ( nsvfl < 0  .OR.  ncsfl < 0 )  THEN
    93249327                 WRITE( message_string, * ) 'Wrong number of SVF or CSF'
     
    93309333                 IF ( debug_output )  CALL debug_message( debug_string, 'info' )
    93319334             ENDIF
    9332              
     9335
    93339336             IF ( nsurfl_from_file /= nsurfl )  THEN
    93349337                 WRITE( message_string, * ) 'nsurfl from SVF file does not ',  &
     
    93379340                 CALL message( 'radiation_read_svf', 'PA0490', 1, 2, 0, 6, 0 )
    93389341             ENDIF
    9339              
     9342
    93409343             IF ( npcbl_from_file /= npcbl )  THEN
    93419344                 WRITE( message_string, * ) 'npcbl from SVF file does not ',   &
     
    93449347                 CALL message( 'radiation_read_svf', 'PA0493', 1, 2, 0, 6, 0 )
    93459348             ENDIF
    9346              
     9349
    93479350             IF ( ndsidir_from_file /= ndsidir )  THEN
    93489351                 WRITE( message_string, * ) 'ndsidir from SVF file does not ', &
     
    93609363                 IF ( debug_output )  CALL debug_message( debug_string, 'info' )
    93619364             ENDIF
    9362              
    9363 !
    9364 !--          Arrays skyvf, skyvft, dsitrans and dsitransc are allready 
    9365 !--          allocated in radiation_interaction_init and 
     9365
     9366!
     9367!--          Arrays skyvf, skyvft, dsitrans and dsitransc are allready
     9368!--          allocated in radiation_interaction_init and
    93669369!--          radiation_presimulate_solar_pos
    93679370             IF ( nsurfl > 0 )  THEN
    93689371                READ(88) skyvf
    93699372                READ(88) skyvft
    9370                 READ(88) dsitrans 
     9373                READ(88) dsitrans
    93719374             ENDIF
    9372              
     9375
    93739376             IF ( plant_canopy  .AND.  npcbl > 0 ) THEN
    93749377                READ ( 88 )  dsitransc
    93759378             ENDIF
    9376              
    9377 !
    9378 !--          The allocation of svf, svfsurf, csf, csfsurf, mrtf, mrtft, and 
    9379 !--          mrtfsurf happens in routine radiation_calc_svf which is not 
    9380 !--          called if the program enters radiation_read_svf. Therefore 
     9379
     9380!
     9381!--          The allocation of svf, svfsurf, csf, csfsurf, mrtf, mrtft, and
     9382!--          mrtfsurf happens in routine radiation_calc_svf which is not
     9383!--          called if the program enters radiation_read_svf. Therefore
    93819384!--          these arrays has to allocate in the following
    93829385             IF ( nsvfl > 0 )  THEN
     
    94089411                READ(88) mrtfsurf
    94099412             ENDIF
    9410              
    9411 !
    9412 !--          Close binary file                 
     9413
     9414!
     9415!--          Close binary file
    94139416             CALL close_file( 88 )
    9414                
     9417
    94159418          ENDIF
    94169419#if defined( __parallel )
     
    94369439
    94379440       IMPLICIT NONE
    9438        
     9441
    94399442       INTEGER(iwp)        :: i
    94409443
     
    94459448          IF ( i == io_group )  THEN
    94469449!
    9447 !--          Open binary file 
     9450!--          Open binary file
    94489451             CALL check_open( 89 )
    94499452
     
    94739476             IF ( nmrtf > 0 )  THEN
    94749477                 WRITE ( 89 )  mrtf
    9475                  WRITE ( 89 )  mrtft               
     9478                 WRITE ( 89 )  mrtft
    94769479                 WRITE ( 89 )  mrtfsurf
    94779480             ENDIF
    94789481!
    9479 !--          Close binary file                 
     9482!--          Close binary file
    94809483             CALL close_file( 89 )
    94819484
     
    96279630    END SUBROUTINE quicksort_csf
    96289631
    9629    
     9632
    96309633!------------------------------------------------------------------------------!
    96319634!
     
    97129715    END SUBROUTINE merge_and_grow_csf
    97139716
    9714    
     9717
    97159718!-- quicksort.f -*-f90-*-
    97169719!-- Author: t-nissie, adaptation J.Resler
     
    97469749        IF ( j+1 < last )  CALL quicksort_csf2(kpcsflt, pcsflt, j+1, last)
    97479750    END SUBROUTINE quicksort_csf2
    9748    
     9751
    97499752
    97509753    PURE FUNCTION csf_lt2(item1, item2) result(res)
     
    97819784!------------------------------------------------------------------------------!
    97829785SUBROUTINE radiation_3d_data_averaging( mode, variable )
    9783  
     9786
    97849787
    97859788    USE control_parameters
     
    97919794    IMPLICIT NONE
    97929795
    9793     CHARACTER (LEN=*) ::  mode    !< 
    9794     CHARACTER (LEN=*) :: variable !< 
     9796    CHARACTER (LEN=*) ::  mode    !<
     9797    CHARACTER (LEN=*) :: variable !<
    97959798
    97969799    LOGICAL      ::  match_lsm !< flag indicating natural-type surface
    97979800    LOGICAL      ::  match_usm !< flag indicating urban-type surface
    9798    
    9799     INTEGER(iwp) ::  i !< 
    9800     INTEGER(iwp) ::  j !< 
    9801     INTEGER(iwp) ::  k !< 
     9801
     9802    INTEGER(iwp) ::  i !<
     9803    INTEGER(iwp) ::  j !<
     9804    INTEGER(iwp) ::  k !<
    98029805    INTEGER(iwp) ::  l, m !< index of current surface element
    98039806
     
    98359838                ENDIF
    98369839                rad_net_av = 0.0_wp
    9837              
     9840
    98389841             CASE ( 'rad_lw_in*' )
    98399842                IF ( .NOT. ALLOCATED( rad_lw_in_xy_av ) )  THEN
     
    98419844                ENDIF
    98429845                rad_lw_in_xy_av = 0.0_wp
    9843                
     9846
    98449847             CASE ( 'rad_lw_out*' )
    98459848                IF ( .NOT. ALLOCATED( rad_lw_out_xy_av ) )  THEN
     
    98479850                ENDIF
    98489851                rad_lw_out_xy_av = 0.0_wp
    9849                
     9852
    98509853             CASE ( 'rad_sw_in*' )
    98519854                IF ( .NOT. ALLOCATED( rad_sw_in_xy_av ) )  THEN
     
    98539856                ENDIF
    98549857                rad_sw_in_xy_av = 0.0_wp
    9855                
     9858
    98569859             CASE ( 'rad_sw_out*' )
    98579860                IF ( .NOT. ALLOCATED( rad_sw_out_xy_av ) )  THEN
    98589861                   ALLOCATE( rad_sw_out_xy_av(nysg:nyng,nxlg:nxrg) )
    98599862                ENDIF
    9860                 rad_sw_out_xy_av = 0.0_wp               
     9863                rad_sw_out_xy_av = 0.0_wp
    98619864
    98629865             CASE ( 'rad_lw_in' )
     
    1009710100                ENDDO
    1009810101             ENDIF
    10099              
     10102
    1010010103          CASE ( 'rad_lw_out*' )
    1010110104             IF ( ALLOCATED( rad_lw_out_xy_av ) ) THEN
     
    1011910122                ENDDO
    1012010123             ENDIF
    10121              
     10124
    1012210125          CASE ( 'rad_sw_in*' )
    1012310126             IF ( ALLOCATED( rad_sw_in_xy_av ) ) THEN
     
    1014110144                ENDDO
    1014210145             ENDIF
    10143              
     10146
    1014410147          CASE ( 'rad_sw_out*' )
    1014510148             IF ( ALLOCATED( rad_sw_out_xy_av ) ) THEN
     
    1016310166                ENDDO
    1016410167             ENDIF
    10165              
     10168
    1016610169          CASE ( 'rad_lw_in' )
    1016710170             IF ( ALLOCATED( rad_lw_in_av ) ) THEN
     
    1041810421                ENDDO
    1041910422             ENDIF
    10420              
     10423
    1042110424          CASE ( 'rad_lw_in*' )
    1042210425             IF ( ALLOCATED( rad_lw_in_xy_av ) ) THEN
     
    1042810431                ENDDO
    1042910432             ENDIF
    10430              
     10433
    1043110434          CASE ( 'rad_lw_out*' )
    1043210435             IF ( ALLOCATED( rad_lw_out_xy_av ) ) THEN
     
    1043810441                ENDDO
    1043910442             ENDIF
    10440              
     10443
    1044110444          CASE ( 'rad_sw_in*' )
    1044210445             IF ( ALLOCATED( rad_sw_in_xy_av ) ) THEN
     
    1044810451                ENDDO
    1044910452             ENDIF
    10450              
     10453
    1045110454          CASE ( 'rad_sw_out*' )
    1045210455             IF ( ALLOCATED( rad_sw_out_xy_av ) ) THEN
     
    1070210705!------------------------------------------------------------------------------!
    1070310706SUBROUTINE radiation_define_netcdf_grid( variable, found, grid_x, grid_y, grid_z )
    10704    
     10707
    1070510708    IMPLICIT NONE
    1070610709
    1070710710    CHARACTER (LEN=*), INTENT(IN)  ::  variable    !<
    10708     LOGICAL, INTENT(OUT)           ::  found       !< 
    10709     CHARACTER (LEN=*), INTENT(OUT) ::  grid_x      !< 
    10710     CHARACTER (LEN=*), INTENT(OUT) ::  grid_y      !< 
    10711     CHARACTER (LEN=*), INTENT(OUT) ::  grid_z      !< 
     10711    LOGICAL, INTENT(OUT)           ::  found       !<
     10712    CHARACTER (LEN=*), INTENT(OUT) ::  grid_x      !<
     10713    CHARACTER (LEN=*), INTENT(OUT) ::  grid_y      !<
     10714    CHARACTER (LEN=*), INTENT(OUT) ::  grid_z      !<
    1071210715
    1071310716    CHARACTER (len=varnamelength)  :: var
     
    1077910782 SUBROUTINE radiation_data_output_2d( av, variable, found, grid, mode,         &
    1078010783                                      local_pf, two_d, nzb_do, nzt_do )
    10781  
     10784
    1078210785    USE indices
    1078310786
     
    1078710790    IMPLICIT NONE
    1078810791
    10789     CHARACTER (LEN=*) ::  grid     !< 
    10790     CHARACTER (LEN=*) ::  mode     !< 
    10791     CHARACTER (LEN=*) ::  variable !< 
    10792 
    10793     INTEGER(iwp) ::  av !< 
    10794     INTEGER(iwp) ::  i  !< 
    10795     INTEGER(iwp) ::  j  !< 
    10796     INTEGER(iwp) ::  k  !< 
     10792    CHARACTER (LEN=*) ::  grid     !<
     10793    CHARACTER (LEN=*) ::  mode     !<
     10794    CHARACTER (LEN=*) ::  variable !<
     10795
     10796    INTEGER(iwp) ::  av !<
     10797    INTEGER(iwp) ::  i  !<
     10798    INTEGER(iwp) ::  j  !<
     10799    INTEGER(iwp) ::  k  !<
    1079710800    INTEGER(iwp) ::  m  !< index of surface element at grid point (j,i)
    10798     INTEGER(iwp) ::  nzb_do   !< 
    10799     INTEGER(iwp) ::  nzt_do   !< 
    10800 
    10801     LOGICAL      ::  found !< 
     10801    INTEGER(iwp) ::  nzb_do   !<
     10802    INTEGER(iwp) ::  nzt_do   !<
     10803
     10804    LOGICAL      ::  found !<
    1080210805    LOGICAL      ::  two_d !< flag parameter that indicates 2D variables (horizontal cross sections)
    1080310806
    1080410807    REAL(wp) ::  fill_value = -999.0_wp    !< value for the _FillValue attribute
    1080510808
    10806     REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !< 
     10809    REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !<
    1080710810
    1080810811    found = .TRUE.
     
    1081810821!--                Natural-type surfaces
    1081910822                   DO  m = surf_lsm_h%start_index(j,i),                        &
    10820                            surf_lsm_h%end_index(j,i) 
     10823                           surf_lsm_h%end_index(j,i)
    1082110824                      local_pf(i,j,nzb+1) = surf_lsm_h%rad_net(m)
    1082210825                   ENDDO
     
    1082410827!--                Urban-type surfaces
    1082510828                   DO  m = surf_usm_h%start_index(j,i),                        &
    10826                            surf_usm_h%end_index(j,i) 
     10829                           surf_usm_h%end_index(j,i)
    1082710830                      local_pf(i,j,nzb+1) = surf_usm_h%rad_net(m)
    1082810831                   ENDDO
     
    1083510838             ENDIF
    1083610839             DO  i = nxl, nxr
    10837                 DO  j = nys, nyn 
     10840                DO  j = nys, nyn
    1083810841                   local_pf(i,j,nzb+1) = rad_net_av(j,i)
    1083910842                ENDDO
     
    1084210845          two_d = .TRUE.
    1084310846          grid = 'zu1'
    10844          
     10847
    1084510848       CASE ( 'rad_lw_in*_xy' )        ! 2d-array
    1084610849          IF ( av == 0 ) THEN
     
    1085110854!--                Natural-type surfaces
    1085210855                   DO  m = surf_lsm_h%start_index(j,i),                        &
    10853                            surf_lsm_h%end_index(j,i) 
     10856                           surf_lsm_h%end_index(j,i)
    1085410857                      local_pf(i,j,nzb+1) = surf_lsm_h%rad_lw_in(m)
    1085510858                   ENDDO
     
    1085710860!--                Urban-type surfaces
    1085810861                   DO  m = surf_usm_h%start_index(j,i),                        &
    10859                            surf_usm_h%end_index(j,i) 
     10862                           surf_usm_h%end_index(j,i)
    1086010863                      local_pf(i,j,nzb+1) = surf_usm_h%rad_lw_in(m)
    1086110864                   ENDDO
     
    1086810871             ENDIF
    1086910872             DO  i = nxl, nxr
    10870                 DO  j = nys, nyn 
     10873                DO  j = nys, nyn
    1087110874                   local_pf(i,j,nzb+1) = rad_lw_in_xy_av(j,i)
    1087210875                ENDDO
     
    1087510878          two_d = .TRUE.
    1087610879          grid = 'zu1'
    10877          
     10880
    1087810881       CASE ( 'rad_lw_out*_xy' )        ! 2d-array
    1087910882          IF ( av == 0 ) THEN
     
    1088410887!--                Natural-type surfaces
    1088510888                   DO  m = surf_lsm_h%start_index(j,i),                        &
    10886                            surf_lsm_h%end_index(j,i) 
     10889                           surf_lsm_h%end_index(j,i)
    1088710890                      local_pf(i,j,nzb+1) = surf_lsm_h%rad_lw_out(m)
    1088810891                   ENDDO
     
    1089010893!--                Urban-type surfaces
    1089110894                   DO  m = surf_usm_h%start_index(j,i),                        &
    10892                            surf_usm_h%end_index(j,i) 
     10895                           surf_usm_h%end_index(j,i)
    1089310896                      local_pf(i,j,nzb+1) = surf_usm_h%rad_lw_out(m)
    1089410897                   ENDDO
     
    1090110904             ENDIF
    1090210905             DO  i = nxl, nxr
    10903                 DO  j = nys, nyn 
     10906                DO  j = nys, nyn
    1090410907                   local_pf(i,j,nzb+1) = rad_lw_out_xy_av(j,i)
    1090510908                ENDDO
     
    1090810911          two_d = .TRUE.
    1090910912          grid = 'zu1'
    10910          
     10913
    1091110914       CASE ( 'rad_sw_in*_xy' )        ! 2d-array
    1091210915          IF ( av == 0 ) THEN
     
    1091710920!--                Natural-type surfaces
    1091810921                   DO  m = surf_lsm_h%start_index(j,i),                        &
    10919                            surf_lsm_h%end_index(j,i) 
     10922                           surf_lsm_h%end_index(j,i)
    1092010923                      local_pf(i,j,nzb+1) = surf_lsm_h%rad_sw_in(m)
    1092110924                   ENDDO
     
    1092310926!--                Urban-type surfaces
    1092410927                   DO  m = surf_usm_h%start_index(j,i),                        &
    10925                            surf_usm_h%end_index(j,i) 
     10928                           surf_usm_h%end_index(j,i)
    1092610929                      local_pf(i,j,nzb+1) = surf_usm_h%rad_sw_in(m)
    1092710930                   ENDDO
     
    1093410937             ENDIF
    1093510938             DO  i = nxl, nxr
    10936                 DO  j = nys, nyn 
     10939                DO  j = nys, nyn
    1093710940                   local_pf(i,j,nzb+1) = rad_sw_in_xy_av(j,i)
    1093810941                ENDDO
     
    1094110944          two_d = .TRUE.
    1094210945          grid = 'zu1'
    10943          
     10946
    1094410947       CASE ( 'rad_sw_out*_xy' )        ! 2d-array
    1094510948          IF ( av == 0 ) THEN
     
    1095010953!--                Natural-type surfaces
    1095110954                   DO  m = surf_lsm_h%start_index(j,i),                        &
    10952                            surf_lsm_h%end_index(j,i) 
     10955                           surf_lsm_h%end_index(j,i)
    1095310956                      local_pf(i,j,nzb+1) = surf_lsm_h%rad_sw_out(m)
    1095410957                   ENDDO
     
    1095610959!--                Urban-type surfaces
    1095710960                   DO  m = surf_usm_h%start_index(j,i),                        &
    10958                            surf_usm_h%end_index(j,i) 
     10961                           surf_usm_h%end_index(j,i)
    1095910962                      local_pf(i,j,nzb+1) = surf_usm_h%rad_sw_out(m)
    1096010963                   ENDDO
     
    1096710970             ENDIF
    1096810971             DO  i = nxl, nxr
    10969                 DO  j = nys, nyn 
     10972                DO  j = nys, nyn
    1097010973                   local_pf(i,j,nzb+1) = rad_sw_out_xy_av(j,i)
    1097110974                ENDDO
     
    1097310976          ENDIF
    1097410977          two_d = .TRUE.
    10975           grid = 'zu1'         
    10976          
     10978          grid = 'zu1'
     10979
    1097710980       CASE ( 'rad_lw_in_xy', 'rad_lw_in_xz', 'rad_lw_in_yz' )
    1097810981          IF ( av == 0 ) THEN
     
    1099010993            ENDIF
    1099110994             DO  i = nxl, nxr
    10992                 DO  j = nys, nyn 
     10995                DO  j = nys, nyn
    1099310996                   DO  k = nzb_do, nzt_do
    1099410997                      local_pf(i,j,k) = rad_lw_in_av(k,j,i)
     
    1101411017            ENDIF
    1101511018             DO  i = nxl, nxr
    11016                 DO  j = nys, nyn 
     11019                DO  j = nys, nyn
    1101711020                   DO  k = nzb_do, nzt_do
    1101811021                      local_pf(i,j,k) = rad_lw_out_av(k,j,i)
     
    1102011023                ENDDO
    1102111024             ENDDO
    11022           ENDIF   
     11025          ENDIF
    1102311026          IF ( mode == 'xy' )  grid = 'zu'
    1102411027
     
    1103811041            ENDIF
    1103911042             DO  i = nxl, nxr
    11040                 DO  j = nys, nyn 
     11043                DO  j = nys, nyn
    1104111044                   DO  k = nzb_do, nzt_do
    1104211045                      local_pf(i,j,k) = rad_lw_cs_hr_av(k,j,i)
     
    1106211065            ENDIF
    1106311066             DO  i = nxl, nxr
    11064                 DO  j = nys, nyn 
     11067                DO  j = nys, nyn
    1106511068                   DO  k = nzb_do, nzt_do
    1106611069                      local_pf(i,j,k) = rad_lw_hr_av(k,j,i)
     
    1108611089            ENDIF
    1108711090             DO  i = nxl, nxr
    11088                 DO  j = nys, nyn 
     11091                DO  j = nys, nyn
    1108911092                   DO  k = nzb_do, nzt_do
    1109011093                      local_pf(i,j,k) = rad_sw_in_av(k,j,i)
     
    1111011113            ENDIF
    1111111114             DO  i = nxl, nxr
    11112                 DO  j = nys, nyn 
     11115                DO  j = nys, nyn
    1111311116                   DO  k = nzb, nzt+1
    1111411117                      local_pf(i,j,k) = rad_sw_out_av(k,j,i)
     
    1113411137            ENDIF
    1113511138             DO  i = nxl, nxr
    11136                 DO  j = nys, nyn 
     11139                DO  j = nys, nyn
    1113711140                   DO  k = nzb_do, nzt_do
    1113811141                      local_pf(i,j,k) = rad_sw_cs_hr_av(k,j,i)
     
    1115811161            ENDIF
    1115911162             DO  i = nxl, nxr
    11160                 DO  j = nys, nyn 
     11163                DO  j = nys, nyn
    1116111164                   DO  k = nzb_do, nzt_do
    1116211165                      local_pf(i,j,k) = rad_sw_hr_av(k,j,i)
     
    1117211175
    1117311176    END SELECT
    11174  
     11177
    1117511178 END SUBROUTINE radiation_data_output_2d
    1117611179
     
    1118311186!------------------------------------------------------------------------------!
    1118411187 SUBROUTINE radiation_data_output_3d( av, variable, found, local_pf, nzb_do, nzt_do )
    11185  
     11188
    1118611189
    1118711190    USE indices
     
    1119211195    IMPLICIT NONE
    1119311196
    11194     CHARACTER (LEN=*) ::  variable !< 
     11197    CHARACTER (LEN=*) ::  variable !<
    1119511198
    1119611199    INTEGER(iwp) ::  av          !<
     
    1121711220       RETURN
    1121811221    ENDIF
    11219    
     11222
    1122011223    IF ( var(1:3) /= 'rad'  .AND.  var(1:3) /= 'rtm' )  THEN
    1122111224       found = .FALSE.
     
    1168411687            ENDIF
    1168511688         ENDIF
    11686 !         
     11689!
    1168711690!--   block of RTM output variables
    1168811691!--   variables are intended mainly for debugging and detailed analyse purposes
    1168911692      CASE ( 'rtm_skyvf' )
    11690 !     
     11693!
    1169111694!--      sky view factor
    1169211695         DO isurf = dirstart(ids), dirend(ids)
     
    1175811761!------------------------------------------------------------------------------!
    1175911762 SUBROUTINE radiation_data_output_mask( av, variable, found, local_pf, mid )
    11760  
     11763
    1176111764    USE control_parameters
    11762        
     11765
    1176311766    USE indices
    11764    
     11767
    1176511768    USE kinds
    11766    
     11769
    1176711770
    1176811771    IMPLICIT NONE
    1176911772
    11770     CHARACTER (LEN=*) ::  variable   !< 
     11773    CHARACTER (LEN=*) ::  variable   !<
    1177111774
    1177211775    CHARACTER(LEN=5) ::  grid        !< flag to distinquish between staggered grids
    1177311776
    11774     INTEGER(iwp) ::  av              !< 
    11775     INTEGER(iwp) ::  i               !< 
    11776     INTEGER(iwp) ::  j               !< 
    11777     INTEGER(iwp) ::  k               !< 
     11777    INTEGER(iwp) ::  av              !<
     11778    INTEGER(iwp) ::  i               !<
     11779    INTEGER(iwp) ::  j               !<
     11780    INTEGER(iwp) ::  k               !<
    1177811781    INTEGER(iwp) ::  mid             !< masked output running index
    1177911782    INTEGER(iwp) ::  topo_top_index  !< k index of highest horizontal surface
     
    1178511788    REAL(wp),                                                                  &
    1178611789       DIMENSION(mask_size_l(mid,1),mask_size_l(mid,2),mask_size_l(mid,3)) ::  &
    11787           local_pf   !< 
     11790          local_pf   !<
    1178811791
    1178911792    REAL(wp), DIMENSION(:,:,:), POINTER ::  to_be_resorted  !< points to array which needs to be resorted for output
     
    1191811921       WRITE ( 14 )  rad_net_av
    1191911922    ENDIF
    11920    
     11923
    1192111924    IF ( ALLOCATED( rad_lw_in_xy_av ) )  THEN
    1192211925       CALL wrd_write_string( 'rad_lw_in_xy_av' )
    1192311926       WRITE ( 14 )  rad_lw_in_xy_av
    1192411927    ENDIF
    11925    
     11928
    1192611929    IF ( ALLOCATED( rad_lw_out_xy_av ) )  THEN
    1192711930       CALL wrd_write_string( 'rad_lw_out_xy_av' )
    1192811931       WRITE ( 14 )  rad_lw_out_xy_av
    1192911932    ENDIF
    11930    
     11933
    1193111934    IF ( ALLOCATED( rad_sw_in_xy_av ) )  THEN
    1193211935       CALL wrd_write_string( 'rad_sw_in_xy_av' )
    1193311936       WRITE ( 14 )  rad_sw_in_xy_av
    1193411937    ENDIF
    11935    
     11938
    1193611939    IF ( ALLOCATED( rad_sw_out_xy_av ) )  THEN
    1193711940       CALL wrd_write_string( 'rad_sw_out_xy_av' )
     
    1203012033                                nxr_on_file, nynf, nync, nyn_on_file, nysf,    &
    1203112034                                nysc, nys_on_file, tmp_2d, tmp_3d, found )
    12032  
     12035
    1203312036
    1203412037    USE control_parameters
    12035        
     12038
    1203612039    USE indices
    12037    
     12040
    1203812041    USE kinds
    12039    
     12042
    1204012043    USE pegrid
    1204112044
     
    1204312046    IMPLICIT NONE
    1204412047
    12045     INTEGER(iwp) ::  k               !< 
    12046     INTEGER(iwp) ::  nxlc            !< 
    12047     INTEGER(iwp) ::  nxlf            !< 
    12048     INTEGER(iwp) ::  nxl_on_file     !< 
    12049     INTEGER(iwp) ::  nxrc            !< 
    12050     INTEGER(iwp) ::  nxrf            !< 
    12051     INTEGER(iwp) ::  nxr_on_file     !< 
    12052     INTEGER(iwp) ::  nync            !< 
    12053     INTEGER(iwp) ::  nynf            !< 
    12054     INTEGER(iwp) ::  nyn_on_file     !< 
    12055     INTEGER(iwp) ::  nysc            !< 
    12056     INTEGER(iwp) ::  nysf            !< 
    12057     INTEGER(iwp) ::  nys_on_file     !< 
     12048    INTEGER(iwp) ::  k               !<
     12049    INTEGER(iwp) ::  nxlc            !<
     12050    INTEGER(iwp) ::  nxlf            !<
     12051    INTEGER(iwp) ::  nxl_on_file     !<
     12052    INTEGER(iwp) ::  nxrc            !<
     12053    INTEGER(iwp) ::  nxrf            !<
     12054    INTEGER(iwp) ::  nxr_on_file     !<
     12055    INTEGER(iwp) ::  nync            !<
     12056    INTEGER(iwp) ::  nynf            !<
     12057    INTEGER(iwp) ::  nyn_on_file     !<
     12058    INTEGER(iwp) ::  nysc            !<
     12059    INTEGER(iwp) ::  nysf            !<
     12060    INTEGER(iwp) ::  nys_on_file     !<
    1205812061
    1205912062    LOGICAL, INTENT(OUT)  :: found
    1206012063
    12061     REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d   !< 
    12062 
    12063     REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d   !< 
    12064 
    12065     REAL(wp), DIMENSION(0:0,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d2   !< 
     12064    REAL(wp), DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_2d   !<
     12065
     12066    REAL(wp), DIMENSION(nzb:nzt+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d   !<
     12067
     12068    REAL(wp), DIMENSION(0:0,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) :: tmp_3d2   !<
    1206612069
    1206712070
     
    1207412077          IF ( .NOT. ALLOCATED( rad_net_av ) )  THEN
    1207512078             ALLOCATE( rad_net_av(nysg:nyng,nxlg:nxrg) )
    12076           ENDIF 
     12079          ENDIF
    1207712080          IF ( k == 1 )  READ ( 13 )  tmp_2d
    1207812081          rad_net_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =           &
    1207912082                        tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
    12080                        
     12083
    1208112084       CASE ( 'rad_lw_in_xy_av' )
    1208212085          IF ( .NOT. ALLOCATED( rad_lw_in_xy_av ) )  THEN
    1208312086             ALLOCATE( rad_lw_in_xy_av(nysg:nyng,nxlg:nxrg) )
    12084           ENDIF 
     12087          ENDIF
    1208512088          IF ( k == 1 )  READ ( 13 )  tmp_2d
    1208612089          rad_lw_in_xy_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =      &
    1208712090                        tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
    12088                        
     12091
    1208912092       CASE ( 'rad_lw_out_xy_av' )
    1209012093          IF ( .NOT. ALLOCATED( rad_lw_out_xy_av ) )  THEN
    1209112094             ALLOCATE( rad_lw_out_xy_av(nysg:nyng,nxlg:nxrg) )
    12092           ENDIF 
     12095          ENDIF
    1209312096          IF ( k == 1 )  READ ( 13 )  tmp_2d
    1209412097          rad_lw_out_xy_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =     &
    1209512098                        tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
    12096                        
     12099
    1209712100       CASE ( 'rad_sw_in_xy_av' )
    1209812101          IF ( .NOT. ALLOCATED( rad_sw_in_xy_av ) )  THEN
    1209912102             ALLOCATE( rad_sw_in_xy_av(nysg:nyng,nxlg:nxrg) )
    12100           ENDIF 
     12103          ENDIF
    1210112104          IF ( k == 1 )  READ ( 13 )  tmp_2d
    1210212105          rad_sw_in_xy_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =      &
    1210312106                        tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
    12104                        
     12107
    1210512108       CASE ( 'rad_sw_out_xy_av' )
    1210612109          IF ( .NOT. ALLOCATED( rad_sw_out_xy_av ) )  THEN
    1210712110             ALLOCATE( rad_sw_out_xy_av(nysg:nyng,nxlg:nxrg) )
    12108           ENDIF 
     12111          ENDIF
    1210912112          IF ( k == 1 )  READ ( 13 )  tmp_2d
    1211012113          rad_sw_out_xy_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp)  =     &
    1211112114                        tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp)
    12112                        
     12115
    1211312116       CASE ( 'rad_lw_in' )
    1211412117          IF ( .NOT. ALLOCATED( rad_lw_in ) )  THEN
     
    1212012123                ALLOCATE( rad_lw_in(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1212112124             ENDIF
    12122           ENDIF 
     12125          ENDIF
    1212312126          IF ( k == 1 )  THEN
    1212412127             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1214412147                ALLOCATE( rad_lw_in_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1214512148             ENDIF
    12146           ENDIF 
     12149          ENDIF
    1214712150          IF ( k == 1 )  THEN
    1214812151             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1216812171                ALLOCATE( rad_lw_out(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1216912172             ENDIF
    12170           ENDIF 
     12173          ENDIF
    1217112174          IF ( k == 1 )  THEN
    1217212175             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1219212195                ALLOCATE( rad_lw_out_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1219312196             ENDIF
    12194           ENDIF 
     12197          ENDIF
    1219512198          IF ( k == 1 )  THEN
    1219612199             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1224812251                ALLOCATE( rad_sw_in(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1224912252             ENDIF
    12250           ENDIF 
     12253          ENDIF
    1225112254          IF ( k == 1 )  THEN
    1225212255             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1227212275                ALLOCATE( rad_sw_in_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1227312276             ENDIF
    12274           ENDIF 
     12277          ENDIF
    1227512278          IF ( k == 1 )  THEN
    1227612279             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1229612299                ALLOCATE( rad_sw_out(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1229712300             ENDIF
    12298           ENDIF 
     12301          ENDIF
    1229912302          IF ( k == 1 )  THEN
    1230012303             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
     
    1232012323                ALLOCATE( rad_sw_out_av(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    1232112324             ENDIF
    12322           ENDIF 
     12325          ENDIF
    1232312326          IF ( k == 1 )  THEN
    1232412327             IF ( radiation_scheme == 'clear-sky'  .OR.                    &
Note: See TracChangeset for help on using the changeset viewer.