Context Navigation

readspecies.f90 @ 37

Last change on this file since 37 was 37, checked in by pesei, 9 years ago
Wet dep quick fix and other small changes. Wet depo quick fix not final yet.
File size: 11.8 KB

Line
1	!**********************************************************************
2	! Copyright 1998-2015 *
3	! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa, *
4	! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann *
5	! *
6	! This file is part of FLEXPART. *
7	! *
8	! FLEXPART is free software: you can redistribute it and/or modify *
9	! it under the terms of the GNU General Public License as published by*
10	! the Free Software Foundation, either version 3 of the License, or *
11	! (at your option) any later version. *
12	! *
13	! FLEXPART is distributed in the hope that it will be useful, *
14	! but WITHOUT ANY WARRANTY; without even the implied warranty of *
15	! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
16	! GNU General Public License for more details. *
17	! *
18	! You should have received a copy of the GNU General Public License *
19	! along with FLEXPART. If not, see <http://www.gnu.org/licenses/>. *
20	!**********************************************************************
21
22	subroutine readspecies(id_spec,id_pos)
23
24	!*****************************************************************************
25	! *
26	! This routine reads names and physical constants of chemical species/ *
27	! radionuclides given in the parameter id_pos *
28	! *
29	! Author: A. Stohl *
30	! *
31	! 11 July 1996 *
32	!
33	! Changes: *
34	! *
35	! N. Kristiansen, 31.01.2013: Including parameters for in-cloud scavenging *
36	! *
37	! HSO, 13 August 2013: added optional namelist input
38	! PS, 2/2015: access= -> position=
39	! PS, 4/2015: quick fix for bug associated with end=22
40	! *
41	!*****************************************************************************
42	! *
43	! Variables: *
44	! decaytime(maxtable) half time for radiological decay *
45	! specname(maxtable) names of chemical species, radionuclides *
46	! weta, wetb Parameters for determining below-cloud scavenging *
47	! weta_in Parameter for determining in-cloud scavenging *
48	! wetb_in Parameter for determining in-cloud scavenging *
49	! wetc_in Parameter for determining in-cloud scavenging *
50	! wetd_in Parameter for determining in-cloud scavenging *
51	! ohreact OH reaction rate *
52	! id_spec SPECIES number as referenced in RELEASE file *
53	! id_pos position where SPECIES data shall be stored *
54	! *
55	! Constants: *
56	! *
57	!*****************************************************************************
58
59	use par_mod
60	use com_mod
61
62	implicit none
63
64	integer :: i, id_pos,j
65	integer :: idow,ihour,id_spec
66	character(len=3) :: aspecnumb
67	logical :: spec_found
68
69	character(len=16) :: pspecies
70	real :: pdecay, pweta, pwetb, preldiff, phenry, pf0, pdensity, pdquer
71	real :: pdsigma, pdryvel, pweightmolar, pohreact, pspec_ass, pkao
72	real :: pweta_in, pwetb_in, pwetc_in, pwetd_in
73	integer :: ios
74
75	! declare namelist
76	namelist /species_params/ &
77	pspecies, pdecay, pweta, pwetb, &
78	pweta_in, pwetb_in, pwetc_in, pwetd_in, &
79	preldiff, phenry, pf0, pdensity, pdquer, &
80	pdsigma, pdryvel, pweightmolar, pohreact, pspec_ass, pkao
81
82	pspecies=" "
83	pdecay=-999.9
84	pweta=-9.9E-09
85	pwetb=0.0
86	pweta_in=-9.9E-09
87	pwetb_in=-9.9E-09
88	pwetc_in=-9.9E-09
89	pwetd_in=-9.9E-09
90	preldiff=-9.9
91	phenry=0.0
92	pf0=0.0
93	pdensity=-9.9E09
94	pdquer=0.0
95	pdsigma=0.0
96	pdryvel=-9.99
97	pohreact=-9.9E-09
98	pspec_ass=-9
99	pkao=-99.99
100	pweightmolar=-789.0 ! read failure indicator value
101
102	! Open the SPECIES file and read species names and properties
103	!************************************************************
104	specnum(id_pos)=id_spec
105	write(aspecnumb,'(i3.3)') specnum(id_pos)
106	open(unitspecies,file=path(1)(1:length(1))//'SPECIES/SPECIES_'//aspecnumb,status='old',form='formatted',err=998)
107	!write(,) 'reading SPECIES',specnum(id_pos)
108
109	ASSSPEC=.FALSE.
110
111	! try namelist input
112	read(unitspecies,species_params,iostat=ios)
113	close(unitspecies)
114
115	if ((pweightmolar.eq.-789.0).or.(ios.ne.0)) then ! no namelist found
116
117	open(unitspecies,file=path(1) (1:length(1))//'SPECIES/SPECIES_'//aspecnumb, status='old', err=998)
118
119	do i=1,6
120	read(unitspecies,*)
121	end do
122	read(unitspecies,'(a10)',end=20) species(id_pos)
123	! write(,) species(id_pos)
124	read(unitspecies,'(f18.1)',end=20) decay(id_pos)
125	! write(,) decay(id_pos)
126	read(unitspecies,'(e18.1)',end=20) weta(id_pos)
127	! write(,) weta(id_pos)
128	read(unitspecies,'(f18.2)',end=20) wetb(id_pos)
129	! write(,) wetb(id_pos)
130
131	!*** NIK 31.01.2013: including in-cloud scavening parameters
132	read(unitspecies,'(e18.1)',end=20) weta_in(id_pos)
133	! write(,) weta_in(id_pos)
134	read(unitspecies,'(f18.2)',end=20) wetb_in(id_pos)
135	! write(,) wetb_in(id_pos)
136	read(unitspecies,'(f18.2)',end=20) wetc_in(id_pos)
137	! write(,) wetc_in(id_pos)
138	read(unitspecies,'(f18.2)',end=20) wetd_in(id_pos)
139	! write(,) wetd_in(id_pos)
140
141	read(unitspecies,'(f18.1)',end=20) reldiff(id_pos)
142	! write(,) reldiff(id_pos)
143	read(unitspecies,'(e18.1)',end=20) henry(id_pos)
144	! write(,) henry(id_pos)
145	read(unitspecies,'(f18.1)',end=20) f0(id_pos)
146	! write(,) f0(id_pos)
147	read(unitspecies,'(e18.1)',end=20) density(id_pos)
148	! write(,) density(id_pos)
149	read(unitspecies,'(e18.1)',end=20) dquer(id_pos)
150	! write(,) dquer(id_pos)
151	read(unitspecies,'(e18.1)',end=20) dsigma(id_pos)
152	! write(,) dsigma(id_pos)
153	read(unitspecies,'(f18.2)',end=20) dryvel(id_pos)
154	! write(,) dryvel(id_pos)
155	read(unitspecies,'(f18.2)',end=20) weightmolar(id_pos)
156	! write(,) weightmolar(id_pos)
157	read(unitspecies,'(e18.1)',end=20) ohreact(id_pos)
158	! write(,) ohreact(id_pos)
159	read(unitspecies,'(i18)',end=20) spec_ass(id_pos)
160	! write(,) spec_ass(id_pos)
161	read(unitspecies,'(f18.2)',end=20) kao(id_pos)
162	! write(,) kao(id_pos)
163
164	pspecies=species(id_pos)
165	pdecay=decay(id_pos)
166	pweta=weta(id_pos)
167	pwetb=wetb(id_pos)
168	pweta_in=weta_in(id_pos)
169	pwetb_in=wetb_in(id_pos)
170	pwetc_in=wetc_in(id_pos)
171	pwetd_in=wetd_in(id_pos)
172	preldiff=reldiff(id_pos)
173	phenry=henry(id_pos)
174	pf0=f0(id_pos)
175	pdensity=density(id_pos)
176	pdquer=dquer(id_pos)
177	pdsigma=dsigma(id_pos)
178	pdryvel=dryvel(id_pos)
179	pweightmolar=weightmolar(id_pos)
180	pohreact=ohreact(id_pos)
181	pspec_ass=spec_ass(id_pos)
182	pkao=kao(id_pos)
183
184	else
185
186	species(id_pos)=pspecies
187	decay(id_pos)=pdecay
188	weta(id_pos)=pweta
189	wetb(id_pos)=pwetb
190	weta_in(id_pos)=pweta_in
191	wetb_in(id_pos)=pwetb_in
192	wetc_in(id_pos)=pwetc_in
193	wetd_in(id_pos)=pwetd_in
194	reldiff(id_pos)=preldiff
195	henry(id_pos)=phenry
196	f0(id_pos)=pf0
197	density(id_pos)=pdensity
198	dquer(id_pos)=pdquer
199	dsigma(id_pos)=pdsigma
200	dryvel(id_pos)=pdryvel
201	weightmolar(id_pos)=pweightmolar
202	ohreact(id_pos)=pohreact
203	spec_ass(id_pos)=pspec_ass
204	kao(id_pos)=pkao
205
206	endif
207
208	if ((weta(id_pos).gt.0).and.(henry(id_pos).le.0)) then
209	if (dquer(id_pos).le.0) goto 996 ! no particle, no henry set
210	endif
211
212	if (spec_ass(id_pos).gt.0) then
213	spec_found=.FALSE.
214	do j=1,id_pos-1
215	if (spec_ass(id_pos).eq.specnum(j)) then
216	spec_ass(id_pos)=j
217	spec_found=.TRUE.
218	ASSSPEC=.TRUE.
219	endif
220	end do
221	if (spec_found.eqv..FALSE.) then
222	goto 997
223	endif
224	endif
225
226	if (dsigma(id_pos).eq.1.) dsigma(id_pos)=1.0001 ! avoid floating exception
227	if (dsigma(id_pos).eq.0.) dsigma(id_pos)=1.0001 ! avoid floating exception
228
229	if ((reldiff(id_pos).gt.0.).and.(density(id_pos).gt.0.)) then
230	write(,) '#### FLEXPART MODEL ERROR! FILE "SPECIES" ####'
231	write(,) '#### IS CORRUPT. SPECIES CANNOT BE BOTH ####'
232	write(,) '#### PARTICLE AND GAS. ####'
233	write(,) '#### SPECIES NUMBER',aspecnumb
234	stop
235	endif
236	20 continue
237
238
239	! Read in daily and day-of-week variation of emissions, if available
240	!*******************************************************************
241	! HSO: This is not yet implemented as namelist parameters
242	! default values set to 1
243
244	do j=1,24 ! initialize everything to no variation
245	area_hour(id_pos,j)=1.
246	point_hour(id_pos,j)=1.
247	end do
248	do j=1,7
249	area_dow(id_pos,j)=1.
250	point_dow(id_pos,j)=1.
251	end do
252
253	if (ios.ne.0) then ! text format input
254
255	read(unitspecies,*,iostat=ios)
256	if (ios .ne. 0) goto 22 ! ifort does not like err=
257	do j=1,24 ! 24 hours, starting with 0-1 local time
258	read(unitspecies,*) ihour,area_hour(id_pos,j),point_hour(id_pos,j)
259	end do
260	read(unitspecies,*)
261	do j=1,7 ! 7 days of the week, starting with Monday
262	read(unitspecies,*) idow,area_dow(id_pos,j),point_dow(id_pos,j)
263	end do
264
265	endif
266
267	22 close(unitspecies)
268
269	! namelist output if requested
270	if (nmlout.eqv..true.) then
271	open(unitspecies,file=path(2)(1:length(2))//'SPECIES_'//aspecnumb//'.namelist',position='append',status='new',err=1000)
272	write(unitspecies,nml=species_params)
273	close(unitspecies)
274	endif
275
276	return
277
278	996 write(,) '#####################################################'
279	write(,) '#### FLEXPART MODEL ERROR! #### '
280	write(,) '#### WET DEPOSITION SWITCHED ON, BUT NO HENRYS #### '
281	write(,) '#### CONSTANT IS SET ####'
282	write(,) '#### PLEASE MODIFY SPECIES DESCR. FILE! #### '
283	write(,) '#####################################################'
284	stop
285
286
287	997 write(,) '#####################################################'
288	write(,) '#### FLEXPART MODEL ERROR! #### '
289	write(,) '#### THE ASSSOCIATED SPECIES HAS TO BE DEFINED #### '
290	write(,) '#### BEFORE THE ONE WHICH POINTS AT IT #### '
291	write(,) '#### PLEASE CHANGE ORDER IN RELEASES OR ADD #### '
292	write(,) '#### THE ASSOCIATED SPECIES IN RELEASES #### '
293	write(,) '#####################################################'
294	stop
295
296
297	998 write(,) '#####################################################'
298	write(,) '#### FLEXPART MODEL ERROR! #### '
299	write(,) '#### THE SPECIES FILE FOR SPECIES ', id_spec
300	write(,) '#### CANNOT BE FOUND: CREATE FILE'
301	write(,) '#### ',path(1)(1:length(1)),'SPECIES/SPECIES_',aspecnumb
302	write(,) '#####################################################'
303	stop
304
305	1000 write(,) ' #### FLEXPART MODEL ERROR! FILE "SPECIES_',aspecnumb,'.namelist'
306	write(,) ' #### CANNOT BE OPENED IN THE DIRECTORY #### '
307	write(*,'(a)') path(2)(1:length(2))
308	stop
309
310	end subroutine readspecies

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

Download in other formats:

Original Format