-
Notifications
You must be signed in to change notification settings - Fork 10
/
Copy pathgetgbh.f
180 lines (179 loc) · 7.36 KB
/
getgbh.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
C> @file
C> @brief Find a grib message.
C> @author Mark Iredell @date 1994-04-01
C> Read a grib index file (or optionally the grib file itself)
C> to get the index buffer (i.e. table of contents) for the grib file.
C> (The index buffer is saved for use by future prospective calls.)
C> Find in the index buffer a reference to the grib message requested.
C> The grib message request specifies the number of messages to skip
C> and the unpacked pds and gds parameters. (A requested parameter
C> of -1 means to allow any value of this parameter to be found.)
C> If the requested grib message is found, then its message number is
C> returned along with the unpacked pds and gds parameters. If the
C> grib message is not found, then the return code will be nonzero.
C>
C> @param[in] lugb integer unit of the unblocked grib data file
C> (only used if lugi=0).
C> @param[in] lugi integer unit of the unblocked grib index file
C> (=0 to get index buffer from the grib file).
C> @param[in] j integer number of messages to skip
C> (=0 to search from beginning)
C> (<0 to read index buffer and skip -1-j messages)
C> @param[in] jpds integer (200) pds parameters for which to search
C> (=-1 for wildcard).
C> - 1: id of center.
C> - 2: generating process id number.
C> - 3: grid definition.
C> - 4: gds/bms flag (right adj copy of octet 8).
C> - 5: indicator of parameter.
C> - 6: type of level.
C> - 7: height/pressure , etc of level.
C> - 8: year including (century-1).
C> - 9: month of year.
C> - 10: day of month.
C> - 11: hour of day.
C> - 12: minute of hour.
C> - 13: indicator of forecast time unit.
C> - 14: time range 1.
C> - 15: time range 2.
C> - 16: time range flag.
C> - 17: number included in average.
C> - 18: version nr of grib specification.
C> - 19: version nr of parameter table.
C> - 20: nr missing from average/accumulation.
C> - 21: century of reference time of data.
C> - 22: units decimal scale factor.
C> - 23: subcenter number.
C> - 24: pds byte 29, for nmc ensemble products.
C> - 128 if forecast field error.
C> - 64 if bias corrected fcst field.
C> - 32 if smoothed field.
C> - warning: can be combination of more than 1.
C> - 25: pds byte 30, not used
C> @param[in] jgds integer (200) gds parameters for which to search
C> (only searched if jpds(3)=255)
C> (=-1 for wildcard).
C> - 1: data representation type.
C> - 19: number of vertical coordinate parameters.
C> - 20: octet number of the list of vertical coordinate parameters or
C> octet number of the list of numbers of points in each row or
C> 255 if neither are present.
C> - 21: for grids with pl, number of points in grid.
C> - 22: number of words in each row.
C> - tu: ngitude grids.
C> - 2: n(i) nr points on latitude circle.
C> - 3: n(j) nr points on longitude meridian.
C> - 4: la(1) latitude of origin.
C> - 5: lo(1) longitude of origin.
C> - 6: resolution flag (right adj copy of octet 17).
C> - 7: la(2) latitude of extreme point.
C> - 8: lo(2) longitude of extreme point.
C> - 9: di longitudinal direction of increment.
C> - 10: dj latitudinal direction increment.
C> - 11: scanning mode flag (right adj copy of octet 28).
C> - Gaussian grids.
C> - 2: n(i) nr points on latitude circle.
C> - 3: n(j) nr points on longitude meridian.
C> - 4: la(1) latitude of origin.
C> - 5: lo(1) longitude of origin.
C> - 6: resolution flag (right adj copy of octet 17).
C> - 7: la(2) latitude of extreme point.
C> - 8: lo(2) longitude of extreme point.
C> - 9: di longitudinal direction of increment.
C> - 10: n - nr of circles pole to equator.
C> - 11: scanning mode flag (right adj copy of octet 28).
C> - 12: nv - nr of vert coord parameters.
C> - 13: pv - octet nr of list of vert coord parameters or
C> pl - location of the list of numbers of points in
C> each row (if no vert coord parameters are present) or
C> 255 if neither are present
C> - Polar stereographic grids.
C> - 2: n(i) nr points along lat circle.
C> - 3: n(j) nr points along lon circle.
C> - 4: la(1) latitude of origin.
C> - 5: lo(1) longitude of origin.
C> - 6: resolution flag (right adj copy of octet 17).
C> - 7: lov grid orientation.
C> - 8: dx - x direction increment.
C> - 9: dy - y direction increment.
C> - 10: projection center flag.
C> - 11: scanning mode (right adj copy of octet 28).
C> - Spherical harmonic coefficients.
C> - 2: j pentagonal resolution parameter.
C> - 3: k pentagonal resolution parameter.
C> - 4: m pentagonal resolution parameter.
C> - 5: representation type.
C> - 6: coefficient storage mode.
C> - Mercator grids.
C> - 2: n(i) nr points on latitude circle.
C> - 3: n(j) nr points on longitude meridian.
C> - 4: la(1) latitude of origin.
C> - 5: lo(1) longitude of origin.
C> - 6: resolution flag (right adj copy of octet 17).
C> - 7: la(2) latitude of last grid point.
C> - 8: lo(2) longitude of last grid point.
C> - 9: latit - latitude of projection intersection.
C> - 10: reserved.
C> - 11: scanning mode flag (right adj copy of octet 28).
C> - 12: longitudinal dir grid length.
C> - 13: latitudinal dir grid length.
C> - Lambert conformal grids.
C> - 2: nx nr points along x-axis.
C> - 3: ny nr points along y-axis.
C> - 4: la1 lat of origin (lower left).
C> - 5: lo1 lon of origin (lower left).
C> - 6: resolution (right adj copy of octet 17).
C> - 7: lov - orientation of grid.
C> - 8: dx - x-dir increment.
C> - 9: dy - y-dir increment.
C> - 10: projection center flag.
C> - 11: scanning mode flag (right adj copy of octet 28).
C> - 12: latin 1 - first lat from pole of secant cone inter.
C> - 13: latin 2 - second lat from pole of secant cone inter.
C> @param[out] kg integer number of bytes in the grib message.
C> @param[out] kf integer number of data points in the message.
C> @param[out] k integer message number unpacked
C> (can be same as j in calling program in order to facilitate multiple searches).
C> @param[out] kpds integer (200) unpacked pds parameters.
C> @param[out] kgds integer (200) unpacked gds parameters.
C> @param[out] iret integer return code.
C> - 0: all ok.
C> - 96: error reading index file.
C> - 99: request not found.
C>
C> @note In order to unpack grib from a multiprocessing environment
C> where each processor is attempting to read from its own pair of
C> logical units, one must directly call subprogram getgbmh as below,
C> allocating a private copy of cbuf, nlen and nnum to each processor.
C> Do not engage the same logical unit from more than one processor.
C>
C> @author Mark Iredell @date 1994-04-01
C-----------------------------------------------------------------------
SUBROUTINE GETGBH(LUGB,LUGI,J,JPDS,JGDS,
& KG,KF,K,KPDS,KGDS,IRET)
INTEGER JPDS(200),JGDS(200)
INTEGER KPDS(200),KGDS(200)
PARAMETER(MBUF=256*1024)
CHARACTER CBUF(MBUF)
SAVE CBUF,NLEN,NNUM,MNUM
DATA LUX/0/
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C DETERMINE WHETHER INDEX BUFFER NEEDS TO BE INITIALIZED
IF(LUGI.GT.0.AND.(J.LT.0.OR.LUGI.NE.LUX)) THEN
LUX=LUGI
JJ=MIN(J,-1-J)
ELSEIF(LUGI.LE.0.AND.(J.LT.0.OR.LUGB.NE.LUX)) THEN
LUX=LUGB
JJ=MIN(J,-1-J)
ELSE
JJ=J
ENDIF
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C FIND AND UNPACK GRIB MESSAGE
CALL GETGBMH(LUGB,LUGI,JJ,JPDS,JGDS,
& MBUF,CBUF,NLEN,NNUM,MNUM,
& KG,KF,K,KPDS,KGDS,IRET)
IF(IRET.EQ.96) LUX=0
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
RETURN
END