Home

Context Navigation

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

Last change on this file since 3997 was 3968, checked in by suehring, 6 years ago
Updates from chemistriy branched merged into trunk: code cleaning and formatting, code structure optimizations
Property svn:keywords set to `Id`
File size: 65.0 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |