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

Last change on this file since 3445 was 3373, checked in by kanani, 5 years ago

Fix cpp directives and error messages

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