Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |