Home

Context Navigation

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

Last change on this file since 3428 was 3373, checked in by kanani, 6 years ago
Fix cpp directives and error messages
Property svn:keywords set to `Id`
File size: 78.6 KB

Rev	Line
[3190]	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	!
[3298]	17	! Copyright 2018-2018 Leibniz Universitaet Hannover
	18	! Copyright 2018-2018 Freie Universitaet Berlin
	19	! Copyright 2018-2018 Karlsruhe Institute of Technology
[3190]	20	!--------------------------------------------------------------------------------!
	21	!
	22	! Current revisions:
[3298]	23	! ------------------
[3190]	24	!
[3298]	25	!
[3190]	26	! Former revisions:
	27	! -----------------
	28	! $Id: chem_emissions_mod.f90 3373 2018-10-18 15:25:56Z gronemeier $
[3373]	29	! Fix wrong location of __netcdf directive
	30	!
	31	! 3337 2018-10-12 15:17:09Z kanani
[3337]	32	! (from branch resler)
	33	! Formatting
	34	!
	35	! 3312 2018-10-06 14:15:46Z knoop
[3298]	36	! Initial revision
[3190]	37	!
[3298]	38	! 3286 2018-09-28 07:41:39Z forkel
	39	!
[3190]	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
[3266]	59	USE arrays_3d, &
	60	ONLY: rho_air
	61
[3190]	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
[3373]	73	#if defined ( __netcdf )
[3190]	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
[3228]	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.
[3190]	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'
[3281]	298	CALL message( 'chem_emissions_check_parameters', 'CM0437', 2, 2, 0, 6, 0 )
[3190]	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
[3228]	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
[3190]	325	END SUBROUTINE chem_emissions_check_parameters
	326
	327	!------------------------------------------------------------------------------!
	328	! Description:
	329	! ------------
[3228]	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.
[3190]	334	!------------------------------------------------------------------------------!
	335	SUBROUTINE 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
[3221]	343	INTEGER(iwp) :: nspec_emis_inp !< Variable where to store the number of the emission species in input
[3190]	344
	345	INTEGER(iwp) :: ind_voc !< Indices to check whether a split for voc should be done
	346
[3266]	347	INTEGER(iwp) :: ispec !> index for cycle over effective number of emission species
[3190]	348
[3266]	349
[3190]	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
[3266]	357	!> Check the emission mode: DEFAULT, PRE-PROCESSED or PARAMETERIZED !TBD: Add option for capital or small letters
[3190]	358	SELECT CASE(TRIM(mode_emis))
	359
[3266]	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")
[3190]	362
[3266]	363	CALL location_message( 'chem_emissions PRE-PROCESSED_mode_matching', .FALSE. )
[3190]	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
[3266]	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
[3190]	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' // &
[3266]	449	' DO NOT MATCH' // &
	450	' model chemical species:' // &
	451	' Chemistry Emissions routine is not called'
[3281]	452	CALL message( 'chem_emissions_matching', 'CM0438', 0, 0, 0, 6, 0 )
[3190]	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:' // &
[3266]	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'
[3281]	463	CALL message( 'chem_emissions_matching', 'CM0439', 0, 2, 0, 6, 0 )
[3190]	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' // &
[3266]	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
[3281]	664	CALL message( 'chem_emissions_matching', 'CM0440', 0, 0, 0, 6, 0 )
[3190]	665
	666	ENDIF
	667
	668	ELSE
	669
	670	message_string = 'Array of Emission species not allocated: ' // &
[3266]	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'
[3281]	674	CALL message( 'chem_emissions_matching', 'CM0441', 0, 2, 0, 6, 0 )
[3190]	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
[3266]	687
[3190]	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
[3266]	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
[3337]	726	IF ( emiss_factor_main(match_spec_input(ispec)) .LT. 0 .AND. &
	727	emiss_factor_side(match_spec_input(ispec)) .LT. 0 ) THEN
[3266]	728
[3337]	729	message_string = 'PARAMETERIZED emissions mode selected:' // &
[3266]	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'
[3281]	736	CALL message( 'chem_emissions_matching', 'CM0442', 2, 2, 0, 6, 0 )
[3266]	737	ENDIF
	738	END DO
	739
	740
[3190]	741	ELSE
	742
	743	message_string = 'Given Chemistry Emission Species' // &
[3266]	744	' DO NOT MATCH' // &
	745	' model chemical species' // &
	746	' Chemistry Emissions routine is not called'
[3281]	747	CALL message( 'chem_emissions_matching', 'CM0443', 0, 0, 0, 6, 0 )
[3190]	748	ENDIF
	749
	750	ELSE
[3266]	751
[3190]	752	message_string = 'Array of Emission species not allocated: ' // &
[3266]	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'
[3281]	756	CALL message( 'chem_emissions_matching', 'CM0444', 0, 2, 0, 6, 0 )
[3190]	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' // &
[3266]	765	' either no emission mode specified or incorrectly given :' // &
	766	' please, pass the correct value to the namelist parameter "mode_emis"'
[3281]	767	CALL message( 'chem_emissions_matching', 'CM0445', 2, 2, 0, 6, 0 )
[3190]	768
	769	END SELECT
	770
	771	END 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
[3373]	783
[3190]	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
[3373]	801	#if defined( __netcdf )
[3190]	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
[3266]	867	!> PRE-PROCESSED case
	868	CASE ("PRE-PROCESSED")
[3190]	869
	870	IF ( .NOT. ALLOCATED( emis_distribution) ) ALLOCATE(emis_distribution(nzb:nzt+1,0:ny,0:nx,nspec_out))
	871
[3266]	872	CALL location_message( 'emis_distribution array allocated in PRE-PROCESSED mode', .FALSE. )
[3190]	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
[3266]	906	USE arrays_3d, &
	907	ONLY: dzw
[3190]	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
[3373]	921
[3190]	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
[3266]	932	INTEGER,INTENT(IN) :: nspec_out !> Output of routine chem_emis_matching with number
	933	! of matched species
[3190]	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
[3266]	948	REAL(wp),ALLOCATABLE, DIMENSION(:) :: time_factor !> factor for time scaling of emissions
[3190]	949	REAL(wp),ALLOCATABLE, DIMENSION(:,:) :: emis
	950
[3266]	951	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: conv_to_ratio !> factor used for converting input
[3190]	952	! to adimensional concentration ratio
[3266]	953
[3190]	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
[3221]	962	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: tmp_temp
[3190]	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
[3266]	984	REAL(wp), PARAMETER :: ratio2ppm = 1.0e06_wp
[3190]	985	!------------------------------------------------------
[3373]	986	#if defined( __netcdf )
[3190]	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	!>-----------------------------------------------------
[3286]	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.
[3190]	1008	!>-----------------------------------------------------
	1009
[3266]	1010	IF (TRIM(mode_emis)=="DEFAULT" .OR. TRIM(mode_emis)=="PRE-PROCESSED" ) THEN
[3190]	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 ' // &
[3266]	1085	' chemistry emission module documentation.'
[3312]	1086	CALL message( 'chem_emissions_setup', 'CM0446', 2, 2, 0, 6, 0 )
[3190]	1087
	1088	END SELECT
	1089
	1090
[3266]	1091	!> PRE-PROCESSED and parameterized mode
[3190]	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 '
[3281]	1096	CALL message( 'chem_emissions_setup', 'CM0447', 0, 0, 0, 6, 0 )
[3190]	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.
[3266]	1107	! Input is Nmole/(m^2*s)
[3190]	1108	! To go to ppmdz basically we need to multiply the input by (m2/N)dz
[3266]	1109	! (m*2/N)dz == V/N
[3190]	1110	! V/N=RT/P
	1111
[3266]	1112	!> m*3/Nmole (J/mol)K^-1 K Pa
[3221]	1113	conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) )
[3190]	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	!>-----------------------------------------------
[3266]	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
[3190]	1124
[3266]	1125	!> PRE-PROCESSED MODE
	1126	IF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
[3190]	1127
[3266]	1128	CALL location_message( 'PRE-PROCESSED MODE: No time-factor specification required', .FALSE. )
[3190]	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
[3266]	1154	message_string = 'The "HOUR" time-factors (DEFAULT mode) of the chemistry emission species' // &
[3190]	1155	' are not provided for each hour of the total simulation time'
[3281]	1156	CALL message( 'chem_emissions_setup', 'CM0448', 2, 2, 0, 6, 0 )
[3190]	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
[3266]	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)
[3190]	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'
[3281]	1178	CALL message( 'chem_emissions_setup', 'CM0449', 2, 2, 0, 6, 0 )
[3190]	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'
[3281]	1187	CALL message( 'chem_emissions_setup', 'CM0450', 2, 2, 0, 6, 0 )
[3190]	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
[3266]	1214	!> Values are still micromoles/(m*2s). Units are not given by the user for this mode and are always micromoles (or micrograms for PMs)
[3190]	1215	emis_distribution(1,nys:nyn,nxl:nxr,ispec)=surface_csflux(match_spec_input(ispec))*time_factor(1)
	1216
	1217	ENDDO
	1218
[3266]	1219	!> PRE-PROCESSED CASE
	1220	ELSEIF (TRIM(mode_emis)=="PRE-PROCESSED") THEN
[3190]	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
[3267]	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
[3190]	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
[3337]	1254	delta_emis(nys:nyn,nxl:nxr) = &
	1255	emis(nys:nyn,nxl:nxr)time_factor(icat)emt_att%nox_comp(icat,1)*con_factor
[3190]	1256
[3337]	1257	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1258	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	1259
	1260	ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="NO2") THEN
	1261
[3337]	1262	delta_emis(nys:nyn,nxl:nxr) = &
	1263	emis(nys:nyn,nxl:nxr)time_factor(icat)emt_att%nox_comp(icat,2)*con_factor
[3190]	1264
[3337]	1265	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1266	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	1267
	1268	!> SOX Compositions
	1269
	1270	ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO2") THEN
	1271
[3337]	1272	delta_emis(nys:nyn,nxl:nxr) = &
	1273	emis(nys:nyn,nxl:nxr)time_factor(icat)emt_att%sox_comp(icat,1)*con_factor
[3190]	1274
[3337]	1275	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1276	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	1277
	1278	ELSE IF (TRIM(spc_names(match_spec_model(ispec)))=="SO4") THEN
	1279
[3337]	1280	delta_emis(nys:nyn,nxl:nxr) = &
	1281	emis(nys:nyn,nxl:nxr)time_factor(icat)emt_att%sox_comp(icat,2)*con_factor
[3190]	1282
[3337]	1283	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1284	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	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
[3337]	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
[3190]	1296
	1297
[3337]	1298	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1299	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	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
[3337]	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
[3190]	1311
	1312
[3337]	1313	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1314	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	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
[3337]	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
[3190]	1326
	1327
[3337]	1328	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1329	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	1330
	1331	ENDDO
	1332
	1333	!> VOCs Compositions: for VOCs, the input value is provided in kg/m*2s but the composition is provided in mole/kg: a conversion factor for the input that could be eventually provided in, for example, tons/(m*2s) 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
[3337]	1341	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr)time_factor(icat) &
[3190]	1342	emt_att%voc_comp(icat,match_spec_voc_input(ivoc))*con_factor
	1343
[3337]	1344	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1345	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	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
[3337]	1356	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
	1357	emis_distribution(1,nys:nyn,nxl:nxr,ispec)+delta_emis(nys:nyn,nxl:nxr)
[3190]	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
[3266]	1375	!-------------------------------------------------------------------------------------------------------
[3190]	1376	!--- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to cssws
[3266]	1377	!-------------------------------------------------------------------------------------------------------
[3190]	1378
	1379	!-- PARAMETERIZED mode
	1380	!> For the PARAMETERIZED mode units of inputs are always micromoles/(m*2s). 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
[3266]	1398	!> PMs are already in mass units:micrograms:they have to be converted to kilograms
[3337]	1399	IF (TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
	1400	.OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25" &
[3190]	1401	.OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
	1402
[3266]	1403	! kg/(m^2s) kg/m^3
[3190]	1404	surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec)) * &
[3276]	1405	! kg/(m^2*s)
[3266]	1406	emis_distribution(1,j,i,ispec)* &
	1407	! kg/m^3
	1408	rho_air(k)
[3190]	1409
	1410	ELSE
	1411
[3266]	1412	!> Other Species: inputs are micromoles: they have to be converted in moles
	1413	! ppm/s m kg/m^3
[3337]	1414	surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_main(match_spec_input(ispec))* &
[3266]	1415	! micromoles/(m^2*s)
[3337]	1416	emis_distribution(1,j,i,ispec) * &
[3266]	1417	! m**3/Nmole
[3337]	1418	conv_to_ratio(k,j,i)* &
[3266]	1419	! kg/m^3
	1420	rho_air(k)
[3190]	1421
[3266]	1422
[3190]	1423	ENDIF
	1424
	1425	ENDDO
	1426
	1427
[3337]	1428	ELSEIF ( street_type_f%var(j,i) >= side_street_id .AND. &
	1429	street_type_f%var(j,i) < main_street_id ) THEN
[3190]	1430
	1431	!> Cycle over already matched species
	1432	DO ispec=1,nspec_out
	1433
	1434	!> PMs are already in mass units: micrograms
[3337]	1435	IF ( TRIM(spc_names(match_spec_model(ispec)))=="PM1" &
	1436	.OR. TRIM(spc_names(match_spec_model(ispec)))=="PM25" &
[3190]	1437	.OR. TRIM(spc_names(match_spec_model(ispec)))=="PM10") THEN
	1438
[3266]	1439	! kg/(m^2s) kg/m^3
[3337]	1440	surf_lsm_h%cssws(match_spec_model(ispec),m)= emiss_factor_side(match_spec_input(ispec)) * &
[3276]	1441	! kg/(m^2*s)
[3266]	1442	emis_distribution(1,j,i,ispec)* &
	1443	! kg/m^3
	1444	rho_air(k)
[3190]	1445	ELSE
	1446
	1447
	1448	!>Other Species: inputs are micromoles
[3266]	1449	! ppm/s m kg/m^3
[3337]	1450	surf_lsm_h%cssws(match_spec_model(ispec),m) = emiss_factor_side(match_spec_input(ispec)) * &
[3266]	1451	! micromoles/(m^2*s)
[3337]	1452	emis_distribution(1,j,i,ispec) * &
[3266]	1453	! m**3/Nmole
[3337]	1454	conv_to_ratio(k,j,i)* &
[3266]	1455	! kg/m^3
	1456	rho_air(k)
[3190]	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
[3266]	1474	!> For both DEFAULT and PRE-PROCESSED
[3190]	1475	ELSE
	1476
	1477
[3266]	1478	DO ispec=1,nspec_out
	1479	!TBD: for the PRE-PROCESSED mode in the future, values at different heights should be included!
[3190]	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*2s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
	1490	! other Species (still expressed in Kg/(m*2s) at this point)
	1491
	1492	!> PMs
[3337]	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
[3190]	1496
[3266]	1497	! kg/(m^2s) 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)
[3190]	1501
[3266]	1502
[3190]	1503	ELSE
	1504
	1505	!> VOCs
	1506	IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
[3266]	1507	! ( ppm/s) * m * kg/m^3 mole/(m**2/s)
[3190]	1508	surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) * &
	1509	! m**3/mole ppm
[3266]	1510	conv_to_ratio(k,j,i)ratio2ppm &
	1511	! kg/m^3
	1512	rho_air(k)
[3190]	1513
[3266]	1514
[3190]	1515	!> OTHER SPECIES
	1516	ELSE
	1517
[3266]	1518	! ( ppm/s )m kg/m^3 kg/(m**2/s)
[3190]	1519	surf_def_h(0)%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) * &
[3266]	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
[3190]	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*2s) using conv_factors: they do not need molar masses for their conversion to PPMs ) and
	1547	! other Species (still expressed in Kg/(m*2s) at this point)
	1548
	1549	!> PMs
[3337]	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
[3190]	1553
[3266]	1554	! kg/(m^2s) 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)
[3190]	1558
	1559	ELSE
	1560
	1561	!> VOCs
	1562	IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
[3266]	1563	! ( ppm/s) * m * kg/m^3 mole/(m**2/s)
[3190]	1564	surf_lsm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) * &
[3266]	1565	! m**3/mole ppm
	1566	conv_to_ratio(k,j,i)ratio2ppm &
	1567	! kg/m^3
	1568	rho_air(k)
[3190]	1569
[3266]	1570
[3190]	1571	!> OTHER SPECIES
	1572	ELSE
	1573
[3266]	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
[3190]	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
[3337]	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
[3190]	1605
[3266]	1606	! kg/(m^2s) 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
[3190]	1611
	1612	ELSE
	1613
	1614	!> VOCs
	1615	IF ( len_index_voc .GT. 0 .AND. emt_att%species_name(match_spec_input(ispec))=="VOC" ) THEN
[3266]	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)
[3190]	1622
	1623	!> OTHER SPECIES
	1624	ELSE
	1625
	1626
[3266]	1627	! ( ppm/s) * m * kg/m^3 kg/(m**2/s)
[3190]	1628	surf_usm_h%cssws(match_spec_model(ispec),m) = emis_distribution(1,j,i,ispec) * &
[3266]	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
[3190]	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.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |