Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

getfield.F @ 3567

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

Line
1	subroutine getfield(cgrib,lcgrib,ifldnum,igds,igdstmpl,igdslen,
2	& ideflist,idefnum,ipdsnum,ipdstmpl,ipdslen,
3	& coordlist,numcoord,ndpts,idrsnum,idrstmpl,
4	& idrslen,ibmap,bmap,fld,ierr)
5	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
6	! . . . .
7	! SUBPROGRAM: getfield
8	! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-26
9	!
10	! ABSTRACT: This subroutine returns the Grid Definition, Product Definition,
11	! Bit-map ( if applicable ), and the unpacked data for a given data
12	! field. Since there can be multiple data fields packed into a GRIB2
13	! message, the calling routine indicates which field is being requested
14	! with the ifldnum argument.
15	!
16	! PROGRAM HISTORY LOG:
17	! 2000-05-26 Gilbert
18	!
19	! USAGE: CALL getfield(cgrib,lcgrib,ifldnum,igds,igdstmpl,igdslen,
20	! & ideflist,idefnum,ipdsnum,ipdstmpl,ipdslen,
21	! & coordlist,numcoord,ndpts,idrsnum,idrstmpl,
22	! & idrslen,ibmap,bmap,fld,ierr)
23	! INPUT ARGUMENT LIST:
24	! cgrib - Character array that contains the GRIB2 message
25	! lcgrib - Length (in bytes) of GRIB message array cgrib.
26	! ifldnum - Specifies which field in the GRIB2 message to return.
27	!
28	! OUTPUT ARGUMENT LIST:
29	! igds - Contains information read from the appropriate GRIB Grid
30	! Definition Section 3 for the field being returned.
31	! Must be dimensioned >= 5.
32	! igds(1)=Source of grid definition (see Code Table 3.0)
33	! igds(2)=Number of grid points in the defined grid.
34	! igds(3)=Number of octets needed for each
35	! additional grid points definition.
36	! Used to define number of
37	! points in each row ( or column ) for
38	! non-regular grids.
39	! = 0, if using regular grid.
40	! igds(4)=Interpretation of list for optional points
41	! definition. (Code Table 3.11)
42	! igds(5)=Grid Definition Template Number (Code Table 3.1)
43	! igdstmpl - Contains the data values for the specified Grid Definition
44	! Template ( NN=igds(5) ). Each element of this integer
45	! array contains an entry (in the order specified) of Grid
46	! Defintion Template 3.NN
47	! A safe dimension for this array can be obtained in advance
48	! from maxvals(2), which is returned from subroutine gribinfo.
49	! igdslen - Number of elements in igdstmpl(). i.e. number of entries
50	! in Grid Defintion Template 3.NN ( NN=igds(5) ).
51	! ideflist - (Used if igds(3) .ne. 0) This array contains the
52	! number of grid points contained in each row ( or column ).
53	! (part of Section 3)
54	! A safe dimension for this array can be obtained in advance
55	! from maxvals(3), which is returned from subroutine gribinfo.
56	! idefnum - (Used if igds(3) .ne. 0) The number of entries
57	! in array ideflist. i.e. number of rows ( or columns )
58	! for which optional grid points are defined.
59	! ipdsnum - Product Definition Template Number ( see Code Table 4.0)
60	! ipdstmpl - Contains the data values for the specified Product Definition
61	! Template ( N=ipdsnum ). Each element of this integer
62	! array contains an entry (in the order specified) of Product
63	! Defintion Template 4.N
64	! A safe dimension for this array can be obtained in advance
65	! from maxvals(4), which is returned from subroutine gribinfo.
66	! ipdslen - Number of elements in ipdstmpl(). i.e. number of entries
67	! in Product Defintion Template 4.N ( N=ipdsnum ).
68	! coordlist- Array containg floating point values intended to document
69	! the vertical discretisation associated to model data
70	! on hybrid coordinate vertical levels. (part of Section 4)
71	! The dimension of this array can be obtained in advance
72	! from maxvals(5), which is returned from subroutine gribinfo.
73	! numcoord - number of values in array coordlist.
74	! ndpts - Number of data points unpacked and returned.
75	! idrsnum - Data Representation Template Number ( see Code Table 5.0)
76	! idrstmpl - Contains the data values for the specified Data Representation
77	! Template ( N=idrsnum ). Each element of this integer
78	! array contains an entry (in the order specified) of Product
79	! Defintion Template 5.N
80	! A safe dimension for this array can be obtained in advance
81	! from maxvals(6), which is returned from subroutine gribinfo.
82	! idrslen - Number of elements in idrstmpl(). i.e. number of entries
83	! in Data Representation Template 5.N ( N=idrsnum ).
84	! ibmap - Bitmap indicator ( see Code Table 6.0 )
85	! 0 = bitmap applies and is included in Section 6.
86	! 1-253 = Predefined bitmap applies
87	! 254 = Previously defined bitmap applies to this field
88	! 255 = Bit map does not apply to this product.
89	! bmap() - Logical*1 array containing decoded bitmap. ( if ibmap=0 )
90	! The dimension of this array can be obtained in advance
91	! from maxvals(7), which is returned from subroutine gribinfo.
92	! fld() - Array of ndpts unpacked data points.
93	! A safe dimension for this array can be obtained in advance
94	! from maxvals(7), which is returned from subroutine gribinfo.
95	! ierr - Error return code.
96	! 0 = no error
97	! 1 = Beginning characters "GRIB" not found.
98	! 2 = GRIB message is not Edition 2.
99	! 3 = The data field request number was not positive.
100	! 4 = End string "7777" found, but not where expected.
101	! 6 = GRIB message did not contain the requested number of
102	! data fields.
103	! 7 = End string "7777" not found at end of message.
104	! 9 = Data Representation Template 5.NN not yet implemented.
105	! 10 = Error unpacking Section 3.
106	! 11 = Error unpacking Section 4.
107	! 12 = Error unpacking Section 5.
108	! 13 = Error unpacking Section 6.
109	! 14 = Error unpacking Section 7.
110	!
111	! REMARKS: Note that subroutine gribinfo can be used to first determine
112	! how many data fields exist in a given GRIB message.
113	!
114	! ATTRIBUTES:
115	! LANGUAGE: Fortran 90
116	! MACHINE: IBM SP
117	!
118	!$$$
119
120	character(len=1),intent(in) :: cgrib(lcgrib)
121	integer,intent(in) :: lcgrib,ifldnum
122	integer,intent(out) :: igds(),igdstmpl(),ideflist(*)
123	integer,intent(out) :: ipdsnum,ipdstmpl(*)
124	integer,intent(out) :: idrsnum,idrstmpl(*)
125	integer,intent(out) :: ndpts,ibmap,idefnum,numcoord
126	integer,intent(out) :: ierr
127	logical1,intent(out) :: bmap()
128	real,intent(out) :: fld(),coordlist()
129
130	character(len=4),parameter :: grib='GRIB',c7777='7777'
131	character(len=4) :: ctemp
132	integer:: listsec0(2)
133	integer iofst,ibeg,istart
134	integer(4) :: ieee
135	logical have3,have4,have5,have6,have7
136
137	have3=.false.
138	have4=.false.
139	have5=.false.
140	have6=.false.
141	have7=.false.
142	ierr=0
143	numfld=0
144	!
145	! Check for valid request number
146	!
147	if (ifldnum.le.0) then
148	print *,'getfield: Request for field number must be positive.'
149	ierr=3
150	return
151	endif
152	!
153	! Check for beginning of GRIB message in the first 100 bytes
154	!
155	istart=0
156	do j=1,100
157	ctemp=cgrib(j)//cgrib(j+1)//cgrib(j+2)//cgrib(j+3)
158	if (ctemp.eq.grib ) then
159	istart=j
160	exit
161	endif
162	enddo
163	if (istart.eq.0) then
164	print *,'getfield: Beginning characters GRIB not found.'
165	ierr=1
166	return
167	endif
168	!
169	! Unpack Section 0 - Indicator Section
170	!
171	iofst=8*(istart+5)
172	call g2lib_gbyte(cgrib,listsec0(1),iofst,8) ! Discipline
173	iofst=iofst+8
174	call g2lib_gbyte(cgrib,listsec0(2),iofst,8) ! GRIB edition number
175	iofst=iofst+8
176	iofst=iofst+32
177	call g2lib_gbyte(cgrib,lengrib,iofst,32) ! Length of GRIB message
178	iofst=iofst+32
179	lensec0=16
180	ipos=istart+lensec0
181	!
182	! Currently handles only GRIB Edition 2.
183	!
184	if (listsec0(2).ne.2) then
185	print *,'getfield: can only decode GRIB edition 2.'
186	ierr=2
187	return
188	endif
189	!
190	! Loop through the remaining sections keeping track of the
191	! length of each. Also keep the latest Grid Definition Section info.
192	! Unpack the requested field number.
193	!
194	do
195	! Check to see if we are at end of GRIB message
196	ctemp=cgrib(ipos)//cgrib(ipos+1)//cgrib(ipos+2)//cgrib(ipos+3)
197	if (ctemp.eq.c7777 ) then
198	ipos=ipos+4
199	! If end of GRIB message not where expected, issue error
200	if (ipos.ne.(istart+lengrib)) then
201	print *,'getfield: "7777" found, but not where expected.'
202	ierr=4
203	return
204	endif
205	exit
206	endif
207	! Get length of Section and Section number
208	iofst=(ipos-1)*8
209	call g2lib_gbyte(cgrib,lensec,iofst,32) ! Get Length of Section
210	iofst=iofst+32
211	call g2lib_gbyte(cgrib,isecnum,iofst,8) ! Get Section number
212	iofst=iofst+8
213	!print *,' lensec= ',lensec,' secnum= ',isecnum
214	!
215	! If found Section 3, unpack the GDS info using the
216	! appropriate template. Save in case this is the latest
217	! grid before the requested field.
218	!
219	if (isecnum.eq.3) then
220	iofst=iofst-40 ! reset offset to beginning of section
221	call unpack3(cgrib,lcgrib,iofst,igds,igdstmpl,igdslen,
222	& ideflist,idefnum,jerr)
223	if (jerr.eq.0) then
224	have3=.true.
225	else
226	ierr=10
227	return
228	endif
229	endif
230	!
231	! If found Section 4, check to see if this field is the
232	! one requested.
233	!
234	if (isecnum.eq.4) then
235	numfld=numfld+1
236	if (numfld.eq.ifldnum) then
237	iofst=iofst-40 ! reset offset to beginning of section
238	call unpack4(cgrib,lcgrib,iofst,ipdsnum,ipdstmpl,ipdslen,
239	& coordlist,numcoord,jerr)
240	if (jerr.eq.0) then
241	have4=.true.
242	else
243	ierr=11
244	return
245	endif
246	endif
247	endif
248	!
249	! If found Section 5, check to see if this field is the
250	! one requested.
251	!
252	if ((isecnum.eq.5).and.(numfld.eq.ifldnum)) then
253	iofst=iofst-40 ! reset offset to beginning of section
254	call unpack5(cgrib,lcgrib,iofst,ndpts,idrsnum,idrstmpl,
255	& idrslen,jerr)
256	if (jerr.eq.0) then
257	have5=.true.
258	else
259	ierr=12
260	return
261	endif
262	endif
263	!
264	! If found Section 6, Unpack bitmap.
265	! Save in case this is the latest
266	! bitmap before the requested field.
267	!
268	if (isecnum.eq.6) then
269	iofst=iofst-40 ! reset offset to beginning of section
270	call unpack6(cgrib,lcgrib,iofst,igds(2),ibmap,bmap,jerr)
271	if (jerr.eq.0) then
272	have6=.true.
273	else
274	ierr=13
275	return
276	endif
277	endif
278	!
279	! If found Section 7, check to see if this field is the
280	! one requested.
281	!
282	if ((isecnum.eq.7).and.(numfld.eq.ifldnum)) then
283	if (idrsnum.eq.0) then
284	call simunpack(cgrib(ipos+5),lensec-6,idrstmpl,ndpts,fld)
285	have7=.true.
286	elseif (idrsnum.eq.2.or.idrsnum.eq.3) then
287	call comunpack(cgrib(ipos+5),lensec-6,lensec,idrsnum,
288	& idrstmpl,ndpts,fld,ier)
289	if ( ier .ne. 0 ) then
290	ierr=14
291	return
292	endif
293	have7=.true.
294	elseif (idrsnum.eq.50) then
295	call simunpack(cgrib(ipos+5),lensec-6,idrstmpl,ndpts-1,
296	& fld(2))
297	ieee=idrstmpl(5)
298	call rdieee(ieee,fld(1),1)
299	have7=.true.
300	else
301	print *,'getfield: Data Representation Template ',idrsnum,
302	& ' not yet implemented.'
303	ierr=9
304	return
305	endif
306	endif
307	!
308	! Check to see if we read pass the end of the GRIB
309	! message and missed the terminator string '7777'.
310	!
311	ipos=ipos+lensec ! Update beginning of section pointer
312	if (ipos.gt.(istart+lengrib)) then
313	print *,'getfield: "7777" not found at end of GRIB message.'
314	ierr=7
315	return
316	endif
317
318	if (have3.and.have4.and.have5.and.have6.and.have7) return
319
320	enddo
321
322	!
323	! If exited from above loop, the end of the GRIB message was reached
324	! before the requested field was found.
325	!
326	print *,'getfield: GRIB message contained ',numlocal,
327	& ' different fields.'
328	print *,'getfield: The request was for the ',ifldnum,
329	& ' field.'
330	ierr=6
331
332	return
333	end
334
335
336	subroutine unpack3(cgrib,lcgrib,iofst,igds,igdstmpl,
337	& mapgridlen,ideflist,idefnum,ierr)
338	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
339	! . . . .
340	! SUBPROGRAM: unpack3
341	! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-26
342	!
343	! ABSTRACT: This subroutine unpacks Section 3 (Grid Definition Section)
344	! starting at octet 6 of that Section.
345	!
346	! PROGRAM HISTORY LOG:
347	! 2000-05-26 Gilbert
348	!
349	! USAGE: CALL unpack3(cgrib,lcgrib,lensec,iofst,igds,igdstmpl,
350	! & mapgridlen,ideflist,idefnum,ierr)
351	! INPUT ARGUMENT LIST:
352	! cgrib - Character array that contains the GRIB2 message
353	! lcgrib - Length (in bytes) of GRIB message array cgrib.
354	! iofst - Bit offset of the beginning of Section 3.
355	!
356	! OUTPUT ARGUMENT LIST:
357	! iofst - Bit offset at the end of Section 3, returned.
358	! igds - Contains information read from the appropriate GRIB Grid
359	! Definition Section 3 for the field being returned.
360	! Must be dimensioned >= 5.
361	! igds(1)=Source of grid definition (see Code Table 3.0)
362	! igds(2)=Number of grid points in the defined grid.
363	! igds(3)=Number of octets needed for each
364	! additional grid points definition.
365	! Used to define number of
366	! points in each row ( or column ) for
367	! non-regular grids.
368	! = 0, if using regular grid.
369	! igds(4)=Interpretation of list for optional points
370	! definition. (Code Table 3.11)
371	! igds(5)=Grid Definition Template Number (Code Table 3.1)
372	! igdstmpl - Contains the data values for the specified Grid Definition
373	! Template ( NN=igds(5) ). Each element of this integer
374	! array contains an entry (in the order specified) of Grid
375	! Defintion Template 3.NN
376	! mapgridlen- Number of elements in igdstmpl(). i.e. number of entries
377	! in Grid Defintion Template 3.NN ( NN=igds(5) ).
378	! ideflist - (Used if igds(3) .ne. 0) This array contains the
379	! number of grid points contained in each row ( or column ).
380	! (part of Section 3)
381	! idefnum - (Used if igds(3) .ne. 0) The number of entries
382	! in array ideflist. i.e. number of rows ( or columns )
383	! for which optional grid points are defined.
384	! ierr - Error return code.
385	! 0 = no error
386	! 5 = "GRIB" message contains an undefined Grid Definition
387	! Template.
388	!
389	! REMARKS: Uses Fortran 90 module gridtemplates.
390	!
391	! ATTRIBUTES:
392	! LANGUAGE: Fortran 90
393	! MACHINE: IBM SP
394	!
395	!$$$
396
397	use gridtemplates
398
399	character(len=1),intent(in) :: cgrib(lcgrib)
400	integer,intent(in) :: lcgrib
401	integer,intent(inout) :: iofst
402	integer,intent(out) :: igds(),igdstmpl(),ideflist(*)
403	integer,intent(out) :: ierr,idefnum
404
405	integer,allocatable :: mapgrid(:)
406	integer :: mapgridlen,ibyttem
407	logical needext
408
409	ierr=0
410
411	call g2lib_gbyte(cgrib,lensec,iofst,32) ! Get Length of Section
412	iofst=iofst+32
413	iofst=iofst+8 ! skip section number
414
415	call g2lib_gbyte(cgrib,igds(1),iofst,8) ! Get source of Grid def.
416	iofst=iofst+8
417	call g2lib_gbyte(cgrib,igds(2),iofst,32) ! Get number of grid pts.
418	iofst=iofst+32
419	call g2lib_gbyte(cgrib,igds(3),iofst,8) ! Get num octets for opt. list
420	iofst=iofst+8
421	call g2lib_gbyte(cgrib,igds(4),iofst,8) ! Get interpret. for opt. list
422	iofst=iofst+8
423	call g2lib_gbyte(cgrib,igds(5),iofst,16) ! Get Grid Def Template num.
424	iofst=iofst+16
425	if (igds(1).eq.0) then
426	! if (igds(1).eq.0.OR.igds(1).eq.255) then ! FOR ECMWF TEST ONLY
427	allocate(mapgrid(lensec))
428	! Get Grid Definition Template
429	call getgridtemplate(igds(5),mapgridlen,mapgrid,needext,
430	& iret)
431	if (iret.ne.0) then
432	ierr=5
433	return
434	endif
435	else
436	! igdstmpl=-1
437	mapgridlen=0
438	needext=.false.
439	endif
440	!
441	! Unpack each value into array igdstmpl from the
442	! the appropriate number of octets, which are specified in
443	! corresponding entries in array mapgrid.
444	!
445	ibyttem=0
446	do i=1,mapgridlen
447	nbits=iabs(mapgrid(i))*8
448	if ( mapgrid(i).ge.0 ) then
449	call g2lib_gbyte(cgrib,igdstmpl(i),iofst,nbits)
450	else
451	call g2lib_gbyte(cgrib,isign,iofst,1)
452	call g2lib_gbyte(cgrib,igdstmpl(i),iofst+1,nbits-1)
453	if (isign.eq.1) igdstmpl(i)=-igdstmpl(i)
454	endif
455	iofst=iofst+nbits
456	ibyttem=ibyttem+iabs(mapgrid(i))
457	enddo
458	!
459	! Check to see if the Grid Definition Template needs to be
460	! extended.
461	! The number of values in a specific template may vary
462	! depending on data specified in the "static" part of the
463	! template.
464	!
465	if ( needext ) then
466	call extgridtemplate(igds(5),igdstmpl,newmapgridlen,mapgrid)
467	! Unpack the rest of the Grid Definition Template
468	do i=mapgridlen+1,newmapgridlen
469	nbits=iabs(mapgrid(i))*8
470	if ( mapgrid(i).ge.0 ) then
471	call g2lib_gbyte(cgrib,igdstmpl(i),iofst,nbits)
472	else
473	call g2lib_gbyte(cgrib,isign,iofst,1)
474	call g2lib_gbyte(cgrib,igdstmpl(i),iofst+1,nbits-1)
475	if (isign.eq.1) igdstmpl(i)=-igdstmpl(i)
476	endif
477	iofst=iofst+nbits
478	ibyttem=ibyttem+iabs(mapgrid(i))
479	enddo
480	mapgridlen=newmapgridlen
481	endif
482	!
483	! Unpack optional list of numbers defining number of points
484	! in each row or column, if included. This is used for non regular
485	! grids.
486	!
487	if ( igds(3).ne.0 ) then
488	nbits=igds(3)*8
489	idefnum=(lensec-14-ibyttem)/igds(3)
490	call g2lib_gbytes(cgrib,ideflist,iofst,nbits,0,idefnum)
491	iofst=iofst+(nbits*idefnum)
492	else
493	idefnum=0
494	endif
495	if( allocated(mapgrid) ) deallocate(mapgrid)
496	return ! End of Section 3 processing
497	end
498
499
500	subroutine unpack4(cgrib,lcgrib,iofst,ipdsnum,ipdstmpl,mappdslen,
501	& coordlist,numcoord,ierr)
502	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
503	! . . . .
504	! SUBPROGRAM: unpack4
505	! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-26
506	!
507	! ABSTRACT: This subroutine unpacks Section 4 (Product Definition Section)
508	! starting at octet 6 of that Section.
509	!
510	! PROGRAM HISTORY LOG:
511	! 2000-05-26 Gilbert
512	!
513	! USAGE: CALL unpack4(cgrib,lcgrib,iofst,ipdsnum,ipdstmpl,mappdslen,
514	! & coordlist,numcoord,ierr)
515	! INPUT ARGUMENT LIST:
516	! cgrib - Character array that contains the GRIB2 message
517	! lcgrib - Length (in bytes) of GRIB message array cgrib.
518	! iofst - Bit offset of the beginning of Section 4.
519	!
520	! OUTPUT ARGUMENT LIST:
521	! iofst - Bit offset of the end of Section 4, returned.
522	! ipdsnum - Product Definition Template Number ( see Code Table 4.0)
523	! ipdstmpl - Contains the data values for the specified Product Definition
524	! Template ( N=ipdsnum ). Each element of this integer
525	! array contains an entry (in the order specified) of Product
526	! Defintion Template 4.N
527	! mappdslen- Number of elements in ipdstmpl(). i.e. number of entries
528	! in Product Defintion Template 4.N ( N=ipdsnum ).
529	! coordlist- Array containg floating point values intended to document
530	! the vertical discretisation associated to model data
531	! on hybrid coordinate vertical levels. (part of Section 4)
532	! numcoord - number of values in array coordlist.
533	! ierr - Error return code.
534	! 0 = no error
535	! 5 = "GRIB" message contains an undefined Product Definition
536	! Template.
537	!
538	! REMARKS: Uses Fortran 90 module pdstemplates.
539	!
540	! ATTRIBUTES:
541	! LANGUAGE: Fortran 90
542	! MACHINE: IBM SP
543	!
544	!$$$
545
546	use pdstemplates
547
548	character(len=1),intent(in) :: cgrib(lcgrib)
549	integer,intent(in) :: lcgrib
550	integer,intent(inout) :: iofst
551	real,intent(out) :: coordlist(*)
552	integer,intent(out) :: ipdsnum,ipdstmpl(*)
553	integer,intent(out) :: ierr,numcoord
554
555	real(4),allocatable :: coordieee(:)
556	integer,allocatable :: mappds(:)
557	integer :: mappdslen
558	logical needext
559
560	ierr=0
561
562	call g2lib_gbyte(cgrib,lensec,iofst,32) ! Get Length of Section
563	iofst=iofst+32
564	iofst=iofst+8 ! skip section number
565	allocate(mappds(lensec))
566
567	call g2lib_gbyte(cgrib,numcoord,iofst,16) ! Get num of coordinate values
568	iofst=iofst+16
569	call g2lib_gbyte(cgrib,ipdsnum,iofst,16) ! Get Prod. Def Template num.
570	iofst=iofst+16
571	! Get Product Definition Template
572	call getpdstemplate(ipdsnum,mappdslen,mappds,needext,iret)
573	if (iret.ne.0) then
574	ierr=5
575	return
576	endif
577	!
578	! Unpack each value into array ipdstmpl from the
579	! the appropriate number of octets, which are specified in
580	! corresponding entries in array mappds.
581	!
582	do i=1,mappdslen
583	nbits=iabs(mappds(i))*8
584	if ( mappds(i).ge.0 ) then
585	call g2lib_gbyte(cgrib,ipdstmpl(i),iofst,nbits)
586	else
587	call g2lib_gbyte(cgrib,isign,iofst,1)
588	call g2lib_gbyte(cgrib,ipdstmpl(i),iofst+1,nbits-1)
589	if (isign.eq.1) ipdstmpl(i)=-ipdstmpl(i)
590	endif
591	iofst=iofst+nbits
592	enddo
593	!
594	! Check to see if the Product Definition Template needs to be
595	! extended.
596	! The number of values in a specific template may vary
597	! depending on data specified in the "static" part of the
598	! template.
599	!
600	if ( needext ) then
601	call extpdstemplate(ipdsnum,ipdstmpl,newmappdslen,mappds)
602	! Unpack the rest of the Product Definition Template
603	do i=mappdslen+1,newmappdslen
604	nbits=iabs(mappds(i))*8
605	if ( mappds(i).ge.0 ) then
606	call g2lib_gbyte(cgrib,ipdstmpl(i),iofst,nbits)
607	else
608	call g2lib_gbyte(cgrib,isign,iofst,1)
609	call g2lib_gbyte(cgrib,ipdstmpl(i),iofst+1,nbits-1)
610	if (isign.eq.1) ipdstmpl(i)=-ipdstmpl(i)
611	endif
612	iofst=iofst+nbits
613	enddo
614	mappdslen=newmappdslen
615	endif
616	!
617	! Get Optional list of vertical coordinate values
618	! after the Product Definition Template, if necessary.
619	!
620	if ( numcoord .ne. 0 ) then
621	allocate (coordieee(numcoord))
622	call g2lib_gbytes(cgrib,coordieee,iofst,32,0,numcoord)
623	call rdieee(coordieee,coordlist,numcoord)
624	deallocate (coordieee)
625	iofst=iofst+(32*numcoord)
626	endif
627	if( allocated(mappds) ) deallocate(mappds)
628	return ! End of Section 4 processing
629	end
630
631
632	subroutine unpack5(cgrib,lcgrib,iofst,ndpts,idrsnum,idrstmpl,
633	& mapdrslen,ierr)
634	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
635	! . . . .
636	! SUBPROGRAM: unpack5
637	! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-26
638	!
639	! ABSTRACT: This subroutine unpacks Section 5 (Data Representation Section)
640	! starting at octet 6 of that Section.
641	!
642	! PROGRAM HISTORY LOG:
643	! 2000-05-26 Gilbert
644	!
645	! USAGE: CALL unpack5(cgrib,lcgrib,iofst,ndpts,idrsnum,idrstmpl,
646	! mapdrslen,ierr)
647	! INPUT ARGUMENT LIST:
648	! cgrib - Character array that contains the GRIB2 message
649	! lcgrib - Length (in bytes) of GRIB message array cgrib.
650	! iofst - Bit offset of the beginning of Section 5.
651	!
652	! OUTPUT ARGUMENT LIST:
653	! iofst - Bit offset at the end of Section 5, returned.
654	! ndpts - Number of data points unpacked and returned.
655	! idrsnum - Data Representation Template Number ( see Code Table 5.0)
656	! idrstmpl - Contains the data values for the specified Data Representation
657	! Template ( N=idrsnum ). Each element of this integer
658	! array contains an entry (in the order specified) of Data
659	! Representation Template 5.N
660	! mapdrslen- Number of elements in idrstmpl(). i.e. number of entries
661	! in Data Representation Template 5.N ( N=idrsnum ).
662	! ierr - Error return code.
663	! 0 = no error
664	! 7 = "GRIB" message contains an undefined Data
665	! Representation Template.
666	!
667	! REMARKS: None
668	!
669	! ATTRIBUTES:
670	! LANGUAGE: Fortran 90
671	! MACHINE: IBM SP
672	!
673	!$$$
674
675	use drstemplates
676
677	character(len=1),intent(in) :: cgrib(lcgrib)
678	integer,intent(in) :: lcgrib
679	integer,intent(inout) :: iofst
680	integer,intent(out) :: ndpts,idrsnum,idrstmpl(*)
681	integer,intent(out) :: ierr
682
683	C integer,allocatable :: mapdrs(:)
684	integer,allocatable :: mapdrs(:)
685	integer :: mapdrslen
686	logical needext
687
688	ierr=0
689
690	call g2lib_gbyte(cgrib,lensec,iofst,32) ! Get Length of Section
691	iofst=iofst+32
692	iofst=iofst+8 ! skip section number
693	allocate(mapdrs(lensec))
694
695	call g2lib_gbyte(cgrib,ndpts,iofst,32) ! Get num of data points
696	iofst=iofst+32
697	call g2lib_gbyte(cgrib,idrsnum,iofst,16) ! Get Data Rep Template Num.
698	iofst=iofst+16
699	! Gen Data Representation Template
700	call getdrstemplate(idrsnum,mapdrslen,mapdrs,needext,iret)
701	if (iret.ne.0) then
702	ierr=7
703	return
704	endif
705	!
706	! Unpack each value into array ipdstmpl from the
707	! the appropriate number of octets, which are specified in
708	! corresponding entries in array mappds.
709	!
710	do i=1,mapdrslen
711	nbits=iabs(mapdrs(i))*8
712	if ( mapdrs(i).ge.0 ) then
713	call g2lib_gbyte(cgrib,idrstmpl(i),iofst,nbits)
714	else
715	call g2lib_gbyte(cgrib,isign,iofst,1)
716	call g2lib_gbyte(cgrib,idrstmpl(i),iofst+1,nbits-1)
717	if (isign.eq.1) idrstmpl(i)=-idrstmpl(i)
718	endif
719	iofst=iofst+nbits
720	enddo
721	!
722	! Check to see if the Data Representation Template needs to be
723	! extended.
724	! The number of values in a specific template may vary
725	! depending on data specified in the "static" part of the
726	! template.
727	!
728	if ( needext ) then
729	call extdrstemplate(idrsnum,idrstmpl,newmapdrslen,mapdrs)
730	! Unpack the rest of the Data Representation Template
731	do i=mapdrslen+1,newmapdrslen
732	nbits=iabs(mapdrs(i))*8
733	if ( mapdrs(i).ge.0 ) then
734	call g2lib_gbyte(cgrib,idrstmpl(i),iofst,nbits)
735	else
736	call g2lib_gbyte(cgrib,isign,iofst,1)
737	call g2lib_gbyte(cgrib,idrstmpl(i),iofst+1,nbits-1)
738	if (isign.eq.1) idrstmpl(i)=-idrstmpl(i)
739	endif
740	iofst=iofst+nbits
741	enddo
742	mapdrslen=newmapdrslen
743	endif
744	if( allocated(mapdrs) ) deallocate(mapdrs)
745	return ! End of Section 5 processing
746	end
747
748
749	subroutine unpack6(cgrib,lcgrib,iofst,ngpts,ibmap,bmap,ierr)
750	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
751	! . . . .
752	! SUBPROGRAM: unpack6
753	! PRGMMR: Gilbert ORG: W/NP11 DATE: 2000-05-26
754	!
755	! ABSTRACT: This subroutine unpacks Section 6 (Bit-Map Section)
756	! starting at octet 6 of that Section.
757	!
758	! PROGRAM HISTORY LOG:
759	! 2000-05-26 Gilbert
760	!
761	! USAGE: CALL unpack6(cgrib,lcgrib,iofst,ngpts,ibmap,bmap,ierr)
762	! INPUT ARGUMENT LIST:
763	! cgrib - Character array that contains the GRIB2 message
764	! lcgrib - Length (in bytes) of GRIB message array cgrib.
765	! iofst - Bit offset of the beginning of Section 6.
766	! ngpts - Number of grid points specified in the bit-map
767	!
768	! OUTPUT ARGUMENT LIST:
769	! iofst - Bit offset at the end of Section 6, returned.
770	! ibmap - Bitmap indicator ( see Code Table 6.0 )
771	! 0 = bitmap applies and is included in Section 6.
772	! 1-253 = Predefined bitmap applies
773	! 254 = Previously defined bitmap applies to this field
774	! 255 = Bit map does not apply to this product.
775	! bmap() - Logical*1 array containing decoded bitmap. ( if ibmap=0 )
776	! ierr - Error return code.
777	! 0 = no error
778	! 4 = Unrecognized pre-defined bit-map.
779	!
780	! REMARKS: None
781	!
782	! ATTRIBUTES:
783	! LANGUAGE: Fortran 90
784	! MACHINE: IBM SP
785	!
786	!$$$
787
788	character(len=1),intent(in) :: cgrib(lcgrib)
789	integer,intent(in) :: lcgrib,ngpts
790	integer,intent(inout) :: iofst
791	integer,intent(out) :: ibmap
792	integer,intent(out) :: ierr
793	logical*1,intent(out) :: bmap(ngpts)
794
795	integer :: intbmap(ngpts)
796
797	ierr=0
798
799	iofst=iofst+32 ! skip Length of Section
800	iofst=iofst+8 ! skip section number
801
802	call g2lib_gbyte(cgrib,ibmap,iofst,8) ! Get bit-map indicator
803	iofst=iofst+8
804
805	if (ibmap.eq.0) then ! Unpack bitmap
806	call g2lib_gbytes(cgrib,intbmap,iofst,1,0,ngpts)
807	iofst=iofst+ngpts
808	do j=1,ngpts
809	bmap(j)=.true.
810	if (intbmap(j).eq.0) bmap(j)=.false.
811	enddo
812	elseif (ibmap.eq.254) then ! Use previous bitmap
813	return
814	elseif (ibmap.eq.255) then ! No bitmap in message
815	bmap(1:ngpts)=.true.
816	else
817	print *,'unpack6: Predefined bitmap ',ibmap,' not recognized.'
818	ierr=4
819	endif
820
821	return ! End of Section 6 processing
822	end
823

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

Download in other formats:

Original Format