Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

ungrib.F90 @ 87

Last change on this file since 87 was 11, checked in by aslmd, 14 years ago
spiga@svn-planeto:ajoute le modele meso-echelle martien
File size: 14.7 KB

Line
1	!*****************************************************************************!
2	! program ungrib !
3	! !
4	! Questions, comments, suggestions, even complaints should be directed to: !
5	! wrfhelp@ucar.edu !
6	! Externals: !
7	! Module TABLE !
8	! Module GRIDINFO !
9	! Module FILELIST !
10	! Subroutine READ_NAMELIST !
11	! Subroutine PARSE_TABLE !
12	! Subroutine CLEAR_STORAGE !
13	! Subroutine RD_GRIB !
14	! Subroutine RD_GRIB2 !
15	! Subroutine PUT_STORAGE !
16	! Subroutine OUTPUT !
17	! Subroutine CCLOSE !
18	! Subroutine RRPR !
19	! Subroutine DATINT !
20	! Subroutine FILE_DELETE !
21	! !
22	! Kevin W. Manning, NCAR/MMM - original 'pregrid' code, 1998-2001 !
23	! Jim Bresch, Michael Duda, Dave Gill, NCAR/MMM - adapted for WPS, 2006 !
24	! !
25	!*****************************************************************************!
26	! !
27	!*****************************************************************************!
28	! !
29	! This program reads GRIB-formatted data and puts it into intermediate format !
30	! for passing data to a horizontal interpolation program. The intermediate !
31	! format can be for WPS, SI, or MM5. !
32	! !
33	! The program tries to read from files called "GRIBFILE.AAA", "GRIBFILE.AAB", !
34	! "GRIBFILE.AAC", ... "GRIBFILE.ABA", "GRIBFILE.ABB", ... "GRIBFILE.ZZZ" until!
35	! it cannot find a file. This naming format allows for up to 17576 files, !
36	! which should be enough for most applications. !
37	! !
38	! The program searches through those "GRIBFILE.???" files, and pulls out all !
39	! the requested fields which fall between a starting date and an ending date. !
40	! It writes the fields from a given time period to a file named according to !
41	! the date and hour, i.e., "FILE:YYYY-MM-DD_HH" !
42	! !
43	!*****************************************************************************!
44	program ungrib
45
46	use table
47	use gridinfo
48	use storage_module
49	use filelist
50	use datarray
51	use module_debug
52	use stringutil
53
54	implicit none
55
56	integer :: nunit1 = 12
57	character(LEN=132) :: gribflnm = 'GRIBFILE.AAA ' ! won't work with len=12
58
59	integer :: debug_level
60
61	integer , parameter :: maxlvl = 100
62
63	real , dimension(maxlvl) :: plvl
64	integer :: iplvl
65
66	integer :: nlvl
67
68	real :: startlat, startlon, deltalat, deltalon
69	real :: level
70	character (LEN=9) :: field
71	character (LEN=3) :: out_format
72	character (LEN=256) :: prefix
73
74	logical :: readit
75
76	integer, dimension(255) :: iuarr = 0
77
78	character (LEN=19) :: HSTART, HEND, HDATE
79	character(LEN=19) :: hsave = '0000-00-00_00:00:00'
80	integer :: itime
81	integer :: ntimes
82	integer :: interval
83	integer :: ierr
84	logical :: ordered_by_date
85	integer :: grib_version
86	integer :: vtable_columns
87
88
89	call mprintf(.true.,STDOUT,' * Starting program ungrib.exe * ')
90	call mprintf(.true.,LOGFILE,' * Starting program ungrib.exe * ')
91	! -----------------
92	! Read the namelist, and return the information we want:
93
94	call read_namelist(hstart, hend, interval, ntimes, &
95	ordered_by_date, debug_level, out_format, prefix)
96
97	call mprintf(.true.,INFORM,"Interval value: %i seconds or %f hours", &
98	i1=interval, f1=float(interval)/3600.)
99
100	call mprintf(.true.,STDOUT,'Path to intermediate files is %s',s1=get_path(prefix))
101	call mprintf(.true.,LOGFILE,'Path to intermediate files is %s',s1=get_path(prefix))
102
103	! -----------------
104	! Determine GRIB Edition number
105	grib_version=0
106	call edition_num(nunit1, gribflnm, grib_version, ierr)
107	call mprintf((ierr.eq.3),ERROR,"GRIB file problem")
108	if (grib_version.eq.2) then
109	vtable_columns=11
110	#if defined (USE_PNG) && (USE_JPEG2000)
111	call mprintf(.true.,INFORM, &
112	"Linked in png and jpeg libraries for Grib Edition 2")
113	#else
114	call mprintf(.true.,STDOUT,"WARNING - Grib Edition 2 data detected, and")
115	call mprintf(.true.,STDOUT," - png and jpeg libs were NOT selected")
116	call mprintf(.true.,STDOUT," - during the build.")
117	call mprintf(.true.,STDOUT,"Stopping")
118	call mprintf(.true.,LOGFILE,"WARNING - Grib Edition 2 data detected, and")
119	call mprintf(.true.,LOGFILE," - png and jpeg libs were NOT selected")
120	call mprintf(.true.,LOGFILE," - during the build.")
121	call mprintf(.true.,LOGFILE,"Stopping")
122	call mprintf(.true.,ERROR,"NEED_GRIB2_LIBS")
123	#endif
124	else
125	vtable_columns=7
126	endif
127	call mprintf(.true.,INFORM,"Reading Grib Edition %i", i1=grib_version)
128
129	! -----------------
130	! Read the "Vtable" file, and put the information contained therein into
131	! the module "table".
132
133	call parse_table(debug_level,vtable_columns)
134
135	call mprintf(.true.,DEBUG,"Parsed the vtable.")
136
137	! -----------------
138	! Initialize the input filename to GRIBFILE.AA{character just before A}
139	! That way, when we update the filename below for the first time, it will
140	! have the correct value of GRIBFILE.AAA.
141
142	gribflnm(12:12)=char(ichar(gribflnm(12:12))-1)
143
144	! -----------------
145	! LOOP2 cycles through the list of files named GRIBFILE.???, until it finds
146	! a non-existent file. Then we exit
147
148	LOOP2 : DO
149
150	! At the beginning of LOOP2 update the input filename.
151	if (gribflnm(12:12).eq.'Z') then
152	if (gribflnm(11:11).eq.'Z') then
153	gribflnm(10:10) = char(ichar(gribflnm(10:10))+1)
154	gribflnm(11:11) = 'A'
155	else
156	gribflnm(11:11) = char(ichar(gribflnm(11:11))+1)
157	endif
158	gribflnm(12:12) = 'A'
159	else
160	gribflnm(12:12) = char(ichar(gribflnm(12:12))+1)
161	endif
162
163	! Set READIT to .TRUE., meaning that we have not read any records yet
164	! from the file GRIBFLNM.
165
166	call mprintf(.true.,DEBUG,"Reading from gribflnm %s ",s1=gribflnm)
167
168	readit = .TRUE. ! i.e., "Yes, we want to read a record."
169
170
171	! LOOP1 reads through the file GRIBFLNM, exiting under two conditions:
172	! 1) We have hit the end-of-file
173	! 2) We have read past the ending date HEND.
174	!
175	! Condition 2 assumes that the data in file GRIBFLNM are ordered in time.
176
177	LOOP1 : DO
178	! At the beginning of LOOP1, we are at a new time period.
179	! Clear the storage arrays and associated level information.
180	nlvl = 0
181	plvl = -999.
182	call clear_storage
183
184	! LOOP0 reads through the file GRIBFLNM, looking for data of the current
185	! date. It exits under the following conditions.
186	! 1) We have hit the end-of-file
187	! 2) The GRIBFLNM variable has been updated to a nonexistent file.
188	! 3) We start reading a new date and the data are assumed to be
189	! ordered by date.
190	! 4) We have a valid record and the data are not assumed to be
191	! ordered by date.
192
193	LOOP0 : DO
194
195	! If we need to read a new grib record, then read one.
196	if (READIT) then
197
198	if (grib_version.ne.2) then
199	call mprintf(.true.,DEBUG, &
200	"Calling rd_grib1 with iunit %i", i1=nunit1)
201	call mprintf(.true.,DEBUG, &
202	"flnm = %s",s1=gribflnm)
203	! Read one record at a time from GRIB1 (and older Editions)
204	call rd_grib1(nunit1, gribflnm, level, field, &
205	hdate, ierr, iuarr, debug_level)
206	else
207
208	! Read one file of records from GRIB2.
209	call mprintf(.true.,DEBUG,"Calling rd_grib2")
210	call rd_grib2(nunit1, gribflnm, hdate, &
211	grib_version, ierr, debug_level)
212	FIELD='NULL'
213
214	endif
215
216	call mprintf(.true.,DEBUG,"ierr = %i ",i1=ierr)
217	if (ierr.eq.1) then
218	! We have hit the end of a file. Exit LOOP0 so we can
219	! write output for date HDATE, and then exit LOOP1
220	! to update the GRIBFLNM
221	hsave = hdate
222	exit LOOP0
223	endif
224
225	if (ierr.eq.2) then
226	! We have hit the end of all the files. We can exit LOOP2
227	! because there are no more files to read.
228	exit LOOP2
229	endif
230
231	call mprintf(.true.,DEBUG, &
232	"Read a record %s with date %s", s1=field,s2=hdate(1:13))
233
234	endif
235
236	call mprintf(.true.,DEBUG, &
237	"hdate = %s , hsave = %s ",s1=hdate(1:13), s2=hsave(1:13) )
238
239	! if (hdate < hstart) then
240	! ! The data read has a date HDATE earlier than the starting
241	! ! date HSTART. So cycle LOOP0 to read the the next GRIB record.
242	! READIT = .TRUE.
243	! cycle LOOP0
244	! endif
245
246	if (FIELD.EQ.'NULL') then
247	! The data read does not match any fields requested
248	! in the Vtable. So cycle LOOP0 to read the next GRIB record.
249	READIT = .TRUE.
250	cycle LOOP0
251	endif
252
253	if (ordered_by_date .and. (hdate > hsave)) then
254
255	! Exit LOOP0, because we started to read data from another
256	! date.
257
258	call mprintf(.true.,DEBUG, &
259	"hdate %s > hsave %s so exit LOOP0",s1=hdate,s2=hsave)
260
261	! We set READIT to FALSE because we have not yet processed
262	! the data from this record, and we will want to process this
263	! data on the next pass of LOOP1 (referring to the next time
264	! period of interest.
265
266	READIT = .FALSE.
267
268	exit LOOP0
269
270	endif
271
272	! When we have reached this point, we have a data array ARRAY which has
273	! some data we want to save, with field name FIELD at pressure level
274	! LEVEL (Pa). Dimensions of this data are map%nx and map%ny. Put
275	! this data into storage.
276
277	if (((field == "SST").or.(field == "SKINTEMP")) .and. &
278	(level /= 200100.)) level = 200100.
279	iplvl = int(level)
280	call mprintf((.not.allocated(rdatarray)),ERROR, &
281	"GRIB data slab not allocated in ungrib.F before call to put_storage.")
282	call put_storage(iplvl,field, &
283	reshape(rdatarray(1:map%nx*map%ny),(/map%nx, map%ny/)),&
284	map%nx,map%ny)
285	deallocate(rdatarray)
286
287	! Since we processed the record we just read, we set
288	! READIT to .TRUE. so that LOOP0 will read the next record.
289	READIT = .TRUE.
290
291	if (.not. ordered_by_date) then
292	if (hdate >= hstart) then
293	hsave = hdate
294	endif
295	exit LOOP0
296	endif
297
298	enddo LOOP0
299
300	! When we have reached this point, we have either hit the end of file, or
301	! hit the end of the current date. Either way, we want to output
302	! the data found for this date. This current date is in HSAVE.
303	! However, if (HSAVE == 0000-00-00_00:00:00), no output is necessary,
304	! because that 0000 date is the flag for the very opening of a file.
305
306	if ((hsave(1:4).ne.'0000').and.(hsave.le.hend)) then
307	if (debug_level .gt. 100) print*, 'Calling output: '//hsave(1:13)
308	call output(hsave, nlvl, maxlvl, plvl, interval, 1, out_format, prefix, debug_level)
309	hsave=hdate
310
311	! If the next record we process has a date later than HEND,
312	! then time to exit LOOP1.
313	if ((ordered_by_date) .and. (hdate.gt.hend)) exit LOOP1
314
315	else
316	hsave = hdate
317	endif
318
319	! If we hit the end-of-file, its time to exit LOOP1 so we can
320	! increment the GRIBFLNM to read the next file.
321	if (ierr.eq.1) exit LOOP1
322
323	enddo LOOP1
324
325	! When we have reached this point, we read all the data we want to from
326	! file GRIBFLNM (either because we reached the end-of-file, or we
327	! read past the date HEND). So we close the file and cycle LOOP2 to read
328	! the next file.
329
330	if (grib_version.ne.2) then
331	call cclose(iuarr(nunit1), debug_level, ierr)
332	iuarr(nunit1) = 0
333	endif
334	hsave = '0000-00-00_00:00:00'
335
336	enddo LOOP2
337
338	! Now Reread, process, and reoutput.
339
340	call mprintf(.true.,INFORM,"First pass done, doing a reprocess")
341
342	call rrpr(hstart, ntimes, interval, nlvl, maxlvl, plvl, debug_level, out_format, prefix)
343
344	! Make sure the filedates are in order, with an inefficient sort:
345
346	call sort_filedates
347
348	! Interpolate temporally to fill in missing times:
349
350	call datint(filedates, nfiles, hstart, ntimes, interval, out_format, prefix)
351
352	! Now delete the temporary files:
353
354	call file_delete(filedates, nfiles, trim(get_path(prefix))//'PFILE:', interval)
355
356	! And Now we are done:
357
358	call mprintf(.true.,STDOUT,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!')
359	call mprintf(.true.,STDOUT,'! Successful completion of ungrib. !')
360	! call mprintf(.true.,STDOUT,"! We're hauling gear at Bandimere. !")
361	call mprintf(.true.,STDOUT,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!')
362
363	call mprintf(.true.,LOGFILE,' * Successful completion of program ungrib.exe * ')
364
365
366
367	contains
368	subroutine sort_filedates
369	implicit none
370
371	integer :: n
372	logical :: done
373	if (nfiles > 1) then
374	done = .FALSE.
375	do while ( .not. done)
376	done = .TRUE.
377	do n = 1, nfiles-1
378	if (filedates(n) > filedates(n+1)) then
379	filedates(size(filedates)) = filedates(n)
380	filedates(n) = filedates(n+1)
381	filedates(n+1) = filedates(size(filedates))
382	filedates(size(filedates)) = '0000-00-00 00:00:00.0000'
383	done = .FALSE.
384	endif
385	enddo
386	enddo
387	endif
388	end subroutine sort_filedates
389
390	end program ungrib

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

Download in other formats:

Original Format