Context Navigation

readavailable.f90 @ 5f2e8f6

flexpart-noresm

Last change on this file since 5f2e8f6 was 5f2e8f6, checked in by Ignacio Pisso <Ignacio.Pisso@…>, 8 years ago
new flexpart noresm code in src
Property mode set to `100755`
File size: 12.6 KB

Line
1	!**********************************************************************
2	! Copyright 2016 *
3	! Andreas Stohl, Massimo Cassiani, Petra Seibert, A. Frank, *
4	! Gerhard Wotawa, Caroline Forster, Sabine Eckhardt, John Burkhart, *
5	! Harald Sodemann, Ignacio Pisso *
6	! *
7	! This file is part of FLEXPART-NorESM *
8	! *
9	! FLEXPART-NorESM is free software: you can redistribute it *
10	! and/or modify *
11	! it under the terms of the GNU General Public License as published by*
12	! the Free Software Foundation, either version 3 of the License, or *
13	! (at your option) any later version. *
14	! *
15	! FLEXPART-NorESM is distributed in the hope that it will be useful, *
16	! but WITHOUT ANY WARRANTY; without even the implied warranty of *
17	! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
18	! GNU General Public License for more details. *
19	! *
20	! You should have received a copy of the GNU General Public License *
21	! along with FLEXPART-NorESM. *
22	! If not, see <http://www.gnu.org/licenses/>. *
23	!**********************************************************************
24
25	subroutine readavailable
26
27	!*****************************************************************************
28	! *
29	! This routine reads the dates and times for which windfields are *
30	! available. *
31	! *
32	! Authors: A. Stohl *
33	! *
34	! 6 February 1994 *
35	! 8 February 1999, Use of nested fields, A. Stohl *
36	! *
37	! modified by M. Cassiani, 2016 *
38	! adapted for NorESM output with no-leap calendar *
39	! *
40	!*****************************************************************************
41	! *
42	! Variables: *
43	! bdate beginning date as Julian date *
44	! beg beginning date for windfields *
45	! end ending date for windfields *
46	! fname filename of wind field, help variable *
47	! ideltas [s] duration of modelling period *
48	! idiff time difference between 2 wind fields *
49	! idiffnorm normal time difference between 2 wind fields *
50	! idiffmax [s] maximum allowable time between 2 wind fields *
51	! jul julian date, help variable *
52	! numbwf actual number of wind fields *
53	! wfname(maxwf) file names of needed wind fields *
54	! wfspec(maxwf) file specifications of wind fields (e.g., if on disc) *
55	! wftime(maxwf) [s]times of wind fields relative to beginning time *
56	! wfname1,wfspec1,wftime1 = same as above, but only local (help variables) *
57	! *
58	! Constants: *
59	! maxwf maximum number of wind fields *
60	! unitavailab unit connected to file AVAILABLE *
61	! *
62	!*****************************************************************************
63
64	use par_mod
65	use com_mod
66
67	implicit none
68	!integer YYYYMMDD,HHMISS !added by mc mc for test reason
69	integer :: i,idiff,ldat,ltim,wftime1(maxwf),numbwfn(maxnests),k
70	integer :: wftime1n(maxnests,maxwf),wftimen(maxnests,maxwf)
71	real(kind=dp) :: juldate,jul,beg,end
72	character(len=255) :: fname,spec,wfname1(maxwf),wfspec1(maxwf)
73	character(len=255) :: wfname1n(maxnests,maxwf)
74	character(len=40) :: wfspec1n(maxnests,maxwf)
75
76
77	! Windfields are only used, if they are within the modelling period.
78	! However, 1 additional day at the beginning and at the end is used for
79	! interpolation. -> Compute beginning and ending date for the windfields.
80	!************************************************************************
81
82	if (ideltas.gt.0) then ! forward trajectories
83	beg=bdate-1._dp
84	end=bdate+real(ideltas,kind=dp)/86400._dp+real(idiffmax,kind=dp)/ &
85	86400._dp
86	else ! backward trajectories
87	beg=bdate+real(ideltas,kind=dp)/86400._dp-real(idiffmax,kind=dp)/ &
88	86400._dp
89	end=bdate+1._dp
90	endif
91
92	! Open the wind field availability file and read available wind fields
93	! within the modelling period.
94	!*********************************************************************
95
96	open(unitavailab,file=path(4)(1:length(4)),status='old', &
97	err=999)
98
99	do i=1,3
100	read(unitavailab,*)
101	end do
102
103	numbwf=0
104	100 read(unitavailab,'(i8,1x,i6,2(6x,a255))',end=99) &
105	ldat,ltim,fname,spec
106	if (ldat.lt.yearorigin) then !modified by mc
107	write(,) ' #### FLEXPART MODEL ERROR! date ####' !modified by mc
108	write(,) ' #### must be gt yearorign: Difference in time ####' !modified by mc
109	write(,) ' #### computed using no-leap calendar from an ####' !modified by mc
110	write(,) ' #### arbitrary year zero #### ' !modified by mc
111	write(,) ' #### "AVAILABLE FILE". #### ' !modified by mc
112	stop
113	end if
114	jul=juldate(ldat,ltim)
115	!call CALDATE(JUL,YYYYMMDD,HHMISS) !added by mc for test reason
116	if ((jul.ge.beg).and.(jul.le.end)) then
117	numbwf=numbwf+1
118	if (numbwf.gt.maxwf) then ! check exceedance of dimension
119	write(,) 'Number of wind fields needed is too great.'
120	write(,) 'Reduce modelling period (file "COMMAND") or'
121	write(,) 'reduce number of wind fields (file "AVAILABLE").'
122	stop
123	endif
124
125	wfname1(numbwf)=fname(1:index(fname,' '))
126	wfspec1(numbwf)=spec
127	wftime1(numbwf)=nint((jul-bdate)*86400._dp)
128	endif
129	goto 100 ! next wind field
130
131	99 continue
132
133	close(unitavailab)
134
135	! Open the wind field availability file and read available wind fields
136	! within the modelling period (nested grids)
137	!*********************************************************************
138
139	do k=1,numbnests
140	open(unitavailab,file=path(numpath+2*(k-1)+2) &
141	(1:length(numpath+2*(k-1)+2)),status='old',err=998)
142
143	do i=1,3
144	read(unitavailab,*)
145	end do
146
147	numbwfn(k)=0
148	700 read(unitavailab,'(i8,1x,i6,2(6x,a255))',end=699) ldat, &
149	ltim,fname,spec
150	jul=juldate(ldat,ltim)
151	if ((jul.ge.beg).and.(jul.le.end)) then
152	numbwfn(k)=numbwfn(k)+1
153	if (numbwfn(k).gt.maxwf) then ! check exceedance of dimension
154	write(,) 'Number of nested wind fields is too great.'
155	write(,) 'Reduce modelling period (file "COMMAND") or'
156	write(,) 'reduce number of wind fields (file "AVAILABLE").'
157	stop
158	endif
159
160	wfname1n(k,numbwfn(k))=fname
161	wfspec1n(k,numbwfn(k))=spec
162	wftime1n(k,numbwfn(k))=nint((jul-bdate)*86400._dp)
163	endif
164	goto 700 ! next wind field
165
166	699 continue
167
168	close(unitavailab)
169	end do
170
171
172	! Check wind field times of file AVAILABLE (expected to be in temporal order)
173	!****************************************************************************
174
175	if (numbwf.eq.0) then
176	write(,) ' #### FLEXPART MODEL ERROR! NO WIND FIELDS #### '
177	write(,) ' #### AVAILABLE FOR SELECTED TIME PERIOD. #### '
178	stop
179	endif
180
181	do i=2,numbwf
182	if (wftime1(i).le.wftime1(i-1)) then
183	write(,) 'FLEXPART ERROR: FILE AVAILABLE IS CORRUPT.'
184	write(,) 'THE WIND FIELDS ARE NOT IN TEMPORAL ORDER.'
185	write(,) 'PLEASE CHECK FIELD ',wfname1(i)
186	stop
187	endif
188	end do
189
190	! Check wind field times of file AVAILABLE for the nested fields
191	! (expected to be in temporal order)
192	!***************************************************************
193
194	do k=1,numbnests
195	if (numbwfn(k).eq.0) then
196	write(,) '#### FLEXPART MODEL ERROR! NO WIND FIELDS ####'
197	write(,) '#### AVAILABLE FOR SELECTED TIME PERIOD. ####'
198	stop
199	endif
200
201	do i=2,numbwfn(k)
202	if (wftime1n(k,i).le.wftime1n(k,i-1)) then
203	write(,) 'FLEXPART ERROR: FILE AVAILABLE IS CORRUPT. '
204	write(,) 'THE NESTED WIND FIELDS ARE NOT IN TEMPORAL ORDER.'
205	write(,) 'PLEASE CHECK FIELD ',wfname1n(k,i)
206	write(,) 'AT NESTING LEVEL ',k
207	stop
208	endif
209	end do
210
211	end do
212
213
214	! For backward trajectories, reverse the order of the windfields
215	!***************************************************************
216
217	if (ideltas.ge.0) then
218	do i=1,numbwf
219	wfname(i)=wfname1(i)
220	wfspec(i)=wfspec1(i)
221	wftime(i)=wftime1(i)
222	end do
223	do k=1,numbnests
224	do i=1,numbwfn(k)
225	wfnamen(k,i)=wfname1n(k,i)
226	wfspecn(k,i)=wfspec1n(k,i)
227	wftimen(k,i)=wftime1n(k,i)
228	end do
229	end do
230	else
231	do i=1,numbwf
232	wfname(numbwf-i+1)=wfname1(i)
233	wfspec(numbwf-i+1)=wfspec1(i)
234	wftime(numbwf-i+1)=wftime1(i)
235	end do
236	do k=1,numbnests
237	do i=1,numbwfn(k)
238	wfnamen(k,numbwfn(k)-i+1)=wfname1n(k,i)
239	wfspecn(k,numbwfn(k)-i+1)=wfspec1n(k,i)
240	wftimen(k,numbwfn(k)-i+1)=wftime1n(k,i)
241	end do
242	end do
243	endif
244
245	! Check the time difference between the wind fields. If it is big,
246	! write a warning message. If it is too big, terminate the trajectory.
247	!*********************************************************************
248
249	do i=2,numbwf
250	idiff=abs(wftime(i)-wftime(i-1))
251	if (idiff.gt.idiffmax) then
252	write(,) 'FLEXPART WARNING: TIME DIFFERENCE BETWEEN TWO'
253	write(,) 'WIND FIELDS IS TOO BIG FOR TRANSPORT CALCULATION.&
254	&'
255	write(,) 'THEREFORE, TRAJECTORIES HAVE TO BE SKIPPED.'
256	else if (idiff.gt.idiffnorm) then
257	print *,idiff,idiffnorm
258	write(,) 'FLEXPART WARNING: TIME DIFFERENCE BETWEEN TWO'
259	write(,) 'WIND FIELDS IS BIG. THIS MAY CAUSE A DEGRADATION'
260	write(,) 'OF SIMULATION QUALITY.'
261	endif
262	end do
263
264	do k=1,numbnests
265	if (numbwfn(k).ne.numbwf) then
266	write(,) 'FLEXPART ERROR: THE AVAILABLE FILES FOR THE'
267	write(,) 'NESTED WIND FIELDS ARE NOT CONSISTENT WITH'
268	write(,) 'THE AVAILABLE FILE OF THE MOTHER DOMAIN. '
269	write(,) 'ERROR AT NEST LEVEL: ',k
270	stop
271	endif
272	do i=1,numbwf
273	if (wftimen(k,i).ne.wftime(i)) then
274	write(,) 'FLEXPART ERROR: THE AVAILABLE FILES FOR THE'
275	write(,) 'NESTED WIND FIELDS ARE NOT CONSISTENT WITH'
276	write(,) 'THE AVAILABLE FILE OF THE MOTHER DOMAIN. '
277	write(,) 'ERROR AT NEST LEVEL: ',k
278	stop
279	endif
280	end do
281	end do
282
283	! Reset the times of the wind fields that are kept in memory to no time
284	!**********************************************************************
285
286	do i=1,2
287	memind(i)=i
288	memtime(i)=999999999
289	end do
290
291	return
292
293	998 write(,) ' #### FLEXPART MODEL ERROR! FILE #### '
294	write(,'(a)') ' '//path(numpath+2(k-1)+2) &
295	(1:length(numpath+2*(k-1)+2))
296	write(,) ' #### CANNOT BE OPENED #### '
297	stop
298
299	999 write(,) ' #### FLEXPART MODEL ERROR! FILE #### '
300	write(*,'(a)') ' '//path(4)(1:length(4))
301	write(,) ' #### CANNOT BE OPENED #### '
302	stop
303
304	end subroutine readavailable

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

Download in other formats:

Original Format