Home

Context Navigation

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

Last change on this file since 3814 was 3788, checked in by banzhafs, 6 years ago
Removed unused variables from chem_emissions_mod
Property svn:keywords set to `Id`
File size: 63.6 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |