source: palm/trunk/SOURCE/chem_emissions_mod.f90 @ 3582

Last change on this file since 3582 was 3582, checked in by suehring, 3 years ago

Merge branch salsa with trunk

  • Property svn:keywords set to Id
File size: 77.9 KB
Line 
1!> @file chem_emissions_mod.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 2018-2018 Leibniz Universitaet Hannover
18! Copyright 2018-2018 Freie Universitaet Berlin
19! Copyright 2018-2018 Karlsruhe Institute of Technology
20!--------------------------------------------------------------------------------!
21!
22! Current revisions:
23! ------------------
24! - Removed salsa dependency.
25! - Enabled PARAMETRIZED mode for default surfaces when LSM is not applied but
26!   salsa is (M. Kurppa)
27!
28! Former revisions:
29! -----------------
30! $Id: chem_emissions_mod.f90 3582 2018-11-29 19:16:36Z suehring $
31! resler:
32! Break lines at 132 characters
33!
34! 3483 2018-11-02 14:19:26Z raasch
35! bugfix: wrong locations of netCDF directives fixed
36!
37! 3467 2018-10-30 19:05:21Z suehring
38! Enabled PARAMETRIZED mode for default surfaces when LSM is not applied but
39! salsa is used
40!
41! 3458 2018-10-30 14:51:23Z kanani
42! from chemistry branch r3443, banzhafs, Russo:
43! Additional correction for index of input file of pre-processed mode
44! Removed atomic and molecular weights as now available in chem_modules and
45! added par_emis_time_factor (formerly in netcdf_data_input_mod)
46! Added correction for index of input file of pre-processed mode
47! Added a condition for default mode necessary for information read from ncdf file
48! in pre-processed and default mode
49! Correction of multiplying factor necessary for scaling emission values in time
50! Introduced correction for scaling units in the case of DEFAULT emission mode
51!
52! 3373 2018-10-18 15:25:56Z kanani
53! Fix wrong location of __netcdf directive
54!
55! 3337 2018-10-12 15:17:09Z kanani
56! (from branch resler)
57! Formatting
58!
59! 3312 2018-10-06 14:15:46Z knoop
60! Initial revision
61!
62! 3286 2018-09-28 07:41:39Z forkel
63!
64! Authors:
65! --------
66! @author Emmanuele Russo (Fu-Berlin)
67! @author Sabine Banzhaf  (FU-Berlin)
68! @author Martijn Schaap  (FU-Berlin, TNO Utrecht)
69!
70! Description:
71! ------------
72!>  MODULE for reading-in Chemistry Emissions data
73!>
74!> @todo Check_parameters routine should be developed: for now it includes just one condition
75!> @todo Use of Restart files not contempled at the moment
76!> @todo revise indices of files read from the netcdf: preproc_emission_data and expert_emission_data
77!> @todo for now emission data may be passed on a singular vertical level: need to be more flexible
78!> @note <Enter notes on the module>
79!> @bug  <Enter known bugs here>
80!------------------------------------------------------------------------------!
81
82 MODULE chem_emissions_mod
83
84    USE arrays_3d,                                                             &
85       ONLY: rho_air
86
87    USE control_parameters,                                                    &
88       ONLY:  initializing_actions, end_time, message_string,                  &
89              intermediate_timestep_count, dt_3d
90 
91    USE indices
92
93    USE kinds
94
95#if defined ( __netcdf )
96    USE NETCDF
97#endif
98
99    USE netcdf_data_input_mod,                                                  &
100       ONLY: chem_emis_att_type, chem_emis_val_type
101
102    USE date_and_time_mod,                                                      &
103       ONLY: time_default_indices, time_preprocessed_indices,                  &
104             month_of_year, day_of_month, hour_of_day,                          &
105             index_mm, index_dd, index_hh
106
107    USE chem_gasphase_mod,                                                      &
108       ONLY: spc_names
109 
110    USE chem_modules
111
112    USE statistics,                                                             &   
113           ONLY:  weight_pres
114
115    IMPLICIT NONE
116
117!-- Declare all global variables within the module
118   
119    CHARACTER (LEN=80)                               :: filename_emis                   !> Variable for the name of the netcdf input file
120
121    INTEGER(iwp)                                     :: i                               !> index 1st selected dimension (some dims are not spatial)
122    INTEGER(iwp)                                     :: j                               !> index 2nd selected dimension
123    INTEGER(iwp)                                     :: i_start                         !> Index to start read variable from netcdf along one dims
124    INTEGER(iwp)                                     :: i_end                           !> Index to end read variable from netcdf in one dims
125    INTEGER(iwp)                                     :: j_start                         !> Index to start read variable from netcdf in additional dims
126    INTEGER(iwp)                                     :: j_end                           !> Index to end read variable from netcdf in additional dims
127    INTEGER(iwp)                                     :: z_start                         !> Index to start read variable from netcdf in additional dims
128    INTEGER(iwp)                                     :: z_end                           !> Index to end read variable from netcdf in additional     dims
129    INTEGER(iwp)                                     :: dt_emis                         !> Time Step Emissions
130    INTEGER(iwp)                                     :: len_index                       !> length of index (used for several indices)
131    INTEGER(iwp)                                     :: len_index_voc                   !> length of voc index
132    INTEGER(iwp)                                     :: len_index_pm                    !> length of PMs index
133    REAL(wp)                                         :: con_factor                      !> Units Conversion Factor
134 
135
136    ! ---------------------------------------------------------------
137    ! Set Values of molecules, mols, etc
138    ! ---------------------------------------------------------------
139 
140    !> Avogadro number
141    REAL, PARAMETER        ::  Avog = 6.02205e23    ! mlc/mol
142 
143    !> Dobson units:
144    REAL, PARAMETER        ::  Dobs = 2.68668e16    ! (mlc/cm2) / DU
145
146    !> sesalt composition:
147    ! (Seinfeld and Pandis, "Atmospheric Chemistry and Physics",
148    !  table 7.8 "Composition of Sea-Salt", p. 444)
149    REAL, PARAMETER        ::  massfrac_Cl_in_seasalt  = 0.5504       ! (kg Cl )/(kg seasalt)
150    REAL, PARAMETER        ::  massfrac_Na_in_seasalt  = 0.3061       ! (kg Na )/(kg seasalt)
151    REAL, PARAMETER        ::  massfrac_SO4_in_seasalt = 0.0768       ! (kg SO4)/(kg seasalt)
152 
153    !> other numbers
154    REAL, PARAMETER        ::  xm_seasalt =  74.947e-3                ! kg/mol : NaCl, SO4, ..
155    REAL, PARAMETER        ::  rho_seasalt = 2.2e3                    ! kg/m3
156
157    !> * amonium sulphate
158 
159    REAL, PARAMETER        ::  xm_NH4HSO4  =  xm_NH4 + xm_H + xm_SO4  ! kg/mol
160    REAL, PARAMETER        ::  rho_NH4HSO4a = 1.8e3                   ! kg/m3
161
162    ! ---------------------------------------------------------------
163    ! gas
164    ! ---------------------------------------------------------------
165 
166    !> gas constant                       
167    REAL, PARAMETER        ::  Rgas = 8.3144     ! J/mol/K
168 
169    !> gas constant for dry air
170    REAL, PARAMETER        ::  Rgas_air = Rgas / xm_air   ! J/kg/K
171 
172    !> water vapour
173    REAL, PARAMETER        ::  Rgas_h2o = Rgas / xm_h2o   ! J/kg/K
174 
175    !> standard pressure
176    REAL, PARAMETER        ::  p0 = 1.0e5    ! Pa
177
178    !> sea-level pressure
179    REAL, PARAMETER        ::  p_sealevel = 1.01325e5    ! Pa  <-- suggestion Bram Bregman
180
181    !> global mean pressure
182    REAL, PARAMETER        ::  p_global = 98500.0   ! Pa
183
184    !> specific heat of dry air at constant pressure
185    REAL, PARAMETER        ::  cp_air = 1004.0           ! J/kg/K
186
187    !> Latent heat of evaporation
188    REAL, PARAMETER        ::  lvap = 2.5e6     !  [J kg-1]
189
190    !> Latent heat of condensation at 0 deg Celcius
191    !  (heat (J) necesarry to evaporate 1 kg of water)
192    REAL, PARAMETER        ::  Lc = 22.6e5           ! J/kg
193 
194    !> kappa = R/cp = 2/7
195    REAL, PARAMETER        ::  kappa = 2.0/7.0
196
197    !> Von Karman constant (dry_dep)
198    REAL, PARAMETER        ::  vkarman = 0.4
199
200    !> Boltzmann constant:
201    REAL(wp), PARAMETER     ::  kbolz = 1.38066e-23_wp    ! J/K
202
203    !> Inverse Reference Pressure (1/Pa)   
204    REAL(wp), PARAMETER     ::  pref_i  = 1.0_wp / 100000.0_wp       
205 
206    !>
207    REAL(wp), PARAMETER     ::  r_cp    = 0.286_wp                  !< R / cp (exponent for potential temperature)
208
209
210    SAVE
211
212!-- Interfaces Part
213
214!-- Checks Input parameters
215    INTERFACE chem_emissions_check_parameters
216       MODULE PROCEDURE chem_emissions_check_parameters
217    END INTERFACE chem_emissions_check_parameters
218
219!-- Reading of NAMELIST parameters
220!    INTERFACE chem_emissions_parin
221!       MODULE PROCEDURE chem_emissions_parin
222!    END INTERFACE chem_emissions_parin
223
224!-- Output of information to the header file
225!    INTERFACE chem_emissions_header
226!       MODULE PROCEDURE chem_emissions_header
227!    END INTERFACE chem_emissions_header
228
229!-- Matching Emissions actions 
230    INTERFACE chem_emissions_match
231       MODULE PROCEDURE chem_emissions_match
232    END INTERFACE chem_emissions_match
233
234!-- Initialization actions 
235    INTERFACE chem_emissions_init
236       MODULE PROCEDURE chem_emissions_init
237    END INTERFACE chem_emissions_init
238
239!-- Setup of Emissions
240    INTERFACE chem_emissions_setup
241       MODULE PROCEDURE chem_emissions_setup
242    END INTERFACE chem_emissions_setup
243
244
245!-- Public Interfaces
246    PUBLIC chem_emissions_init,chem_emissions_match,chem_emissions_setup
247
248!-- Public Variables
249
250    PUBLIC con_factor, len_index,len_index_voc,len_index_pm
251
252 CONTAINS
253
254!------------------------------------------------------------------------------!
255! Description:
256! ------------
257!> Routine for checking input parameters
258!------------------------------------------------------------------------------!
259 SUBROUTINE chem_emissions_check_parameters
260
261
262    !TBD: Where Should we put the call to chem_emissions_check_parameters? In chem_init or in check_parameters?
263
264    IMPLICIT NONE
265
266    INTEGER(iwp)           ::  tmp
267
268    TYPE(chem_emis_att_type)   ::  emt
269!
270!--    Call routine for controlling chem_emissions namelist
271!    CALL chem_emissions_parin
272
273!TBD: In case the given emission values do not coincide with the passed names we should think of a solution:
274!  I would personally do that if the name of the species differ from the number of emission values, I would
275!  print a warning that says that in correspondance of that particular species the value is zero.
276!  An alternative would be to initialize the cs_surface_flux values to a negative number.
277
278!-- Check Emission Species Number equal to number of passed names for the chemistry species:
279   IF ( SIZE(emt%species_name) /= emt%nspec  )  THEN
280
281       message_string = 'Numbers of input emission species names and number of species'             //          &
282                         'for which emission values are given do not match'                 
283       CALL message( 'chem_emissions_check_parameters', 'CM0437', 2, 2, 0, 6, 0 ) 
284
285
286    ENDIF
287
288    !-- Check Emission Species Number equals to number of passed names for the species
289    !IF ( SIZE(emt%species_name) /= SIZE(emt%species_index)  )  THEN
290    !      message_string = 'Number of input emission species names and '             //          &
291    !                       ' number of passed species indices do not match'                 
292    !      CALL message( 'chem_emissions_check_parameters', 'CM0101', 2, 2, 0, 6, 0 )
293
294    !ENDIF
295
296    !-- Check Emission Categories
297!    IF ( SIZE(chem_emis%cat_name) /= SIZE(chem_emis%cat_index)  )  THEN
298!       WRITE( message_string, * )                                                        &
299!       'Number of input emission categories name and categories index does not match\\' 
300!       CALL message( 'chem_emissions_check_parameters', 'CM0101', 1, 2, 0, 6, 0 )
301!    ENDIF
302
303
304
305    !TBD: Check which other consistency controls should be included
306
307   !TBD: Include check for spatial dimension of netcdf file. If they do not match with the ones
308   !     of the simulation, the model should print an error. 
309
310END SUBROUTINE chem_emissions_check_parameters
311
312!------------------------------------------------------------------------------!
313! Description:
314! ------------
315!> Matching the chemical species indices. The routine checks what are the indices of the emission input species
316!  and the corresponding ones of the model species. The routine gives as output a vector containing the number
317!  of common species: it is important to note that while the model species are distinct, their values could be
318!  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.
319!------------------------------------------------------------------------------!
320SUBROUTINE chem_emissions_match(emt_att,len_index)   
321
322
323    INTEGER(iwp), INTENT(INOUT)              ::  len_index        !< Variable where to store the number of common species between the input dataset and the model species   
324
325    TYPE(chem_emis_att_type), INTENT(INOUT)  ::  emt_att          !< Chemistry Emission Array containing information for all the input chemical emission species
326   
327    INTEGER(iwp)                             ::  ind_mod, ind_inp !< Parameters for cycling through chemical model and input species
328    INTEGER(iwp)                             ::  nspec_emis_inp   !< Variable where to store the number of the emission species in input
329
330    INTEGER(iwp)                             ::  ind_voc          !< Indices to check whether a split for voc should be done
331
332    INTEGER(iwp)                             ::  ispec            !> index for cycle over effective number of emission species
333
334
335    !> Tell the user where we are!!
336    CALL location_message( 'Matching input emissions and model chemistry species', .FALSE. )
337
338    !> Number of input emission species.
339
340    nspec_emis_inp=emt_att%nspec
341
342    !> Check the emission mode: DEFAULT, PRE-PROCESSED or PARAMETERIZED  !TBD: Add option for capital or small letters
343    SELECT CASE(TRIM(mode_emis)) 
344
345       !> 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.
346       CASE ("PRE-PROCESSED")
347
348          CALL location_message( 'chem_emissions PRE-PROCESSED_mode_matching', .FALSE. )
349
350          len_index=0
351          len_index_voc=0 
352
353          !> The first condition is that both the number of model and input emissions species are not null
354          IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0)) THEN
355
356             !> Cycle over model species
357             DO ind_mod = 1, SIZE(spc_names)
358                !> Cycle over Input species 
359                DO ind_inp = 1, nspec_emis_inp
360
361           !> 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
362                   ! > Check for VOC Species: In input in this case we only have one species (VOC) 
363                   IF (TRIM(emt_att%species_name(ind_inp))=="VOC") THEN
364                      !> 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
365                      DO ind_voc= 1,emt_att%nvoc
366                         !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
367                         IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
368                            len_index=len_index+1
369                            len_index_voc=len_index_voc+1
370                         ENDIF
371                      END DO
372                   ENDIF
373                   !> Other Species
374                   IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
375                      len_index=len_index+1
376                   ENDIF
377                ENDDO
378             ENDDO
379
380             !> Allocate array for storing the indices of the matched species
381             IF (len_index>0) THEN
382 
383                ALLOCATE(match_spec_input(len_index)) 
384 
385                ALLOCATE(match_spec_model(len_index))
386
387                IF (len_index_voc>0) THEN
388
389                   ALLOCATE(match_spec_voc_model(len_index_voc))   !> Contains indices of the VOC model species
390
391                   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
392
393                ENDIF
394
395                !> Repeat the same cycle of above but passing the species indices to the newly declared arrays
396                len_index=0
397
398                !> Cycle over model species
399                DO ind_mod = 1, SIZE(spc_names) 
400                   !> Cycle over Input species 
401                   DO ind_inp = 1, nspec_emis_inp
402
403                   !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
404
405                      ! > VOCs
406                      IF ( TRIM(emt_att%species_name(ind_inp))=="VOC" .AND. ALLOCATED(match_spec_voc_input) ) THEN     
407                         DO ind_voc= 1,emt_att%nvoc
408                            IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
409                               len_index=len_index+1
410                               len_index_voc=len_index_voc+1
411                       
412                               match_spec_input(len_index)=ind_inp
413                               match_spec_model(len_index)=ind_mod
414
415                               match_spec_voc_input(len_index_voc)=ind_voc
416                               match_spec_voc_model(len_index_voc)=ind_mod                         
417                            ENDIF
418                         END DO
419                      ENDIF
420
421                      !>Other Species
422                      IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
423                         len_index=len_index+1
424                         match_spec_input(len_index)=ind_inp
425                         match_spec_model(len_index)=ind_mod
426                      ENDIF
427                   END DO
428                END DO
429
430             !> In the case there are no species matching, the emission module should not be called
431             ELSE
432
433                message_string = 'Given Chemistry Emission Species'            //          &
434                                 ' DO NOT MATCH'                                //          &
435                                 ' model chemical species:'                      //          &
436                                 ' Chemistry Emissions routine is not called' 
437                CALL message( 'chem_emissions_matching', 'CM0438', 0, 0, 0, 6, 0 )
438             ENDIF !> IF (len_index>0)
439 
440          ELSE
441
442             !In this case either spc_names is null or nspec_emis_inp is not allocated
443
444             message_string = 'Array of Emission species not allocated:'             //          &
445                              ' Either no emission species are provided as input or'  //          &
446                              ' no chemical species are used by PALM:'                //          &
447                              ' Chemistry Emissions routine is not called'                 
448             CALL message( 'chem_emissions_matching', 'CM0439', 0, 2, 0, 6, 0 ) 
449
450          ENDIF !> IF ( (SIZE(spc_names) .GT. 0) .AND. ALLOCATED(nspec_emis_inp))
451
452       !> ------------------------------------------------------------------
453       !> DEFAULT case
454
455       CASE ("DEFAULT")
456
457          CALL location_message( 'chem_emissions DEFAULT_mode_matching', .FALSE. )
458
459          len_index=0      !>  index for TOTAL number of species   
460          len_index_voc=0  !>  index for TOTAL number of VOCs
461          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.
462
463          !> The first condition is that both the number of model and input emissions species are not null
464          IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0) ) THEN
465
466             !> Cycle over model species
467             DO ind_mod = 1, SIZE(spc_names)
468                !> Cycle over input species
469                DO ind_inp = 1, nspec_emis_inp
470
471                   ! > Check for VOC Species: In input in this case we only have one species (VOC) 
472                   IF (TRIM(emt_att%species_name(ind_inp))=="VOC") THEN
473                      !> 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
474                      DO ind_voc= 1,emt_att%nvoc
475                         !> Determine the indices of the coinciding species in the two cases and assign them to matching arrays
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
479                         ENDIF
480                      END DO
481                   ENDIF
482
483                   !> 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.
484                   IF (TRIM(emt_att%species_name(ind_inp))=="PM") THEN
485                      !>pm1
486                      IF (TRIM(spc_names(ind_mod))=="PM1") THEN
487                         len_index=len_index+1
488                      !>pm2.5
489                      ELSE IF (TRIM(spc_names(ind_mod))=="PM25") THEN
490                         len_index=len_index+1
491                      !>pm10
492                      ELSE IF (TRIM(spc_names(ind_mod))=="PM10") THEN
493                         len_index=len_index+1   
494                      ENDIF
495                   ENDIF
496
497                   !> 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   
498                   IF (TRIM(emt_att%species_name(ind_inp))=="NOx") THEN 
499                      !>no
500                      IF (TRIM(spc_names(ind_mod))=="NO") THEN
501                         len_index=len_index+1
502                      !>no2
503                      ELSE IF (TRIM(spc_names(ind_mod))=="NO2") THEN
504                         len_index=len_index+1                       
505                      ENDIF
506                   ENDIF
507
508                   !>SOX: same as for NOx, but with SO2 and SO4
509                   IF (TRIM(emt_att%species_name(ind_inp))=="SOx") THEN
510                      !>no
511                      IF (TRIM(spc_names(ind_mod))=="SO2") THEN
512                         len_index=len_index+1
513                      !>no2
514                      ELSE IF (TRIM(spc_names(ind_mod))=="SO4") THEN
515                         len_index=len_index+1                       
516                      ENDIF
517                   ENDIF
518
519                   !> Other Species
520                   IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
521                      len_index=len_index+1
522                   ENDIF
523                END DO   !> number of emission input species
524             END DO   !> number of model species
525
526
527             !> Allocate Arrays
528             IF (len_index>0) THEN
529
530                ALLOCATE(match_spec_input(len_index))
531                ALLOCATE(match_spec_model(len_index))
532
533                IF (len_index_voc>0) THEN
534                   ALLOCATE(match_spec_voc_model(len_index_voc))   !> contains indices of the VOC model species
535                   ALLOCATE(match_spec_voc_input(len_index_voc))   !> In input there is only one value for VOCs in the DEFAULT mode.
536                                                                   !  This array contains the indices of the different values of VOC compositions of the input variable VOC_composition
537                ENDIF
538
539                !> ASSIGN INDICES
540                !> Repeat the same cycles of above, but passing the species indices to the newly declared arrays
541                len_index=0
542                len_index_voc=0
543               
544                DO ind_mod = 1, SIZE(spc_names)
545                   DO ind_inp = 1, nspec_emis_inp
546
547                      ! > VOCs
548                      IF ( TRIM(emt_att%species_name(ind_inp))=="VOC" .AND. ALLOCATED(match_spec_voc_input) ) THEN     
549                         DO ind_voc= 1,emt_att%nvoc
550                            IF (TRIM(emt_att%voc_name(ind_voc))==TRIM(spc_names(ind_mod))) THEN
551                               len_index=len_index+1
552                               len_index_voc=len_index_voc+1
553                       
554                               match_spec_input(len_index)=ind_inp
555                               match_spec_model(len_index)=ind_mod
556
557                               match_spec_voc_input(len_index_voc)=ind_voc
558                               match_spec_voc_model(len_index_voc)=ind_mod                         
559                            ENDIF
560                         END DO
561                      ENDIF
562
563                      !> 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.
564                      IF (TRIM(emt_att%species_name(ind_inp))=="PM") THEN
565                         !>pm1
566                         IF (TRIM(spc_names(ind_mod))=="PM1") THEN
567                            len_index=len_index+1
568
569                            match_spec_input(len_index)=ind_inp
570                            match_spec_model(len_index)=ind_mod
571 
572                            !match_spec_pm(1)=ind_mod
573                         !>pm2.5
574                         ELSE IF (TRIM(spc_names(ind_mod))=="PM25") THEN
575                            len_index=len_index+1
576
577                            match_spec_input(len_index)=ind_inp
578                            match_spec_model(len_index)=ind_mod
579 
580                            !match_spec_pm(2)=ind_mod
581                         !>pm10
582                         ELSE IF (TRIM(spc_names(ind_mod))=="PM10") THEN
583                            len_index=len_index+1
584   
585                            match_spec_input(len_index)=ind_inp
586                            match_spec_model(len_index)=ind_mod
587 
588                            !match_spec_pm(3)=ind_mod
589                         ENDIF
590                      ENDIF
591
592                      !> NOx - The same as for PMs but here the model species are only 2: NO and NO2
593                      IF (TRIM(emt_att%species_name(ind_inp))=="NOx") THEN
594                         !>no
595                         IF (TRIM(spc_names(ind_mod))=="NO") THEN
596                            len_index=len_index+1
597
598                            match_spec_input(len_index)=ind_inp
599                            match_spec_model(len_index)=ind_mod
600 
601                            !match_spec_nox(1)=ind_mod
602                         !>no2
603                         ELSE IF (TRIM(spc_names(ind_mod))=="NO2") THEN
604                            len_index=len_index+1
605
606                            match_spec_input(len_index)=ind_inp
607                            match_spec_model(len_index)=ind_mod
608 
609                           ! match_spec_nox(2)=ind_mod
610                         ENDIF
611                      ENDIF
612
613                      !> SOx
614                      IF (TRIM(emt_att%species_name(ind_inp))=="SOx") THEN
615                         !>so2
616                         IF (TRIM(spc_names(ind_mod))=="SO2") THEN
617                            len_index=len_index+1
618
619                            match_spec_input(len_index)=ind_inp
620                            match_spec_model(len_index)=ind_mod
621 
622                           ! match_spec_sox(1)=ind_mod
623                         !>so4
624                         ELSE IF (TRIM(spc_names(ind_mod))=="SO4") THEN
625                            len_index=len_index+1
626
627                            match_spec_input(len_index)=ind_inp
628                            match_spec_model(len_index)=ind_mod
629 
630                           ! match_spec_sox(2)=ind_mod
631                         ENDIF
632                      ENDIF
633
634                      !> Other Species
635                      IF (TRIM(emt_att%species_name(ind_inp))==TRIM(spc_names(ind_mod))) THEN
636                         len_index=len_index+1
637                           
638                         match_spec_input(len_index)=ind_inp
639                         match_spec_model(len_index)=ind_mod
640                      ENDIF
641                   END DO
642                END DO
643
644             ELSE
645
646                message_string = 'Given Chemistry Emission Species'            //          &
647                                 ' DO NOT MATCH'                                //          &
648                                 ' model chemical species'                      //          &
649                                 ' Chemistry Emissions routine is not called'         
650                CALL message( 'chem_emissions_matching', 'CM0440', 0, 0, 0, 6, 0 ) 
651
652             ENDIF
653
654          ELSE
655
656             message_string = 'Array of Emission species not allocated: '            //          &
657                              ' Either no emission species are provided as input or'  //          &
658                              ' no chemical species are used by PALM:'                //          &
659                              ' Chemistry Emissions routine is not called'                                   
660             CALL message( 'chem_emissions_matching', 'CM0441', 0, 2, 0, 6, 0 ) 
661 
662          ENDIF
663
664!-- 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.
665
666       CASE ("PARAMETERIZED")
667
668          len_index=0
669
670         !spc_names have to be greater than zero for proceeding
671          IF ( (SIZE(spc_names) .GT. 0) .AND. (nspec_emis_inp .GT. 0) ) THEN
672
673
674             !cycle over model species
675             DO ind_mod = 1, SIZE(spc_names)
676                ind_inp=1
677                DO  WHILE ( TRIM( surface_csflux_name( ind_inp ) ) /= 'novalue' )   !<'novalue' is the default 
678                   IF (TRIM(surface_csflux_name( ind_inp ))==TRIM(spc_names(ind_mod))) THEN
679                      len_index=len_index+1
680                   ENDIF
681                   ind_inp=ind_inp+1
682
683                ENDDO
684             ENDDO
685
686             !Proceed only if there are matched species
687             IF ( len_index .GT. 0 ) THEN
688
689                !Allocation of Arrays of the matched species
690                ALLOCATE(match_spec_input(len_index)) 
691 
692                ALLOCATE(match_spec_model(len_index))
693
694                !Assign corresponding indices of input and model matched species
695                len_index=0
696
697                DO ind_mod = 1, SIZE(spc_names) 
698                   ind_inp = 1
699                   DO  WHILE ( TRIM( surface_csflux_name( ind_inp ) ) /= 'novalue' )   !<'novalue' is the default 
700                      IF (TRIM( surface_csflux_name( ind_inp ) ) == TRIM(spc_names(ind_mod)))           THEN
701                         len_index=len_index+1
702                         match_spec_input(len_index)=ind_inp
703                         match_spec_model(len_index)=ind_mod
704                      ENDIF
705                   ind_inp=ind_inp+1
706                   END DO
707                END DO
708
709                ! also, Add AN if to check that when the surface_csflux are passed to the namelist, also the street factor values are provided
710                DO  ispec = 1 , len_index
711
712                   IF ( emiss_factor_main(match_spec_input(ispec)) .LT. 0 .AND. &
713                        emiss_factor_side(match_spec_input(ispec)) .LT. 0 ) THEN
714
715                      message_string = 'PARAMETERIZED emissions mode selected:'            //          &
716                                       ' EMISSIONS POSSIBLE ONLY ON STREET SURFACES'        //          &
717                                       ' but values of scaling factors for street types'    //          &
718                                       ' emiss_factor_main AND emiss_factor_side'           //          &
719                                       ' not provided for each of the emissions species'    //          &
720                                       ' or not provided at all: PLEASE set a finite value' //          &
721                                       ' for these parameters in the chemistry namelist'         
722                      CALL message( 'chem_emissions_matching', 'CM0442', 2, 2, 0, 6, 0 ) 
723                   ENDIF
724                END DO
725
726
727             ELSE
728               
729                message_string = 'Given Chemistry Emission Species'            //          &
730                                 ' DO NOT MATCH'                                //          &
731                                 ' model chemical species'                      //          &
732                                 ' Chemistry Emissions routine is not called'         
733                CALL message( 'chem_emissions_matching', 'CM0443', 0, 0, 0, 6, 0 ) 
734             ENDIF
735
736          ELSE
737     
738             message_string = 'Array of Emission species not allocated: '            //          &
739                              ' Either no emission species are provided as input or'  //          &
740                              ' no chemical species are used by PALM.'                //          &
741                              ' Chemistry Emissions routine is not called'                                   
742             CALL message( 'chem_emissions_matching', 'CM0444', 0, 2, 0, 6, 0 ) 
743 
744          ENDIF             
745
746
747!-- CASE when emission module is switched on but mode_emis is not specified or it is given the wrong name
748       CASE DEFAULT
749
750          message_string = 'Chemistry Emissions Module switched ON, but'            //          &
751                           ' either no emission mode specified or incorrectly given :'  //          &
752                           ' please, pass the correct value to the namelist parameter "mode_emis"'                 
753          CALL message( 'chem_emissions_matching', 'CM0445', 2, 2, 0, 6, 0 )             
754
755       END SELECT
756
757END SUBROUTINE chem_emissions_match
758
759!------------------------------------------------------------------------------!
760! Description:
761! ------------
762!> Initialization:
763!> Netcdf reading, arrays allocation and first calculation of cssws fluxes at timestep 0
764!> 
765!------------------------------------------------------------------------------!
766 SUBROUTINE chem_emissions_init(emt_att,emt,nspec_out)
767
768    USE surface_mod,                                                           &
769       ONLY:  surf_lsm_h,surf_def_h,surf_usm_h
770   
771    IMPLICIT NONE
772 
773    TYPE(chem_emis_att_type),INTENT(INOUT)                            :: emt_att   !> variable where to store all the information of
774                                                                                   !  emission inputs whose values do not depend
775                                                                                   !  on the considered species
776
777    TYPE(chem_emis_val_type),INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt      !> variable where to store emission inputs values,
778                                                                                   !  depending on the considered species
779   
780    INTEGER(iwp),INTENT(INOUT)                                        :: nspec_out !> number of outputs of chemistry emission routines
781
782    CHARACTER (LEN=80)                                                :: units     !> units of chemistry inputs
783 
784    INTEGER(iwp)                                                      :: ispec     !> Index to go through the emission chemical species
785
786
787!-- Actions for initial runs : TBD: needs to be updated
788  IF (  TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
789!--    ...
790   
791!
792!-- Actions for restart runs
793  ELSE
794!--    ...
795
796  ENDIF
797
798
799  CALL location_message( 'Starting initialization of chemistry emissions arrays', .FALSE. )
800
801
802  !-- Match Input and Model emissions
803  CALL  chem_emissions_match(emt_att,nspec_out)
804
805  !> If nspec_out==0, then do not use emission module: The corresponding message is already printed in the matching routine
806  IF ( nspec_out == 0 ) THEN
807 
808     emission_output_required=.FALSE.
809
810  ELSE
811
812     emission_output_required=.TRUE.
813
814
815!-----------------------------------------------------------------
816     !> Set molecule masses'
817     ALLOCATE(emt_att%xm(nspec_out))
818     ! loop over emisisons:
819
820     DO ispec = 1, nspec_out
821       ! switch:
822        SELECT CASE ( TRIM(spc_names(match_spec_model(ispec))) )
823           CASE ( 'SO2','SO4' ) ; emt_att%xm(ispec) = xm_S + xm_O * 2      ! kg/mole
824           CASE ( 'NO','NO2' )  ; emt_att%xm(ispec) = xm_N + xm_O * 2      ! kg/mole  NO2 equivalent
825           CASE ( 'NH3' ) ; emt_att%xm(ispec) = xm_N + xm_H * 3  ! kg/mole
826           CASE ( 'CO'  ) ; emt_att%xm(ispec) = xm_C + xm_O      ! kg/mole
827           CASE ( 'CO2' ) ; emt_att%xm(ispec) = xm_C + xm_O * 2  ! kg/mole
828           CASE ( 'CH4' ) ; emt_att%xm(ispec) = xm_C + xm_H * 4  ! kg/mole     
829           CASE ( 'HNO3' ); emt_att%xm(ispec) = xm_H + xm_N + xm_O*3 !kg/mole 
830           CASE DEFAULT
831             !--  TBD: check if this hase to be removed
832              !emt_att%xm(ispec) = 1.0_wp
833        END SELECT
834     ENDDO
835
836   
837     !> ASSIGN EMISSION VALUES for the different emission modes
838     SELECT CASE(TRIM(mode_emis))   !TBD: Add the option for CApital or small letters
839
840
841        !> PRE-PROCESSED case
842        CASE ("PRE-PROCESSED")
843
844           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE(emis_distribution(nzb:nzt+1,0:ny,0:nx,nspec_out)) 
845
846           CALL location_message( 'emis_distribution array allocated in PRE-PROCESSED mode', .FALSE. )
847 
848           !> Calculate the values of the emissions at the first time step
849           CALL chem_emissions_setup(emt_att,emt,nspec_out)
850
851        !> Default case
852        CASE ("DEFAULT")
853
854           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE( emis_distribution( 1, 0:ny, 0:nx, nspec_out ) ) 
855       
856           CALL location_message( 'emis_distribution array allocated in DEFAULT mode', .FALSE. )
857
858           !> Calculate the values of the emissions at the first time step
859           CALL chem_emissions_setup(emt_att,emt,nspec_out)
860
861        !> PARAMETERIZED case
862        CASE ("PARAMETERIZED")
863
864           CALL location_message( 'emis_distribution array allocated in PARAMETERIZED mode', .FALSE. )
865
866           ! For now for PAR and DEF values only, first vertical level of emis_distribution is allocated, while for PRE-PROCESSED all.
867           IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE(emis_distribution(1,0:ny,0:nx,nspec_out))
868
869           !> Calculate the values of the emissions at the first time step
870           CALL chem_emissions_setup(emt_att,emt,nspec_out)
871
872     END SELECT
873
874  ENDIF   
875
876 END SUBROUTINE chem_emissions_init
877
878
879
880!------------------------------------------------------------------------------!
881! Description:
882! ------------
883!> Routine for Update of Emission values at each timestep
884!-------------------------------------------------------------------------------!
885
886 SUBROUTINE chem_emissions_setup(emt_att,emt,nspec_out)
887
888   USE arrays_3d,                                                    &
889       ONLY:  dzw
890   USE grid_variables,                                                        &
891       ONLY: dx, dy
892   USE indices,                                                               &
893       ONLY: nnx,nny,nnz
894   USE surface_mod,                                                           &
895       ONLY:  surf_lsm_h,surf_def_h,surf_usm_h
896   USE netcdf_data_input_mod,                                                 &
897       ONLY:  street_type_f
898   USE arrays_3d,                                                             &       
899       ONLY: hyp, pt
900
901 IMPLICIT NONE
902 
903    !--- IN/OUT
904
905    TYPE(chem_emis_att_type),INTENT(INOUT)                            ::  emt_att    !> variable where to store all the information of
906                                                                                     !  emission inputs whose values do not depend
907                                                                                     !  on the considered species
908
909    TYPE(chem_emis_val_type),INTENT(INOUT), ALLOCATABLE, DIMENSION(:) ::  emt        !> variable where to store emission inputs values,
910                                                                                     !  depending on the considered species
911
912    INTEGER,INTENT(IN)                                                ::  nspec_out  !> Output of routine chem_emis_matching with number
913                                                                                     !  of matched species
914    !---
915
916    REAL(wp),ALLOCATABLE,DIMENSION(:,:)                               ::  delta_emis !> Term to add to the emission_distribution array
917                                                                                     !  for each of the categories at each time step
918                                                                                     !  in the default mode
919
920    CHARACTER(LEN=80)                                                 ::  units      !> Units of the emissions
921
922    INTEGER(iwp)                                                      ::  icat       !> Index for number of categories
923    INTEGER(iwp)                                                      ::  ispec      !> index for number of species
924    INTEGER(iwp)                                                      ::  i_pm_comp  !> index for number of PM components
925    INTEGER(iwp)                                                      ::  ivoc       !> Index for number of components of  VOCs
926    INTEGER(iwp)                                                      ::  time_index !> Index for time
927   
928    REAL(wp),ALLOCATABLE, DIMENSION(:)                                ::  time_factor !> factor for time scaling of emissions
929    REAL(wp),ALLOCATABLE, DIMENSION(:,:)                              ::  emis
930
931    REAL(wp), DIMENSION(24)                                           :: par_emis_time_factor      !< time factors
932                                                                                      !  for the parameterized mode: these are fixed for each hour
933                                                                                      !  of a single day.
934    REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr)                    ::  conv_to_ratio !> factor used for converting input
935                                                                                        !  to adimensional concentration ratio
936
937   ! ------------------------------------------
938    ! variables for the conversion of indices i and j according to surface_mod
939
940    INTEGER(iwp)                                                      ::  i          !> running index for grid in x-direction
941    INTEGER(iwp)                                                      ::  j          !> running index for grid in y-direction
942    INTEGER(iwp)                                                      ::  k          !> running index for grid in z-direction
943    INTEGER(iwp)                                                      ::  m          !> running index for horizontal surfaces
944
945    REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr)                    ::  tmp_temp
946
947    ! --- const -------------------------------
948    !-CONVERSION FACTORS: TIME
949    ! number of sec per hour = 3600   
950    REAL, PARAMETER   ::  s_per_hour = 3600.0  !  (s)/(hour)
951    ! number of sec per day = 86400   
952    REAL, PARAMETER   ::  s_per_day = 86400.0  !  (s)/(day)
953    ! number of hours in a year of 365 days:
954    REAL, PARAMETER   ::  hour_per_year = 8760.0 !> TBD: What about leap years?
955    ! number of hours in a day:
956    REAL, PARAMETER   ::  hour_per_day = 24.0
957
958    ! conversion from hours to seconds (s/hour) = 1/3600.0 ~ 0.2777778   
959    REAL, PARAMETER   ::  hour_to_s = 1/s_per_hour  !  (hour)/(s)
960    ! conversion from day to seconds (s/day) = 1/86400.0 ~ 1.157407e-05   
961    REAL, PARAMETER   ::  day_to_s = 1/s_per_day   !  (day)/(s)
962    ! conversion from year to sec (s/year) = 1/31536000.0 ~ 3.170979e-08   
963    REAL, PARAMETER   ::  year_to_s = 1/(s_per_hour*hour_per_year)  !  (year)/(s)
964
965    !-CONVERSION FACTORS: WEIGHT
966    !  Conversion from tons to kg (kg/tons) = 100.0/1 ~ 100
967    REAL, PARAMETER   ::  tons_to_kg = 100  !  (tons)/(kg)
968    !  Conversion from g to kg (kg/g) = 1/1000.0 ~ 0.001
969    REAL, PARAMETER   ::  g_to_kg = 0.001  !  (g)/(kg)
970    !  Conversion from g to kg (kg/g) = 1/1000.0 ~ 0.001
971    REAL, PARAMETER   ::  miug_to_kg = 0.000000001  !  (microng)/(kg)
972
973    !< CONVERSION FACTORS: fraction to ppm
974    REAL(wp), PARAMETER   ::  ratio2ppm  = 1.0e06_wp 
975    !------------------------------------------------------   
976
977    IF ( emission_output_required ) THEN
978
979    !>  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.
980
981       IF ( call_chem_at_all_substeps )  THEN
982
983          dt_emis = dt_3d * weight_pres(intermediate_timestep_count)
984
985       ELSE
986
987          dt_emis = dt_3d
988
989       ENDIF
990
991
992    ! --- CHECK UNITS
993    !>-----------------------------------------------------
994    !> Conversion of the units to the ones employed in PALM.
995    !> Possible temporal units of input emissions data are: year, hour and second;
996    !> the mass could be expressed in: tons, kilograms or grams.
997    !> IN the PARAMETERIZED mode no conversion is possible: in this case INPUTS have to have fixed units.
998    !>-----------------------------------------------------
999
1000        IF (TRIM(mode_emis)=="DEFAULT" .OR. TRIM(mode_emis)=="PRE-PROCESSED" ) THEN
1001
1002          SELECT CASE(TRIM(emt_att%units))
1003          !> kilograms
1004             CASE ( 'kg/m2/s','KG/M2/S') 
1005                CALL location_message( 'Units of the emissions are consistent with'  //          &
1006                                       ' the ones employed in the PALM-4U model ', .FALSE. )
1007                con_factor=1
1008
1009             CASE ('kg/m2/hour','KG/M2/HOUR')
1010                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1011
1012                con_factor=hour_to_s
1013
1014             CASE ( 'kg/m2/day','KG/M2/DAY' ) 
1015                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1016
1017                con_factor=day_to_s
1018
1019             CASE ( 'kg/m2/year','KG/M2/YEAR' )
1020                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1021
1022                con_factor=year_to_s
1023
1024          !> Tons
1025             CASE ( 'ton/m2/s','TON/M2/S' ) 
1026                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1027
1028                con_factor=tons_to_kg
1029
1030             CASE ( 'ton/m2/hour','TON/M2/HOUR' ) 
1031                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1032
1033                con_factor=tons_to_kg*hour_to_s
1034
1035             CASE ( 'ton/m2/year','TON/M2/YEAR' ) 
1036                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1037
1038                con_factor=tons_to_kg*year_to_s
1039
1040          !> Grams
1041             CASE ( 'g/m2/s','G/M2/S' )
1042                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1043
1044                con_factor=g_to_kg
1045
1046             CASE ( 'g/m2/hour','G/M2/HOUR' ) 
1047                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1048
1049                con_factor=g_to_kg*hour_to_s
1050
1051             CASE ( 'g/m2/year','G/M2/YEAR' )
1052                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1053
1054                con_factor=g_to_kg*year_to_s
1055
1056          !> Micro Grams
1057             CASE ( 'micrograms/m2/s','MICROGRAMS/M2/S' )
1058                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1059
1060                con_factor=miug_to_kg
1061
1062             CASE ( 'micrograms/m2/hour','MICROGRAMS/M2/HOUR' ) 
1063                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1064
1065                con_factor=miug_to_kg*hour_to_s
1066
1067             CASE ( 'micrograms/m2/year','MICROGRAMS/M2/YEAR' )
1068                CALL location_message( 'Units of emission inputs need conversion', .FALSE. )
1069
1070                con_factor=miug_to_kg*year_to_s
1071
1072             CASE DEFAULT 
1073                message_string = 'The Units of the provided input chemistry emission species' // &
1074                                 ' are not the ones required by PALM-4U: please check '      // &
1075                                 ' chemistry emission module documentation.'                                 
1076                CALL message( 'chem_emissions_setup', 'CM0446', 2, 2, 0, 6, 0 )
1077
1078          END SELECT
1079
1080         
1081       !> PRE-PROCESSED and parameterized mode
1082       ELSE
1083 
1084          message_string = 'No Units conversion required for units of chemistry emissions' // &
1085                           ' of the PARAMETERIZED mode: units have to be provided in'     //  & 
1086                           ' micromole/m**2/day for GASES and'                            //  &
1087                           ' kg/m**2/day for PMs'                     
1088          CALL message( 'chem_emissions_setup', 'CM0447', 0, 0, 0, 6, 0 ) 
1089
1090       ENDIF
1091
1092    !> Conversion factors (if the units are kg/m**2/s) we have to convert them to ppm/s
1093     DO i=nxl,nxr
1094          DO j=nys,nyn
1095          !> Derive Temperature from Potential Temperature
1096             tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * ( hyp(nzb:nzt+1) * pref_i )**r_cp
1097
1098             !> We need to pass to cssws<- (ppm/s) * dz.
1099             !  Input is Nmole/(m^2*s)
1100             !  To go to ppm*dz basically we need to multiply the input by (m**2/N)*dz
1101             !  (m**2/N)*dz == V/N
1102             !  V/N=RT/P
1103
1104             !>    m**3/Nmole              (J/mol)*K^-1           K                      Pa         
1105             conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) ) 
1106          ENDDO
1107       ENDDO
1108
1109    !>------------------------------------------------
1110    !> Start The Processing of the input data
1111
1112        emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp
1113
1114    !>-----------------------------------------------
1115    !> 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
1116 
1117    !> PRE-PROCESSED MODE
1118       IF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
1119
1120          !> Update time indices
1121          CALL time_preprocessed_indices(index_hh)
1122
1123          CALL location_message( 'PRE-PROCESSED MODE: No time-factor specification required', .FALSE. )
1124
1125       ELSEIF (TRIM(mode_emis)=="DEFAULT") THEN
1126
1127          CALL location_message( 'DEFAULT MODE: time-factor specification required', .FALSE. )
1128 
1129       !> Allocate array where to store temporary emission values     
1130          IF(.NOT. ALLOCATED(emis)) ALLOCATE(emis(nys:nyn,nxl:nxr))
1131
1132       !> Allocate time factor per emitted component category
1133          ALLOCATE(time_factor(emt_att%ncat)) 
1134
1135       !> Read-in HOURLY emission time factor
1136          IF (TRIM(time_fac_type)=="HOUR") THEN
1137
1138          !> Update time indices
1139             CALL time_default_indices(month_of_year,day_of_month,hour_of_day,index_hh)
1140
1141          !> Check if the index is less or equal to the temporal dimension of HOURLY emission files             
1142             IF (index_hh .LE. SIZE(emt_att%hourly_emis_time_factor(1,:))) THEN
1143
1144             !> Read-in the correspondant time factor             
1145                time_factor(:)= emt_att%hourly_emis_time_factor(:,index_hh)     
1146
1147             ELSE
1148
1149                message_string = 'The "HOUR" time-factors (DEFAULT mode) of the chemistry emission species'           // &
1150                              ' are not provided for each hour of the total simulation time'     
1151                CALL message( 'chem_emissions_setup', 'CM0448', 2, 2, 0, 6, 0 ) 
1152
1153             ENDIF
1154
1155       !> Read-in MDH emissions time factors
1156          ELSEIF (TRIM(time_fac_type)=="MDH") THEN
1157
1158          !> Update time indices     
1159             CALL time_default_indices(daytype_mdh,month_of_year,day_of_month,hour_of_day,index_mm,index_dd,index_hh)
1160
1161
1162          !> Check if the index is less or equal to the temporal dimension of MDH emission files             
1163             IF ((index_hh+index_dd+index_mm) .LE. SIZE(emt_att%mdh_emis_time_factor(1,:))) THEN
1164
1165             !> Read-in the correspondant time factor             
1166                time_factor(:)=emt_att%mdh_emis_time_factor(:,index_mm)*emt_att%mdh_emis_time_factor(:,index_dd)*   &
1167                                                                         emt_att%mdh_emis_time_factor(:,index_hh)
1168     
1169             ELSE
1170
1171                message_string = 'The "MDH" time-factors (DEFAULT mode) of the chemistry emission species'           // &
1172                              ' are not provided for each hour/day/month  of the total simulation time'     
1173                CALL message( 'chem_emissions_setup', 'CM0449', 2, 2, 0, 6, 0 )
1174
1175             ENDIF
1176
1177          ELSE
1178          !> condition when someone used the DEFAULT mode but forgets to indicate the time-factor in the namelist
1179             
1180             message_string = 'The time-factor (DEFAULT mode) of the chemistry emission species'           // &
1181                              ' is not provided in the NAMELIST'     
1182             CALL message( 'chem_emissions_setup', 'CM0450', 2, 2, 0, 6, 0 )
1183         
1184          ENDIF
1185
1186       !> PARAMETERIZED MODE
1187       ELSEIF (TRIM(mode_emis)=="PARAMETERIZED") THEN
1188          CALL location_message( 'PARAMETERIZED MODE: time-factor specification is fixed: '                         // &
1189                                 ' 24 values for every day of the year ', .FALSE. )
1190       
1191          !Assign Constant Values of time factors, diurnal time profile for traffic sector:
1192          par_emis_time_factor( : ) = &
1193            (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, 0.056, 0.053, 0.051, 0.051, 0.052, 0.055, &
1194               0.059, 0.061, 0.064, 0.067, 0.069, 0.069, 0.049, 0.039, 0.039, 0.029, 0.024, 0.019 /)
1195         
1196          !> in this case allocate time factor each hour in a day
1197          IF (.NOT. ALLOCATED(time_factor)) ALLOCATE(time_factor(1))
1198
1199          !>Pass the values of time factors in the parameterized mode to the time_factor variable. in this case index_hh==hour_of_day
1200          index_hh=hour_of_day
1201
1202          time_factor(1) = par_emis_time_factor(index_hh)
1203
1204       ENDIF
1205
1206!--------------------------------
1207!--  EMISSION DISTRIBUTION Calculation
1208
1209       !> PARAMETERIZED CASE
1210       IF ( mode_emis == "PARAMETERIZED" ) THEN
1211
1212          DO ispec=1,nspec_out
1213
1214             !> Values are still micromoles/(m**2*s). Units in this case are always micromoles/m**2*day (or kilograms/m**2*day for PMs)
1215             emis_distribution(1,nys:nyn,nxl:nxr,ispec)=surface_csflux(match_spec_input(ispec))*time_factor(1)*hour_to_s
1216
1217          ENDDO
1218
1219       !> PRE-PROCESSED CASE
1220       ELSEIF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
1221
1222          !> Start Cycle over Species
1223          DO ispec=1,nspec_out !> nspec_out represents the number of species in common between
1224                               !  the emission input data and the chemistry mechanism used
1225   
1226             emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emt(match_spec_input(ispec))%                               &
1227                                                                   preproc_emission_data(index_hh,1,nys+1:nyn+1,nxl+1:nxr+1)* &
1228                                                                      con_factor
1229         
1230          ENDDO
1231
1232      !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.
1233
1234       ELSE IF (TRIM(mode_emis)=="DEFAULT") THEN
1235
1236       !> Allocate array for the emission value corresponding to a specific category and time factor
1237          ALLOCATE(delta_emis(nys:nyn,nxl:nxr)) 
1238
1239
1240       !> Start Cycle over categories
1241          DO icat = 1, emt_att%ncat
1242 
1243          !> Start Cycle over Species
1244             DO ispec=1,nspec_out !> nspec_out represents the number of species in common between
1245                                  !  the emission input data and the chemistry mechanism used
1246
1247                emis(nys:nyn,nxl:nxr) = emt(match_spec_input(ispec))%default_emission_data(icat,nys+1:nyn+1,nxl+1:nxr+1)
1248
1249!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.
1250
1251             !> 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.
1252             !  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 )
1253             !                                                                                 ((value/sec)*(24*3600)*time_factor)/3600=24*(value/sec)*time_factor                         
1254
1255             !> NOX Compositions
1256                IF (TRIM(spc_names(match_spec_model(ispec)))=="NO") THEN
1257                !>             Kg/m2*s                   kg/m2*s
1258                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1259                                                 emt_att%nox_comp(icat,1)*con_factor*hour_per_day
1260
1261                   emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
1262
1263                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="NO2") THEN
1264
1265                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1266                                                 emt_att%nox_comp(icat,2)*con_factor*hour_per_day
1267
1268                   emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
1269 
1270             !> SOX Compositions
1271
1272                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO2") THEN
1273                   !>             Kg/m2*s                   kg/m2*s
1274                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1275                                                 emt_att%sox_comp(icat,1)*con_factor*hour_per_day
1276
1277                   emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
1278
1279                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO4") THEN
1280                   !>             Kg/m2*s                   kg/m2*s
1281                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1282                                                 emt_att%sox_comp(icat,2)*con_factor*hour_per_day
1283
1284                   emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
1285 
1286
1287             !> PMs should be in kg/m**2/s, so conversion factor is here still required
1288             !> PM1 Compositions
1289                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1") THEN
1290
1291                !> Cycle over the different pm components for PM1 type
1292                   DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,1))
1293
1294                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1295                                                    emt_att%pm_comp(icat,i_pm_comp,1)*con_factor*hour_per_day
1296
1297                      emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
1298                                                                 delta_emis(nys:nyn,nxl:nxr)
1299 
1300                   ENDDO
1301
1302             !> PM2.5 Compositions
1303                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM25") THEN
1304
1305                !> Cycle over the different pm components for PM2.5 type
1306                   DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,2))
1307
1308                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1309                                                    emt_att%pm_comp(icat,i_pm_comp,2)*con_factor*hour_per_day
1310
1311                      emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
1312                                                                 delta_emis(nys:nyn,nxl:nxr)
1313 
1314                   ENDDO
1315
1316             !> PM10 Compositions
1317                ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1318
1319                !> Cycle over the different pm components for PM10 type
1320                   DO i_pm_comp= 1,SIZE(emt_att%pm_comp(1,:,3)) 
1321
1322                      delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)* &
1323                                                    emt_att%pm_comp(icat,i_pm_comp,3)*con_factor*hour_per_day
1324
1325                      emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
1326                                                                 delta_emis(nys:nyn,nxl:nxr) 
1327
1328                   ENDDO
1329
1330             !> 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
1331                ELSE IF  (SIZE(match_spec_voc_input) .GT. 0) THEN
1332
1333                  !TBD: maybe we can avoid the cycle
1334                   DO ivoc= 1,SIZE(match_spec_voc_input)
1335
1336                      IF (TRIM(spc_names(match_spec_model(ispec)))==TRIM(emt_att%voc_name(ivoc))) THEN   
1337
1338                         delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)*                                    &
1339                                                       emt_att%voc_comp(icat,match_spec_voc_input(ivoc))*con_factor*hour_per_day
1340
1341                         emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+ &
1342                                                                    delta_emis(nys:nyn,nxl:nxr) 
1343
1344                      ENDIF                       
1345
1346                   ENDDO
1347               
1348             !> General case (other species)
1349                ELSE
1350
1351                   delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)*time_factor(icat)*con_factor*hour_per_day
1352 
1353                   emis_distribution(1,nys:nyn,nxl:nxr,ispec)=emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr) 
1354
1355                ENDIF  ! IF (spc_names==)
1356
1357                !> for every species and category set emis to 0 so to avoid overwriting. TBD: discuss whether necessary
1358
1359                emis(:,:)= 0 
1360
1361             ENDDO
1362
1363          !> for every category set delta_emis to 0 so to avoid overwriting. TBD: discuss whether necessary
1364
1365          delta_emis(:,:)=0 
1366 
1367          ENDDO
1368
1369       ENDIF  !> mode_emis
1370
1371!-------------------------------------------------------------------------------------------------------
1372!--- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to cssws
1373!-------------------------------------------------------------------------------------------------------
1374
1375!-- PARAMETERIZED mode
1376       !> 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
1377       IF (TRIM(mode_emis)=="PARAMETERIZED") THEN
1378
1379          IF ( street_type_f%from_file )  THEN
1380
1381       !> Streets are lsm surfaces, hence, no usm surface treatment required
1382             IF (surf_lsm_h%ns .GT. 0 ) THEN
1383                DO  m = 1, surf_lsm_h%ns
1384                   i = surf_lsm_h%i(m)
1385                   j = surf_lsm_h%j(m)
1386                   k = surf_lsm_h%k(m)
1387
1388
1389                   IF ( street_type_f%var(j,i) >= main_street_id  .AND. street_type_f%var(j,i) < max_street_id )  THEN
1390
1391                   !> Cycle over already matched species
1392                      DO  ispec=1,nspec_out
1393
1394                         !> PMs are already in mass units: kilograms
1395                         IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
1396                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
1397                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1398
1399                            !              kg/(m^2*s) *kg/m^3
1400                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec)) * &
1401                            !                                                       kg/(m^2*s)
1402                                                                                emis_distribution(1,j,i,ispec)*        &
1403                            !                                                    kg/m^3
1404                                                                                rho_air(k)   
1405
1406                         ELSE
1407
1408                         !> Other Species: inputs are micromoles: they have to be converted in moles
1409                         !                 ppm/s *m *kg/m^3               
1410                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec))*   &
1411                         !                                                    micromoles/(m^2*s)
1412                                                                          emis_distribution(1,j,i,ispec) *              &
1413                         !                                                    m^3/Nmole
1414                                                                          conv_to_ratio(k,j,i)*                         &       
1415                         !                                                    kg/m^3
1416                                                                          rho_air(k)   
1417 
1418
1419                         ENDIF
1420
1421                      ENDDO
1422
1423
1424                   ELSEIF ( street_type_f%var(j,i) >= side_street_id  .AND. street_type_f%var(j,i) < main_street_id )  THEN
1425
1426                   !> Cycle over already matched species
1427                      DO  ispec=1,nspec_out
1428
1429                         !> PMs are already in mass units: micrograms
1430                         IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
1431                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
1432                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1433
1434                            !              kg/(m^2*s) *kg/m^3                               
1435                            surf_lsm_h%cssws(match_spec_model(ispec),m)= emiss_factor_side(match_spec_input(ispec)) *   &
1436                            !                                                       kg/(m^2*s)
1437                                                                                emis_distribution(1,j,i,ispec)*        &
1438                            !                                                    kg/m^3
1439                                                                                rho_air(k)   
1440                         ELSE
1441                       
1442
1443                         !>Other Species: inputs are micromoles
1444                         !                 ppm/s *m *kg/m^3               
1445                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_side(match_spec_input(ispec)) *   &
1446                         !                                                    micromoles/(m^2*s)
1447                                                                          emis_distribution(1,j,i,ispec) *              &
1448                         !                                                    m^3/Nmole
1449                                                                          conv_to_ratio(k,j,i)*                         &       
1450                         !                                                    kg/m^3
1451                                                                          rho_air(k)   
1452                         ENDIF
1453
1454                      ENDDO
1455
1456                   ELSE
1457
1458                   !> If no street type is defined, then assign null emissions to all the species
1459                      surf_lsm_h%cssws(:,m) = 0.0_wp
1460
1461                   ENDIF
1462
1463                ENDDO
1464             ENDIF
1465
1466          ENDIF
1467
1468       !> For both DEFAULT and PRE-PROCESSED
1469       ELSE   
1470       
1471
1472          DO ispec=1,nspec_out
1473                                !TBD: for the PRE-PROCESSED mode in the future, values at different heights should be included!           
1474             !> Default surface type
1475             IF (surf_def_h(0)%ns .GT. 0) THEN
1476
1477                DO  m = 1, surf_def_h(0)%ns
1478
1479                   i = surf_def_h(0)%i(m)
1480                   j = surf_def_h(0)%j(m)
1481
1482                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
1483
1484
1485                      !> Distinguish between PMs (no needing conversion in ppms),
1486                      !  VOC (already converted to moles/(m**2*s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
1487                      ! other Species (still expressed in Kg/(m**2*s) at this point)
1488
1489                      !> PMs
1490                      IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
1491                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1492                   
1493                         !            kg/(m^2*s) *kg/m^3                         kg/(m^2*s)                 
1494                         surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*   &
1495                         !                                                  kg/m^3
1496                                                                          rho_air(nzb)   
1497 
1498 
1499                      ELSE
1500
1501                         !> VOCs
1502                         IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
1503                            !          ( ppm/s) * m * kg/m^3                         mole/(m^2/s)   
1504                            surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
1505                                                                            !    m^3/mole          ppm
1506                                                                             conv_to_ratio(nzb,j,i)*ratio2ppm *      &
1507                         !                                                    kg/m^3
1508                                                                             rho_air(nzb)   
1509
1510
1511                         !> OTHER SPECIES
1512                         ELSE
1513
1514                            !               ( ppm/s )*m  * kg/m^3                      kg/(m^2/s)                     
1515                            surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
1516                                                                               !  mole/Kg       
1517                                                                             (1/emt_att%xm(ispec))*                &
1518                            !                                                    m^3/mole          ppm
1519                                                                             conv_to_ratio(nzb,j,i)*ratio2ppm*       &
1520                            !                                                  kg/m^3
1521                                                                             rho_air(nzb)   
1522 
1523
1524                         ENDIF
1525
1526                      ENDIF
1527
1528                   ENDIF
1529
1530                ENDDO
1531
1532             END IF
1533         
1534             !-- Land Surface Mode surface type
1535             IF (surf_lsm_h%ns .GT. 0) THEN
1536
1537                DO m = 1, surf_lsm_h%ns
1538
1539                   i = surf_lsm_h%i(m)
1540                   j = surf_lsm_h%j(m)
1541                   k = surf_lsm_h%k(m)
1542
1543                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
1544
1545                      !> Distinguish between PMs (no needing conversion in ppms),
1546                      !  VOC (already converted to moles/(m**2*s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
1547                      ! other Species (still expressed in Kg/(m**2*s) at this point)
1548
1549                      !> PMs
1550                      IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
1551                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1552   
1553                         !         kg/(m^2*s) * kg/m^3                           kg/(m^2*s)           
1554                         surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *              &
1555                         !                                                  kg/m^3
1556                                                                          rho_air(k)   
1557 
1558                      ELSE
1559
1560                         !> VOCs
1561                         IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
1562                            !          ( ppm/s) * m * kg/m^3                        mole/(m^2/s)   
1563                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
1564                                                                          !    m^3/mole          ppm
1565                                                                          conv_to_ratio(k,j,i)*ratio2ppm*    &
1566                         !                                                 kg/m^3
1567                                                                          rho_air(k)   
1568
1569                         !> OTHER SPECIES
1570                         ELSE
1571   
1572                            !         ( ppm/s) * m * kg/m^3                        kg/(m^2/s)                     
1573                            surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *               &
1574                                                                               !  mole/Kg   
1575                                                                         (1/emt_att%xm(ispec))*                          &
1576                            !                                                m^3/mole           ppm
1577                                                                         conv_to_ratio(k,j,i)*ratio2ppm*                 &
1578                            !                                            kg/m^3
1579                                                                         rho_air(k)   
1580                                                     
1581                         ENDIF
1582
1583                      ENDIF
1584
1585                   ENDIF
1586
1587                ENDDO
1588
1589             END IF
1590
1591             !-- Urban Surface Mode surface type
1592             IF (surf_usm_h%ns .GT. 0) THEN
1593
1594
1595                DO m = 1, surf_usm_h%ns
1596
1597                   i = surf_usm_h%i(m)
1598                   j = surf_usm_h%j(m)
1599                   k = surf_usm_h%k(m)
1600
1601                   IF ( emis_distribution(1,j,i,ispec) > 0.0_wp )  THEN
1602
1603                      !> PMs
1604                      IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25"  &
1605                             .OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
1606                     
1607                         !          kg/(m^2*s) *kg/m^3                             kg/(m^2*s)                     
1608                         surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec)*        & 
1609                         !                                              kg/m^3
1610                                                                       rho_air(k)   
1611
1612 
1613                      ELSE
1614
1615                         !> VOCs
1616                         IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
1617                            !          ( ppm/s) * m * kg/m^3                        mole/(m^2/s)   
1618                            surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *   &
1619                                                                          !    m^3/mole          ppm
1620                                                                          conv_to_ratio(k,j,i)*ratio2ppm*    &
1621                         !                                                kg/m^3
1622                                                                          rho_air(k)   
1623
1624                         !> OTHER SPECIES
1625                         ELSE
1626
1627
1628                         !            ( ppm/s) * m * kg/m^3                        kg/(m^2/s)                     
1629                            surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) *      &
1630                                                                             !  mole/Kg   
1631                                                                         (1/emt_att%xm(ispec))*                 &
1632                            !                                                m^3/mole           ppm
1633                                                                         conv_to_ratio(k,j,i)*ratio2ppm*        &
1634                            !                                            kg/m^3
1635                                                                         rho_air(k)   
1636
1637
1638                         ENDIF
1639
1640                      ENDIF
1641
1642                   ENDIF
1643 
1644                ENDDO
1645             END IF
1646          ENDDO
1647
1648       ENDIF
1649
1650
1651    !> At the end of each call to chem_emissions setup, the time_factor is deallocated (ALLOCATED ONLY in the DEFAULT mode)
1652
1653       IF ( ALLOCATED ( time_factor ) ) DEALLOCATE( time_factor )
1654
1655   ENDIF !> emis_output_required
1656
1657
1658 END SUBROUTINE chem_emissions_setup
1659
1660 END MODULE chem_emissions_mod
Note: See TracBrowser for help on using the repository browser.