Home

Context Navigation

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

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |