Ignore:
Timestamp:
Dec 7, 2018 2:14:11 PM (3 years ago)
Author:
banzhafs
Message:

chem_emissions_mod and chem_modules update to comply PALM coding rules

File:
1 edited

Legend:

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

    r3589 r3611  
    2222! Current revisions:
    2323! ------------------
    24 ! 
     24!
    2525!
    2626! Former revisions:
    2727! -----------------
    2828! $Id$
     29! Code update to comply PALM coding rules
     30! removed unnecessary informative messsages/location
     31! messages and corrected comments on PM units from to kg 
     32! bug fix: spcs_name replaced by nvar in DO loops
     33!
     34! 3591 2018-11-30 17:37:32Z suehring
    2935! - Removed salsa dependency.
    3036! - Enabled PARAMETRIZED mode for default surfaces when LSM is not applied but
    3137!   salsa is (M. Kurppa)
    32 ! 
     38!
    3339! 3582 2018-11-29 19:16:36Z suehring
    3440! resler:
     
    6773! Authors:
    6874! --------
    69 ! @author Emmanuele Russo (Fu-Berlin)
     75! @author Emmanuele Russo (FU-Berlin)
    7076! @author Sabine Banzhaf  (FU-Berlin)
    7177! @author Martijn Schaap  (FU-Berlin, TNO Utrecht)
     
    7985!> @todo revise indices of files read from the netcdf: preproc_emission_data and expert_emission_data
    8086!> @todo for now emission data may be passed on a singular vertical level: need to be more flexible
     87!> @todo fill/activate restart option in chem_emissions_init
     88!> @todo discuss dt_emis
    8189!> @note <Enter notes on the module>
    8290!> @bug  <Enter known bugs here>
     
    8694
    8795    USE arrays_3d,                                                             &
    88        ONLY: rho_air
     96        ONLY: rho_air
    8997
    9098    USE control_parameters,                                                    &
    91        ONLY:  initializing_actions, end_time, message_string,                  &
    92               intermediate_timestep_count, dt_3d
     99        ONLY:  end_time, message_string, initializing_actions,                 &
     100               intermediate_timestep_count, dt_3d
    93101 
    94102    USE indices
     
    96104    USE kinds
    97105
    98 #if defined ( __netcdf )
    99     USE NETCDF
     106#if defined( __netcdf )
     107    USE netcdf
    100108#endif
    101109
    102110    USE netcdf_data_input_mod,                                                  &
    103        ONLY: chem_emis_att_type, chem_emis_val_type
     111        ONLY: chem_emis_att_type, chem_emis_val_type
    104112
    105113    USE date_and_time_mod,                                                      &
    106        ONLY: time_default_indices, time_preprocessed_indices,                  &
    107              month_of_year, day_of_month, hour_of_day,                          &
    108              index_mm, index_dd, index_hh
    109 
     114        ONLY: day_of_month, hour_of_day,                                        &
     115             index_mm, index_dd, index_hh,                                      &
     116             month_of_year, hour_of_day,                                        &
     117             time_default_indices, time_preprocessed_indices
     118   
    110119    USE chem_gasphase_mod,                                                      &
    111        ONLY: spc_names
     120        ONLY: spc_names
    112121 
    113122    USE chem_modules
    114123
    115124    USE statistics,                                                             &   
    116            ONLY:  weight_pres
    117 
     125        ONLY:  weight_pres
     126
     127   
    118128    IMPLICIT NONE
    119129
     130!
    120131!-- Declare all global variables within the module
    121132   
    122     CHARACTER (LEN=80)                               :: filename_emis                   !> Variable for the name of the netcdf input file
    123 
    124     INTEGER(iwp)                                     :: i                               !> index 1st selected dimension (some dims are not spatial)
    125     INTEGER(iwp)                                     :: j                               !> index 2nd selected dimension
    126     INTEGER(iwp)                                     :: i_start                         !> Index to start read variable from netcdf along one dims
    127     INTEGER(iwp)                                     :: i_end                           !> Index to end read variable from netcdf in one dims
    128     INTEGER(iwp)                                     :: j_start                         !> Index to start read variable from netcdf in additional dims
    129     INTEGER(iwp)                                     :: j_end                           !> Index to end read variable from netcdf in additional dims
    130     INTEGER(iwp)                                     :: z_start                         !> Index to start read variable from netcdf in additional dims
    131     INTEGER(iwp)                                     :: z_end                           !> Index to end read variable from netcdf in additional     dims
    132     INTEGER(iwp)                                     :: dt_emis                         !> Time Step Emissions
    133     INTEGER(iwp)                                     :: len_index                       !> length of index (used for several indices)
    134     INTEGER(iwp)                                     :: len_index_voc                   !> length of voc index
    135     INTEGER(iwp)                                     :: len_index_pm                    !> length of PMs index
    136     REAL(wp)                                         :: con_factor                      !> Units Conversion Factor
    137  
    138 
    139     ! ---------------------------------------------------------------
    140     ! Set Values of molecules, mols, etc
    141     ! ---------------------------------------------------------------
    142  
    143     !> Avogadro number
    144     REAL, PARAMETER        ::  Avog = 6.02205e23    ! mlc/mol
    145  
    146     !> Dobson units:
    147     REAL, PARAMETER        ::  Dobs = 2.68668e16    ! (mlc/cm2) / DU
    148 
    149     !> sesalt composition:
    150     ! (Seinfeld and Pandis, "Atmospheric Chemistry and Physics",
    151     !  table 7.8 "Composition of Sea-Salt", p. 444)
    152     REAL, PARAMETER        ::  massfrac_Cl_in_seasalt  = 0.5504       ! (kg Cl )/(kg seasalt)
    153     REAL, PARAMETER        ::  massfrac_Na_in_seasalt  = 0.3061       ! (kg Na )/(kg seasalt)
    154     REAL, PARAMETER        ::  massfrac_SO4_in_seasalt = 0.0768       ! (kg SO4)/(kg seasalt)
    155  
    156     !> other numbers
    157     REAL, PARAMETER        ::  xm_seasalt =  74.947e-3                ! kg/mol : NaCl, SO4, ..
    158     REAL, PARAMETER        ::  rho_seasalt = 2.2e3                    ! kg/m3
    159 
    160     !> * amonium sulphate
    161  
    162     REAL, PARAMETER        ::  xm_NH4HSO4  =  xm_NH4 + xm_H + xm_SO4  ! kg/mol
    163     REAL, PARAMETER        ::  rho_NH4HSO4a = 1.8e3                   ! kg/m3
    164 
    165     ! ---------------------------------------------------------------
    166     ! gas
    167     ! ---------------------------------------------------------------
    168  
    169     !> gas constant                       
    170     REAL, PARAMETER        ::  Rgas = 8.3144     ! J/mol/K
    171  
    172     !> gas constant for dry air
    173     REAL, PARAMETER        ::  Rgas_air = Rgas / xm_air   ! J/kg/K
    174  
    175     !> water vapour
    176     REAL, PARAMETER        ::  Rgas_h2o = Rgas / xm_h2o   ! J/kg/K
    177  
    178     !> standard pressure
    179     REAL, PARAMETER        ::  p0 = 1.0e5    ! Pa
    180 
    181     !> sea-level pressure
    182     REAL, PARAMETER        ::  p_sealevel = 1.01325e5    ! Pa  <-- suggestion Bram Bregman
    183 
    184     !> global mean pressure
    185     REAL, PARAMETER        ::  p_global = 98500.0   ! Pa
    186 
    187     !> specific heat of dry air at constant pressure
    188     REAL, PARAMETER        ::  cp_air = 1004.0           ! J/kg/K
    189 
    190     !> Latent heat of evaporation
    191     REAL, PARAMETER        ::  lvap = 2.5e6     !  [J kg-1]
    192 
    193     !> Latent heat of condensation at 0 deg Celcius
    194     !  (heat (J) necesarry to evaporate 1 kg of water)
    195     REAL, PARAMETER        ::  Lc = 22.6e5           ! J/kg
    196  
    197     !> kappa = R/cp = 2/7
    198     REAL, PARAMETER        ::  kappa = 2.0/7.0
    199 
    200     !> Von Karman constant (dry_dep)
    201     REAL, PARAMETER        ::  vkarman = 0.4
    202 
    203     !> Boltzmann constant:
    204     REAL(wp), PARAMETER     ::  kbolz = 1.38066e-23_wp    ! J/K
    205 
    206     !> Inverse Reference Pressure (1/Pa)   
    207     REAL(wp), PARAMETER     ::  pref_i  = 1.0_wp / 100000.0_wp       
    208  
    209     !>
    210     REAL(wp), PARAMETER     ::  r_cp    = 0.286_wp                  !< R / cp (exponent for potential temperature)
     133    CHARACTER (LEN=80) ::  filename_emis             !< Variable for the name of the netcdf input file
     134
     135    INTEGER(iwp) ::  i                               !< index 1st selected dimension (some dims are not spatial)
     136    INTEGER(iwp) ::  j                               !< index 2nd selected dimension
     137    INTEGER(iwp) ::  i_start                         !< Index to start read variable from netcdf along one dims
     138    INTEGER(iwp) ::  i_end                           !< Index to end read variable from netcdf in one dims
     139    INTEGER(iwp) ::  j_start                         !< Index to start read variable from netcdf in additional dims
     140    INTEGER(iwp) ::  j_end                           !< Index to end read variable from netcdf in additional dims
     141    INTEGER(iwp) ::  z_start                         !< Index to start read variable from netcdf in additional dims
     142    INTEGER(iwp) ::  z_end                           !< Index to end read variable from netcdf in additional dims
     143    INTEGER(iwp) ::  dt_emis                         !< Time Step Emissions
     144    INTEGER(iwp) ::  len_index                       !< length of index (used for several indices)
     145    INTEGER(iwp) ::  len_index_voc                   !< length of voc index
     146    INTEGER(iwp) ::  len_index_pm                    !< length of PMs index
     147
     148    REAL(wp) ::  con_factor                          !< Units Conversion Factor
     149                           
     150    REAL(wp), PARAMETER ::  Rgas = 8.3144                    !< gas constant in J/mol/K           
     151    REAL(wp), PARAMETER ::  pref_i  = 1.0_wp / 100000.0_wp   !< Inverse Reference Pressure (1/Pa)
     152    REAL(wp), PARAMETER ::  r_cp    = 0.286_wp               !< R / cp (exponent for potential temperature)
    211153
    212154
    213155    SAVE
    214156
    215 !-- Interfaces Part
    216 
     157
     158!   
    217159!-- Checks Input parameters
    218160    INTERFACE chem_emissions_check_parameters
    219161       MODULE PROCEDURE chem_emissions_check_parameters
    220162    END INTERFACE chem_emissions_check_parameters
    221 
    222 !-- Reading of NAMELIST parameters
    223 !    INTERFACE chem_emissions_parin
    224 !       MODULE PROCEDURE chem_emissions_parin
    225 !    END INTERFACE chem_emissions_parin
    226 
    227 !-- Output of information to the header file
    228 !    INTERFACE chem_emissions_header
    229 !       MODULE PROCEDURE chem_emissions_header
    230 !    END INTERFACE chem_emissions_header
    231 
     163!
    232164!-- Matching Emissions actions 
    233165    INTERFACE chem_emissions_match
    234166       MODULE PROCEDURE chem_emissions_match
    235167    END INTERFACE chem_emissions_match
    236 
     168!
    237169!-- Initialization actions 
    238170    INTERFACE chem_emissions_init
    239171       MODULE PROCEDURE chem_emissions_init
    240172    END INTERFACE chem_emissions_init
    241 
     173!
    242174!-- Setup of Emissions
    243175    INTERFACE chem_emissions_setup
     
    246178
    247179
    248 !-- Public Interfaces
    249     PUBLIC chem_emissions_init,chem_emissions_match,chem_emissions_setup
    250 
     180   
     181    PUBLIC chem_emissions_init, chem_emissions_match, chem_emissions_setup
     182!
    251183!-- Public Variables
    252 
    253     PUBLIC con_factor, len_index,len_index_voc,len_index_pm
     184    PUBLIC con_factor, len_index, len_index_pm, len_index_voc
    254185
    255186 CONTAINS
     
    263194
    264195
    265     !TBD: Where Should we put the call to chem_emissions_check_parameters? In chem_init or in check_parameters?
    266 
    267196    IMPLICIT NONE
    268197
    269     INTEGER(iwp)           ::  tmp
    270 
    271     TYPE(chem_emis_att_type)   ::  emt
    272 !
    273 !--    Call routine for controlling chem_emissions namelist
    274 !    CALL chem_emissions_parin
    275 
    276 !TBD: In case the given emission values do not coincide with the passed names we should think of a solution:
    277 !  I would personally do that if the name of the species differ from the number of emission values, I would
    278 !  print a warning that says that in correspondance of that particular species the value is zero.
    279 !  An alternative would be to initialize the cs_surface_flux values to a negative number.
    280 
    281 !-- Check Emission Species Number equal to number of passed names for the chemistry species:
    282    IF ( SIZE(emt%species_name) /= emt%nspec  )  THEN
    283 
    284        message_string = 'Numbers of input emission species names and number of species'             //          &
     198    INTEGER(iwp) ::  tmp
     199
     200    TYPE(chem_emis_att_type) ::  emt
     201
     202    !
     203    !-- Check Emission Species Number equal to number of passed names for the chemistry species:
     204    IF ( SIZE(emt%species_name) /= emt%nspec  )  THEN
     205
     206       message_string = 'Numbers of input emission species names and number of species'     //      &
    285207                         'for which emission values are given do not match'                 
    286208       CALL message( 'chem_emissions_check_parameters', 'CM0437', 2, 2, 0, 6, 0 )
    287 
    288 
     209           
    289210    ENDIF
    290 
    291     !-- Check Emission Species Number equals to number of passed names for the species
    292     !IF ( SIZE(emt%species_name) /= SIZE(emt%species_index)  )  THEN
    293     !      message_string = 'Number of input emission species names and '             //          &
    294     !                       ' number of passed species indices do not match'                 
    295     !      CALL message( 'chem_emissions_check_parameters', 'CM0101', 2, 2, 0, 6, 0 )
    296 
    297     !ENDIF
    298 
    299     !-- Check Emission Categories
    300 !    IF ( SIZE(chem_emis%cat_name) /= SIZE(chem_emis%cat_index)  )  THEN
    301 !       WRITE( message_string, * )                                                        &
    302 !       'Number of input emission categories name and categories index does not match\\' 
    303 !       CALL message( 'chem_emissions_check_parameters', 'CM0101', 1, 2, 0, 6, 0 )
    304 !    ENDIF
    305 
    306 
    307 
    308     !TBD: Check which other consistency controls should be included
    309 
    310    !TBD: Include check for spatial dimension of netcdf file. If they do not match with the ones
    311    !     of the simulation, the model should print an error. 
    312 
    313 END SUBROUTINE chem_emissions_check_parameters
     211   
     212
     213 END SUBROUTINE chem_emissions_check_parameters
    314214
    315215!------------------------------------------------------------------------------!
     
    317217! ------------
    318218!> Matching the chemical species indices. The routine checks what are the indices of the emission input species
    319 !  and the corresponding ones of the model species. The routine gives as output a vector containing the number
    320 !  of common species: it is important to note that while the model species are distinct, their values could be
    321 !  given to a single species in input: for example, in the case of NO2 and NO, values may be passed in input as NOx values.
     219!> and the corresponding ones of the model species. The routine gives as output a vector containing the number
     220!> of common species: it is important to note that while the model species are distinct, their values could be
     221!> given to a single species in input: for example, in the case of NO2 and NO, values may be passed in input as
     222!> NOx values.
    322223!------------------------------------------------------------------------------!
    323 SUBROUTINE chem_emissions_match(emt_att,len_index)   
    324 
    325 
    326     INTEGER(iwp), INTENT(INOUT)              ::  len_index        !< Variable where to store the number of common species between the input dataset and the model species   
    327 
    328     TYPE(chem_emis_att_type), INTENT(INOUT)  ::  emt_att          !< Chemistry Emission Array containing information for all the input chemical emission species
     224SUBROUTINE chem_emissions_match( emt_att,len_index )   
     225
     226
     227    INTEGER(iwp), INTENT(INOUT)             ::  len_index   !< Variable where to store the number of common species between the input dataset and the model species   
     228
     229    TYPE(chem_emis_att_type), INTENT(INOUT) ::  emt_att     !< Chemistry Emission Array containing information for all the input chemical emission species
    329230   
    330     INTEGER(iwp)                             ::  ind_mod, ind_inp !< Parameters for cycling through chemical model and input species
    331     INTEGER(iwp)                             ::  nspec_emis_inp   !< Variable where to store the number of the emission species in input
    332 
    333     INTEGER(iwp)                             ::  ind_voc          !< Indices to check whether a split for voc should be done
    334 
    335     INTEGER(iwp)                             ::  ispec            !> index for cycle over effective number of emission species
    336 
    337 
    338     !> Tell the user where we are!!
     231    INTEGER(iwp) ::  ind_mod, ind_inp                       !< Parameters for cycling through chemical model and input species
     232    INTEGER(iwp) ::  nspec_emis_inp                         !< Variable where to store the number of the emission species in input
     233    INTEGER(iwp) ::  ind_voc                                !< Indices to check whether a split for voc should be done
     234    INTEGER(iwp) ::  ispec                                  !< index for cycle over effective number of emission species
     235
     236
    339237    CALL location_message( 'Matching input emissions and model chemistry species', .FALSE. )
    340238
    341     !> Number of input emission species.
    342 
     239    !
     240    !-- Number of input emission species.
    343241    nspec_emis_inp=emt_att%nspec
    344242
    345     !> Check the emission mode: DEFAULT, PRE-PROCESSED or PARAMETERIZED  !TBD: Add option for capital or small letters
    346     SELECT CASE(TRIM(mode_emis))
    347 
    348        !> PRE-PROCESSED case: in this case the input species have to coincide with the ones of the model (except VOCs for which we have the VOC split): NO and NO2 in input and not NOx.
    349        CASE ("PRE-PROCESSED")
    350 
    351           CALL location_message( 'chem_emissions PRE-PROCESSED_mode_matching', .FALSE. )
    352 
    353           len_index=0
    354           len_index_voc=0
    355 
    356           !> The first condition is that both the number of model and input emissions species are not null
    357           IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0)) THEN
    358 
    359              !> Cycle over model species
    360              DO ind_mod = 1, SIZE(spc_names)
    361                 !> Cycle over Input species 
     243    !
     244    !-- Check the emission mode: DEFAULT, PRE-PROCESSED or PARAMETERIZED
     245    SELECT CASE( TRIM( mode_emis ) )
     246
     247       !
     248       !-- PRE-PROCESSED mode
     249       CASE ( "PRE-PROCESSED" )
     250
     251          len_index = 0
     252          len_index_voc = 0
     253
     254          IF ( nvar > 0 .AND. (nspec_emis_inp > 0) )  THEN
     255             !
     256             !-- Cycle over model species
     257             DO ind_mod = 1,  nvar
     258                !
     259                !-- Cycle over input species 
    362260                DO ind_inp = 1, nspec_emis_inp
    363261
    364            !> In the PRE-PROCESSED mode each emission species is given input values, made exception for the VOCs, having the total number of VOCs in input: the model then executes a split through the different VOC species
    365                    ! > Check for VOC Species: In input in this case we only have one species (VOC) 
    366                    IF (TRIM(emt_att%species_name(ind_inp))=="VOC") THEN
    367                       !> Cycle over the input voc species: they have to be one of the VOCs of the mechanism used. A list of VOC species for different mechanisms is provided in the module documentation
    368                       DO ind_voc= 1,emt_att%nvoc
    369                          !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
    370                          IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
    371                             len_index=len_index+1
    372                             len_index_voc=len_index_voc+1
     262                   !
     263                   !-- Check for VOC Species 
     264                   IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" )  THEN       
     265                      DO ind_voc = 1, emt_att%nvoc
     266             
     267                         IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) )  THEN
     268                            len_index = len_index + 1
     269                            len_index_voc = len_index_voc + 1
    373270                         ENDIF
    374271                      END DO
    375272                   ENDIF
    376                    !> Other Species
    377                    IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
    378                       len_index=len_index+1
     273                   !
     274                   !-- Other Species
     275                   IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     276                      len_index = len_index + 1
    379277                   ENDIF
    380278                ENDDO
    381279             ENDDO
    382280
    383              !> Allocate array for storing the indices of the matched species
    384              IF (len_index>0) THEN
    385  
    386                 ALLOCATE(match_spec_input(len_index))
    387  
    388                 ALLOCATE(match_spec_model(len_index))
    389 
    390                 IF (len_index_voc>0) THEN
    391 
    392                    ALLOCATE(match_spec_voc_model(len_index_voc))   !> Contains indices of the VOC model species
    393 
    394                    ALLOCATE(match_spec_voc_input(len_index_voc))   !> In input there is only one value for VOCs in the DEFAULT mode. This array contains the indices of the different values of VOC compositions of the input variable VOC_composition
     281             !
     282             !-- Allocate array for storing the indices of the matched species
     283             IF ( len_index > 0 )  THEN
     284 
     285                ALLOCATE( match_spec_input(len_index) )
     286 
     287                ALLOCATE( match_spec_model(len_index) )
     288
     289                IF ( len_index_voc > 0 )  THEN
     290                   !
     291                   !-- contains indices of the VOC model species
     292                   ALLOCATE( match_spec_voc_model(len_index_voc) )
     293                   !
     294                   !-- contains the indices of different values of VOC composition of input variable VOC_composition
     295                   ALLOCATE( match_spec_voc_input(len_index_voc) )
    395296
    396297                ENDIF
    397298
    398                 !> Repeat the same cycle of above but passing the species indices to the newly declared arrays
    399                 len_index=0
    400 
    401                 !> Cycle over model species
    402                 DO ind_mod = 1, SIZE(spc_names) 
    403                    !> Cycle over Input species 
     299                !
     300                !-- pass the species indices to declared arrays
     301                len_index = 0
     302
     303                !
     304                !-- Cycle over model species
     305                DO ind_mod = 1, nvar
     306                   !
     307                   !-- Cycle over Input species 
    404308                   DO ind_inp = 1, nspec_emis_inp
    405 
    406                    !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
    407 
    408                       ! > VOCs
    409                       IF ( TRIM(emt_att%species_name(ind_inp))=="VOC" .AND. ALLOCATED(match_spec_voc_input) ) THEN      
    410                          DO ind_voc= 1,emt_att%nvoc
    411                             IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
    412                                len_index=len_index+1
    413                                len_index_voc=len_index_voc+1
     309                      !
     310                      !-- VOCs
     311                      IF ( TRIM(emt_att%species_name(ind_inp) ) == "VOC" .AND.    &
     312                           ALLOCATED(match_spec_voc_input) )  THEN
     313                         
     314                         DO ind_voc= 1, emt_att%nvoc
     315                            IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN
     316                               len_index = len_index + 1
     317                               len_index_voc = len_index_voc + 1
    414318                       
    415                                match_spec_input(len_index)=ind_inp
    416                                match_spec_model(len_index)=ind_mod
    417 
    418                                match_spec_voc_input(len_index_voc)=ind_voc
    419                                match_spec_voc_model(len_index_voc)=ind_mod                         
     319                               match_spec_input(len_index) = ind_inp
     320                               match_spec_model(len_index) = ind_mod
     321
     322                               match_spec_voc_input(len_index_voc) = ind_voc
     323                               match_spec_voc_model(len_index_voc) = ind_mod                         
    420324                            ENDIF
    421325                         END DO
    422326                      ENDIF
    423327
    424                       !>Other Species
    425                       IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
    426                          len_index=len_index+1
    427                          match_spec_input(len_index)=ind_inp
    428                          match_spec_model(len_index)=ind_mod
     328                      !
     329                      !-- Other Species
     330                      IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     331                         len_index = len_index + 1
     332                         match_spec_input(len_index) = ind_inp
     333                         match_spec_model(len_index) = ind_mod
    429334                      ENDIF
    430335                   END DO
    431336                END DO
    432337
    433              !> In the case there are no species matching, the emission module should not be called
    434338             ELSE
    435 
    436                 message_string = 'Given Chemistry Emission Species'            //          &
    437                                  ' DO NOT MATCH'                                //          &
    438                                  ' model chemical species:'                      //          &
    439                                  ' Chemistry Emissions routine is not called' 
     339                !
     340                !-- in case there are no species matching
     341                message_string = 'Non of given emission species'            //         &
     342                                 ' matches'                                //          &
     343                                 ' model chemical species:'                //          &
     344                                 ' Emission routine is not called' 
    440345                CALL message( 'chem_emissions_matching', 'CM0438', 0, 0, 0, 6, 0 )
    441              ENDIF !> IF (len_index>0)
     346             ENDIF
    442347 
    443348          ELSE
    444349
    445              !In this case either spc_names is null or nspec_emis_inp is not allocated
    446 
     350             !
     351             !-- either spc_names is zero or nspec_emis_inp is not allocated
    447352             message_string = 'Array of Emission species not allocated:'             //          &
    448                               ' Either no emission species are provided as input or'  //          &
    449                               ' no chemical species are used by PALM:'                //          &
    450                               ' Chemistry Emissions routine is not called'                 
     353                              ' Either no emission species are provided as input or'  //         &
     354                              ' no chemical species are used by PALM:'                //         &
     355                              ' Emission routine is not called'                 
    451356             CALL message( 'chem_emissions_matching', 'CM0439', 0, 2, 0, 6, 0 )
    452357
    453           ENDIF !> IF ( (SIZE(spc_names) .GT. 0) .AND. ALLOCATED(nspec_emis_inp))
    454 
    455        !> ------------------------------------------------------------------
    456        !> DEFAULT case
    457 
     358          ENDIF
     359
     360       !
     361       !-- DEFAULT mode
    458362       CASE ("DEFAULT")
    459363
    460           CALL location_message( 'chem_emissions DEFAULT_mode_matching', .FALSE. )
    461 
    462           len_index=0      !>  index for TOTAL number of species   
    463           len_index_voc=0  !>  index for TOTAL number of VOCs
    464           len_index_pm=3   !>  index for TOTAL number of PM Types: PM1, PM2.5, PM10. It is needed because the number of emission inputs and the one of PMs in the model may not be the same.
    465 
    466           !> The first condition is that both the number of model and input emissions species are not null
    467           IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0) ) THEN
    468 
    469              !> Cycle over model species
    470              DO ind_mod = 1, SIZE(spc_names)
    471                 !> Cycle over input species
     364          len_index = 0          !<  index for TOTAL number of species   
     365          len_index_voc = 0      !<  index for TOTAL number of VOCs
     366          len_index_pm = 3       !<  index for TOTAL number of PMs: PM1, PM2.5, PM10.
     367 
     368          IF ( nvar > 0 .AND. nspec_emis_inp > 0 )  THEN
     369
     370             !
     371             !-- Cycle over model species
     372             DO ind_mod = 1, nvar
     373                !
     374                !-- Cycle over input species
    472375                DO ind_inp = 1, nspec_emis_inp
    473376
    474                    ! > Check for VOC Species: In input in this case we only have one species (VOC) 
    475                    IF (TRIM(emt_att%species_name(ind_inp))=="VOC") THEN
    476                       !> Cycle over the voc species given in input: they have to be one of the VOCs of the mechanism used. A list of VOC species for different mechanisms is provided in the module documentation
    477                       DO ind_voc= 1,emt_att%nvoc
    478                          !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
    479                          IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
    480                             len_index=len_index+1
    481                             len_index_voc=len_index_voc+1
     377                   !
     378                   !-- Check for VOC Species 
     379                   IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" )  THEN
     380                      DO ind_voc= 1, emt_att%nvoc
     381                       
     382                         IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN
     383                            len_index = len_index + 1
     384                            len_index_voc = len_index_voc + 1
    482385                         ENDIF
     386                         
    483387                      END DO
    484388                   ENDIF
    485389
    486                    !> PMs: For PMs there is only one input species name for all the PM types. This variable has 3 dimensions, one for each PM type.
    487                    IF (TRIM(emt_att%species_name(ind_inp))=="PM") THEN
    488                       !>pm1
    489                       IF (TRIM(spc_names(ind_mod))=="PM1") THEN
    490                          len_index=len_index+1
    491                       !>pm2.5
    492                       ELSE IF (TRIM(spc_names(ind_mod))=="PM25") THEN
    493                          len_index=len_index+1
    494                       !>pm10
    495                       ELSE IF (TRIM(spc_names(ind_mod))=="PM10") THEN
    496                          len_index=len_index+1   
     390                   !
     391                   !-- PMs: There is one input species name for all PM
     392                   !-- This variable has 3 dimensions, one for PM1, PM2.5 and PM10
     393                   IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" )  THEN
     394                      !
     395                      !-- PM1
     396                      IF ( TRIM( spc_names(ind_mod) ) == "PM1" )  THEN
     397                         len_index = len_index + 1
     398                      !
     399                      !-- PM2.5
     400                      ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" )  THEN
     401                         len_index = len_index + 1
     402                      !
     403                      !-- PM10
     404                      ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" )  THEN
     405                         len_index = len_index + 1
    497406                      ENDIF
    498407                   ENDIF
    499408
    500                    !> NOx: for NOx by definition we have 2 splits. The Emission Input Name in this case is only one: NOx, while in the model we can have NO2 and NO   
    501                    IF (TRIM(emt_att%species_name(ind_inp))=="NOx") THEN
    502                       !>no
    503                       IF (TRIM(spc_names(ind_mod))=="NO") THEN
    504                          len_index=len_index+1
    505                       !>no2
    506                       ELSE IF (TRIM(spc_names(ind_mod))=="NO2") THEN
    507                          len_index=len_index+1                       
     409                   !
     410                   !-- NOx: NO2 and NO   
     411                   IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOx" )  THEN
     412                      !
     413                      !-- NO
     414                      IF ( TRIM( spc_names(ind_mod) ) == "NO" )  THEN
     415                         len_index = len_index + 1
     416                      !
     417                      !-- NO2
     418                      ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" )  THEN
     419                         len_index = len_index + 1
    508420                      ENDIF
    509421                   ENDIF
    510422
    511                    !>SOX: same as for NOx, but with SO2 and SO4
    512                    IF (TRIM(emt_att%species_name(ind_inp))=="SOx") THEN
    513                       !>no
    514                       IF (TRIM(spc_names(ind_mod))=="SO2") THEN
    515                          len_index=len_index+1
    516                       !>no2
    517                       ELSE IF (TRIM(spc_names(ind_mod))=="SO4") THEN
    518                          len_index=len_index+1                       
     423                   !
     424                   !-- SOx: SO2 and SO4
     425                   IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOx" )  THEN
     426                      !
     427                      !-- SO2
     428                      IF ( TRIM( spc_names(ind_mod) ) == "SO2" )  THEN
     429                         len_index = len_index + 1
     430                      !
     431                      !-- SO4
     432                      ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" )  THEN
     433                         len_index = len_index + 1
    519434                      ENDIF
    520435                   ENDIF
    521436
    522                    !> Other Species
    523                    IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
    524                       len_index=len_index+1
     437                   !
     438                   !-- Other Species
     439                   IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     440                      len_index = len_index + 1
    525441                   ENDIF
    526                 END DO   !> number of emission input species
    527              END DO   !> number of model species
    528 
    529 
    530              !> Allocate Arrays
    531              IF (len_index>0) THEN
    532 
    533                 ALLOCATE(match_spec_input(len_index))
    534                 ALLOCATE(match_spec_model(len_index))
    535 
    536                 IF (len_index_voc>0) THEN
    537                    ALLOCATE(match_spec_voc_model(len_index_voc))   !> contains indices of the VOC model species
    538                    ALLOCATE(match_spec_voc_input(len_index_voc))   !> In input there is only one value for VOCs in the DEFAULT mode.
    539                                                                    !  This array contains the indices of the different values of VOC compositions of the input variable VOC_composition
     442                END DO
     443             END DO
     444
     445
     446             !
     447             !-- Allocate arrays
     448             IF ( len_index > 0 )  THEN
     449
     450                ALLOCATE( match_spec_input(len_index) )
     451                ALLOCATE( match_spec_model(len_index) )
     452
     453                IF ( len_index_voc > 0 )  THEN
     454                   !
     455                   !-- Contains indices of the VOC model species
     456                   ALLOCATE( match_spec_voc_model(len_index_voc) )
     457                   !
     458                   !-- Contains the indices of different values of VOC composition
     459                   !-- of input variable VOC_composition
     460                   ALLOCATE( match_spec_voc_input(len_index_voc) )                                                 
    540461                ENDIF
    541462
    542                 !> ASSIGN INDICES
    543                 !> Repeat the same cycles of above, but passing the species indices to the newly declared arrays
    544                 len_index=0
    545                 len_index_voc=0
     463                !
     464                !-- Pass the species indices to declared arrays
     465                len_index = 0
     466                len_index_voc = 0
    546467               
    547                 DO ind_mod = 1, SIZE(spc_names)
     468                DO ind_mod = 1, nvar
    548469                   DO ind_inp = 1, nspec_emis_inp
    549470
    550                       ! > VOCs
    551                       IF ( TRIM(emt_att%species_name(ind_inp))=="VOC" .AND. ALLOCATED(match_spec_voc_input) ) THEN     
    552                          DO ind_voc= 1,emt_att%nvoc
    553                             IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
    554                                len_index=len_index+1
    555                                len_index_voc=len_index_voc+1
     471                      !
     472                      !-- VOCs
     473                      IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND.   &
     474                         ALLOCATED(match_spec_voc_input) )  THEN     
     475                         DO ind_voc= 1, emt_att%nvoc
     476                            IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) )  THEN
     477                               len_index = len_index + 1
     478                               len_index_voc = len_index_voc + 1
    556479                       
    557                                match_spec_input(len_index)=ind_inp
    558                                match_spec_model(len_index)=ind_mod
    559 
    560                                match_spec_voc_input(len_index_voc)=ind_voc
    561                                match_spec_voc_model(len_index_voc)=ind_mod                         
     480                               match_spec_input(len_index) = ind_inp
     481                               match_spec_model(len_index) = ind_mod
     482
     483                               match_spec_voc_input(len_index_voc) = ind_voc
     484                               match_spec_voc_model(len_index_voc) = ind_mod                         
    562485                            ENDIF
    563486                         END DO
    564487                      ENDIF
    565488
    566                       !> PMs:In this case the Inputs have only PM while the model has three different possible types: PM1, PM2.5, PM10. We need an additional index for matching each PM index in the model.
    567                       IF (TRIM(emt_att%species_name(ind_inp))=="PM") THEN
    568                          !>pm1
    569                          IF (TRIM(spc_names(ind_mod))=="PM1") THEN
    570                             len_index=len_index+1
    571 
    572                             match_spec_input(len_index)=ind_inp
    573                             match_spec_model(len_index)=ind_mod
     489                      !
     490                      !-- PMs
     491                      IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" )  THEN
     492                         !
     493                         !-- PM1
     494                         IF ( TRIM( spc_names(ind_mod) ) == "PM1" )  THEN
     495                            len_index = len_index + 1
     496
     497                            match_spec_input(len_index) = ind_inp
     498                            match_spec_model(len_index) = ind_mod
     499                         !
     500                         !-- PM2.5
     501                         ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" )  THEN
     502                            len_index = len_index + 1
     503
     504                            match_spec_input(len_index) = ind_inp
     505                            match_spec_model(len_index) = ind_mod
     506                         !
     507                         !-- PM10
     508                         ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" )  THEN
     509                            len_index = len_index + 1
     510   
     511                            match_spec_input(len_index) = ind_inp
     512                            match_spec_model(len_index) = ind_mod
    574513 
    575                             !match_spec_pm(1)=ind_mod
    576                          !>pm2.5
    577                          ELSE IF (TRIM(spc_names(ind_mod))=="PM25") THEN
    578                             len_index=len_index+1
    579 
    580                             match_spec_input(len_index)=ind_inp
    581                             match_spec_model(len_index)=ind_mod
    582  
    583                             !match_spec_pm(2)=ind_mod
    584                          !>pm10
    585                          ELSE IF (TRIM(spc_names(ind_mod))=="PM10") THEN
    586                             len_index=len_index+1
    587    
    588                             match_spec_input(len_index)=ind_inp
    589                             match_spec_model(len_index)=ind_mod
    590  
    591                             !match_spec_pm(3)=ind_mod
    592514                         ENDIF
    593515                      ENDIF
    594516
    595                       !> NOx - The same as for PMs but here the model species are only 2: NO and NO2
    596                       IF (TRIM(emt_att%species_name(ind_inp))=="NOx") THEN
    597                          !>no
    598                          IF (TRIM(spc_names(ind_mod))=="NO") THEN
    599                             len_index=len_index+1
    600 
    601                             match_spec_input(len_index)=ind_inp
    602                             match_spec_model(len_index)=ind_mod
     517                      !
     518                      !-- NOx
     519                      IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOx" )  THEN
     520                         !
     521                         !-- NO
     522                         IF ( TRIM( spc_names(ind_mod) ) == "NO" )  THEN
     523                            len_index = len_index + 1
     524
     525                            match_spec_input(len_index) = ind_inp
     526                            match_spec_model(len_index) = ind_mod
     527                         !
     528                         !-- NO2
     529                         ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" )  THEN
     530                            len_index = len_index + 1
     531
     532                            match_spec_input(len_index) = ind_inp
     533                            match_spec_model(len_index) = ind_mod
    603534 
    604                             !match_spec_nox(1)=ind_mod
    605                          !>no2
    606                          ELSE IF (TRIM(spc_names(ind_mod))=="NO2") THEN
    607                             len_index=len_index+1
    608 
    609                             match_spec_input(len_index)=ind_inp
    610                             match_spec_model(len_index)=ind_mod
    611  
    612                            ! match_spec_nox(2)=ind_mod
    613535                         ENDIF
    614536                      ENDIF
    615537
    616                       !> SOx
    617                       IF (TRIM(emt_att%species_name(ind_inp))=="SOx") THEN
    618                          !>so2
    619                          IF (TRIM(spc_names(ind_mod))=="SO2") THEN
    620                             len_index=len_index+1
    621 
    622                             match_spec_input(len_index)=ind_inp
    623                             match_spec_model(len_index)=ind_mod
     538                      !
     539                      !-- SOx
     540                      IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOx" ) THEN
     541                         !
     542                         !-- SO2
     543                         IF ( TRIM( spc_names(ind_mod) ) == "SO2" )  THEN
     544                            len_index = len_index + 1
     545
     546                            match_spec_input(len_index) = ind_inp
     547                            match_spec_model(len_index) = ind_mod
    624548 
    625                            ! match_spec_sox(1)=ind_mod
    626                          !>so4
    627                          ELSE IF (TRIM(spc_names(ind_mod))=="SO4") THEN
    628                             len_index=len_index+1
    629 
    630                             match_spec_input(len_index)=ind_inp
    631                             match_spec_model(len_index)=ind_mod
     549                         !
     550                         !-- SO4
     551                         ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" ) THEN
     552                            len_index = len_index + 1
     553
     554                            match_spec_input(len_index) = ind_inp
     555                            match_spec_model(len_index) = ind_mod
    632556 
    633                            ! match_spec_sox(2)=ind_mod
    634557                         ENDIF
    635558                      ENDIF
    636559
    637                       !> Other Species
    638                       IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
    639                          len_index=len_index+1
     560                      !
     561                      !-- Other Species
     562                      IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     563                         len_index = len_index + 1
    640564                           
    641                          match_spec_input(len_index)=ind_inp
    642                          match_spec_model(len_index)=ind_mod
     565                         match_spec_input(len_index) = ind_inp
     566                         match_spec_model(len_index) = ind_mod
    643567                      ENDIF
    644568                   END DO
     
    647571             ELSE
    648572
    649                 message_string = 'Given Chemistry Emission Species'            //          &
    650                                  ' DO NOT MATCH'                                //          &
    651                                  ' model chemical species'                      //          &
    652                                  ' Chemistry Emissions routine is not called'         
     573                message_string = 'Non of given Emission Species'            //         &
     574                                 ' matches'                                //          &
     575                                 ' model chemical species'                 //          &
     576                                 ' Emission routine is not called'         
    653577                CALL message( 'chem_emissions_matching', 'CM0440', 0, 0, 0, 6, 0 )
    654578
     
    658582
    659583             message_string = 'Array of Emission species not allocated: '            //          &
    660                               ' Either no emission species are provided as input or'  //          &
    661                               ' no chemical species are used by PALM:'                //          &
    662                               ' Chemistry Emissions routine is not called'                                   
     584                              ' Either no emission species are provided as input or'  //         &
     585                              ' no chemical species are used by PALM:'                //         &
     586                              ' Emission routine is not called'                                   
    663587             CALL message( 'chem_emissions_matching', 'CM0441', 0, 2, 0, 6, 0 )
    664588 
    665589          ENDIF
    666 
    667 !-- CASE parameterized: In the parameterized case the user gives in input the exact species names of the chemical mechanism: no split is required for NOx, SOx, PMs and VOCs, but units of emissions inputs must be in (mole/m**)/s, made exception for PMs.
    668 
     590 
     591       !
     592       !-- PARAMETERIZED mode
    669593       CASE ("PARAMETERIZED")
    670594
    671           len_index=0
    672 
    673          !spc_names have to be greater than zero for proceeding
    674           IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0) ) THEN
    675 
    676 
    677              !cycle over model species
    678              DO ind_mod = 1, SIZE(spc_names)
    679                 ind_inp=1
    680                 DO  WHILE ( TRIM( surface_csflux_name( ind_inp ) ) /= 'novalue' )   !<'novalue' is the default 
    681                    IF (TRIM(surface_csflux_name( ind_inp ))==TRIM(spc_names(ind_mod))) THEN
    682                       len_index=len_index+1
     595          len_index = 0
     596
     597          IF ( nvar > 0 .AND. nspec_emis_inp > 0 )  THEN
     598
     599             !
     600             !-- Cycle over model species
     601             DO ind_mod = 1, nvar
     602                ind_inp = 1
     603                DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' )    !< 'novalue' is the default 
     604                   IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     605                      len_index = len_index + 1
    683606                   ENDIF
    684                    ind_inp=ind_inp+1
    685 
     607                   ind_inp = ind_inp + 1
    686608                ENDDO
    687609             ENDDO
    688610
    689              !Proceed only if there are matched species
    690              IF ( len_index .GT. 0 ) THEN
    691 
    692                 !Allocation of Arrays of the matched species
    693                 ALLOCATE(match_spec_input(len_index))
    694  
    695                 ALLOCATE(match_spec_model(len_index))
    696 
    697                 !Assign corresponding indices of input and model matched species
    698                 len_index=0
    699 
    700                 DO ind_mod = 1, SIZE(spc_names) 
     611             IF ( len_index > 0 )  THEN
     612
     613                !
     614                !-- Allocation of Arrays of the matched species
     615                ALLOCATE( match_spec_input(len_index) )
     616 
     617                ALLOCATE( match_spec_model(len_index) )
     618
     619                !
     620                !-- Pass the species indices to declared arrays   
     621                len_index = 0
     622
     623                DO ind_mod = 1, nvar 
    701624                   ind_inp = 1
    702                    DO  WHILE ( TRIM( surface_csflux_name( ind_inp ) ) /= 'novalue' )   !<'novalue' is the default 
    703                       IF (TRIM( surface_csflux_name( ind_inp ) ) == TRIM(spc_names(ind_mod)))           THEN
    704                          len_index=len_index+1
    705                          match_spec_input(len_index)=ind_inp
    706                          match_spec_model(len_index)=ind_mod
     625                   DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' )    
     626                      IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) )  THEN
     627                         len_index = len_index + 1
     628                         match_spec_input(len_index) = ind_inp
     629                         match_spec_model(len_index) = ind_mod
    707630                      ENDIF
    708                    ind_inp=ind_inp+1
     631                   ind_inp = ind_inp + 1
    709632                   END DO
    710633                END DO
    711634
    712                 ! also, Add AN if to check that when the surface_csflux are passed to the namelist, also the street factor values are provided
    713                 DO  ispec = 1 , len_index
    714 
    715                    IF ( emiss_factor_main(match_spec_input(ispec)) .LT. 0 .AND. &
    716                         emiss_factor_side(match_spec_input(ispec)) .LT. 0 ) THEN
    717 
    718                       message_string = 'PARAMETERIZED emissions mode selected:'            //          &
     635                !
     636                !-- Check
     637                DO ispec = 1, len_index
     638
     639                   IF ( emiss_factor_main(match_spec_input(ispec) ) < 0 .AND.    &
     640                        emiss_factor_side(match_spec_input(ispec) ) < 0 )  THEN
     641
     642                      message_string = 'PARAMETERIZED emissions mode selected:'            //           &
    719643                                       ' EMISSIONS POSSIBLE ONLY ON STREET SURFACES'        //          &
    720644                                       ' but values of scaling factors for street types'    //          &
     
    730654             ELSE
    731655               
    732                 message_string = 'Given Chemistry Emission Species'            //          &
    733                                  ' DO NOT MATCH'                                //          &
    734                                  ' model chemical species'                      //          &
    735                                  ' Chemistry Emissions routine is not called'         
     656                message_string = 'Non of given Emission Species'            //          &
     657                                 ' matches'                                //           &
     658                                 ' model chemical species'                 //           &
     659                                 ' Emission routine is not called'         
    736660                CALL message( 'chem_emissions_matching', 'CM0443', 0, 0, 0, 6, 0 )
    737661             ENDIF
     
    739663          ELSE
    740664     
    741              message_string = 'Array of Emission species not allocated: '            //          &
     665             message_string = 'Array of Emission species not allocated: '            //           &
    742666                              ' Either no emission species are provided as input or'  //          &
    743667                              ' no chemical species are used by PALM.'                //          &
    744                               ' Chemistry Emissions routine is not called'                                   
     668                              ' Emission routine is not called'                                   
    745669             CALL message( 'chem_emissions_matching', 'CM0444', 0, 2, 0, 6, 0 )
    746670 
     
    748672
    749673
    750 !-- CASE when emission module is switched on but mode_emis is not specified or it is given the wrong name
     674       !
     675       !-- If emission module is switched on but mode_emis is not specified or it is given the wrong name
    751676       CASE DEFAULT
    752677
    753           message_string = 'Chemistry Emissions Module switched ON, but'            //          &
     678          message_string = 'Emission Module switched ON, but'            //                         &
    754679                           ' either no emission mode specified or incorrectly given :'  //          &
    755680                           ' please, pass the correct value to the namelist parameter "mode_emis"'                 
     
    758683       END SELECT
    759684
    760 END SUBROUTINE chem_emissions_match
    761 
     685 END SUBROUTINE chem_emissions_match
     686
     687 
    762688!------------------------------------------------------------------------------!
    763689! Description:
    764690! ------------
    765691!> Initialization:
    766 !> Netcdf reading, arrays allocation and first calculation of cssws fluxes at timestep 0
    767 !>
     692!> Netcdf reading, arrays allocation and first calculation of cssws
     693!> fluxes at timestep 0
    768694!------------------------------------------------------------------------------!
    769  SUBROUTINE chem_emissions_init(emt_att,emt,nspec_out)
     695 SUBROUTINE chem_emissions_init( emt_att, emt, nspec_out )
    770696
    771697    USE surface_mod,                                                           &
    772        ONLY:  surf_lsm_h,surf_def_h,surf_usm_h
     698        ONLY:  surf_def_h, surf_lsm_h, surf_usm_h
    773699   
    774700    IMPLICIT NONE
    775701 
    776     TYPE(chem_emis_att_type),INTENT(INOUT)                            :: emt_att   !> variable where to store all the information of
    777                                                                                    !  emission inputs whose values do not depend
    778                                                                                    !  on the considered species
    779 
    780     TYPE(chem_emis_val_type),INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt      !> variable where to store emission inputs values,
    781                                                                                    !  depending on the considered species
     702    CHARACTER (LEN=80)          ::  units                                           !< units of inputs
    782703   
    783     INTEGER(iwp),INTENT(INOUT)                                        :: nspec_out !> number of outputs of chemistry emission routines
    784 
    785     CHARACTER (LEN=80)                                                :: units     !> units of chemistry inputs
    786  
    787     INTEGER(iwp)                                                      :: ispec     !> Index to go through the emission chemical species
    788 
    789 
    790 !-- Actions for initial runs : TBD: needs to be updated
    791   IF (  TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
     704    INTEGER(iwp), INTENT(INOUT) ::  nspec_out                                       !< number of outputs
     705
     706    INTEGER(iwp)                :: ispec                                            !< running index
     707
     708    TYPE(chem_emis_att_type), INTENT(INOUT) ::  emt_att                             !< variable where to store all the information
     709                                                                                    !< of emission inputs whose values do not
     710                                                                                    !< depend on the considered species
     711
     712    TYPE(chem_emis_val_type), INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt      !< variable where to store emission input 
     713                                                                                    !< values depending on the considered species
     714
     715!   
     716!-- Actions for initial runs
     717!  IF (  TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
    792718!--    ...
    793    
     719!   
    794720!
    795721!-- Actions for restart runs
     722!  ELSE
     723!--    ...
     724!
     725!  ENDIF
     726
     727
     728  CALL location_message( 'Starting initialization of emission arrays', .FALSE. )
     729
     730  !
     731  !-- Matching
     732  CALL  chem_emissions_match( emt_att, nspec_out )
     733 
     734  IF ( nspec_out == 0 )  THEN
     735 
     736     emission_output_required = .FALSE.
     737
    796738  ELSE
    797 !--    ...
    798 
    799   ENDIF
    800 
    801 
    802   CALL location_message( 'Starting initialization of chemistry emissions arrays', .FALSE. )
    803 
    804 
    805   !-- Match Input and Model emissions
    806   CALL  chem_emissions_match(emt_att,nspec_out)
    807 
    808   !> If nspec_out==0, then do not use emission module: The corresponding message is already printed in the matching routine
    809   IF ( nspec_out == 0 ) THEN
    810  
    811      emission_output_required=.FALSE.
    812 
    813   ELSE
    814 
    815      emission_output_required=.TRUE.
    816 
    817 
    818 !-----------------------------------------------------------------
    819      !> Set molecule masses'
    820      ALLOCATE(emt_att%xm(nspec_out))
    821      ! loop over emisisons:
     739
     740     emission_output_required = .TRUE.
     741
     742
     743     !
     744     !-- Set molecule masses'
     745     ALLOCATE( emt_att%xm(nspec_out) )
    822746
    823747     DO ispec = 1, nspec_out
    824        ! switch:
    825         SELECT CASE ( TRIM(spc_names(match_spec_model(ispec))) )
    826            CASE ( 'SO2','SO4' ) ; emt_att%xm(ispec) = xm_S + xm_O * 2      ! kg/mole
    827            CASE ( 'NO','NO2' )  ; emt_att%xm(ispec) = xm_N + xm_O * 2      ! kg/mole  NO2 equivalent
    828            CASE ( 'NH3' ) ; emt_att%xm(ispec) = xm_N + xm_H * 3  ! kg/mole
    829            CASE ( 'CO'  ) ; emt_att%xm(ispec) = xm_C + xm_O      ! kg/mole
    830            CASE ( 'CO2' ) ; emt_att%xm(ispec) = xm_C + xm_O * 2  ! kg/mole
    831            CASE ( 'CH4' ) ; emt_att%xm(ispec) = xm_C + xm_H * 4  ! kg/mole     
    832            CASE ( 'HNO3' ); emt_att%xm(ispec) = xm_H + xm_N + xm_O*3 !kg/mole 
     748        SELECT CASE ( TRIM( spc_names(match_spec_model(ispec)) ) )
     749           CASE ( 'SO2' ); emt_att%xm(ispec) = xm_S + xm_O * 2        !< kg/mole
     750           CASE ( 'SO4' ); emt_att%xm(ispec) = xm_S + xm_O * 4        !< kg/mole
     751           CASE ( 'NO' ); emt_att%xm(ispec) = xm_N + xm_O             !< kg/mole
     752           CASE ( 'NO2' ); emt_att%xm(ispec) = xm_N + xm_O * 2        !< kg/mole   
     753           CASE ( 'NH3' ); emt_att%xm(ispec) = xm_N + xm_H * 3        !< kg/mole
     754           CASE ( 'CO'  ); emt_att%xm(ispec) = xm_C + xm_O            !< kg/mole
     755           CASE ( 'CO2' ); emt_att%xm(ispec) = xm_C + xm_O * 2        !< kg/mole
     756           CASE ( 'CH4' ); emt_att%xm(ispec) = xm_C + xm_H * 4        !< kg/mole     
     757           CASE ( 'HNO3' ); emt_att%xm(ispec) = xm_H + xm_N + xm_O*3  !< kg/mole 
    833758           CASE DEFAULT
    834              !--  TBD: check if this hase to be removed
    835               !emt_att%xm(ispec) = 1.0_wp
     759              emt_att%xm(ispec) = 1.0_wp
    836760        END SELECT
    837761     ENDDO
    838762
    839763   
    840      !> ASSIGN EMISSION VALUES for the different emission modes
    841      SELECT CASE(TRIM(mode_emis))   !TBD: Add the option for CApital or small letters
    842 
    843 
    844         !> PRE-PROCESSED case
    845         CASE ("PRE-PROCESSED")
    846 
    847            IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE(emis_distribution(nzb:nzt+1,0:ny,0:nx,nspec_out)) 
    848 
    849            CALL location_message( 'emis_distribution array allocated in PRE-PROCESSED mode', .FALSE. )
    850  
    851            !> Calculate the values of the emissions at the first time step
    852            CALL chem_emissions_setup(emt_att,emt,nspec_out)
    853 
    854         !> Default case
    855         CASE ("DEFAULT")
    856 
    857            IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE( emis_distribution( 1, 0:ny, 0:nx, nspec_out ) )
    858        
    859            CALL location_message( 'emis_distribution array allocated in DEFAULT mode', .FALSE. )
    860 
    861            !> Calculate the values of the emissions at the first time step
    862            CALL chem_emissions_setup(emt_att,emt,nspec_out)
    863 
    864         !> PARAMETERIZED case
    865         CASE ("PARAMETERIZED")
    866 
    867            CALL location_message( 'emis_distribution array allocated in PARAMETERIZED mode', .FALSE. )
    868 
    869            ! For now for PAR and DEF values only, first vertical level of emis_distribution is allocated, while for PRE-PROCESSED all.
    870            IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE(emis_distribution(1,0:ny,0:nx,nspec_out))
    871 
    872            !> Calculate the values of the emissions at the first time step
    873            CALL chem_emissions_setup(emt_att,emt,nspec_out)
     764     !
     765     !-- assign emission values
     766     SELECT CASE ( TRIM( mode_emis ) )   
     767
     768
     769        !
     770        !-- PRE-PROCESSED case
     771        CASE ( "PRE-PROCESSED" )
     772
     773           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE( emis_distribution(nzb:nzt+1,0:ny,0:nx,nspec_out) ) 
     774 
     775           !
     776           !-- Get emissions at the first time step
     777           CALL chem_emissions_setup( emt_att, emt, nspec_out )
     778
     779        !
     780        !-- Default case
     781        CASE ( "DEFAULT" )
     782
     783           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE( emis_distribution(1,0:ny,0:nx,nspec_out) )
     784
     785           !
     786           !-- Get emissions at the first time step
     787           CALL chem_emissions_setup( emt_att, emt, nspec_out )
     788
     789        !
     790        !-- PARAMETERIZED case
     791        CASE ( "PARAMETERIZED" )
     792
     793           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE( emis_distribution(1,0:ny,0:nx,nspec_out) )
     794
     795           !
     796           !-- Get emissions at the first time step
     797           CALL chem_emissions_setup( emt_att, emt, nspec_out)
    874798
    875799     END SELECT
     
    887811!-------------------------------------------------------------------------------!
    888812
    889  SUBROUTINE chem_emissions_setup(emt_att,emt,nspec_out)
     813 SUBROUTINE chem_emissions_setup( emt_att, emt, nspec_out )
    890814
    891815   USE arrays_3d,                                                    &
    892816       ONLY:  dzw
    893    USE grid_variables,                                                        &
     817   USE grid_variables,                                               &
    894818       ONLY: dx, dy
    895    USE indices,                                                               &
    896        ONLY: nnx,nny,nnz
    897    USE surface_mod,                                                           &
    898        ONLY:  surf_lsm_h,surf_def_h,surf_usm_h
    899    USE netcdf_data_input_mod,                                                 &
     819   USE indices,                                                      &
     820       ONLY: nnx, nny, nnz
     821   USE surface_mod,                                                  &
     822       ONLY:  surf_def_h, surf_lsm_h, surf_usm_h
     823   USE netcdf_data_input_mod,                                        &
    900824       ONLY:  street_type_f
    901    USE arrays_3d,                                                             &       
     825   USE arrays_3d,                                                    &       
    902826       ONLY: hyp, pt
    903827
     828   
    904829 IMPLICIT NONE
    905830 
    906     !--- IN/OUT
    907 
    908     TYPE(chem_emis_att_type),INTENT(INOUT)                            ::  emt_att    !> variable where to store all the information of
    909                                                                                      !  emission inputs whose values do not depend
    910                                                                                      !  on the considered species
    911 
    912     TYPE(chem_emis_val_type),INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt        !> variable where to store emission inputs values,
    913                                                                                      !  depending on the considered species
    914 
    915     INTEGER,INTENT(IN)                                                ::  nspec_out  !> Output of routine chem_emis_matching with number
    916                                                                                      !  of matched species
    917     !---
    918 
    919     REAL(wp),ALLOCATABLE,DIMENSION(:,:)                               ::  delta_emis !> Term to add to the emission_distribution array
    920                                                                                      !  for each of the categories at each time step
    921                                                                                      !  in the default mode
    922 
    923     CHARACTER(LEN=80)                                                 ::  units      !> Units of the emissions
    924 
    925     INTEGER(iwp)                                                      ::  icat       !> Index for number of categories
    926     INTEGER(iwp)                                                      ::  ispec      !> index for number of species
    927     INTEGER(iwp)                                                      ::  i_pm_comp  !> index for number of PM components
    928     INTEGER(iwp)                                                      ::  ivoc       !> Index for number of components of  VOCs
    929     INTEGER(iwp)                                                      ::  time_index !> Index for time
     831
     832    TYPE(chem_emis_att_type), INTENT(INOUT) ::  emt_att                         !< variable to store emission information                                                                           
     833
     834    TYPE(chem_emis_val_type), INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt  !< variable to store emission input values,
     835                                                                                !< depending on the considered species
     836
     837    INTEGER,INTENT(IN) ::  nspec_out                                            !< Output of matching routine with number
     838                                                                                !< of matched species
     839
     840    CHARACTER(LEN=80) ::  units                                                 !< Units of the emission data
     841
     842    INTEGER(iwp) ::  i                                                          !< running index for grid in x-direction
     843    INTEGER(iwp) ::  j                                                          !< running index for grid in y-direction
     844    INTEGER(iwp) ::  k                                                          !< running index for grid in z-direction
     845    INTEGER(iwp) ::  m                                                          !< running index for horizontal surfaces
    930846   
    931     REAL(wp),ALLOCATABLE, DIMENSION(:)                                ::  time_factor !> factor for time scaling of emissions
    932     REAL(wp),ALLOCATABLE, DIMENSION(:,:)                              ::  emis
    933 
    934     REAL(wp), DIMENSION(24)                                           :: par_emis_time_factor      !< time factors
    935                                                                                       !  for the parameterized mode: these are fixed for each hour
    936                                                                                       !  of a single day.
    937     REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr)                    ::  conv_to_ratio !> factor used for converting input
    938                                                                                         !  to adimensional concentration ratio
    939 
    940    ! ------------------------------------------
    941     ! variables for the conversion of indices i and j according to surface_mod
    942 
    943     INTEGER(iwp)                                                      ::  i          !> running index for grid in x-direction
    944     INTEGER(iwp)                                                      ::  j          !> running index for grid in y-direction
    945     INTEGER(iwp)                                                      ::  k          !> running index for grid in z-direction
    946     INTEGER(iwp)                                                      ::  m          !> running index for horizontal surfaces
    947 
    948     REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr)                    ::  tmp_temp
    949 
    950     ! --- const -------------------------------
    951     !-CONVERSION FACTORS: TIME
    952     ! number of sec per hour = 3600   
    953     REAL, PARAMETER   ::  s_per_hour = 3600.0  !  (s)/(hour)
    954     ! number of sec per day = 86400   
    955     REAL, PARAMETER   ::  s_per_day = 86400.0  !  (s)/(day)
    956     ! number of hours in a year of 365 days:
    957     REAL, PARAMETER   ::  hour_per_year = 8760.0 !> TBD: What about leap years?
    958     ! number of hours in a day:
    959     REAL, PARAMETER   ::  hour_per_day = 24.0
    960 
    961     ! conversion from hours to seconds (s/hour) = 1/3600.0 ~ 0.2777778   
    962     REAL, PARAMETER   ::  hour_to_s = 1/s_per_hour  !  (hour)/(s)
    963     ! conversion from day to seconds (s/day) = 1/86400.0 ~ 1.157407e-05   
    964     REAL, PARAMETER   ::  day_to_s = 1/s_per_day   !  (day)/(s)
    965     ! conversion from year to sec (s/year) = 1/31536000.0 ~ 3.170979e-08   
    966     REAL, PARAMETER   ::  year_to_s = 1/(s_per_hour*hour_per_year)  !  (year)/(s)
    967 
    968     !-CONVERSION FACTORS: WEIGHT
    969     !  Conversion from tons to kg (kg/tons) = 100.0/1 ~ 100
    970     REAL, PARAMETER   ::  tons_to_kg = 100  !  (tons)/(kg)
    971     !  Conversion from g to kg (kg/g) = 1/1000.0 ~ 0.001
    972     REAL, PARAMETER   ::  g_to_kg = 0.001  !  (g)/(kg)
    973     !  Conversion from g to kg (kg/g) = 1/1000.0 ~ 0.001
    974     REAL, PARAMETER   ::  miug_to_kg = 0.000000001  !  (microng)/(kg)
    975 
    976     !< CONVERSION FACTORS: fraction to ppm
    977     REAL(wp), PARAMETER   ::  ratio2ppm  = 1.0e06_wp 
     847    INTEGER(iwp) ::  icat                                                       !< Index for number of categories
     848    INTEGER(iwp) ::  ispec                                                      !< index for number of species
     849    INTEGER(iwp) ::  i_pm_comp                                                  !< index for number of PM components
     850    INTEGER(iwp) ::  ivoc                                                       !< Index for number of VOCs
     851    INTEGER(iwp) ::  time_index                                                 !< Index for time
     852
     853    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::  delta_emis                       
     854    REAL(wp), ALLOCATABLE, DIMENSION(:)   ::  time_factor                       !< factor for time scaling of emissions
     855    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::  emis
     856
     857    REAL(wp), DIMENSION(24) :: par_emis_time_factor                             !< time factors for the parameterized mode:
     858                                                                               !< fixed houlry profile for example day
     859    REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  conv_to_ratio          !< factor used for converting input
     860                                                                                !< to concentration ratio
     861    REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) ::  tmp_temp
     862
     863    !
     864    !-- CONVERSION FACTORS: TIME 
     865    REAL(wp), PARAMETER ::  s_per_hour = 3600.0                       !< number of sec per hour (s)/(hour)   
     866    REAL(wp), PARAMETER ::  s_per_day = 86400.0                       !< number of sec per day (s)/(day) 
     867    REAL(wp), PARAMETER ::  hour_per_year = 8760.0                    !< number of hours in a year of 365 days 
     868    REAL(wp), PARAMETER ::  hour_per_day = 24.0                       !< number of hours in a day
     869   
     870    REAL(wp), PARAMETER ::  hour_to_s = 1/s_per_hour                  !< conversion from hours to seconds (s/hour) ~ 0.2777778     
     871    REAL(wp), PARAMETER ::  day_to_s = 1/s_per_day                    !< conversion from day to seconds (s/day) ~ 1.157407e-05
     872    REAL(wp), PARAMETER ::  year_to_s = 1/(s_per_hour*hour_per_year)  !< conversion from year to sec (s/year) ~ 3.170979e-08
     873    !
     874    !-- CONVERSION FACTORS: WEIGHT 
     875    REAL(wp), PARAMETER ::  tons_to_kg = 100                          !< Conversion from tons to kg (kg/tons)   
     876    REAL(wp), PARAMETER ::  g_to_kg = 0.001                           !< Conversion from g to kg (kg/g)
     877    REAL(wp), PARAMETER ::  miug_to_kg = 0.000000001                  !< Conversion from g to kg (kg/g)
     878    !
     879    !-- CONVERSION FACTORS: fraction to ppm
     880    REAL(wp), PARAMETER ::  ratio2ppm  = 1.0e06 
    978881    !------------------------------------------------------   
    979882
    980     IF ( emission_output_required ) THEN
    981 
    982     !>  Set emis_dt  !TBD: this is the same as dt_chem. We should consider the fact that dt_emis should be the timestep of input emissions or better defined, the timestep at which the emission routines are called: for now one hour. It should be made changeable.
    983 
     883    IF ( emission_output_required )  THEN
     884
     885       !
     886       !-- Set emis_dt 
    984887       IF ( call_chem_at_all_substeps )  THEN
    985888
     
    993896
    994897
    995     ! --- CHECK UNITS
    996     !>-----------------------------------------------------
    997     !> Conversion of the units to the ones employed in PALM.
    998     !> Possible temporal units of input emissions data are: year, hour and second;
    999     !> the mass could be expressed in: tons, kilograms or grams.
    1000     !> IN the PARAMETERIZED mode no conversion is possible: in this case INPUTS have to have fixed units.
    1001     !>-----------------------------------------------------
    1002 
    1003         IF (TRIM(mode_emis)=="DEFAULT" .OR. TRIM(mode_emis)=="PRE-PROCESSED" ) THEN
    1004 
    1005           SELECT CASE(TRIM(emt_att%units))
    1006           !> kilograms
    1007              CASE ( 'kg/m2/s','KG/M2/S')
    1008                 CALL location_message( 'Units of the emissions are consistent with'  //          &
    1009                                        ' the ones employed in the PALM-4U model ', .FALSE. )
     898       !
     899       !-- Conversion of units to the ones employed in PALM
     900       !-- In PARAMETERIZED mode no conversion is performed: in this case input units are fixed
     901
     902        IF ( TRIM( mode_emis ) == "DEFAULT" .OR. TRIM( mode_emis ) == "PRE-PROCESSED" )  THEN
     903
     904          SELECT CASE ( TRIM( emt_att%units ) )
     905             !
     906             !-- kilograms
     907             CASE ( 'kg/m2/s', 'KG/M2/S' )
     908               
    1010909                con_factor=1
    1011910
    1012              CASE ('kg/m2/hour','KG/M2/HOUR')
    1013                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
     911             CASE ('kg/m2/hour', 'KG/M2/HOUR' )
    1014912
    1015913                con_factor=hour_to_s
    1016914
    1017              CASE ( 'kg/m2/day','KG/M2/DAY' )
    1018                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1019 
     915             CASE ( 'kg/m2/day', 'KG/M2/DAY' )
     916               
    1020917                con_factor=day_to_s
    1021918
    1022              CASE ( 'kg/m2/year','KG/M2/YEAR' )
    1023                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1024 
     919             CASE ( 'kg/m2/year', 'KG/M2/YEAR' )
     920             
    1025921                con_factor=year_to_s
    1026922
    1027           !> Tons
    1028              CASE ( 'ton/m2/s','TON/M2/S' )
    1029                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1030 
     923             !
     924             !-- Tons
     925             CASE ( 'ton/m2/s', 'TON/M2/S' )
     926               
    1031927                con_factor=tons_to_kg
    1032928
    1033              CASE ( 'ton/m2/hour','TON/M2/HOUR' )
    1034                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1035 
     929             CASE ( 'ton/m2/hour', 'TON/M2/HOUR' )
     930               
    1036931                con_factor=tons_to_kg*hour_to_s
    1037932
    1038              CASE ( 'ton/m2/year','TON/M2/YEAR' )
    1039                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1040 
     933             CASE ( 'ton/m2/year', 'TON/M2/YEAR' )
     934               
    1041935                con_factor=tons_to_kg*year_to_s
    1042936
    1043           !> Grams
    1044              CASE ( 'g/m2/s','G/M2/S' )
    1045                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1046 
     937             !
     938             !-- Grams
     939             CASE ( 'g/m2/s', 'G/M2/S' )
     940         
    1047941                con_factor=g_to_kg
    1048942
    1049              CASE ( 'g/m2/hour','G/M2/HOUR' )
    1050                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
     943             CASE ( 'g/m2/hour', 'G/M2/HOUR' )
    1051944
    1052945                con_factor=g_to_kg*hour_to_s
    1053946
    1054              CASE ( 'g/m2/year','G/M2/YEAR' )
    1055                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1056 
     947             CASE ( 'g/m2/year', 'G/M2/YEAR' )
     948 
    1057949                con_factor=g_to_kg*year_to_s
    1058950
    1059           !> Micro Grams
    1060              CASE ( 'micrograms/m2/s','MICROGRAMS/M2/S' )
    1061                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1062 
     951             !
     952             !-- Micrograms
     953             CASE ( 'micrograms/m2/s', 'MICROGRAMS/M2/S' )
     954 
    1063955                con_factor=miug_to_kg
    1064956
    1065              CASE ( 'micrograms/m2/hour','MICROGRAMS/M2/HOUR' )
    1066                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1067 
     957             CASE ( 'micrograms/m2/hour', 'MICROGRAMS/M2/HOUR' )
     958 
    1068959                con_factor=miug_to_kg*hour_to_s
    1069960
    1070              CASE ( 'micrograms/m2/year','MICROGRAMS/M2/YEAR' )
    1071                 CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
    1072 
     961             CASE ( 'micrograms/m2/year', 'MICROGRAMS/M2/YEAR' )
     962 
    1073963                con_factor=miug_to_kg*year_to_s
    1074964
    1075965             CASE DEFAULT 
    1076                 message_string = 'The Units of the provided input chemistry emission species' // &
     966                message_string = 'The Units of the provided emission input' // &
    1077967                                 ' are not the ones required by PALM-4U: please check '      // &
    1078                                  ' chemistry emission module documentation.'                                 
     968                                 ' emission module documentation.'                                 
    1079969                CALL message( 'chem_emissions_setup', 'CM0446', 2, 2, 0, 6, 0 )
    1080970
     
    1082972
    1083973         
    1084        !> PRE-PROCESSED and parameterized mode
    1085        ELSE
    1086  
    1087           message_string = 'No Units conversion required for units of chemistry emissions' // &
    1088                            ' of the PARAMETERIZED mode: units have to be provided in'     //  &
    1089                            ' micromole/m**2/day for GASES and'                            //  &
    1090                            ' kg/m**2/day for PMs'                     
    1091           CALL message( 'chem_emissions_setup', 'CM0447', 0, 0, 0, 6, 0 )
    1092 
    1093974       ENDIF
    1094975
    1095     !> Conversion factors (if the units are kg/m**2/s) we have to convert them to ppm/s
    1096      DO i=nxl,nxr
    1097           DO j=nys,nyn
    1098           !> Derive Temperature from Potential Temperature
     976       !
     977       !-- Conversion factor to convert  kg/m**2/s to ppm/s
     978       DO i = nxl, nxr
     979          DO j = nys, nyn
     980             !
     981             !-- Derive Temperature from Potential Temperature
    1099982             tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * ( hyp(nzb:nzt+1) * pref_i )**r_cp
    1100 
    1101              !> We need to pass to cssws<- (ppm/s) * dz.
    1102              !  Input is Nmole/(m^2*s)
    1103              !  To go to ppm*dz basically we need to multiply the input by (m**2/N)*dz
    1104              !  (m**2/N)*dz == V/N
    1105              !  V/N=RT/P
    1106 
    1107              !>    m**3/Nmole              (J/mol)*K^-1           K                      Pa         
     983             
     984             
     985             !>  We need to pass to cssws <- (ppm/s) * dz
     986             !>  Input is Nmole/(m^2*s)
     987             !>  To go to ppm*dz multiply the input by (m**2/N)*dz
     988             !>  (m**2/N)*dz == V/N
     989             !>  V/N = RT/P
     990             !>  m**3/Nmole               (J/mol)*K^-1           K                      Pa         
    1108991             conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) ) 
    1109992          ENDDO
    1110993       ENDDO
    1111994
    1112     !>------------------------------------------------
    1113     !> Start The Processing of the input data
    1114 
    1115         emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp
    1116 
    1117     !>-----------------------------------------------
    1118     !> Distinguish between DEFAULT, PRE-PROCESSED and PARAMETERIZED mode when calculating time_factors: only done for DEFAULT mode. For the PARAMETERIZED mode there is a time factor but it is fixed in the model
     995
     996       !
     997       !-- Initialize
     998       emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp
     999
    11191000 
    1120     !> PRE-PROCESSED MODE
    1121        IF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
    1122 
    1123           !> Update time indices
    1124           CALL time_preprocessed_indices(index_hh)
    1125 
    1126           CALL location_message( 'PRE-PROCESSED MODE: No time-factor specification required', .FALSE. )
    1127 
    1128        ELSEIF (TRIM(mode_emis)=="DEFAULT") THEN
    1129 
    1130           CALL location_message( 'DEFAULT MODE: time-factor specification required', .FALSE. )
    1131  
    1132        !> Allocate array where to store temporary emission values     
    1133           IF(.NOT. ALLOCATED(emis)) ALLOCATE(emis(nys:nyn,nxl:nxr))
    1134 
    1135        !> Allocate time factor per emitted component category
    1136           ALLOCATE(time_factor(emt_att%ncat)) 
    1137 
    1138        !> Read-in HOURLY emission time factor
    1139           IF (TRIM(time_fac_type)=="HOUR") THEN
    1140 
    1141           !> Update time indices
    1142              CALL time_default_indices(month_of_year,day_of_month,hour_of_day,index_hh)
    1143 
    1144           !> Check if the index is less or equal to the temporal dimension of HOURLY emission files             
    1145              IF (index_hh .LE. SIZE(emt_att%hourly_emis_time_factor(1,:))) THEN
    1146 
    1147              !> Read-in the correspondant time factor             
    1148                 time_factor(:)= emt_att%hourly_emis_time_factor(:,index_hh)     
     1001       !
     1002       !-- PRE-PROCESSED MODE
     1003       IF ( TRIM( mode_emis ) == "PRE-PROCESSED" )  THEN
     1004
     1005          !
     1006          !-- Update time indices
     1007          CALL time_preprocessed_indices( index_hh )
     1008
     1009       ELSEIF ( TRIM( mode_emis ) == "DEFAULT" ) THEN
     1010
     1011          !
     1012          !-- Allocate array where to store temporary emission values     
     1013          IF ( .NOT. ALLOCATED(emis) ) ALLOCATE( emis(nys:nyn,nxl:nxr) )
     1014          !
     1015          !-- Allocate time factor per category
     1016          ALLOCATE( time_factor(emt_att%ncat) ) 
     1017          !
     1018          !-- Read-in hourly emission time factor
     1019          IF ( TRIM( time_fac_type ) == "HOUR" )  THEN
     1020
     1021             !
     1022             !-- Update time indices
     1023             CALL time_default_indices( month_of_year, day_of_month, hour_of_day, index_hh )
     1024             !
     1025             !-- Check if the index is less or equal to the temporal dimension of HOURLY emission files             
     1026             IF ( index_hh <= SIZE( emt_att%hourly_emis_time_factor(1,:) ) ) THEN
     1027                !
     1028                !-- Read-in the correspondant time factor             
     1029                time_factor(:) = emt_att%hourly_emis_time_factor(:,index_hh)     
    11491030
    11501031             ELSE
    11511032
    1152                 message_string = 'The "HOUR" time-factors (DEFAULT mode) of the chemistry emission species'           // &
     1033                message_string = 'The "HOUR" time-factors in the DEFAULT mode '           // &
    11531034                              ' are not provided for each hour of the total simulation time'     
    11541035                CALL message( 'chem_emissions_setup', 'CM0448', 2, 2, 0, 6, 0 )
    11551036
    11561037             ENDIF
    1157 
    1158        !> Read-in MDH emissions time factors
    1159           ELSEIF (TRIM(time_fac_type)=="MDH") THEN
    1160 
    1161           !> Update time indices     
    1162              CALL time_default_indices(daytype_mdh,month_of_year,day_of_month,hour_of_day,index_mm,index_dd,index_hh)
    1163 
    1164 
    1165           !> Check if the index is less or equal to the temporal dimension of MDH emission files             
    1166              IF ((index_hh+index_dd+index_mm) .LE. SIZE(emt_att%mdh_emis_time_factor(1,:))) THEN
    1167 
    1168              !> Read-in the correspondant time factor             
    1169                 time_factor(:)=emt_att%mdh_emis_time_factor(:,index_mm)*emt_att%mdh_emis_time_factor(:,index_dd)*   &
     1038          !
     1039          !-- Read-in MDH emissions time factors
     1040          ELSEIF ( TRIM( time_fac_type ) == "MDH" )  THEN
     1041
     1042             !
     1043             !-- Update time indices     
     1044             CALL time_default_indices( daytype_mdh, month_of_year, day_of_month,     &
     1045                  hour_of_day, index_mm, index_dd,index_hh )
     1046
     1047             !
     1048             !-- Check if the index is less or equal to the temporal dimension of MDH emission files             
     1049             IF ( ( index_hh + index_dd + index_mm) <= SIZE( emt_att%mdh_emis_time_factor(1,:) ) )  THEN
     1050                !
     1051                !-- Read-in the correspondant time factor             
     1052                time_factor(:) = emt_att%mdh_emis_time_factor(:,index_mm) * emt_att%mdh_emis_time_factor(:,index_dd) *   &
    11701053                                                                         emt_att%mdh_emis_time_factor(:,index_hh)
    11711054     
    11721055             ELSE
    11731056
    1174                 message_string = 'The "MDH" time-factors (DEFAULT mode) of the chemistry emission species'           // &
     1057                message_string = 'The "MDH" time-factors in the DEFAULT mode '           // &
    11751058                              ' are not provided for each hour/day/month  of the total simulation time'     
    11761059                CALL message( 'chem_emissions_setup', 'CM0449', 2, 2, 0, 6, 0 )
     
    11791062
    11801063          ELSE
    1181           !> condition when someone used the DEFAULT mode but forgets to indicate the time-factor in the namelist
    1182              
    1183              message_string = 'The time-factor (DEFAULT mode) of the chemistry emission species'           // &
    1184                               ' is not provided in the NAMELIST'     
     1064                 
     1065             message_string = 'In the DEFAULT mode the time factor'           //  &
     1066                              ' has to be defined in the NAMELIST'     
    11851067             CALL message( 'chem_emissions_setup', 'CM0450', 2, 2, 0, 6, 0 )
    11861068         
    11871069          ENDIF
    11881070
    1189        !> PARAMETERIZED MODE
    1190        ELSEIF (TRIM(mode_emis)=="PARAMETERIZED") THEN
    1191           CALL location_message( 'PARAMETERIZED MODE: time-factor specification is fixed: '                         // &
    1192                                  ' 24 values for every day of the year ', .FALSE. )
     1071       !
     1072       !-- PARAMETERIZED MODE
     1073       ELSEIF ( TRIM( mode_emis ) == "PARAMETERIZED" )  THEN
    11931074       
    1194           !Assign Constant Values of time factors, diurnal time profile for traffic sector:
    1195           par_emis_time_factor( : ) = &
    1196             (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, 0.056, 0.053, 0.051, 0.051, 0.052, 0.055, &
     1075          !
     1076          !-- assign constant values of time factors, diurnal time profile for traffic sector
     1077          par_emis_time_factor( : ) =  &
     1078            (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, 0.056, 0.053, 0.051, 0.051, 0.052, 0.055,  &
    11971079               0.059, 0.061, 0.064, 0.067, 0.069, 0.069, 0.049, 0.039, 0.039, 0.029, 0.024, 0.019 /)
    11981080         
    1199           !> in this case allocate time factor each hour in a day
    1200           IF (.NOT. ALLOCATED(time_factor)) ALLOCATE(time_factor(1))
    1201 
    1202           !>Pass the values of time factors in the parameterized mode to the time_factor variable. in this case index_hh==hour_of_day
    1203           index_hh=hour_of_day
     1081          IF ( .NOT. ALLOCATED( time_factor ) ) ALLOCATE( time_factor(1) )
     1082
     1083          !
     1084          !-- Get time-factor for specific hour
     1085          index_hh = hour_of_day
    12041086
    12051087          time_factor(1) = par_emis_time_factor(index_hh)
     
    12071089       ENDIF
    12081090
    1209 !--------------------------------
    1210 !--  EMISSION DISTRIBUTION Calculation
    1211 
    1212        !> PARAMETERIZED CASE
    1213        IF ( mode_emis == "PARAMETERIZED" ) THEN
    1214 
    1215           DO ispec=1,nspec_out
    1216 
    1217              !> Values are still micromoles/(m**2*s). Units in this case are always micromoles/m**2*day (or kilograms/m**2*day for PMs)
    1218              emis_distribution(1,nys:nyn,nxl:nxr,ispec)=surface_csflux(match_spec_input(ispec))*time_factor(1)*hour_to_s
     1091       
     1092       !
     1093       !--  Emission distribution calculation
     1094
     1095       !
     1096       !-- PARAMETERIZED case
     1097       IF ( TRIM( mode_emis ) == "PARAMETERIZED" )  THEN
     1098
     1099          DO ispec = 1, nspec_out
     1100
     1101             !
     1102             !-- Units are micromoles/m**2*day (or kilograms/m**2*day for PMs)
     1103             emis_distribution(1,nys:nyn,nxl:nxr,ispec) = surface_csflux(match_spec_input(ispec)) *  &
     1104                                                           time_factor(1) * hour_to_s
    12191105
    12201106          ENDDO
    12211107
    1222        !> PRE-PROCESSED CASE
    1223        ELSEIF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
    1224 
    1225           !> Start Cycle over Species
    1226           DO ispec=1,nspec_out !> nspec_out represents the number of species in common between
    1227                                !  the emission input data and the chemistry mechanism used
     1108       !
     1109       !-- PRE-PROCESSED case
     1110       ELSEIF ( TRIM( mode_emis ) == "PRE-PROCESSED" )  THEN
     1111
     1112          !
     1113          !-- Cycle over species:
     1114          !-- nspec_out represents the number of species in common between the emission input data
     1115          !-- and the chemistry mechanism used
     1116          DO ispec=1,nspec_out 
    12281117   
    1229              emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emt(match_spec_input(ispec))%                               &
    1230                                                                    preproc_emission_data(index_hh,1,nys+1:nyn+1,nxl+1:nxr+1)* &
    1231                                                                       con_factor
     1118             emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emt(match_spec_input(ispec))%                                 &
     1119                                                            preproc_emission_data(index_hh,1,nys+1:nyn+1,nxl+1:nxr+1) * &
     1120                                                                con_factor
    12321121         
    12331122          ENDDO
    1234 
    1235       !TBD: At the moment the default mode considers a single vertical level (the surface). So we need to change it accordingly or eventually include the variable vertical_profile to be used in the default mode i we want to consider additional levels.
    1236 
    1237        ELSE IF (TRIM(mode_emis)=="DEFAULT") THEN
    1238 
    1239        !> Allocate array for the emission value corresponding to a specific category and time factor
    1240           ALLOCATE(delta_emis(nys:nyn,nxl:nxr)) 
    1241 
    1242 
    1243        !> Start Cycle over categories
     1123         
     1124       !
     1125       !-- DEFAULT case
     1126       ELSEIF ( TRIM( mode_emis ) == "DEFAULT" )  THEN
     1127
     1128          !
     1129          !-- Allocate array for the emission value corresponding to a specific category and time factor
     1130          ALLOCATE( delta_emis(nys:nyn,nxl:nxr) ) 
     1131
     1132          !
     1133          !-- Cycle over categories
    12441134          DO icat = 1, emt_att%ncat
    12451135 
    1246           !> Start Cycle over Species
    1247              DO ispec=1,nspec_out !> nspec_out represents the number of species in common between
    1248                                   !  the emission input data and the chemistry mechanism used
     1136             !
     1137             !-- Cycle over Species
     1138             !-- nspec_out represents the number of species in common between the emission input data
     1139             !-- and the chemistry mechanism used
     1140             DO ispec = 1, nspec_out
    12491141
    12501142                emis(nys:nyn,nxl:nxr) = emt(match_spec_input(ispec))%default_emission_data(icat,nys+1:nyn+1,nxl+1:nxr+1)
    12511143
    1252 !TBD: The consideration of dt_emis of the input data is still missing. Basically the emissions could be provided every 10, 30 minutes and not always at one hour. This should be eventually solved in the date_and_time mode routine.
    1253 
    1254              !> the time factors are 24 for each day. When multiplied by a daily value, they allow to have an hourly value. Then to convert it to seconds, we still have to divide this value by 3600.
    1255              !  So given any units, we convert them to seconds and finally multiply them by 24 ((value/sec)*(24*3600)=value/day ---- (value/day)*time_factor=value/hour ---(value/hour)/(3600)=value/sec )
    1256              !                                                                                 ((value/sec)*(24*3600)*time_factor)/3600=24*(value/sec)*time_factor                         
    1257 
    1258              !> NOX Compositions
    1259                 IF (TRIM(spc_names(match_spec_model(ispec)))=="NO") THEN
    1260                 !>             Kg/m2*s                   kg/m2*s
    1261                    delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1262                                                  emt_att%nox_comp(icat,1)*con_factor*hour_per_day
    1263 
    1264                    emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
    1265 
    1266                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="NO2") THEN
    1267 
    1268                    delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1269                                                  emt_att%nox_comp(icat,2)*con_factor*hour_per_day
    1270 
    1271                    emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
    1272  
    1273              !> SOX Compositions
    1274 
    1275                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO2") THEN
    1276                    !>             Kg/m2*s                   kg/m2*s
    1277                    delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1278                                                  emt_att%sox_comp(icat,1)*con_factor*hour_per_day
    1279 
    1280                    emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
    1281 
    1282                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO4") THEN
    1283                    !>             Kg/m2*s                   kg/m2*s
    1284                    delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1285                                                  emt_att%sox_comp(icat,2)*con_factor*hour_per_day
    1286 
    1287                    emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
    1288  
    1289 
    1290              !> PMs should be in kg/m**2/s, so conversion factor is here still required
    1291              !> PM1 Compositions
    1292                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1") THEN
    1293 
    1294                 !> Cycle over the different pm components for PM1 type
    1295                    DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,1))
    1296 
    1297                       delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1298                                                     emt_att%pm_comp(icat,i_pm_comp,1)*con_factor*hour_per_day
    1299 
    1300                       emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
     1144
     1145                !
     1146                !-- NOx
     1147                IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "NO" )  THEN
     1148               
     1149                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &                 !<  kg/m2*s
     1150                                                  emt_att%nox_comp(icat,1) * con_factor * hour_per_day
     1151
     1152                   emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
    13011153                                                                 delta_emis(nys:nyn,nxl:nxr)
    1302  
     1154
     1155                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "NO2" )  THEN
     1156
     1157                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &                 !<  kg/m2*s
     1158                                                  emt_att%nox_comp(icat,2) * con_factor * hour_per_day
     1159
     1160                   emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1161                                                                 delta_emis(nys:nyn,nxl:nxr)
     1162 
     1163                !
     1164                !-- SOx
     1165                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "SO2" )  THEN
     1166                 
     1167                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &                 !<  kg/m2*s
     1168                                                 emt_att%sox_comp(icat,1) * con_factor * hour_per_day
     1169
     1170                   emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1171                                                                 delta_emis(nys:nyn,nxl:nxr)
     1172
     1173                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "SO4" )  THEN
     1174                 
     1175                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &                 !<  kg/m2*s
     1176                                                 emt_att%sox_comp(icat,2) * con_factor * hour_per_day
     1177
     1178                   emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1179                                                                 delta_emis(nys:nyn,nxl:nxr)
     1180 
     1181
     1182                !
     1183                !-- PMs
     1184                !-- PM1
     1185                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1" )  THEN
     1186
     1187                   !
     1188                   !-- Cycle over PM1 components
     1189                   DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,1) )
     1190
     1191                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &               !<  kg/m2*s
     1192                                                     emt_att%pm_comp(icat,i_pm_comp,1) * con_factor * hour_per_day
     1193
     1194                      emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1195                                                                    delta_emis(nys:nyn,nxl:nxr)
    13031196                   ENDDO
    13041197
    1305              !> PM2.5 Compositions
    1306                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM25") THEN
    1307 
    1308                 !> Cycle over the different pm components for PM2.5 type
    1309                    DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,2))
    1310 
    1311                       delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1312                                                     emt_att%pm_comp(icat,i_pm_comp,2)*con_factor*hour_per_day
    1313 
    1314                       emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
    1315                                                                  delta_emis(nys:nyn,nxl:nxr)
     1198                !
     1199                !-- PM2.5
     1200                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM25" )  THEN
     1201
     1202                   !
     1203                   !-- Cycle over PM2.5 components
     1204                   DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,2) )
     1205
     1206                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &               !<  kg/m2*s
     1207                                                     emt_att%pm_comp(icat,i_pm_comp,2) * con_factor * hour_per_day
     1208
     1209                      emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1210                                                                    delta_emis(nys:nyn,nxl:nxr)
    13161211 
    13171212                   ENDDO
    13181213
    1319              !> PM10 Compositions
    1320                 ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
    1321 
    1322                 !> Cycle over the different pm components for PM10 type
    1323                    DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,3)) 
    1324 
    1325                       delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
    1326                                                     emt_att%pm_comp(icat,i_pm_comp,3)*con_factor*hour_per_day
     1214                !
     1215                !-- PM10
     1216                ELSEIF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM10" )  THEN
     1217
     1218                   !
     1219                   !-- Cycle over PM10 components
     1220                   DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,3) ) 
     1221
     1222                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &               !<  kg/m2*s
     1223                                                     emt_att%pm_comp(icat,i_pm_comp,3) * con_factor * hour_per_day
    13271224
    13281225                      emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
     
    13311228                   ENDDO
    13321229
    1333              !> VOCs Compositions: for VOCs, the input value is provided in kg/m**2*s but the composition is provided in mole/kg: a conversion factor for the input that could be eventually provided in, for example, tons/(m**2*s) is still required
    1334                 ELSE IF  (SIZE(match_spec_voc_input) .GT. 0) THEN
    1335 
    1336                   !TBD: maybe we can avoid the cycle
    1337                    DO ivoc= 1,SIZE(match_spec_voc_input)
    1338 
    1339                       IF (TRIM(spc_names(match_spec_model(ispec)))==TRIM(emt_att%voc_name(ivoc))) THEN   
    1340 
    1341                          delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)*                                    &
    1342                                                        emt_att%voc_comp(icat,match_spec_voc_input(ivoc))*con_factor*hour_per_day
    1343 
    1344                          emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
    1345                                                                     delta_emis(nys:nyn,nxl:nxr)
     1230                !
     1231                !-- VOCs
     1232                ELSEIF ( SIZE( match_spec_voc_input ) > 0 )  THEN
     1233
     1234                   DO ivoc = 1, SIZE( match_spec_voc_input )
     1235
     1236                      IF ( TRIM( spc_names(match_spec_model(ispec)) ) == TRIM( emt_att%voc_name(ivoc) ) )  THEN   
     1237
     1238                         delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &
     1239                                                        emt_att%voc_comp(icat,match_spec_voc_input(ivoc)) * &
     1240                                                         con_factor * hour_per_day
     1241
     1242                         emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1243                                                                       delta_emis(nys:nyn,nxl:nxr)
    13461244
    13471245                      ENDIF                       
     
    13491247                   ENDDO
    13501248               
    1351              !> General case (other species)
     1249                !
     1250                !-- any other species
    13521251                ELSE
    13531252
    1354                    delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)*con_factor*hour_per_day
    1355  
    1356                    emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
    1357 
    1358                 ENDIF  ! IF (spc_names==)
    1359 
    1360                 !> for every species and category set emis to 0 so to avoid overwriting. TBD: discuss whether necessary
    1361 
     1253                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * &
     1254                                                  con_factor * hour_per_day
     1255 
     1256                   emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
     1257                                                                 delta_emis(nys:nyn,nxl:nxr)
     1258
     1259                ENDIF
     1260               
    13621261                emis(:,:)= 0
    1363 
     1262               
    13641263             ENDDO
    1365 
    1366           !> for every category set delta_emis to 0 so to avoid overwriting. TBD: discuss whether necessary
    1367 
    1368           delta_emis(:,:)=0
    1369  
     1264             
     1265             delta_emis(:,:)=0
     1266         
    13701267          ENDDO
    13711268
    1372        ENDIF  !> mode_emis
    1373 
    1374 !-------------------------------------------------------------------------------------------------------
    1375 !--- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to cssws
    1376 !-------------------------------------------------------------------------------------------------------
    1377 
     1269       ENDIF
     1270
     1271       
     1272!
     1273!-- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to cssws
     1274!
    13781275!-- PARAMETERIZED mode
    1379        !> For the PARAMETERIZED mode units of inputs are always micromoles/(m**2*s). In this case we do not need the molar mass for conversion into ppms
    1380        IF (TRIM(mode_emis)=="PARAMETERIZED") THEN
     1276       !
     1277       !-- Units of inputs are micromoles/(m**2*s)
     1278       IF ( TRIM( mode_emis ) == "PARAMETERIZED" )  THEN
    13811279
    13821280          IF ( street_type_f%from_file )  THEN
    13831281
    1384        !> Streets are lsm surfaces, hence, no usm surface treatment required
    1385              IF (surf_lsm_h%ns .GT. 0 ) THEN
     1282             !
     1283             !-- Streets are lsm surfaces, hence, no usm surface treatment required
     1284             IF ( surf_lsm_h%ns > 0 )  THEN
    13861285                DO  m = 1, surf_lsm_h%ns
    13871286                   i = surf_lsm_h%i(m)
     
    13891288                   k = surf_lsm_h%k(m)
    13901289
    1391 
    13921290                   IF ( street_type_f%var(j,i) >= main_street_id  .AND. street_type_f%var(j,i) < max_street_id )  THEN
    13931291
    1394                    !> Cycle over already matched species
    1395                       DO  ispec=1,nspec_out
    1396 
    1397                          !> PMs are already in mass units: kilograms
    1398                          IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
    1399                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
    1400                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
    1401 
    1402                             !              kg/(m^2*s) *kg/m^3
     1292                      !
     1293                      !-- Cycle over matched species
     1294                      DO  ispec = 1, nspec_out
     1295
     1296                         !
     1297                         !-- PMs are already in kilograms
     1298                         IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1 "  &
     1299                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM25"  &
     1300                             .OR. TRIM( spc_names(match_spec_model(ispec)) )=="PM10")  THEN
     1301
     1302                            !
     1303                            !--           kg/(m^2*s) * kg/m^3
    14031304                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec)) * &
    1404                             !                                                       kg/(m^2*s)
    1405                                                                                 emis_distribution(1,j,i,ispec)*        &
    1406                             !                                                    kg/m^3
    1407                                                                                 rho_air(k)   
    1408 
     1305                                                                           emis_distribution(1,j,i,ispec) *            &  !< in kg/(m^2*s)
     1306                                                                            rho_air(k)                                    !< in kg/m^3
     1307                           
     1308                         !
     1309                         !-- Other Species
     1310                         !-- Inputs are micromoles
    14091311                         ELSE
    14101312
    1411                          !> Other Species: inputs are micromoles: they have to be converted in moles
    1412                          !                 ppm/s *m *kg/m^3               
    1413                             surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec))*   &
    1414                          !                                                    micromoles/(m^2*s)
    1415                                                                           emis_distribution(1,j,i,ispec) *              &
    1416                          !                                                    m^3/Nmole
    1417                                                                           conv_to_ratio(k,j,i)*                         &       
    1418                          !                                                    kg/m^3
    1419                                                                           rho_air(k)   
    1420  
     1313                            !   
     1314                            !--             ppm/s *m *kg/m^3               
     1315                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec)) * &
     1316                                                                           emis_distribution(1,j,i,ispec) *            &  !< in micromoles/(m^2*s)
     1317                                                                            conv_to_ratio(k,j,i) *                     &  !< in m^3/Nmole     
     1318                                                                             rho_air(k)                                   !< in kg/m^3
    14211319
    14221320                         ENDIF
     
    14271325                   ELSEIF ( street_type_f%var(j,i) >= side_street_id  .AND. street_type_f%var(j,i) < main_street_id )  THEN
    14281326
    1429                    !> Cycle over already matched species
    1430                       DO  ispec=1,nspec_out
    1431 
    1432                          !> PMs are already in mass units: micrograms
    1433                          IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
    1434                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
    1435                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
    1436 
    1437                             !              kg/(m^2*s) *kg/m^3                               
    1438                             surf_lsm_h%cssws(match_spec_model(ispec),m)= emiss_factor_side(match_spec_input(ispec)) *   &
    1439                             !                                                       kg/(m^2*s)
    1440                                                                                 emis_distribution(1,j,i,ispec)*        &
    1441                             !                                                    kg/m^3
    1442                                                                                 rho_air(k)   
     1327                      !
     1328                      !-- Cycle over matched species
     1329                      DO  ispec = 1, nspec_out
     1330
     1331                         !
     1332                         !-- PMs are already in kilograms
     1333                         IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1"  &
     1334                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM25"  &
     1335                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM10" )  THEN
     1336
     1337                            !
     1338                            !--           kg/(m^2*s) * kg/m^3
     1339                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_side(match_spec_input(ispec)) * &
     1340                                                                           emis_distribution(1,j,i,ispec) *            &  !< in kg/(m^2*s)
     1341                                                                            rho_air(k)                                    !< in kg/m^3   
     1342                         !
     1343                         !-- Other species
     1344                         !-- Inputs are micromoles
    14431345                         ELSE
    14441346                       
    1445 
    1446                          !>Other Species: inputs are micromoles
    1447                          !                 ppm/s *m *kg/m^3               
    1448                             surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_side(match_spec_input(ispec)) *   &
    1449                          !                                                    micromoles/(m^2*s)
    1450                                                                           emis_distribution(1,j,i,ispec) *              &
    1451                          !                                                    m^3/Nmole
    1452                                                                           conv_to_ratio(k,j,i)*                         &       
    1453                          !                                                    kg/m^3
    1454                                                                           rho_air(k)   
     1347                            !   
     1348                            !--             ppm/s *m *kg/m^3     
     1349                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_side(match_spec_input(ispec)) * &
     1350                                                                           emis_distribution(1,j,i,ispec) *            &  !< in micromoles/(m^2*s)
     1351                                                                            conv_to_ratio(k,j,i) *                     &  !< in m^3/Nmole       
     1352                                                                             rho_air(k)                                   !< in kg/m^3   
    14551353                         ENDIF
    14561354
     
    14591357                   ELSE
    14601358
    1461                    !> If no street type is defined, then assign null emissions to all the species
     1359                      !
     1360                      !-- If no street type is defined, then assign zero emission to all the species
    14621361                      surf_lsm_h%cssws(:,m) = 0.0_wp
    14631362
     
    14651364
    14661365                ENDDO
     1366
    14671367             ENDIF
    14681368
    14691369          ENDIF
    14701370
    1471        !> For both DEFAULT and PRE-PROCESSED
     1371       !
     1372       !-- For both DEFAULT and PRE-PROCESSED mode
    14721373       ELSE   
    14731374       
    14741375
    1475           DO ispec=1,nspec_out
    1476                                 !TBD: for the PRE-PROCESSED mode in the future, values at different heights should be included!           
    1477              !> Default surface type
    1478              IF (surf_def_h(0)%ns .GT. 0) THEN
     1376          DO ispec = 1, nspec_out
     1377                   
     1378             !
     1379             !-- Default surfaces
     1380             IF ( surf_def_h(0)%ns > 0 )  THEN
    14791381
    14801382                DO  m = 1, surf_def_h(0)%ns
     
    14851387                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
    14861388
    1487 
    1488                       !> Distinguish between PMs (no needing conversion in ppms),
    1489                       !  VOC (already converted to moles/(m**2*s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
    1490                       ! other Species (still expressed in Kg/(m**2*s) at this point)
    1491 
    1492                       !> PMs
    1493                       IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
    1494                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
     1389                      !
     1390                      !-- PMs
     1391                      IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1"  &
     1392                           .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM25"  &
     1393                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM10" )  THEN
    14951394                   
    1496                          !            kg/(m^2*s) *kg/m^3                         kg/(m^2*s)                 
    1497                          surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*   &
    1498                          !                                                  kg/m^3
    1499                                                                           rho_air(nzb)  
     1395                         !
     1396                         !--            kg/(m^2*s) *kg/m^3                         kg/(m^2*s)                 
     1397                         surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*  &
     1398                                                                           rho_air(nzb)                           !< in kg/m^3 
    15001399 
    15011400 
    15021401                      ELSE
    15031402
    1504                          !> VOCs
    1505                          IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
    1506                             !          ( ppm/s) * m * kg/m^3                         mole/(m^2/s)   
    1507                             surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
    1508                                                                             !    m^3/mole          ppm
    1509                                                                              conv_to_ratio(nzb,j,i)*ratio2ppm *      &
    1510                          !                                                    kg/m^3
    1511                                                                              rho_air(nzb)   
    1512 
    1513 
    1514                          !> OTHER SPECIES
     1403                         !
     1404                         !-- VOCs
     1405                         IF ( len_index_voc > 0 .AND. emt_att%species_name(match_spec_input(ispec)) == "VOC" )  THEN
     1406                            !
     1407                            !--           (ppm/s) * m * kg/m^3                        mole/(m^2/s)   
     1408                            surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &
     1409                                                                              conv_to_ratio(nzb,j,i) *         &  !< in m^3/mole
     1410                                                                               ratio2ppm *                     &  !< in ppm
     1411                                                                                rho_air(nzb)                      !< in kg/m^3 
     1412
     1413
     1414                         !
     1415                         !-- Other species
    15151416                         ELSE
    15161417
    1517                             !               ( ppm/s )*m  * kg/m^3                      kg/(m^2/s)                     
    1518                             surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
    1519                                                                                !  mole/Kg       
    1520                                                                              (1/emt_att%xm(ispec))*                &
    1521                             !                                                    m^3/mole          ppm
    1522                                                                              conv_to_ratio(nzb,j,i)*ratio2ppm*       &
    1523                             !                                                  kg/m^3
    1524                                                                              rho_air(nzb)   
     1418                            !
     1419                            !--           (ppm/s) * m  * kg/m^3                    kg/(m^2/s)
     1420                            surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &       
     1421                                                                              ( 1.0_wp / emt_att%xm(ispec) ) * &  !< in mole/kg
     1422                                                                                conv_to_ratio(nzb,j,i) *       &  !< in m^3/mole 
     1423                                                                                 ratio2ppm *                   &  !< in ppm
     1424                                                                                  rho_air(nzb)                    !< in  kg/m^3 
    15251425 
    15261426
     
    15331433                ENDDO
    15341434
    1535              END IF
     1435             ENDIF
    15361436         
    1537              !-- Land Surface Mode surface type
    1538              IF (surf_lsm_h%ns .GT. 0) THEN
     1437             !
     1438             !-- LSM surfaces
     1439             IF ( surf_lsm_h%ns > 0 )  THEN
    15391440
    15401441                DO m = 1, surf_lsm_h%ns
     
    15461447                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
    15471448
    1548                       !> Distinguish between PMs (no needing conversion in ppms),
    1549                       !  VOC (already converted to moles/(m**2*s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
    1550                       ! other Species (still expressed in Kg/(m**2*s) at this point)
    1551 
    1552                       !> PMs
    1553                       IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
    1554                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
     1449                      !
     1450                      !-- PMs
     1451                      IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1"  &
     1452                           .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM25"  &
     1453                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM10" )  THEN
    15551454   
    1556                          !         kg/(m^2*s) * kg/m^3                           kg/(m^2*s)           
    1557                          surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *              &
    1558                          !                                                  kg/m^3
    1559                                                                           rho_air(k)   
     1455                         !
     1456                         !--         kg/(m^2*s) * kg/m^3                           kg/(m^2*s)           
     1457                         surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &
     1458                                                                        rho_air(k)                                !< in kg/m^3
    15601459 
    15611460                      ELSE
    15621461
    1563                          !> VOCs
    1564                          IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
    1565                             !          ( ppm/s) * m * kg/m^3                        mole/(m^2/s)   
    1566                             surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
    1567                                                                           !    m^3/mole          ppm
    1568                                                                           conv_to_ratio(k,j,i)*ratio2ppm*    &
    1569                          !                                                 kg/m^3
    1570                                                                           rho_air(k)   
    1571 
    1572                          !> OTHER SPECIES
     1462                         !
     1463                         !-- VOCs
     1464                         IF ( len_index_voc > 0 .AND. emt_att%species_name(match_spec_input(ispec)) == "VOC" )  THEN
     1465                            !
     1466                            !--          (ppm/s) * m * kg/m^3                        mole/(m^2/s)   
     1467                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &     
     1468                                                                           conv_to_ratio(k,j,i) *           &     !< in m^3/mole
     1469                                                                            ratio2ppm *                     &     !< in ppm
     1470                                                                             rho_air(k)                           !< in kg/m^3 
     1471
     1472                         !
     1473                         !-- Other species
    15731474                         ELSE
    1574    
    1575                             !         ( ppm/s) * m * kg/m^3                        kg/(m^2/s)                     
    1576                             surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *               &
    1577                                                                                !  mole/Kg   
    1578                                                                          (1/emt_att%xm(ispec))*                          &
    1579                             !                                                m^3/mole           ppm
    1580                                                                          conv_to_ratio(k,j,i)*ratio2ppm*                 &
    1581                             !                                            kg/m^3
    1582                                                                          rho_air(k)   
     1475                            !
     1476                            !--           (ppm/s) * m  * kg/m^3                    kg/(m^2/s)
     1477                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &
     1478                                                                           ( 1.0_wp / emt_att%xm(ispec) ) * &     !< in mole/kg
     1479                                                                            conv_to_ratio(k,j,i) *          &     !< in m^3/mole
     1480                                                                             ratio2ppm *                    &     !< in ppm
     1481                                                                              rho_air(k)                          !< in kg/m^3     
    15831482                                                     
    15841483                         ENDIF
     
    15901489                ENDDO
    15911490
    1592              END IF
    1593 
    1594              !-- Urban Surface Mode surface type
    1595              IF (surf_usm_h%ns .GT. 0) THEN
     1491             ENDIF
     1492
     1493             !
     1494             !-- USM surfaces
     1495             IF ( surf_usm_h%ns > 0 )  THEN
    15961496
    15971497
     
    16041504                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
    16051505
    1606                       !> PMs
    1607                       IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
    1608                              .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
     1506                      !
     1507                      !-- PMs
     1508                      IF ( TRIM( spc_names(match_spec_model(ispec)) ) == "PM1" &
     1509                           .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM25"  &
     1510                             .OR. TRIM( spc_names(match_spec_model(ispec)) ) == "PM10" )  THEN
    16091511                     
    1610                          !          kg/(m^2*s) *kg/m^3                             kg/(m^2*s)                     
    1611                          surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*        &
    1612                          !                                              kg/m^3
    1613                                                                        rho_air(k)   
     1512                         !
     1513                         !--          kg/(m^2*s) *kg/m^3                             kg/(m^2*s)                   
     1514                         surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*  &
     1515                                                                        rho_air(k)                                !< in kg/m^3
    16141516
    16151517 
    16161518                      ELSE
    1617 
    1618                          !> VOCs
    1619                          IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
    1620                             !          ( ppm/s) * m * kg/m^3                        mole/(m^2/s)   
    1621                             surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *   &
    1622                                                                           !    m^3/mole          ppm
    1623                                                                           conv_to_ratio(k,j,i)*ratio2ppm*    &
    1624                          !                                                kg/m^3
    1625                                                                           rho_air(k)   
    1626 
    1627                          !> OTHER SPECIES
     1519                         
     1520                         !
     1521                         !-- VOCs
     1522                         IF ( len_index_voc > 0 .AND. emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN
     1523                            !
     1524                            !--          (ppm/s) * m * kg/m^3                        mole/(m^2/s)   
     1525                            surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &
     1526                                                                           conv_to_ratio(k,j,i) *           &     !< in m^3/mole
     1527                                                                            ratio2ppm *                     &     !< in ppm
     1528                                                                             rho_air(k)                           !< in kg/m^3   
     1529
     1530                         !
     1531                         !-- Other species
    16281532                         ELSE
    16291533
    1630 
    1631                          !            ( ppm/s) * m * kg/m^3                        kg/(m^2/s)                     
    1632                             surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
    1633                                                                              !  mole/Kg   
    1634                                                                          (1/emt_att%xm(ispec))*                 &
    1635                             !                                                m^3/mole           ppm
    1636                                                                          conv_to_ratio(k,j,i)*ratio2ppm*        &
    1637                             !                                            kg/m^3
    1638                                                                          rho_air(k)   
     1534                            !
     1535                            !--          (ppm/s) * m * kg/m^3                        kg/(m^2/s)                     
     1536                            surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *  &   
     1537                                                                           ( 1.0_wp / emt_att%xm(ispec) ) * &     !< in mole/kg
     1538                                                                            conv_to_ratio(k,j,i) *          &     !< in m^3/mole
     1539                                                                             ratio2ppm*                     &     !< in ppm
     1540                                                                              rho_air(k)                          !< in kg/m^3   
    16391541
    16401542
     
    16461548 
    16471549                ENDDO
    1648              END IF
     1550
     1551             ENDIF
     1552
    16491553          ENDDO
    16501554
    1651        ENDIF
    1652 
    1653 
    1654     !> At the end of each call to chem_emissions setup, the time_factor is deallocated (ALLOCATED ONLY in the DEFAULT mode)
    1655 
     1555       ENDIF
     1556
     1557       !
     1558       !-- Deallocate time_factor in case of DEFAULT mode)
    16561559       IF ( ALLOCATED ( time_factor ) ) DEALLOCATE( time_factor )
    16571560
    1658    ENDIF !> emis_output_required
    1659 
     1561   ENDIF
    16601562
    16611563 END SUBROUTINE chem_emissions_setup
Note: See TracChangeset for help on using the changeset viewer.