grib2io

Introduction

grib2io is a Python package that provides an interface to the NCEP GRIB2 C (g2c) library for the purpose of reading and writing WMO GRIdded Binary, Edition 2 (GRIB2) messages. A physical file can contain one or more GRIB2 messages.

GRIB2 file IO is performed directly in Python. The unpacking/packing of GRIB2 integer, coded metadata and data sections is performed by the g2c library functions via the g2clib Cython wrapper module. The decoding/encoding of GRIB2 metadata is translated into more descriptive, plain language metadata by looking up the integer code values against the appropriate GRIB2 code tables. These code tables are a part of the grib2io module.

Tutorials

The following Jupyter Notebooks are available as tutorials:

 1from ._grib2io import *
 2from ._grib2io import __doc__
 3from ._grib2io import _Grib2Message
 4
 5__all__ = ['open', 'show_config', 'interpolate', 'interpolate_to_stations',
 6           'tables', 'templates', 'utils',
 7           'Grib2Message', '_Grib2Message', 'Grib2GridDef']
 8
 9try:
10    from . import __config__
11    __version__ = __config__.grib2io_version
12except(ImportError):
13    pass
14
15from .g2clib import __version__ as __g2clib_version__
16from .g2clib import _has_jpeg
17from .g2clib import _has_png
18from .g2clib import _has_aec
19
20has_jpeg_support = bool(_has_jpeg)
21has_png_support  = bool(_has_png)
22has_aec_support = bool(_has_aec)
23
24from .tables.originating_centers import _ncep_grib2_table_version
25ncep_grib2_table_version = _ncep_grib2_table_version
26g2c_version = __g2clib_version__
27
28def show_config():
29    """Print grib2io build configuration information."""
30    print(f'grib2io version {__version__} Configuration:\n')
31    print(f'\tg2c library version: {__g2clib_version__}')
32    print(f'\tJPEG compression support: {has_jpeg_support}')
33    print(f'\tPNG compression support: {has_png_support}')
34    print(f'\tAEC compression support: {has_aec_support}')
35    print(f'')
36    print(f'\tNCEP GRIB2 Table Version: {_ncep_grib2_table_version}')
class open:
 61class open():
 62    """
 63    GRIB2 File Object.
 64
 65    A physical file can contain one or more GRIB2 messages.  When instantiated,
 66    class `grib2io.open`, the file named `filename` is opened for reading (`mode
 67    = 'r'`) and is automatically indexed.  The indexing procedure reads some of
 68    the GRIB2 metadata for all GRIB2 Messages.  A GRIB2 Message may contain
 69    submessages whereby Section 2-7 can be repeated.  grib2io accommodates for
 70    this by flattening any GRIB2 submessages into multiple individual messages.
 71
 72    It is important to note that GRIB2 files from some Meteorological agencies
 73    contain other data than GRIB2 messages.  GRIB2 files from ECMWF can contain
 74    GRIB1 and GRIB2 messages.  grib2io checks for these and safely ignores them.
 75
 76    Attributes
 77    ----------
 78    closed : bool
 79        `True` is file handle is close; `False` otherwise.
 80    current_message : int
 81        Current position of the file in units of GRIB2 Messages.
 82    levels : tuple
 83        Tuple containing a unique list of wgrib2-formatted level/layer strings.
 84    messages : int
 85        Count of GRIB2 Messages contained in the file.
 86    mode : str
 87        File IO mode of opening the file.
 88    name : str
 89        Full path name of the GRIB2 file.
 90    size : int
 91        Size of the file in units of bytes.
 92    variables : tuple
 93        Tuple containing a unique list of variable short names (i.e. GRIB2
 94        abbreviation names).
 95    """
 96
 97    __slots__ = ('_fileid', '_filehandle', '_hasindex', '_index', '_nodata',
 98                 '_pos', 'closed', 'current_message', 'messages', 'mode',
 99                 'name', 'size')
100
101    def __init__(self, filename: str, mode: Literal["r", "w", "x"] = "r", **kwargs):
102        """
103        Initialize GRIB2 File object instance.
104
105        Parameters
106        ----------
107        filename
108            File name containing GRIB2 messages.
109        mode: default="r"
110            File access mode where "r" opens the files for reading only; "w"
111            opens the file for overwriting and "x" for writing to a new file.
112        """
113        # Manage keywords
114        if "_xarray_backend" not in kwargs:
115            kwargs["_xarray_backend"] = False
116            self._nodata = False
117        else:
118            self._nodata = kwargs["_xarray_backend"]
119
120        # All write modes are read/write.
121        # All modes are binary.
122        if mode in ("a", "x", "w"):
123            mode += "+"
124        mode = mode + "b"
125
126        # Some GRIB2 files are gzipped, so check for that here, but
127        # raise error when using xarray backend.
128        if 'r' in mode:
129            self._filehandle = builtins.open(filename, mode=mode)
130            # Gzip files contain a 2-byte header b'\x1f\x8b'.
131            if self._filehandle.read(2) == b'\x1f\x8b':
132                self._filehandle.close()
133                if kwargs["_xarray_backend"]:
134                    raise RuntimeError('Gzip GRIB2 files are not supported by the Xarray backend.')
135                import gzip
136                self._filehandle = gzip.open(filename, mode=mode)
137            else:
138                self._filehandle = builtins.open(filename, mode=mode)
139        else:
140            self._filehandle = builtins.open(filename, mode=mode)
141
142        self.name = os.path.abspath(filename)
143        fstat = os.stat(self.name)
144        self._hasindex = False
145        self._index = {}
146        self.mode = mode
147        self.messages = 0
148        self.current_message = 0
149        self.size = fstat.st_size
150        self.closed = self._filehandle.closed
151        self._fileid = hashlib.sha1((self.name+str(fstat.st_ino)+
152                                     str(self.size)).encode('ASCII')).hexdigest()
153        if 'r' in self.mode:
154            self._build_index()
155        # FIX: Cannot perform reads on mode='a'
156        #if 'a' in self.mode and self.size > 0: self._build_index()
157
158
159    def __delete__(self, instance):
160        self.close()
161        del self._index
162
163
164    def __enter__(self):
165        return self
166
167
168    def __exit__(self, atype, value, traceback):
169        self.close()
170
171
172    def __iter__(self):
173        yield from self._index['msg']
174
175
176    def __len__(self):
177        return self.messages
178
179
180    def __repr__(self):
181        strings = []
182        for k in self.__slots__:
183            if k.startswith('_'): continue
184            strings.append('%s = %s\n'%(k,eval('self.'+k)))
185        return ''.join(strings)
186
187
188    def __getitem__(self, key):
189        if isinstance(key,int):
190            if abs(key) >= len(self._index['msg']):
191                raise IndexError("index out of range")
192            else:
193                return self._index['msg'][key]
194        elif isinstance(key,str):
195            return self.select(shortName=key)
196        elif isinstance(key,slice):
197            return self._index['msg'][key]
198        else:
199            raise KeyError('Key must be an integer, slice, or GRIB2 variable shortName.')
200
201
202    def _build_index(self):
203        """Perform indexing of GRIB2 Messages."""
204        # Initialize index dictionary
205        if not self._hasindex:
206            self._index['sectionOffset'] = []
207            self._index['sectionSize'] = []
208            self._index['msgSize'] = []
209            self._index['msgNumber'] = []
210            self._index['msg'] = []
211            self._index['isSubmessage'] = []
212            self._hasindex = True
213
214        # Iterate
215        while True:
216            try:
217                # Initialize
218                bmapflag = None
219                pos = self._filehandle.tell()
220                section2 = b''
221                trailer = b''
222                _secpos = dict.fromkeys(range(8))
223                _secsize = dict.fromkeys(range(8))
224                _isSubmessage = False
225
226                # Ignore headers (usually text) that are not part of the GRIB2
227                # file.  For example, NAVGEM files have a http header at the
228                # beginning that needs to be ignored.
229
230                # Read a byte at a time until "GRIB" is found.  Using
231                # "wgrib2" on a NAVGEM file, the header was 421 bytes and
232                # decided to go to 2048 bytes to be safe. For normal GRIB2
233                # files this should be quick and break out of the first
234                # loop when "test_offset" is 0.
235                for test_offset in range(2048):
236                    self._filehandle.seek(pos + test_offset)
237                    header = struct.unpack(">I", self._filehandle.read(4))[0]
238                    if header.to_bytes(4, "big") == b"GRIB":
239                        pos = pos + test_offset
240                        break
241                else:
242                    # NOTE: Coming here means that no "GRIB" message identifier
243                    # was found in the previous 2048 bytes. So here we continue
244                    # the while True loop.
245                    continue
246
247                # Read the rest of Section 0 using struct.
248                _secpos[0] = self._filehandle.tell()-4
249                _secsize[0] = 16
250                secmsg = self._filehandle.read(12)
251                section0 = np.concatenate(([header],list(struct.unpack('>HBBQ',secmsg))),dtype=np.int64)
252
253                # Make sure message is GRIB2.
254                if section0[3] != 2:
255                    # Check for GRIB1 and ignore.
256                    if secmsg[3] == 1:
257                        warnings.warn("GRIB version 1 message detected.  Ignoring...")
258                        grib1_size = int.from_bytes(secmsg[:3],"big")
259                        self._filehandle.seek(self._filehandle.tell()+grib1_size-16)
260                        continue
261                    else:
262                        raise ValueError("Bad GRIB version number.")
263
264                # Read and unpack sections 1 through 8 which all follow a
265                # pattern of section size, number, and content.
266                while 1:
267                    # Read first 5 bytes of the section which contains the size
268                    # of the section (4 bytes) and the section number (1 byte).
269                    secmsg = self._filehandle.read(5)
270                    secsize, secnum = struct.unpack('>iB',secmsg)
271
272                    # Record the offset of the section number and "append" the
273                    # rest of the section to secmsg.
274                    _secpos[secnum] = self._filehandle.tell()-5
275                    _secsize[secnum] = secsize
276                    if secnum in {1,3,4,5}:
277                        secmsg += self._filehandle.read(secsize-5)
278                    grbpos = 0
279
280                    # Unpack section
281                    if secnum == 1:
282                        # Unpack Section 1
283                        section1, grbpos = g2clib.unpack1(secmsg,grbpos,np.empty)
284                    elif secnum == 2:
285                        # Unpack Section 2
286                        section2 = self._filehandle.read(secsize-5)
287                    elif secnum == 3:
288                        # Unpack Section 3
289                        gds, gdt, deflist, grbpos = g2clib.unpack3(secmsg,grbpos,np.empty)
290                        gds = gds.tolist()
291                        gdt = gdt.tolist()
292                        section3 = np.concatenate((gds,gdt))
293                        section3 = np.where(section3==4294967295,-1,section3)
294                    elif secnum == 4:
295                        # Unpack Section 4
296                        numcoord, pdt, pdtnum, coordlist, grbpos = g2clib.unpack4(secmsg,grbpos,np.empty)
297                        pdt = pdt.tolist()
298                        section4 = np.concatenate((np.array((numcoord,pdtnum)),pdt))
299                    elif secnum == 5:
300                        # Unpack Section 5
301                        drt, drtn, npts, self._pos = g2clib.unpack5(secmsg,grbpos,np.empty)
302                        section5 = np.concatenate((np.array((npts,drtn)),drt))
303                        section5 = np.where(section5==4294967295,-1,section5)
304                    elif secnum == 6:
305                        # Unpack Section 6 - Just the bitmap flag
306                        bmapflag = struct.unpack('>B',self._filehandle.read(1))[0]
307                        if bmapflag == 0:
308                            bmappos = self._filehandle.tell()-6
309                        elif bmapflag == 254:
310                            # Do this to keep the previous position value
311                            pass
312                        else:
313                            bmappos = None
314                        self._filehandle.seek(self._filehandle.tell()+secsize-6)
315                    elif secnum == 7:
316                        # Do not unpack section 7 here, but move the file pointer
317                        # to after section 7.
318                        self._filehandle.seek(self._filehandle.tell()+secsize-5)
319
320                        # Update the file index.
321                        self.messages += 1
322                        self._index['sectionOffset'].append(copy.deepcopy(_secpos))
323                        self._index['sectionSize'].append(copy.deepcopy(_secsize))
324                        self._index['msgSize'].append(section0[-1])
325                        self._index['msgNumber'].append(self.messages)
326                        self._index['isSubmessage'].append(_isSubmessage)
327
328                        # Create Grib2Message with data.
329                        msg = Grib2Message(section0,section1,section2,section3,section4,section5,bmapflag)
330                        msg._msgnum = self.messages-1
331                        msg._deflist = deflist
332                        msg._coordlist = coordlist
333                        if not self._nodata:
334                            msg._ondiskarray = Grib2MessageOnDiskArray((msg.ny,msg.nx), 2,
335                                                                TYPE_OF_VALUES_DTYPE[msg.typeOfValues],
336                                                                self._filehandle,
337                                                                msg, pos, _secpos[6], _secpos[7])
338                        self._index['msg'].append(msg)
339
340                        # If here, then we have moved through GRIB2 section 1-7.
341                        # Now we need to check the next 4 bytes after section 7.
342                        trailer = struct.unpack('>i',self._filehandle.read(4))[0]
343
344                        # If we reach the GRIB2 trailer string ('7777'), then we
345                        # can break begin processing the next GRIB2 message.  If
346                        # not, then we continue within the same iteration to
347                        # process a GRIB2 submessage.
348                        if trailer.to_bytes(4, "big") == b'7777':
349                            break
350                        else:
351                            # If here, trailer should be the size of the first
352                            # section of the next submessage, then the next byte
353                            # is the section number.  Check this value.
354                            nextsec = struct.unpack('>B',self._filehandle.read(1))[0]
355                            if nextsec not in {2,3,4}:
356                                raise ValueError("Bad GRIB2 message structure.")
357                            self._filehandle.seek(self._filehandle.tell()-5)
358                            _isSubmessage = True
359                            continue
360                    else:
361                        raise ValueError("Bad GRIB2 section number.")
362
363            except(struct.error):
364                if 'r' in self.mode:
365                    self._filehandle.seek(0)
366                break
367
368
369    @property
370    def levels(self):
371        if self._hasindex and not self._nodata:
372            return tuple(sorted(set([msg.level for msg in self._index['msg']])))
373        else:
374            return None
375
376
377    @property
378    def variables(self):
379        if self._hasindex and not self._nodata:
380            return tuple(sorted(set([msg.shortName for msg in self._index['msg']])))
381        else:
382            return None
383
384
385    def close(self):
386        """Close the file handle."""
387        if not self._filehandle.closed:
388            self.messages = 0
389            self.current_message = 0
390            self._filehandle.close()
391            self.closed = self._filehandle.closed
392
393
394    def read(self, size: Optional[int]=None):
395        """
396        Read size amount of GRIB2 messages from the current position.
397
398        If no argument is given, then size is None and all messages are returned
399        from the current position in the file. This read method follows the
400        behavior of Python's builtin open() function, but whereas that operates
401        on units of bytes, we operate on units of GRIB2 messages.
402
403        Parameters
404        ----------
405        size: default=None
406            The number of GRIB2 messages to read from the current position. If
407            no argument is give, the default value is None and remainder of
408            the file is read.
409
410        Returns
411        -------
412        read
413            ``Grib2Message`` object when size = 1 or a list of Grib2Messages
414            when size > 1.
415        """
416        if size is not None and size < 0:
417            size = None
418        if size is None or size > 1:
419            start = self.tell()
420            stop = self.messages if size is None else start+size
421            if size is None:
422                self.current_message = self.messages-1
423            else:
424                self.current_message += size
425            return self._index['msg'][slice(start,stop,1)]
426        elif size == 1:
427            self.current_message += 1
428            return self._index['msg'][self.current_message]
429        else:
430            None
431
432
433    def seek(self, pos: int):
434        """
435        Set the position within the file in units of GRIB2 messages.
436
437        Parameters
438        ----------
439        pos
440            The GRIB2 Message number to set the file pointer to.
441        """
442        if self._hasindex:
443            self._filehandle.seek(self._index['sectionOffset'][0][pos])
444            self.current_message = pos
445
446
447    def tell(self):
448        """Returns the position of the file in units of GRIB2 Messages."""
449        return self.current_message
450
451
452    def select(self, **kwargs):
453        """Select GRIB2 messages by `Grib2Message` attributes."""
454        # TODO: Added ability to process multiple values for each keyword (attribute)
455        idxs = []
456        nkeys = len(kwargs.keys())
457        for k,v in kwargs.items():
458            for m in self._index['msg']:
459                if hasattr(m,k) and getattr(m,k) == v: idxs.append(m._msgnum)
460        idxs = np.array(idxs,dtype=np.int32)
461        return [self._index['msg'][i] for i in [ii[0] for ii in collections.Counter(idxs).most_common() if ii[1] == nkeys]]
462
463
464    def write(self, msg):
465        """
466        Writes GRIB2 message object to file.
467
468        Parameters
469        ----------
470        msg
471            GRIB2 message objects to write to file.
472        """
473        if isinstance(msg,list):
474            for m in msg:
475                self.write(m)
476            return
477
478        if issubclass(msg.__class__,_Grib2Message):
479            if hasattr(msg,'_msg'):
480                self._filehandle.write(msg._msg)
481            else:
482                if msg._signature != msg._generate_signature():
483                    msg.pack()
484                    self._filehandle.write(msg._msg)
485                else:
486                    if hasattr(msg._data,'filehandle'):
487                        msg._data.filehandle.seek(msg._data.offset)
488                        self._filehandle.write(msg._data.filehandle.read(msg.section0[-1]))
489                    else:
490                        msg.pack()
491                        self._filehandle.write(msg._msg)
492            self.flush()
493            self.size = os.path.getsize(self.name)
494            self._filehandle.seek(self.size-msg.section0[-1])
495            self._build_index()
496        else:
497            raise TypeError("msg must be a Grib2Message object.")
498        return
499
500
501    def flush(self):
502        """Flush the file object buffer."""
503        self._filehandle.flush()
504
505
506    def levels_by_var(self, name: str):
507        """
508        Return a list of level strings given a variable shortName.
509
510        Parameters
511        ----------
512        name
513            Grib2Message variable shortName
514
515        Returns
516        -------
517        levels_by_var
518            A list of unique level strings.
519        """
520        return list(sorted(set([msg.level for msg in self.select(shortName=name)])))
521
522
523    def vars_by_level(self, level: str):
524        """
525        Return a list of variable shortName strings given a level.
526
527        Parameters
528        ----------
529        level
530            Grib2Message variable level
531
532        Returns
533        -------
534        vars_by_level
535            A list of unique variable shortName strings.
536        """
537        return list(sorted(set([msg.shortName for msg in self.select(level=level)])))

GRIB2 File Object.

A physical file can contain one or more GRIB2 messages. When instantiated, class grib2io.open, the file named filename is opened for reading (mode = 'r') and is automatically indexed. The indexing procedure reads some of the GRIB2 metadata for all GRIB2 Messages. A GRIB2 Message may contain submessages whereby Section 2-7 can be repeated. grib2io accommodates for this by flattening any GRIB2 submessages into multiple individual messages.

It is important to note that GRIB2 files from some Meteorological agencies contain other data than GRIB2 messages. GRIB2 files from ECMWF can contain GRIB1 and GRIB2 messages. grib2io checks for these and safely ignores them.

Attributes
  • closed (bool): True is file handle is close; False otherwise.
  • current_message (int): Current position of the file in units of GRIB2 Messages.
  • levels (tuple): Tuple containing a unique list of wgrib2-formatted level/layer strings.
  • messages (int): Count of GRIB2 Messages contained in the file.
  • mode (str): File IO mode of opening the file.
  • name (str): Full path name of the GRIB2 file.
  • size (int): Size of the file in units of bytes.
  • variables (tuple): Tuple containing a unique list of variable short names (i.e. GRIB2 abbreviation names).
open(filename: str, mode: Literal['r', 'w', 'x'] = 'r', **kwargs)
101    def __init__(self, filename: str, mode: Literal["r", "w", "x"] = "r", **kwargs):
102        """
103        Initialize GRIB2 File object instance.
104
105        Parameters
106        ----------
107        filename
108            File name containing GRIB2 messages.
109        mode: default="r"
110            File access mode where "r" opens the files for reading only; "w"
111            opens the file for overwriting and "x" for writing to a new file.
112        """
113        # Manage keywords
114        if "_xarray_backend" not in kwargs:
115            kwargs["_xarray_backend"] = False
116            self._nodata = False
117        else:
118            self._nodata = kwargs["_xarray_backend"]
119
120        # All write modes are read/write.
121        # All modes are binary.
122        if mode in ("a", "x", "w"):
123            mode += "+"
124        mode = mode + "b"
125
126        # Some GRIB2 files are gzipped, so check for that here, but
127        # raise error when using xarray backend.
128        if 'r' in mode:
129            self._filehandle = builtins.open(filename, mode=mode)
130            # Gzip files contain a 2-byte header b'\x1f\x8b'.
131            if self._filehandle.read(2) == b'\x1f\x8b':
132                self._filehandle.close()
133                if kwargs["_xarray_backend"]:
134                    raise RuntimeError('Gzip GRIB2 files are not supported by the Xarray backend.')
135                import gzip
136                self._filehandle = gzip.open(filename, mode=mode)
137            else:
138                self._filehandle = builtins.open(filename, mode=mode)
139        else:
140            self._filehandle = builtins.open(filename, mode=mode)
141
142        self.name = os.path.abspath(filename)
143        fstat = os.stat(self.name)
144        self._hasindex = False
145        self._index = {}
146        self.mode = mode
147        self.messages = 0
148        self.current_message = 0
149        self.size = fstat.st_size
150        self.closed = self._filehandle.closed
151        self._fileid = hashlib.sha1((self.name+str(fstat.st_ino)+
152                                     str(self.size)).encode('ASCII')).hexdigest()
153        if 'r' in self.mode:
154            self._build_index()
155        # FIX: Cannot perform reads on mode='a'
156        #if 'a' in self.mode and self.size > 0: self._build_index()

Initialize GRIB2 File object instance.

Parameters
  • filename: File name containing GRIB2 messages.
  • mode (default="r"): File access mode where "r" opens the files for reading only; "w" opens the file for overwriting and "x" for writing to a new file.
name
mode
messages
current_message
size
closed
levels
369    @property
370    def levels(self):
371        if self._hasindex and not self._nodata:
372            return tuple(sorted(set([msg.level for msg in self._index['msg']])))
373        else:
374            return None
variables
377    @property
378    def variables(self):
379        if self._hasindex and not self._nodata:
380            return tuple(sorted(set([msg.shortName for msg in self._index['msg']])))
381        else:
382            return None
def close(self):
385    def close(self):
386        """Close the file handle."""
387        if not self._filehandle.closed:
388            self.messages = 0
389            self.current_message = 0
390            self._filehandle.close()
391            self.closed = self._filehandle.closed

Close the file handle.

def read(self, size: Optional[int] = None):
394    def read(self, size: Optional[int]=None):
395        """
396        Read size amount of GRIB2 messages from the current position.
397
398        If no argument is given, then size is None and all messages are returned
399        from the current position in the file. This read method follows the
400        behavior of Python's builtin open() function, but whereas that operates
401        on units of bytes, we operate on units of GRIB2 messages.
402
403        Parameters
404        ----------
405        size: default=None
406            The number of GRIB2 messages to read from the current position. If
407            no argument is give, the default value is None and remainder of
408            the file is read.
409
410        Returns
411        -------
412        read
413            ``Grib2Message`` object when size = 1 or a list of Grib2Messages
414            when size > 1.
415        """
416        if size is not None and size < 0:
417            size = None
418        if size is None or size > 1:
419            start = self.tell()
420            stop = self.messages if size is None else start+size
421            if size is None:
422                self.current_message = self.messages-1
423            else:
424                self.current_message += size
425            return self._index['msg'][slice(start,stop,1)]
426        elif size == 1:
427            self.current_message += 1
428            return self._index['msg'][self.current_message]
429        else:
430            None

Read size amount of GRIB2 messages from the current position.

If no argument is given, then size is None and all messages are returned from the current position in the file. This read method follows the behavior of Python's builtin open() function, but whereas that operates on units of bytes, we operate on units of GRIB2 messages.

Parameters
  • size (default=None): The number of GRIB2 messages to read from the current position. If no argument is give, the default value is None and remainder of the file is read.
Returns
  • read: Grib2Message object when size = 1 or a list of Grib2Messages when size > 1.
def seek(self, pos: int):
433    def seek(self, pos: int):
434        """
435        Set the position within the file in units of GRIB2 messages.
436
437        Parameters
438        ----------
439        pos
440            The GRIB2 Message number to set the file pointer to.
441        """
442        if self._hasindex:
443            self._filehandle.seek(self._index['sectionOffset'][0][pos])
444            self.current_message = pos

Set the position within the file in units of GRIB2 messages.

Parameters
  • pos: The GRIB2 Message number to set the file pointer to.
def tell(self):
447    def tell(self):
448        """Returns the position of the file in units of GRIB2 Messages."""
449        return self.current_message

Returns the position of the file in units of GRIB2 Messages.

def select(self, **kwargs):
452    def select(self, **kwargs):
453        """Select GRIB2 messages by `Grib2Message` attributes."""
454        # TODO: Added ability to process multiple values for each keyword (attribute)
455        idxs = []
456        nkeys = len(kwargs.keys())
457        for k,v in kwargs.items():
458            for m in self._index['msg']:
459                if hasattr(m,k) and getattr(m,k) == v: idxs.append(m._msgnum)
460        idxs = np.array(idxs,dtype=np.int32)
461        return [self._index['msg'][i] for i in [ii[0] for ii in collections.Counter(idxs).most_common() if ii[1] == nkeys]]

Select GRIB2 messages by Grib2Message attributes.

def write(self, msg):
464    def write(self, msg):
465        """
466        Writes GRIB2 message object to file.
467
468        Parameters
469        ----------
470        msg
471            GRIB2 message objects to write to file.
472        """
473        if isinstance(msg,list):
474            for m in msg:
475                self.write(m)
476            return
477
478        if issubclass(msg.__class__,_Grib2Message):
479            if hasattr(msg,'_msg'):
480                self._filehandle.write(msg._msg)
481            else:
482                if msg._signature != msg._generate_signature():
483                    msg.pack()
484                    self._filehandle.write(msg._msg)
485                else:
486                    if hasattr(msg._data,'filehandle'):
487                        msg._data.filehandle.seek(msg._data.offset)
488                        self._filehandle.write(msg._data.filehandle.read(msg.section0[-1]))
489                    else:
490                        msg.pack()
491                        self._filehandle.write(msg._msg)
492            self.flush()
493            self.size = os.path.getsize(self.name)
494            self._filehandle.seek(self.size-msg.section0[-1])
495            self._build_index()
496        else:
497            raise TypeError("msg must be a Grib2Message object.")
498        return

Writes GRIB2 message object to file.

Parameters
  • msg: GRIB2 message objects to write to file.
def flush(self):
501    def flush(self):
502        """Flush the file object buffer."""
503        self._filehandle.flush()

Flush the file object buffer.

def levels_by_var(self, name: str):
506    def levels_by_var(self, name: str):
507        """
508        Return a list of level strings given a variable shortName.
509
510        Parameters
511        ----------
512        name
513            Grib2Message variable shortName
514
515        Returns
516        -------
517        levels_by_var
518            A list of unique level strings.
519        """
520        return list(sorted(set([msg.level for msg in self.select(shortName=name)])))

Return a list of level strings given a variable shortName.

Parameters
  • name: Grib2Message variable shortName
Returns
  • levels_by_var: A list of unique level strings.
def vars_by_level(self, level: str):
523    def vars_by_level(self, level: str):
524        """
525        Return a list of variable shortName strings given a level.
526
527        Parameters
528        ----------
529        level
530            Grib2Message variable level
531
532        Returns
533        -------
534        vars_by_level
535            A list of unique variable shortName strings.
536        """
537        return list(sorted(set([msg.shortName for msg in self.select(level=level)])))

Return a list of variable shortName strings given a level.

Parameters
  • level: Grib2Message variable level
Returns
  • vars_by_level: A list of unique variable shortName strings.
def show_config():
29def show_config():
30    """Print grib2io build configuration information."""
31    print(f'grib2io version {__version__} Configuration:\n')
32    print(f'\tg2c library version: {__g2clib_version__}')
33    print(f'\tJPEG compression support: {has_jpeg_support}')
34    print(f'\tPNG compression support: {has_png_support}')
35    print(f'\tAEC compression support: {has_aec_support}')
36    print(f'')
37    print(f'\tNCEP GRIB2 Table Version: {_ncep_grib2_table_version}')

Print grib2io build configuration information.

def interpolate( a, method: Union[int, str], grid_def_in, grid_def_out, method_options=None, num_threads=1):
1507def interpolate(a, method: Union[int, str], grid_def_in, grid_def_out,
1508                method_options=None, num_threads=1):
1509    """
1510    This is the module-level interpolation function.
1511
1512    This interfaces with the grib2io_interp component package that interfaces to
1513    the [NCEPLIBS-ip library](https://github.com/NOAA-EMC/NCEPLIBS-ip).
1514
1515    Parameters
1516    ----------
1517    a : numpy.ndarray or tuple
1518        Input data.  If `a` is a `numpy.ndarray`, scalar interpolation will be
1519        performed.  If `a` is a `tuple`, then vector interpolation will be
1520        performed with the assumption that u = a[0] and v = a[1] and are both
1521        `numpy.ndarray`.
1522
1523        These data are expected to be in 2-dimensional form with shape (ny, nx)
1524        or 3-dimensional (:, ny, nx) where the 1st dimension represents another
1525        spatial, temporal, or classification (i.e. ensemble members) dimension.
1526        The function will properly flatten the (ny,nx) dimensions into (nx * ny)
1527        acceptable for input into the interpolation subroutines.
1528    method
1529        Interpolate method to use. This can either be an integer or string using
1530        the following mapping:
1531
1532        | Interpolate Scheme | Integer Value |
1533        | :---:              | :---:         |
1534        | 'bilinear'         | 0             |
1535        | 'bicubic'          | 1             |
1536        | 'neighbor'         | 2             |
1537        | 'budget'           | 3             |
1538        | 'spectral'         | 4             |
1539        | 'neighbor-budget'  | 6             |
1540
1541    grid_def_in : grib2io.Grib2GridDef
1542        Grib2GridDef object for the input grid.
1543    grid_def_out : grib2io.Grib2GridDef
1544        Grib2GridDef object for the output grid or station points.
1545    method_options : list of ints, optional
1546        Interpolation options. See the NCEPLIBS-ip documentation for
1547        more information on how these are used.
1548    num_threads : int, optional
1549        Number of OpenMP threads to use for interpolation. The default
1550        value is 1. If grib2io_interp was not built with OpenMP, then
1551        this keyword argument and value will have no impact.
1552
1553    Returns
1554    -------
1555    interpolate
1556        Returns a `numpy.ndarray` when scalar interpolation is performed or a
1557        `tuple` of `numpy.ndarray`s when vector interpolation is performed with
1558        the assumptions that 0-index is the interpolated u and 1-index is the
1559        interpolated v.
1560    """
1561    import grib2io_interp
1562    from grib2io_interp import interpolate
1563
1564    prev_num_threads = 1
1565    try:
1566        import grib2io_interp
1567        if grib2io_interp.has_openmp_support:
1568            prev_num_threads = grib2io_interp.get_openmp_threads()
1569            grib2io_interp.set_openmp_threads(num_threads)
1570    except(AttributeError):
1571        pass
1572
1573    if isinstance(method,int) and method not in _interp_schemes.values():
1574        raise ValueError('Invalid interpolation method.')
1575    elif isinstance(method,str):
1576        if method in _interp_schemes.keys():
1577            method = _interp_schemes[method]
1578        else:
1579            raise ValueError('Invalid interpolation method.')
1580
1581    if method_options is None:
1582        method_options = np.zeros((20),dtype=np.int32)
1583        if method in {3,6}:
1584            method_options[0:2] = -1
1585
1586    ni = grid_def_in.npoints
1587    no = grid_def_out.npoints
1588
1589    # Adjust shape of input array(s)
1590    a,newshp = _adjust_array_shape_for_interp(a,grid_def_in,grid_def_out)
1591
1592    # Set lats and lons if stations, else create array for grids.
1593    if grid_def_out.gdtn == -1:
1594        rlat = np.array(grid_def_out.lats,dtype=np.float32)
1595        rlon = np.array(grid_def_out.lons,dtype=np.float32)
1596    else:
1597        rlat = np.zeros((no),dtype=np.float32)
1598        rlon = np.zeros((no),dtype=np.float32)
1599
1600    # Call interpolation subroutines according to type of a.
1601    if isinstance(a,np.ndarray):
1602        # Scalar
1603        if np.any(np.isnan(a)):
1604            ibi = np.zeros((a.shape[0]),dtype=np.int32)+1
1605            li = np.where(np.isnan(a),0,1).astype(np.int8)
1606        else:
1607            ibi = np.zeros((a.shape[0]),dtype=np.int32)
1608            li = np.zeros(a.shape,dtype=np.int8)
1609        go = np.zeros((a.shape[0],no),dtype=np.float32)
1610        no,ibo,lo,iret = interpolate.interpolate_scalar(method,method_options,
1611                                                 grid_def_in.gdtn,grid_def_in.gdt,
1612                                                 grid_def_out.gdtn,grid_def_out.gdt,
1613                                                 ibi,li.T,a.T,go.T,rlat,rlon)
1614        lo = lo.T.reshape(newshp)
1615        out = go.reshape(newshp)
1616        out = np.where(lo==0,np.nan,out)
1617    elif isinstance(a,tuple):
1618        # Vector
1619        if np.any(np.isnan(a)):
1620            ibi = np.zeros((a.shape[0]),dtype=np.int32)+1
1621            li = np.where(np.isnan(a),0,1).astype(np.int8)
1622        else:
1623            ibi = np.zeros((a.shape[0]),dtype=np.int32)
1624            li = np.zeros(a.shape,dtype=np.int8)
1625        uo = np.zeros((a[0].shape[0],no),dtype=np.float32)
1626        vo = np.zeros((a[1].shape[0],no),dtype=np.float32)
1627        crot = np.ones((no),dtype=np.float32)
1628        srot = np.zeros((no),dtype=np.float32)
1629        no,ibo,lo,iret = interpolate.interpolate_vector(method,method_options,
1630                                                 grid_def_in.gdtn,grid_def_in.gdt,
1631                                                 grid_def_out.gdtn,grid_def_out.gdt,
1632                                                 ibi,li.T,a[0].T,a[1].T,uo.T,vo.T,
1633                                                 rlat,rlon,crot,srot)
1634        del crot
1635        del srot
1636        lo = lo[:,0].reshape(newshp)
1637        uo = uo.reshape(new)
1638        vo = vo.reshape(new)
1639        uo = np.where(lo==0,np.nan,uo)
1640        vo = np.where(lo==0,np.nan,vo)
1641        out = (uo,vo)
1642
1643    del rlat
1644    del rlon
1645
1646    try:
1647        if grib2io_interp.has_openmp_support:
1648            grib2io_interp.set_openmp_threads(prev_num_threads)
1649    except(AttributeError):
1650        pass
1651
1652    return out

This is the module-level interpolation function.

This interfaces with the grib2io_interp component package that interfaces to the NCEPLIBS-ip library.

Parameters
  • a (numpy.ndarray or tuple): Input data. If a is a numpy.ndarray, scalar interpolation will be performed. If a is a tuple, then vector interpolation will be performed with the assumption that u = a[0] and v = a[1] and are both numpy.ndarray.

    These data are expected to be in 2-dimensional form with shape (ny, nx) or 3-dimensional (:, ny, nx) where the 1st dimension represents another spatial, temporal, or classification (i.e. ensemble members) dimension. The function will properly flatten the (ny,nx) dimensions into (nx * ny) acceptable for input into the interpolation subroutines.

  • method: Interpolate method to use. This can either be an integer or string using the following mapping:
Interpolate Scheme Integer Value
'bilinear' 0
'bicubic' 1
'neighbor' 2
'budget' 3
'spectral' 4
'neighbor-budget' 6

  • grid_def_in (grib2io.Grib2GridDef): Grib2GridDef object for the input grid.
  • grid_def_out (grib2io.Grib2GridDef): Grib2GridDef object for the output grid or station points.
  • method_options (list of ints, optional): Interpolation options. See the NCEPLIBS-ip documentation for more information on how these are used.
  • num_threads (int, optional): Number of OpenMP threads to use for interpolation. The default value is 1. If grib2io_interp was not built with OpenMP, then this keyword argument and value will have no impact.
Returns
  • interpolate: Returns a numpy.ndarray when scalar interpolation is performed or a tuple of numpy.ndarrays when vector interpolation is performed with the assumptions that 0-index is the interpolated u and 1-index is the interpolated v.
def interpolate_to_stations( a, method, grid_def_in, lats, lons, method_options=None, num_threads=1):
1655def interpolate_to_stations(a, method, grid_def_in, lats, lons,
1656                            method_options=None, num_threads=1):
1657    """
1658    Module-level interpolation function for interpolation to stations.
1659
1660    Interfaces with the grib2io_interp component package that interfaces to the
1661    [NCEPLIBS-ip
1662    library](https://github.com/NOAA-EMC/NCEPLIBS-ip). It supports scalar and
1663    vector interpolation according to the type of object a.
1664
1665    Parameters
1666    ----------
1667    a : numpy.ndarray or tuple
1668        Input data.  If `a` is a `numpy.ndarray`, scalar interpolation will be
1669        performed.  If `a` is a `tuple`, then vector interpolation will be
1670        performed with the assumption that u = a[0] and v = a[1] and are both
1671        `numpy.ndarray`.
1672
1673        These data are expected to be in 2-dimensional form with shape (ny, nx)
1674        or 3-dimensional (:, ny, nx) where the 1st dimension represents another
1675        spatial, temporal, or classification (i.e. ensemble members) dimension.
1676        The function will properly flatten the (ny,nx) dimensions into (nx * ny)
1677        acceptable for input into the interpolation subroutines.
1678    method
1679        Interpolate method to use. This can either be an integer or string using
1680        the following mapping:
1681
1682        | Interpolate Scheme | Integer Value |
1683        | :---:              | :---:         |
1684        | 'bilinear'         | 0             |
1685        | 'bicubic'          | 1             |
1686        | 'neighbor'         | 2             |
1687        | 'budget'           | 3             |
1688        | 'spectral'         | 4             |
1689        | 'neighbor-budget'  | 6             |
1690
1691    grid_def_in : grib2io.Grib2GridDef
1692        Grib2GridDef object for the input grid.
1693    lats : numpy.ndarray or list
1694        Latitudes for station points
1695    lons : numpy.ndarray or list
1696        Longitudes for station points
1697    method_options : list of ints, optional
1698        Interpolation options. See the NCEPLIBS-ip documentation for
1699        more information on how these are used.
1700    num_threads : int, optional
1701        Number of OpenMP threads to use for interpolation. The default
1702        value is 1. If grib2io_interp was not built with OpenMP, then
1703        this keyword argument and value will have no impact.
1704
1705    Returns
1706    -------
1707    interpolate_to_stations
1708        Returns a `numpy.ndarray` when scalar interpolation is performed or a
1709        `tuple` of `numpy.ndarray`s when vector interpolation is performed with
1710        the assumptions that 0-index is the interpolated u and 1-index is the
1711        interpolated v.
1712    """
1713    import grib2io_interp
1714    from grib2io_interp import interpolate
1715
1716    prev_num_threads = 1
1717    try:
1718        import grib2io_interp
1719        if grib2io_interp.has_openmp_support:
1720            prev_num_threads = grib2io_interp.get_openmp_threads()
1721            grib2io_interp.set_openmp_threads(num_threads)
1722    except(AttributeError):
1723        pass
1724
1725    if isinstance(method,int) and method not in _interp_schemes.values():
1726        raise ValueError('Invalid interpolation method.')
1727    elif isinstance(method,str):
1728        if method in _interp_schemes.keys():
1729            method = _interp_schemes[method]
1730        else:
1731            raise ValueError('Invalid interpolation method.')
1732
1733    if method_options is None:
1734        method_options = np.zeros((20),dtype=np.int32)
1735        if method in {3,6}:
1736            method_options[0:2] = -1
1737
1738    # Check lats and lons
1739    if isinstance(lats,list):
1740        nlats = len(lats)
1741    elif isinstance(lats,np.ndarray) and len(lats.shape) == 1:
1742        nlats = lats.shape[0]
1743    else:
1744        raise ValueError('Station latitudes must be a list or 1-D NumPy array.')
1745    if isinstance(lons,list):
1746        nlons = len(lons)
1747    elif isinstance(lons,np.ndarray) and len(lons.shape) == 1:
1748        nlons = lons.shape[0]
1749    else:
1750        raise ValueError('Station longitudes must be a list or 1-D NumPy array.')
1751    if nlats != nlons:
1752        raise ValueError('Station lats and lons must be same size.')
1753
1754    ni = grid_def_in.npoints
1755    no = nlats
1756
1757    # Adjust shape of input array(s)
1758    a,newshp = _adjust_array_shape_for_interp_stations(a,grid_def_in,no)
1759
1760    # Set lats and lons if stations
1761    rlat = np.array(lats,dtype=np.float32)
1762    rlon = np.array(lons,dtype=np.float32)
1763
1764    # Use gdtn = -1 for stations and an empty template array
1765    gdtn = -1
1766    gdt = np.zeros((200),dtype=np.int32)
1767
1768    # Call interpolation subroutines according to type of a.
1769    if isinstance(a,np.ndarray):
1770        # Scalar
1771        ibi = np.zeros((a.shape[0]),dtype=np.int32)
1772        li = np.zeros(a.shape,dtype=np.int32)
1773        go = np.zeros((a.shape[0],no),dtype=np.float32)
1774        no,ibo,lo,iret = interpolate.interpolate_scalar(method,method_options,
1775                                                 grid_def_in.gdtn,grid_def_in.gdt,
1776                                                 gdtn,gdt,
1777                                                 ibi,li.T,a.T,go.T,rlat,rlon)
1778        out = go.reshape(newshp)
1779    elif isinstance(a,tuple):
1780        # Vector
1781        ibi = np.zeros((a[0].shape[0]),dtype=np.int32)
1782        li = np.zeros(a[0].shape,dtype=np.int32)
1783        uo = np.zeros((a[0].shape[0],no),dtype=np.float32)
1784        vo = np.zeros((a[1].shape[0],no),dtype=np.float32)
1785        crot = np.ones((no),dtype=np.float32)
1786        srot = np.zeros((no),dtype=np.float32)
1787        no,ibo,lo,iret = interpolate.interpolate_vector(method,method_options,
1788                                                 grid_def_in.gdtn,grid_def_in.gdt,
1789                                                 gdtn,gdt,
1790                                                 ibi,li.T,a[0].T,a[1].T,uo.T,vo.T,
1791                                                 rlat,rlon,crot,srot)
1792        del crot
1793        del srot
1794        out = (uo.reshape(newshp),vo.reshape(newshp))
1795
1796    del rlat
1797    del rlon
1798
1799    try:
1800        if grib2io_interp.has_openmp_support:
1801            grib2io_interp.set_openmp_threads(prev_num_threads)
1802    except(AttributeError):
1803        pass
1804
1805    return out

Module-level interpolation function for interpolation to stations.

Interfaces with the grib2io_interp component package that interfaces to the NCEPLIBS-ip library. It supports scalar and vector interpolation according to the type of object a.

Parameters
  • a (numpy.ndarray or tuple): Input data. If a is a numpy.ndarray, scalar interpolation will be performed. If a is a tuple, then vector interpolation will be performed with the assumption that u = a[0] and v = a[1] and are both numpy.ndarray.

    These data are expected to be in 2-dimensional form with shape (ny, nx) or 3-dimensional (:, ny, nx) where the 1st dimension represents another spatial, temporal, or classification (i.e. ensemble members) dimension. The function will properly flatten the (ny,nx) dimensions into (nx * ny) acceptable for input into the interpolation subroutines.

  • method: Interpolate method to use. This can either be an integer or string using the following mapping:
Interpolate Scheme Integer Value
'bilinear' 0
'bicubic' 1
'neighbor' 2
'budget' 3
'spectral' 4
'neighbor-budget' 6

  • grid_def_in (grib2io.Grib2GridDef): Grib2GridDef object for the input grid.
  • lats (numpy.ndarray or list): Latitudes for station points
  • lons (numpy.ndarray or list): Longitudes for station points
  • method_options (list of ints, optional): Interpolation options. See the NCEPLIBS-ip documentation for more information on how these are used.
  • num_threads (int, optional): Number of OpenMP threads to use for interpolation. The default value is 1. If grib2io_interp was not built with OpenMP, then this keyword argument and value will have no impact.
Returns
  • interpolate_to_stations: Returns a numpy.ndarray when scalar interpolation is performed or a tuple of numpy.ndarrays when vector interpolation is performed with the assumptions that 0-index is the interpolated u and 1-index is the interpolated v.
class Grib2Message:
540class Grib2Message:
541    """
542    Creation class for a GRIB2 message.
543
544    This class returns a dynamically-created Grib2Message object that
545    inherits from `_Grib2Message` and grid, product, data representation
546    template classes according to the template numbers for the respective
547    sections. If `section3`, `section4`, or `section5` are omitted, then
548    the appropriate keyword arguments for the template number `gdtn=`,
549    `pdtn=`, or `drtn=` must be provided.
550
551    Parameters
552    ----------
553    section0
554        GRIB2 section 0 array.
555    section1
556        GRIB2 section 1 array.
557    section2
558        Local Use section data.
559    section3
560        GRIB2 section 3 array.
561    section4
562        GRIB2 section 4 array.
563    section5
564        GRIB2 section 5 array.
565
566    Returns
567    -------
568    Msg
569        A dynamically-create Grib2Message object that inherits from
570        _Grib2Message, a grid definition template class, product
571        definition template class, and a data representation template
572        class.
573    """
574    def __new__(self, section0: NDArray = np.array([struct.unpack('>I',b'GRIB')[0],0,0,2,0]),
575                      section1: NDArray = np.zeros((13),dtype=np.int64),
576                      section2: Optional[bytes] = None,
577                      section3: Optional[NDArray] = None,
578                      section4: Optional[NDArray] = None,
579                      section5: Optional[NDArray] = None, *args, **kwargs):
580
581        if np.all(section1==0):
582            try:
583                # Python >= 3.10
584                section1[5:11] = datetime.datetime.fromtimestamp(0, datetime.UTC).timetuple()[:6]
585            except(AttributeError):
586                # Python < 3.10
587                section1[5:11] = datetime.datetime.utcfromtimestamp(0).timetuple()[:6]
588
589        bases = list()
590        if section3 is None:
591            if 'gdtn' in kwargs.keys():
592                gdtn = kwargs['gdtn']
593                Gdt = templates.gdt_class_by_gdtn(gdtn)
594                bases.append(Gdt)
595                section3 = np.zeros((Gdt._len+5),dtype=np.int64)
596                section3[4] = gdtn
597            else:
598                raise ValueError("Must provide GRIB2 Grid Definition Template Number or section 3 array")
599        else:
600            gdtn = section3[4]
601            Gdt = templates.gdt_class_by_gdtn(gdtn)
602            bases.append(Gdt)
603
604        if section4 is None:
605            if 'pdtn' in kwargs.keys():
606                pdtn = kwargs['pdtn']
607                Pdt = templates.pdt_class_by_pdtn(pdtn)
608                bases.append(Pdt)
609                section4 = np.zeros((Pdt._len+2),dtype=np.int64)
610                section4[1] = pdtn
611            else:
612                raise ValueError("Must provide GRIB2 Production Definition Template Number or section 4 array")
613        else:
614            pdtn = section4[1]
615            Pdt = templates.pdt_class_by_pdtn(pdtn)
616            bases.append(Pdt)
617
618        if section5 is None:
619            if 'drtn' in kwargs.keys():
620                drtn = kwargs['drtn']
621                Drt = templates.drt_class_by_drtn(drtn)
622                bases.append(Drt)
623                section5 = np.zeros((Drt._len+2),dtype=np.int64)
624                section5[1] = drtn
625            else:
626                raise ValueError("Must provide GRIB2 Data Representation Template Number or section 5 array")
627        else:
628            drtn = section5[1]
629            Drt = templates.drt_class_by_drtn(drtn)
630            bases.append(Drt)
631
632        # attempt to use existing Msg class if it has already been made with gdtn,pdtn,drtn combo
633        try:
634            Msg = _msg_class_store[f"{gdtn}:{pdtn}:{drtn}"]
635        except KeyError:
636            @dataclass(init=False, repr=False)
637            class Msg(_Grib2Message, *bases):
638                pass
639            _msg_class_store[f"{gdtn}:{pdtn}:{drtn}"] = Msg
640
641        return Msg(section0, section1, section2, section3, section4, section5, *args)

Creation class for a GRIB2 message.

This class returns a dynamically-created Grib2Message object that inherits from _Grib2Message and grid, product, data representation template classes according to the template numbers for the respective sections. If section3, section4, or section5 are omitted, then the appropriate keyword arguments for the template number gdtn=, pdtn=, or drtn= must be provided.

Parameters
  • section0: GRIB2 section 0 array.
  • section1: GRIB2 section 1 array.
  • section2: Local Use section data.
  • section3: GRIB2 section 3 array.
  • section4: GRIB2 section 4 array.
  • section5: GRIB2 section 5 array.
Returns
  • Msg: A dynamically-create Grib2Message object that inherits from _Grib2Message, a grid definition template class, product definition template class, and a data representation template class.
@dataclass
class _Grib2Message:
 644@dataclass
 645class _Grib2Message:
 646    """
 647    GRIB2 Message base class.
 648    """
 649    # GRIB2 Sections
 650    section0: NDArray = field(init=True,repr=False)
 651    section1: NDArray = field(init=True,repr=False)
 652    section2: bytes = field(init=True,repr=False)
 653    section3: NDArray = field(init=True,repr=False)
 654    section4: NDArray = field(init=True,repr=False)
 655    section5: NDArray = field(init=True,repr=False)
 656    bitMapFlag: templates.Grib2Metadata = field(init=True,repr=False,default=255)
 657
 658    # Section 0 looked up attributes
 659    indicatorSection: NDArray = field(init=False,repr=False,default=templates.IndicatorSection())
 660    discipline: templates.Grib2Metadata = field(init=False,repr=False,default=templates.Discipline())
 661
 662    # Section 1 looked up attributes
 663    identificationSection: NDArray = field(init=False,repr=False,default=templates.IdentificationSection())
 664    originatingCenter: templates.Grib2Metadata = field(init=False,repr=False,default=templates.OriginatingCenter())
 665    originatingSubCenter: templates.Grib2Metadata = field(init=False,repr=False,default=templates.OriginatingSubCenter())
 666    masterTableInfo: templates.Grib2Metadata = field(init=False,repr=False,default=templates.MasterTableInfo())
 667    localTableInfo: templates.Grib2Metadata = field(init=False,repr=False,default=templates.LocalTableInfo())
 668    significanceOfReferenceTime: templates.Grib2Metadata = field(init=False,repr=False,default=templates.SignificanceOfReferenceTime())
 669    year: int = field(init=False,repr=False,default=templates.Year())
 670    month: int = field(init=False,repr=False,default=templates.Month())
 671    day: int = field(init=False,repr=False,default=templates.Day())
 672    hour: int = field(init=False,repr=False,default=templates.Hour())
 673    minute: int = field(init=False,repr=False,default=templates.Minute())
 674    second: int = field(init=False,repr=False,default=templates.Second())
 675    refDate: datetime.datetime = field(init=False,repr=False,default=templates.RefDate())
 676    productionStatus: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ProductionStatus())
 677    typeOfData: templates.Grib2Metadata = field(init=False,repr=False,default=templates.TypeOfData())
 678
 679    # Section 3 looked up common attributes.  Other looked up attributes are available according
 680    # to the Grid Definition Template.
 681    gridDefinitionSection: NDArray = field(init=False,repr=False,default=templates.GridDefinitionSection())
 682    sourceOfGridDefinition: int = field(init=False,repr=False,default=templates.SourceOfGridDefinition())
 683    numberOfDataPoints: int = field(init=False,repr=False,default=templates.NumberOfDataPoints())
 684    interpretationOfListOfNumbers: templates.Grib2Metadata = field(init=False,repr=False,default=templates.InterpretationOfListOfNumbers())
 685    gridDefinitionTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.GridDefinitionTemplateNumber())
 686    gridDefinitionTemplate: list = field(init=False,repr=False,default=templates.GridDefinitionTemplate())
 687    _earthparams: dict = field(init=False,repr=False,default=templates.EarthParams())
 688    _dxsign: float = field(init=False,repr=False,default=templates.DxSign())
 689    _dysign: float = field(init=False,repr=False,default=templates.DySign())
 690    _llscalefactor: float = field(init=False,repr=False,default=templates.LLScaleFactor())
 691    _lldivisor: float = field(init=False,repr=False,default=templates.LLDivisor())
 692    _xydivisor: float = field(init=False,repr=False,default=templates.XYDivisor())
 693    shapeOfEarth: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ShapeOfEarth())
 694    earthShape: str = field(init=False,repr=False,default=templates.EarthShape())
 695    earthRadius: float = field(init=False,repr=False,default=templates.EarthRadius())
 696    earthMajorAxis: float = field(init=False,repr=False,default=templates.EarthMajorAxis())
 697    earthMinorAxis: float = field(init=False,repr=False,default=templates.EarthMinorAxis())
 698    resolutionAndComponentFlags: list = field(init=False,repr=False,default=templates.ResolutionAndComponentFlags())
 699    ny: int = field(init=False,repr=False,default=templates.Ny())
 700    nx: int = field(init=False,repr=False,default=templates.Nx())
 701    scanModeFlags: list = field(init=False,repr=False,default=templates.ScanModeFlags())
 702    projParameters: dict = field(init=False,repr=False,default=templates.ProjParameters())
 703
 704    # Section 4
 705    productDefinitionTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ProductDefinitionTemplateNumber())
 706    productDefinitionTemplate: NDArray = field(init=False,repr=False,default=templates.ProductDefinitionTemplate())
 707
 708    # Section 5 looked up common attributes.  Other looked up attributes are
 709    # available according to the Data Representation Template.
 710    numberOfPackedValues: int = field(init=False,repr=False,default=templates.NumberOfPackedValues())
 711    dataRepresentationTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.DataRepresentationTemplateNumber())
 712    dataRepresentationTemplate: list = field(init=False,repr=False,default=templates.DataRepresentationTemplate())
 713    typeOfValues: templates.Grib2Metadata = field(init=False,repr=False,default=templates.TypeOfValues())
 714
 715    def __post_init__(self):
 716        """Set some attributes after init."""
 717        self._auto_nans = _AUTO_NANS
 718        self._coordlist = None
 719        self._data = None
 720        self._deflist = None
 721        self._msgnum = -1
 722        self._ondiskarray = None
 723        self._orig_section5 = np.copy(self.section5)
 724        self._signature = self._generate_signature()
 725        try:
 726            self._sha1_section3 = hashlib.sha1(self.section3).hexdigest()
 727        except(TypeError):
 728            pass
 729        self.bitMapFlag = templates.Grib2Metadata(self.bitMapFlag,table='6.0')
 730        self.bitmap = None
 731
 732
 733    @property
 734    def _isNDFD(self):
 735        """Check if GRIB2 message is from NWS NDFD"""
 736        return np.all(self.section1[0:2]==[8,65535])
 737
 738
 739    @property
 740    def gdtn(self):
 741        """Return Grid Definition Template Number"""
 742        return self.section3[4]
 743
 744
 745    @property
 746    def gdt(self):
 747        """Return Grid Definition Template."""
 748        return self.gridDefinitionTemplate
 749
 750
 751    @property
 752    def pdtn(self):
 753        """Return Product Definition Template Number."""
 754        return self.section4[1]
 755
 756
 757    @property
 758    def pdt(self):
 759        """Return Product Definition Template."""
 760        return self.productDefinitionTemplate
 761
 762
 763    @property
 764    def drtn(self):
 765        """Return Data Representation Template Number."""
 766        return self.section5[1]
 767
 768
 769    @property
 770    def drt(self):
 771        """Return Data Representation Template."""
 772        return self.dataRepresentationTemplate
 773
 774
 775    @property
 776    def pdy(self):
 777        """Return the PDY ('YYYYMMDD')."""
 778        return ''.join([str(i) for i in self.section1[5:8]])
 779
 780
 781    @property
 782    def griddef(self):
 783        """Return a Grib2GridDef instance for a GRIB2 message."""
 784        return Grib2GridDef.from_section3(self.section3)
 785
 786
 787    @property
 788    def lats(self):
 789        """Return grid latitudes."""
 790        return self.latlons()[0]
 791
 792
 793    @property
 794    def lons(self):
 795        """Return grid longitudes."""
 796        return self.latlons()[1]
 797
 798    @property
 799    def min(self):
 800        """Return minimum value of data."""
 801        return np.nanmin(self.data)
 802
 803
 804    @property
 805    def max(self):
 806        """Return maximum value of data."""
 807        return np.nanmax(self.data)
 808
 809
 810    @property
 811    def mean(self):
 812        """Return mean value of data."""
 813        return np.nanmean(self.data)
 814
 815
 816    @property
 817    def median(self):
 818        """Return median value of data."""
 819        return np.nanmedian(self.data)
 820
 821
 822    def __repr__(self):
 823        """
 824        Return an unambiguous string representation of the object.
 825
 826        Returns
 827        -------
 828        repr
 829            A string representation of the object, including information from
 830            sections 0, 1, 3, 4, 5, and 6.
 831        """
 832        info = ''
 833        for sect in [0,1,3,4,5,6]:
 834            for k,v in self.attrs_by_section(sect,values=True).items():
 835                info += f'Section {sect}: {k} = {v}\n'
 836        return info
 837
 838    def __str__(self):
 839        """
 840        Return a readable string representation of the object.
 841
 842        Returns
 843        -------
 844        str
 845            A formatted string representation of the object, including
 846            selected attributes.
 847        """
 848        return (f'{self._msgnum}:d={self.refDate}:{self.shortName}:'
 849                f'{self.fullName} ({self.units}):{self.level}:'
 850                f'{self.leadTime}')
 851
 852
 853    def _generate_signature(self):
 854        """Generature SHA-1 hash string from GRIB2 integer sections."""
 855        return hashlib.sha1(np.concatenate((self.section0,self.section1,
 856                                            self.section3,self.section4,
 857                                            self.section5))).hexdigest()
 858
 859
 860    def attrs_by_section(self, sect: int, values: bool=False):
 861        """
 862        Provide a tuple of attribute names for the given GRIB2 section.
 863
 864        Parameters
 865        ----------
 866        sect
 867            The GRIB2 section number.
 868        values
 869            Optional (default is `False`) argument to return attributes values.
 870
 871        Returns
 872        -------
 873        attrs_by_section
 874            A list of attribute names or dict of name:value pairs if `values =
 875            True`.
 876        """
 877        if sect in {0,1,6}:
 878            attrs = templates._section_attrs[sect]
 879        elif sect in {3,4,5}:
 880            def _find_class_index(n):
 881                _key = {3:'Grid', 4:'Product', 5:'Data'}
 882                for i,c in enumerate(self.__class__.__mro__):
 883                    if _key[n] in c.__name__:
 884                        return i
 885                else:
 886                    return []
 887            if sys.version_info.minor <= 8:
 888                attrs = templates._section_attrs[sect]+\
 889                        [a for a in dir(self.__class__.__mro__[_find_class_index(sect)]) if not a.startswith('_')]
 890            else:
 891                attrs = templates._section_attrs[sect]+\
 892                        self.__class__.__mro__[_find_class_index(sect)]._attrs
 893        else:
 894            attrs = []
 895        if values:
 896            return {k:getattr(self,k) for k in attrs}
 897        else:
 898            return attrs
 899
 900
 901    def pack(self):
 902        """
 903        Pack GRIB2 section data into a binary message.
 904
 905        It is the user's responsibility to populate the GRIB2 section
 906        information with appropriate metadata.
 907        """
 908        # Create beginning of packed binary message with section 0 and 1 data.
 909        self._sections = []
 910        self._msg,self._pos = g2clib.grib2_create(self.indicatorSection[2:4],self.identificationSection)
 911        self._sections += [0,1]
 912
 913        # Add section 2 if present.
 914        if isinstance(self.section2,bytes) and len(self.section2) > 0:
 915            self._msg,self._pos = g2clib.grib2_addlocal(self._msg,self.section2)
 916            self._sections.append(2)
 917
 918        # Add section 3.
 919        self.section3[1] = self.nx * self.ny
 920        self._msg,self._pos = g2clib.grib2_addgrid(self._msg,self.gridDefinitionSection,
 921                                                   self.gridDefinitionTemplate,
 922                                                   self._deflist)
 923        self._sections.append(3)
 924
 925        # Prepare data.
 926        if self._data is None:
 927            if self._ondiskarray is None:
 928                raise ValueError("Grib2Message object has no data, thus it cannot be packed.")
 929        field = np.copy(self.data)
 930        if self.scanModeFlags is not None:
 931            if self.scanModeFlags[3]:
 932                fieldsave = field.astype('f') # Casting makes a copy
 933                field[1::2,:] = fieldsave[1::2,::-1]
 934        fld = field.astype('f')
 935
 936        # Prepare bitmap, if necessary
 937        bitmapflag = self.bitMapFlag.value
 938        if bitmapflag == 0:
 939            if self.bitmap is not None:
 940                bmap = np.ravel(self.bitmap).astype(DEFAULT_NUMPY_INT)
 941            else:
 942                bmap = np.ravel(np.where(np.isnan(fld),0,1)).astype(DEFAULT_NUMPY_INT)
 943        else:
 944            bmap = None
 945
 946        # Prepare optional coordinate list
 947        if self._coordlist is not None:
 948            crdlist = np.array(self._coordlist,'f')
 949        else:
 950            crdlist = None
 951
 952        # Prepare data for packing if nans are present
 953        fld = np.ravel(fld)
 954        if bitmapflag in {0,254}:
 955            fld = np.where(np.isnan(fld),0,fld)
 956        else:
 957            if np.isnan(fld).any():
 958                if hasattr(self,'priMissingValue'):
 959                    fld = np.where(np.isnan(fld),self.priMissingValue,fld)
 960            if hasattr(self,'_missvalmap'):
 961                if hasattr(self,'priMissingValue'):
 962                    fld = np.where(self._missvalmap==1,self.priMissingValue,fld)
 963                if hasattr(self,'secMissingValue'):
 964                    fld = np.where(self._missvalmap==2,self.secMissingValue,fld)
 965
 966        # Add sections 4, 5, 6, and 7.
 967        self._msg,self._pos = g2clib.grib2_addfield(self._msg,self.pdtn,
 968                                                    self.productDefinitionTemplate,
 969                                                    crdlist,
 970                                                    self.drtn,
 971                                                    self.dataRepresentationTemplate,
 972                                                    fld,
 973                                                    bitmapflag,
 974                                                    bmap)
 975        self._sections.append(4)
 976        self._sections.append(5)
 977        self._sections.append(6)
 978        self._sections.append(7)
 979
 980        # Finalize GRIB2 message with section 8.
 981        self._msg, self._pos = g2clib.grib2_end(self._msg)
 982        self._sections.append(8)
 983        self.section0[-1] = len(self._msg)
 984
 985
 986    @property
 987    def data(self) -> np.array:
 988        """Access the unpacked data values."""
 989        if self._data is None:
 990            if self._auto_nans != _AUTO_NANS:
 991                self._data = self._ondiskarray
 992            self._data = np.asarray(self._ondiskarray)
 993        return self._data
 994
 995
 996    @data.setter
 997    def data(self, data):
 998        if not isinstance(data, np.ndarray):
 999            raise ValueError('Grib2Message data only supports numpy arrays')
1000        self._data = data
1001
1002
1003    def flush_data(self):
1004        """
1005        Flush the unpacked data values from the Grib2Message object.
1006
1007        Note: If the Grib2Message object was constructed from "scratch" (i.e.
1008        not read from file), this method will remove the data array from
1009        the object and it cannot be recovered.
1010        """
1011        self._data = None
1012        self.bitmap = None
1013
1014
1015    def __getitem__(self, item):
1016        return self.data[item]
1017
1018
1019    def __setitem__(self, item):
1020        raise NotImplementedError('assignment of data not supported via setitem')
1021
1022
1023    def latlons(self, *args, **kwrgs):
1024        """Alias for `grib2io.Grib2Message.grid` method."""
1025        return self.grid(*args, **kwrgs)
1026
1027
1028    def grid(self, unrotate: bool=True):
1029        """
1030        Return lats,lons (in degrees) of grid.
1031
1032        Currently can handle reg. lat/lon,cglobal Gaussian, mercator,
1033        stereographic, lambert conformal, albers equal-area, space-view and
1034        azimuthal equidistant grids.
1035
1036        Parameters
1037        ----------
1038        unrotate
1039            If `True` [DEFAULT], and grid is rotated lat/lon, then unrotate the
1040            grid, otherwise `False`, do not.
1041
1042        Returns
1043        -------
1044        lats, lons : numpy.ndarray
1045            Returns two numpy.ndarrays with dtype=numpy.float32 of grid
1046            latitudes and longitudes in units of degrees.
1047        """
1048        if self._sha1_section3 in _latlon_datastore.keys():
1049            return (_latlon_datastore[self._sha1_section3]['latitude'],
1050                    _latlon_datastore[self._sha1_section3]['longitude'])
1051        gdtn = self.gridDefinitionTemplateNumber.value
1052        gdtmpl = self.gridDefinitionTemplate
1053        reggrid = self.gridDefinitionSection[2] == 0 # This means regular 2-d grid
1054        if gdtn == 0:
1055            # Regular lat/lon grid
1056            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1057            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1058            dlon = self.gridlengthXDirection
1059            dlat = self.gridlengthYDirection
1060            if lon2 < lon1 and dlon < 0: lon1 = -lon1
1061            lats = np.linspace(lat1,lat2,self.ny)
1062            if reggrid:
1063                lons = np.linspace(lon1,lon2,self.nx)
1064            else:
1065                lons = np.linspace(lon1,lon2,self.ny*2)
1066            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1067        elif gdtn == 1: # Rotated Lat/Lon grid
1068            pj = pyproj.Proj(self.projParameters)
1069            lat1,lon1 = self.latitudeFirstGridpoint,self.longitudeFirstGridpoint
1070            lat2,lon2 = self.latitudeLastGridpoint,self.longitudeLastGridpoint
1071            if lon1 > 180.0: lon1 -= 360.0
1072            if lon2 > 180.0: lon2 -= 360.0
1073            lats = np.linspace(lat1,lat2,self.ny)
1074            lons = np.linspace(lon1,lon2,self.nx)
1075            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1076            if unrotate:
1077                from grib2io.utils import rotated_grid
1078                lats,lons = rotated_grid.unrotate(lats,lons,self.anglePoleRotation,
1079                                                  self.latitudeSouthernPole,
1080                                                  self.longitudeSouthernPole)
1081        elif gdtn == 40: # Gaussian grid (only works for global!)
1082            from grib2io.utils.gauss_grid import gaussian_latitudes
1083            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1084            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1085            nlats = self.ny
1086            if not reggrid: # Reduced Gaussian grid.
1087                nlons = 2*nlats
1088                dlon = 360./nlons
1089            else:
1090                nlons = self.nx
1091                dlon = self.gridlengthXDirection
1092            lons = np.linspace(lon1,lon2,nlons)
1093            # Compute Gaussian lats (north to south)
1094            lats = gaussian_latitudes(nlats)
1095            if lat1 < lat2:  # reverse them if necessary
1096                lats = lats[::-1]
1097            lons,lats = np.meshgrid(lons,lats)
1098        elif gdtn in {10,20,30,31,110}:
1099            # Mercator, Lambert Conformal, Stereographic, Albers Equal Area,
1100            # Azimuthal Equidistant
1101            dx,dy = self.gridlengthXDirection, self.gridlengthYDirection
1102            lon1,lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1103            pj = pyproj.Proj(self.projParameters)
1104            llcrnrx, llcrnry = pj(lon1,lat1)
1105            x = llcrnrx+dx*np.arange(self.nx)
1106            y = llcrnry+dy*np.arange(self.ny)
1107            x,y = np.meshgrid(x, y)
1108            lons,lats = pj(x, y, inverse=True)
1109        elif gdtn == 90:
1110            # Satellite Projection
1111            dx = self.gridlengthXDirection
1112            dy = self.gridlengthYDirection
1113            pj = pyproj.Proj(self.projParameters)
1114            x = dx*np.indices((self.ny,self.nx),'f')[1,:,:]
1115            x -= 0.5*x.max()
1116            y = dy*np.indices((self.ny,self.nx),'f')[0,:,:]
1117            y -= 0.5*y.max()
1118            lons,lats = pj(x,y,inverse=True)
1119            # Set lons,lats to 1.e30 where undefined
1120            abslons = np.fabs(lons)
1121            abslats = np.fabs(lats)
1122            lons = np.where(abslons < 1.e20, lons, 1.e30)
1123            lats = np.where(abslats < 1.e20, lats, 1.e30)
1124        elif gdtn == 32769:
1125            # Special NCEP Grid, Rotated Lat/Lon, Arakawa E-Grid (Non-Staggered)
1126            from grib2io.utils import arakawa_rotated_grid
1127            from grib2io.utils.rotated_grid import DEG2RAD
1128            di, dj = 0.0, 0.0
1129            do_180 = False
1130            idir = 1 if self.scanModeFlags[0] == 0 else -1
1131            jdir = -1 if self.scanModeFlags[1] == 0 else 1
1132            grid_oriented = 0 if self.resolutionAndComponentFlags[4] == 0 else 1
1133            do_rot = 1 if grid_oriented == 1 else 0
1134            la1 = self.latitudeFirstGridpoint
1135            lo1 = self.longitudeFirstGridpoint
1136            clon = self.longitudeCenterGridpoint
1137            clat = self.latitudeCenterGridpoint
1138            lasp = clat - 90.0
1139            losp = clon
1140            llat, llon = arakawa_rotated_grid.ll2rot(la1,lo1,lasp,losp)
1141            la2, lo2 = arakawa_rotated_grid.rot2ll(-llat,-llon,lasp,losp)
1142            rlat = -llat
1143            rlon = -llon
1144            if self.nx == 1:
1145                di = 0.0
1146            elif idir == 1:
1147                ti = rlon
1148                while ti < llon:
1149                    ti += 360.0
1150                di = (ti - llon)/float(self.nx-1)
1151            else:
1152                ti = llon
1153                while ti < rlon:
1154                    ti += 360.0
1155                di = (ti - rlon)/float(self.nx-1)
1156            if self.ny == 1:
1157               dj = 0.0
1158            else:
1159                dj = (rlat - llat)/float(self.ny-1)
1160                if dj < 0.0:
1161                    dj = -dj
1162            if idir == 1:
1163                if llon > rlon:
1164                    llon -= 360.0
1165                if llon < 0 and rlon > 0:
1166                    do_180 = True
1167            else:
1168                if rlon > llon:
1169                    rlon -= 360.0
1170                if rlon < 0 and llon > 0:
1171                    do_180 = True
1172            xlat1d = llat + (np.arange(self.ny)*jdir*dj)
1173            xlon1d = llon + (np.arange(self.nx)*idir*di)
1174            xlons, xlats = np.meshgrid(xlon1d,xlat1d)
1175            rot2ll_vectorized = np.vectorize(arakawa_rotated_grid.rot2ll)
1176            lats, lons = rot2ll_vectorized(xlats,xlons,lasp,losp)
1177            if do_180:
1178                lons = np.where(lons>180.0,lons-360.0,lons)
1179            vector_rotation_angles_vectorized = np.vectorize(arakawa_rotated_grid.vector_rotation_angles)
1180            rots = vector_rotation_angles_vectorized(lats, lons, clat, losp, xlats)
1181            del xlat1d, xlon1d, xlats, xlons
1182        else:
1183            raise ValueError('Unsupported grid')
1184
1185        _latlon_datastore[self._sha1_section3] = dict(latitude=lats,longitude=lons)
1186        try:
1187            _latlon_datastore[self._sha1_section3]['vector_rotation_angles'] = rots
1188        except(NameError):
1189            pass
1190
1191        return lats, lons
1192
1193
1194    def map_keys(self):
1195        """
1196        Unpack data grid replacing integer values with strings.
1197
1198        These types of fields are categorical or classifications where data
1199        values do not represent an observable or predictable physical quantity.
1200        An example of such a field would be [Dominant Precipitation Type -
1201        DPTYPE](https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-201.shtml)
1202
1203        Returns
1204        -------
1205        map_keys
1206            numpy.ndarray of string values per element.
1207        """
1208        hold_auto_nans = _AUTO_NANS
1209        set_auto_nans(False)
1210        if (np.all(self.section1[0:2]==[7,14]) and self.shortName == 'PWTHER') or \
1211        (np.all(self.section1[0:2]==[8,65535]) and self.shortName == 'WX'):
1212            keys = utils.decode_wx_strings(self.section2)
1213            if hasattr(self,'priMissingValue') and self.priMissingValue not in [None,0]:
1214                keys[int(self.priMissingValue)] = 'Missing'
1215            if hasattr(self,'secMissingValue') and self.secMissingValue not in [None,0]:
1216                keys[int(self.secMissingValue)] = 'Missing'
1217            u,inv = np.unique(self.data,return_inverse=True)
1218            fld = np.array([keys[x] for x in u])[inv].reshape(self.data.shape)
1219        else:
1220            # For data whose units are defined in a code table (i.e. classification or mask)
1221            tblname = re.findall(r'\d\.\d+',self.units,re.IGNORECASE)[0]
1222            fld = self.data.astype(np.int32).astype(str)
1223            tbl = tables.get_table(tblname,expand=True)
1224            for val in np.unique(fld):
1225                fld = np.where(fld==val,tbl[val],fld)
1226        set_auto_nans(hold_auto_nans)
1227        return fld
1228
1229
1230    def to_bytes(self, validate: bool=True):
1231        """
1232        Return packed GRIB2 message in bytes format.
1233
1234        This will be useful for exporting data in non-file formats. For example,
1235        can be used to output grib data directly to S3 using the boto3 client
1236        without the need to write a temporary file to upload first.
1237
1238        Parameters
1239        ----------
1240        validate: default=True
1241            If `True` (DEFAULT), validates first/last four bytes for proper
1242            formatting, else returns None. If `False`, message is output as is.
1243
1244        Returns
1245        -------
1246        to_bytes
1247            Returns GRIB2 formatted message as bytes.
1248        """
1249        if hasattr(self,'_msg'):
1250            if validate:
1251                if self.validate():
1252                    return self._msg
1253                else:
1254                    return None
1255            else:
1256                return self._msg
1257        else:
1258            return None
1259
1260
1261    def interpolate(self, method, grid_def_out, method_options=None, drtn=None,
1262                    num_threads=1):
1263        """
1264        Grib2Message Interpolator
1265
1266        Performs spatial interpolation via the [NCEPLIBS-ip
1267        library](https://github.com/NOAA-EMC/NCEPLIBS-ip). This interpolate
1268        method only supports scalar interpolation. If you need to perform
1269        vector interpolation, use the module-level `grib2io.interpolate`
1270        function.
1271
1272        Parameters
1273        ----------
1274        method
1275            Interpolate method to use. This can either be an integer or string
1276            using the following mapping:
1277
1278            | Interpolate Scheme | Integer Value |
1279            | :---:              | :---:         |
1280            | 'bilinear'         | 0             |
1281            | 'bicubic'          | 1             |
1282            | 'neighbor'         | 2             |
1283            | 'budget'           | 3             |
1284            | 'spectral'         | 4             |
1285            | 'neighbor-budget'  | 6             |
1286
1287        grid_def_out : grib2io.Grib2GridDef
1288            Grib2GridDef object of the output grid.
1289        method_options : list of ints, optional
1290            Interpolation options. See the NCEPLIBS-ip documentation for
1291            more information on how these are used.
1292        drtn
1293            Data Representation Template to be used for the returned
1294            interpolated GRIB2 message. When `None`, the data representation
1295            template of the source GRIB2 message is used. Once again, it is the
1296            user's responsibility to properly set the Data Representation
1297            Template attributes.
1298        num_threads : int, optional
1299            Number of OpenMP threads to use for interpolation. The default
1300            value is 1. If grib2io_interp was not built with OpenMP, then
1301            this keyword argument and value will have no impact.
1302
1303        Returns
1304        -------
1305        interpolate
1306            If interpolating to a grid, a new Grib2Message object is returned.
1307            The GRIB2 metadata of the new Grib2Message object is identical to
1308            the input except where required to be different because of the new
1309            grid specs and possibly a new data representation template.
1310
1311            If interpolating to station points, the interpolated data values are
1312            returned as a numpy.ndarray.
1313        """
1314        section0 = self.section0
1315        section0[-1] = 0
1316        gds = [0, grid_def_out.npoints, 0, 255, grid_def_out.gdtn]
1317        section3 = np.concatenate((gds,grid_def_out.gdt))
1318        drtn = self.drtn if drtn is None else drtn
1319
1320        msg = Grib2Message(section0,self.section1,self.section2,section3,
1321                           self.section4,None,self.bitMapFlag.value,drtn=drtn)
1322
1323        msg._msgnum = -1
1324        msg._deflist = self._deflist
1325        msg._coordlist = self._coordlist
1326        shape = (msg.ny,msg.nx)
1327        ndim = 2
1328        if msg.typeOfValues == 0:
1329            dtype = 'float32'
1330        elif msg.typeOfValues == 1:
1331            dtype = 'int32'
1332        msg._data = interpolate(self.data,method,Grib2GridDef.from_section3(self.section3),grid_def_out,
1333                                method_options=method_options,num_threads=num_threads).reshape(msg.ny,msg.nx)
1334        msg.section5[0] = grid_def_out.npoints
1335        return msg
1336
1337
1338    def validate(self):
1339        """
1340        Validate a complete GRIB2 message.
1341
1342        The g2c library does its own internal validation when g2_gribend() is called, but
1343        we will check in grib2io also. The validation checks if the first 4 bytes in
1344        self._msg is 'GRIB' and '7777' as the last 4 bytes and that the message length in
1345        section 0 equals the length of the packed message.
1346
1347        Returns
1348        -------
1349        `True` if the packed GRIB2 message is complete and well-formed, `False` otherwise.
1350        """
1351        valid = False
1352        if hasattr(self,'_msg'):
1353            if self._msg[0:4]+self._msg[-4:] == b'GRIB7777':
1354                if self.section0[-1] == len(self._msg):
1355                    valid = True
1356        return valid

GRIB2 Message base class.

_Grib2Message( section0: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]], section1: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]], section2: bytes, section3: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]], section4: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]], section5: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]], bitMapFlag: grib2io.templates.Grib2Metadata = 255)
section0: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
section1: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
section2: bytes
section3: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
section4: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
section5: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
indicatorSection: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
identificationSection: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]

GRIB2 Section 1, Identification Section

year: int

Year of reference time

month: int

Month of reference time

day: int

Day of reference time

hour: int

Hour of reference time

minute: int

Minute of reference time

second: int

Second of reference time

refDate: datetime.datetime

Reference Date. NOTE: This is a datetime.datetime object.

gridDefinitionSection: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]

GRIB2 Section 3, Grid Definition Section

sourceOfGridDefinition: int
numberOfDataPoints: int

Number of Data Points

interpretationOfListOfNumbers: grib2io.templates.Grib2Metadata

Interpretation of List of Numbers

gridDefinitionTemplate: list

Grid definition template

earthShape: str

Description of the shape of the Earth

earthRadius: float

Radius of the Earth (Assumes "spherical")

earthMajorAxis: float

Major Axis of the Earth (Assumes "oblate spheroid" or "ellipsoid")

earthMinorAxis: float

Minor Axis of the Earth (Assumes "oblate spheroid" or "ellipsoid")

resolutionAndComponentFlags: list
ny: int

Number of grid points in the Y-direction (generally North-South)

nx: int

Number of grid points in the X-direction (generally East-West)

scanModeFlags: list
projParameters: dict

PROJ Parameters to define the reference system

productDefinitionTemplate: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]

Product Definition Template

numberOfPackedValues: int

Number of Packed Values

dataRepresentationTemplate: list

Data Representation Template

gdtn
739    @property
740    def gdtn(self):
741        """Return Grid Definition Template Number"""
742        return self.section3[4]

Return Grid Definition Template Number

gdt
745    @property
746    def gdt(self):
747        """Return Grid Definition Template."""
748        return self.gridDefinitionTemplate

Return Grid Definition Template.

pdtn
751    @property
752    def pdtn(self):
753        """Return Product Definition Template Number."""
754        return self.section4[1]

Return Product Definition Template Number.

pdt
757    @property
758    def pdt(self):
759        """Return Product Definition Template."""
760        return self.productDefinitionTemplate

Return Product Definition Template.

drtn
763    @property
764    def drtn(self):
765        """Return Data Representation Template Number."""
766        return self.section5[1]

Return Data Representation Template Number.

drt
769    @property
770    def drt(self):
771        """Return Data Representation Template."""
772        return self.dataRepresentationTemplate

Return Data Representation Template.

pdy
775    @property
776    def pdy(self):
777        """Return the PDY ('YYYYMMDD')."""
778        return ''.join([str(i) for i in self.section1[5:8]])

Return the PDY ('YYYYMMDD').

griddef
781    @property
782    def griddef(self):
783        """Return a Grib2GridDef instance for a GRIB2 message."""
784        return Grib2GridDef.from_section3(self.section3)

Return a Grib2GridDef instance for a GRIB2 message.

lats
787    @property
788    def lats(self):
789        """Return grid latitudes."""
790        return self.latlons()[0]

Return grid latitudes.

lons
793    @property
794    def lons(self):
795        """Return grid longitudes."""
796        return self.latlons()[1]

Return grid longitudes.

min
798    @property
799    def min(self):
800        """Return minimum value of data."""
801        return np.nanmin(self.data)

Return minimum value of data.

max
804    @property
805    def max(self):
806        """Return maximum value of data."""
807        return np.nanmax(self.data)

Return maximum value of data.

mean
810    @property
811    def mean(self):
812        """Return mean value of data."""
813        return np.nanmean(self.data)

Return mean value of data.

median
816    @property
817    def median(self):
818        """Return median value of data."""
819        return np.nanmedian(self.data)

Return median value of data.

def attrs_by_section(self, sect: int, values: bool = False):
860    def attrs_by_section(self, sect: int, values: bool=False):
861        """
862        Provide a tuple of attribute names for the given GRIB2 section.
863
864        Parameters
865        ----------
866        sect
867            The GRIB2 section number.
868        values
869            Optional (default is `False`) argument to return attributes values.
870
871        Returns
872        -------
873        attrs_by_section
874            A list of attribute names or dict of name:value pairs if `values =
875            True`.
876        """
877        if sect in {0,1,6}:
878            attrs = templates._section_attrs[sect]
879        elif sect in {3,4,5}:
880            def _find_class_index(n):
881                _key = {3:'Grid', 4:'Product', 5:'Data'}
882                for i,c in enumerate(self.__class__.__mro__):
883                    if _key[n] in c.__name__:
884                        return i
885                else:
886                    return []
887            if sys.version_info.minor <= 8:
888                attrs = templates._section_attrs[sect]+\
889                        [a for a in dir(self.__class__.__mro__[_find_class_index(sect)]) if not a.startswith('_')]
890            else:
891                attrs = templates._section_attrs[sect]+\
892                        self.__class__.__mro__[_find_class_index(sect)]._attrs
893        else:
894            attrs = []
895        if values:
896            return {k:getattr(self,k) for k in attrs}
897        else:
898            return attrs

Provide a tuple of attribute names for the given GRIB2 section.

Parameters
  • sect: The GRIB2 section number.
  • values: Optional (default is False) argument to return attributes values.
Returns
  • attrs_by_section: A list of attribute names or dict of name:value pairs if values = True.
def pack(self):
901    def pack(self):
902        """
903        Pack GRIB2 section data into a binary message.
904
905        It is the user's responsibility to populate the GRIB2 section
906        information with appropriate metadata.
907        """
908        # Create beginning of packed binary message with section 0 and 1 data.
909        self._sections = []
910        self._msg,self._pos = g2clib.grib2_create(self.indicatorSection[2:4],self.identificationSection)
911        self._sections += [0,1]
912
913        # Add section 2 if present.
914        if isinstance(self.section2,bytes) and len(self.section2) > 0:
915            self._msg,self._pos = g2clib.grib2_addlocal(self._msg,self.section2)
916            self._sections.append(2)
917
918        # Add section 3.
919        self.section3[1] = self.nx * self.ny
920        self._msg,self._pos = g2clib.grib2_addgrid(self._msg,self.gridDefinitionSection,
921                                                   self.gridDefinitionTemplate,
922                                                   self._deflist)
923        self._sections.append(3)
924
925        # Prepare data.
926        if self._data is None:
927            if self._ondiskarray is None:
928                raise ValueError("Grib2Message object has no data, thus it cannot be packed.")
929        field = np.copy(self.data)
930        if self.scanModeFlags is not None:
931            if self.scanModeFlags[3]:
932                fieldsave = field.astype('f') # Casting makes a copy
933                field[1::2,:] = fieldsave[1::2,::-1]
934        fld = field.astype('f')
935
936        # Prepare bitmap, if necessary
937        bitmapflag = self.bitMapFlag.value
938        if bitmapflag == 0:
939            if self.bitmap is not None:
940                bmap = np.ravel(self.bitmap).astype(DEFAULT_NUMPY_INT)
941            else:
942                bmap = np.ravel(np.where(np.isnan(fld),0,1)).astype(DEFAULT_NUMPY_INT)
943        else:
944            bmap = None
945
946        # Prepare optional coordinate list
947        if self._coordlist is not None:
948            crdlist = np.array(self._coordlist,'f')
949        else:
950            crdlist = None
951
952        # Prepare data for packing if nans are present
953        fld = np.ravel(fld)
954        if bitmapflag in {0,254}:
955            fld = np.where(np.isnan(fld),0,fld)
956        else:
957            if np.isnan(fld).any():
958                if hasattr(self,'priMissingValue'):
959                    fld = np.where(np.isnan(fld),self.priMissingValue,fld)
960            if hasattr(self,'_missvalmap'):
961                if hasattr(self,'priMissingValue'):
962                    fld = np.where(self._missvalmap==1,self.priMissingValue,fld)
963                if hasattr(self,'secMissingValue'):
964                    fld = np.where(self._missvalmap==2,self.secMissingValue,fld)
965
966        # Add sections 4, 5, 6, and 7.
967        self._msg,self._pos = g2clib.grib2_addfield(self._msg,self.pdtn,
968                                                    self.productDefinitionTemplate,
969                                                    crdlist,
970                                                    self.drtn,
971                                                    self.dataRepresentationTemplate,
972                                                    fld,
973                                                    bitmapflag,
974                                                    bmap)
975        self._sections.append(4)
976        self._sections.append(5)
977        self._sections.append(6)
978        self._sections.append(7)
979
980        # Finalize GRIB2 message with section 8.
981        self._msg, self._pos = g2clib.grib2_end(self._msg)
982        self._sections.append(8)
983        self.section0[-1] = len(self._msg)

Pack GRIB2 section data into a binary message.

It is the user's responsibility to populate the GRIB2 section information with appropriate metadata.

data: <built-in function array>
986    @property
987    def data(self) -> np.array:
988        """Access the unpacked data values."""
989        if self._data is None:
990            if self._auto_nans != _AUTO_NANS:
991                self._data = self._ondiskarray
992            self._data = np.asarray(self._ondiskarray)
993        return self._data

Access the unpacked data values.

def flush_data(self):
1003    def flush_data(self):
1004        """
1005        Flush the unpacked data values from the Grib2Message object.
1006
1007        Note: If the Grib2Message object was constructed from "scratch" (i.e.
1008        not read from file), this method will remove the data array from
1009        the object and it cannot be recovered.
1010        """
1011        self._data = None
1012        self.bitmap = None

Flush the unpacked data values from the Grib2Message object.

Note: If the Grib2Message object was constructed from "scratch" (i.e. not read from file), this method will remove the data array from the object and it cannot be recovered.

def latlons(self, *args, **kwrgs):
1023    def latlons(self, *args, **kwrgs):
1024        """Alias for `grib2io.Grib2Message.grid` method."""
1025        return self.grid(*args, **kwrgs)

Alias for grib2io.Grib2Message.grid method.

def grid(self, unrotate: bool = True):
1028    def grid(self, unrotate: bool=True):
1029        """
1030        Return lats,lons (in degrees) of grid.
1031
1032        Currently can handle reg. lat/lon,cglobal Gaussian, mercator,
1033        stereographic, lambert conformal, albers equal-area, space-view and
1034        azimuthal equidistant grids.
1035
1036        Parameters
1037        ----------
1038        unrotate
1039            If `True` [DEFAULT], and grid is rotated lat/lon, then unrotate the
1040            grid, otherwise `False`, do not.
1041
1042        Returns
1043        -------
1044        lats, lons : numpy.ndarray
1045            Returns two numpy.ndarrays with dtype=numpy.float32 of grid
1046            latitudes and longitudes in units of degrees.
1047        """
1048        if self._sha1_section3 in _latlon_datastore.keys():
1049            return (_latlon_datastore[self._sha1_section3]['latitude'],
1050                    _latlon_datastore[self._sha1_section3]['longitude'])
1051        gdtn = self.gridDefinitionTemplateNumber.value
1052        gdtmpl = self.gridDefinitionTemplate
1053        reggrid = self.gridDefinitionSection[2] == 0 # This means regular 2-d grid
1054        if gdtn == 0:
1055            # Regular lat/lon grid
1056            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1057            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1058            dlon = self.gridlengthXDirection
1059            dlat = self.gridlengthYDirection
1060            if lon2 < lon1 and dlon < 0: lon1 = -lon1
1061            lats = np.linspace(lat1,lat2,self.ny)
1062            if reggrid:
1063                lons = np.linspace(lon1,lon2,self.nx)
1064            else:
1065                lons = np.linspace(lon1,lon2,self.ny*2)
1066            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1067        elif gdtn == 1: # Rotated Lat/Lon grid
1068            pj = pyproj.Proj(self.projParameters)
1069            lat1,lon1 = self.latitudeFirstGridpoint,self.longitudeFirstGridpoint
1070            lat2,lon2 = self.latitudeLastGridpoint,self.longitudeLastGridpoint
1071            if lon1 > 180.0: lon1 -= 360.0
1072            if lon2 > 180.0: lon2 -= 360.0
1073            lats = np.linspace(lat1,lat2,self.ny)
1074            lons = np.linspace(lon1,lon2,self.nx)
1075            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1076            if unrotate:
1077                from grib2io.utils import rotated_grid
1078                lats,lons = rotated_grid.unrotate(lats,lons,self.anglePoleRotation,
1079                                                  self.latitudeSouthernPole,
1080                                                  self.longitudeSouthernPole)
1081        elif gdtn == 40: # Gaussian grid (only works for global!)
1082            from grib2io.utils.gauss_grid import gaussian_latitudes
1083            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1084            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1085            nlats = self.ny
1086            if not reggrid: # Reduced Gaussian grid.
1087                nlons = 2*nlats
1088                dlon = 360./nlons
1089            else:
1090                nlons = self.nx
1091                dlon = self.gridlengthXDirection
1092            lons = np.linspace(lon1,lon2,nlons)
1093            # Compute Gaussian lats (north to south)
1094            lats = gaussian_latitudes(nlats)
1095            if lat1 < lat2:  # reverse them if necessary
1096                lats = lats[::-1]
1097            lons,lats = np.meshgrid(lons,lats)
1098        elif gdtn in {10,20,30,31,110}:
1099            # Mercator, Lambert Conformal, Stereographic, Albers Equal Area,
1100            # Azimuthal Equidistant
1101            dx,dy = self.gridlengthXDirection, self.gridlengthYDirection
1102            lon1,lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1103            pj = pyproj.Proj(self.projParameters)
1104            llcrnrx, llcrnry = pj(lon1,lat1)
1105            x = llcrnrx+dx*np.arange(self.nx)
1106            y = llcrnry+dy*np.arange(self.ny)
1107            x,y = np.meshgrid(x, y)
1108            lons,lats = pj(x, y, inverse=True)
1109        elif gdtn == 90:
1110            # Satellite Projection
1111            dx = self.gridlengthXDirection
1112            dy = self.gridlengthYDirection
1113            pj = pyproj.Proj(self.projParameters)
1114            x = dx*np.indices((self.ny,self.nx),'f')[1,:,:]
1115            x -= 0.5*x.max()
1116            y = dy*np.indices((self.ny,self.nx),'f')[0,:,:]
1117            y -= 0.5*y.max()
1118            lons,lats = pj(x,y,inverse=True)
1119            # Set lons,lats to 1.e30 where undefined
1120            abslons = np.fabs(lons)
1121            abslats = np.fabs(lats)
1122            lons = np.where(abslons < 1.e20, lons, 1.e30)
1123            lats = np.where(abslats < 1.e20, lats, 1.e30)
1124        elif gdtn == 32769:
1125            # Special NCEP Grid, Rotated Lat/Lon, Arakawa E-Grid (Non-Staggered)
1126            from grib2io.utils import arakawa_rotated_grid
1127            from grib2io.utils.rotated_grid import DEG2RAD
1128            di, dj = 0.0, 0.0
1129            do_180 = False
1130            idir = 1 if self.scanModeFlags[0] == 0 else -1
1131            jdir = -1 if self.scanModeFlags[1] == 0 else 1
1132            grid_oriented = 0 if self.resolutionAndComponentFlags[4] == 0 else 1
1133            do_rot = 1 if grid_oriented == 1 else 0
1134            la1 = self.latitudeFirstGridpoint
1135            lo1 = self.longitudeFirstGridpoint
1136            clon = self.longitudeCenterGridpoint
1137            clat = self.latitudeCenterGridpoint
1138            lasp = clat - 90.0
1139            losp = clon
1140            llat, llon = arakawa_rotated_grid.ll2rot(la1,lo1,lasp,losp)
1141            la2, lo2 = arakawa_rotated_grid.rot2ll(-llat,-llon,lasp,losp)
1142            rlat = -llat
1143            rlon = -llon
1144            if self.nx == 1:
1145                di = 0.0
1146            elif idir == 1:
1147                ti = rlon
1148                while ti < llon:
1149                    ti += 360.0
1150                di = (ti - llon)/float(self.nx-1)
1151            else:
1152                ti = llon
1153                while ti < rlon:
1154                    ti += 360.0
1155                di = (ti - rlon)/float(self.nx-1)
1156            if self.ny == 1:
1157               dj = 0.0
1158            else:
1159                dj = (rlat - llat)/float(self.ny-1)
1160                if dj < 0.0:
1161                    dj = -dj
1162            if idir == 1:
1163                if llon > rlon:
1164                    llon -= 360.0
1165                if llon < 0 and rlon > 0:
1166                    do_180 = True
1167            else:
1168                if rlon > llon:
1169                    rlon -= 360.0
1170                if rlon < 0 and llon > 0:
1171                    do_180 = True
1172            xlat1d = llat + (np.arange(self.ny)*jdir*dj)
1173            xlon1d = llon + (np.arange(self.nx)*idir*di)
1174            xlons, xlats = np.meshgrid(xlon1d,xlat1d)
1175            rot2ll_vectorized = np.vectorize(arakawa_rotated_grid.rot2ll)
1176            lats, lons = rot2ll_vectorized(xlats,xlons,lasp,losp)
1177            if do_180:
1178                lons = np.where(lons>180.0,lons-360.0,lons)
1179            vector_rotation_angles_vectorized = np.vectorize(arakawa_rotated_grid.vector_rotation_angles)
1180            rots = vector_rotation_angles_vectorized(lats, lons, clat, losp, xlats)
1181            del xlat1d, xlon1d, xlats, xlons
1182        else:
1183            raise ValueError('Unsupported grid')
1184
1185        _latlon_datastore[self._sha1_section3] = dict(latitude=lats,longitude=lons)
1186        try:
1187            _latlon_datastore[self._sha1_section3]['vector_rotation_angles'] = rots
1188        except(NameError):
1189            pass
1190
1191        return lats, lons

Return lats,lons (in degrees) of grid.

Currently can handle reg. lat/lon,cglobal Gaussian, mercator, stereographic, lambert conformal, albers equal-area, space-view and azimuthal equidistant grids.

Parameters
  • unrotate: If True [DEFAULT], and grid is rotated lat/lon, then unrotate the grid, otherwise False, do not.
Returns
  • lats, lons (numpy.ndarray): Returns two numpy.ndarrays with dtype=numpy.float32 of grid latitudes and longitudes in units of degrees.
def map_keys(self):
1194    def map_keys(self):
1195        """
1196        Unpack data grid replacing integer values with strings.
1197
1198        These types of fields are categorical or classifications where data
1199        values do not represent an observable or predictable physical quantity.
1200        An example of such a field would be [Dominant Precipitation Type -
1201        DPTYPE](https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-201.shtml)
1202
1203        Returns
1204        -------
1205        map_keys
1206            numpy.ndarray of string values per element.
1207        """
1208        hold_auto_nans = _AUTO_NANS
1209        set_auto_nans(False)
1210        if (np.all(self.section1[0:2]==[7,14]) and self.shortName == 'PWTHER') or \
1211        (np.all(self.section1[0:2]==[8,65535]) and self.shortName == 'WX'):
1212            keys = utils.decode_wx_strings(self.section2)
1213            if hasattr(self,'priMissingValue') and self.priMissingValue not in [None,0]:
1214                keys[int(self.priMissingValue)] = 'Missing'
1215            if hasattr(self,'secMissingValue') and self.secMissingValue not in [None,0]:
1216                keys[int(self.secMissingValue)] = 'Missing'
1217            u,inv = np.unique(self.data,return_inverse=True)
1218            fld = np.array([keys[x] for x in u])[inv].reshape(self.data.shape)
1219        else:
1220            # For data whose units are defined in a code table (i.e. classification or mask)
1221            tblname = re.findall(r'\d\.\d+',self.units,re.IGNORECASE)[0]
1222            fld = self.data.astype(np.int32).astype(str)
1223            tbl = tables.get_table(tblname,expand=True)
1224            for val in np.unique(fld):
1225                fld = np.where(fld==val,tbl[val],fld)
1226        set_auto_nans(hold_auto_nans)
1227        return fld

Unpack data grid replacing integer values with strings.

These types of fields are categorical or classifications where data values do not represent an observable or predictable physical quantity. An example of such a field would be Dominant Precipitation Type - DPTYPE

Returns
  • map_keys: numpy.ndarray of string values per element.
def to_bytes(self, validate: bool = True):
1230    def to_bytes(self, validate: bool=True):
1231        """
1232        Return packed GRIB2 message in bytes format.
1233
1234        This will be useful for exporting data in non-file formats. For example,
1235        can be used to output grib data directly to S3 using the boto3 client
1236        without the need to write a temporary file to upload first.
1237
1238        Parameters
1239        ----------
1240        validate: default=True
1241            If `True` (DEFAULT), validates first/last four bytes for proper
1242            formatting, else returns None. If `False`, message is output as is.
1243
1244        Returns
1245        -------
1246        to_bytes
1247            Returns GRIB2 formatted message as bytes.
1248        """
1249        if hasattr(self,'_msg'):
1250            if validate:
1251                if self.validate():
1252                    return self._msg
1253                else:
1254                    return None
1255            else:
1256                return self._msg
1257        else:
1258            return None

Return packed GRIB2 message in bytes format.

This will be useful for exporting data in non-file formats. For example, can be used to output grib data directly to S3 using the boto3 client without the need to write a temporary file to upload first.

Parameters
  • validate (default=True): If True (DEFAULT), validates first/last four bytes for proper formatting, else returns None. If False, message is output as is.
Returns
  • to_bytes: Returns GRIB2 formatted message as bytes.
def interpolate( self, method, grid_def_out, method_options=None, drtn=None, num_threads=1):
1261    def interpolate(self, method, grid_def_out, method_options=None, drtn=None,
1262                    num_threads=1):
1263        """
1264        Grib2Message Interpolator
1265
1266        Performs spatial interpolation via the [NCEPLIBS-ip
1267        library](https://github.com/NOAA-EMC/NCEPLIBS-ip). This interpolate
1268        method only supports scalar interpolation. If you need to perform
1269        vector interpolation, use the module-level `grib2io.interpolate`
1270        function.
1271
1272        Parameters
1273        ----------
1274        method
1275            Interpolate method to use. This can either be an integer or string
1276            using the following mapping:
1277
1278            | Interpolate Scheme | Integer Value |
1279            | :---:              | :---:         |
1280            | 'bilinear'         | 0             |
1281            | 'bicubic'          | 1             |
1282            | 'neighbor'         | 2             |
1283            | 'budget'           | 3             |
1284            | 'spectral'         | 4             |
1285            | 'neighbor-budget'  | 6             |
1286
1287        grid_def_out : grib2io.Grib2GridDef
1288            Grib2GridDef object of the output grid.
1289        method_options : list of ints, optional
1290            Interpolation options. See the NCEPLIBS-ip documentation for
1291            more information on how these are used.
1292        drtn
1293            Data Representation Template to be used for the returned
1294            interpolated GRIB2 message. When `None`, the data representation
1295            template of the source GRIB2 message is used. Once again, it is the
1296            user's responsibility to properly set the Data Representation
1297            Template attributes.
1298        num_threads : int, optional
1299            Number of OpenMP threads to use for interpolation. The default
1300            value is 1. If grib2io_interp was not built with OpenMP, then
1301            this keyword argument and value will have no impact.
1302
1303        Returns
1304        -------
1305        interpolate
1306            If interpolating to a grid, a new Grib2Message object is returned.
1307            The GRIB2 metadata of the new Grib2Message object is identical to
1308            the input except where required to be different because of the new
1309            grid specs and possibly a new data representation template.
1310
1311            If interpolating to station points, the interpolated data values are
1312            returned as a numpy.ndarray.
1313        """
1314        section0 = self.section0
1315        section0[-1] = 0
1316        gds = [0, grid_def_out.npoints, 0, 255, grid_def_out.gdtn]
1317        section3 = np.concatenate((gds,grid_def_out.gdt))
1318        drtn = self.drtn if drtn is None else drtn
1319
1320        msg = Grib2Message(section0,self.section1,self.section2,section3,
1321                           self.section4,None,self.bitMapFlag.value,drtn=drtn)
1322
1323        msg._msgnum = -1
1324        msg._deflist = self._deflist
1325        msg._coordlist = self._coordlist
1326        shape = (msg.ny,msg.nx)
1327        ndim = 2
1328        if msg.typeOfValues == 0:
1329            dtype = 'float32'
1330        elif msg.typeOfValues == 1:
1331            dtype = 'int32'
1332        msg._data = interpolate(self.data,method,Grib2GridDef.from_section3(self.section3),grid_def_out,
1333                                method_options=method_options,num_threads=num_threads).reshape(msg.ny,msg.nx)
1334        msg.section5[0] = grid_def_out.npoints
1335        return msg

Grib2Message Interpolator

Performs spatial interpolation via the NCEPLIBS-ip library. This interpolate method only supports scalar interpolation. If you need to perform vector interpolation, use the module-level grib2io.interpolate function.

Parameters
  • method: Interpolate method to use. This can either be an integer or string using the following mapping:
Interpolate Scheme Integer Value
'bilinear' 0
'bicubic' 1
'neighbor' 2
'budget' 3
'spectral' 4
'neighbor-budget' 6

  • grid_def_out (grib2io.Grib2GridDef): Grib2GridDef object of the output grid.
  • method_options (list of ints, optional): Interpolation options. See the NCEPLIBS-ip documentation for more information on how these are used.
  • drtn: Data Representation Template to be used for the returned interpolated GRIB2 message. When None, the data representation template of the source GRIB2 message is used. Once again, it is the user's responsibility to properly set the Data Representation Template attributes.
  • num_threads (int, optional): Number of OpenMP threads to use for interpolation. The default value is 1. If grib2io_interp was not built with OpenMP, then this keyword argument and value will have no impact.
Returns
  • interpolate: If interpolating to a grid, a new Grib2Message object is returned. The GRIB2 metadata of the new Grib2Message object is identical to the input except where required to be different because of the new grid specs and possibly a new data representation template.

If interpolating to station points, the interpolated data values are returned as a numpy.ndarray.

def validate(self):
1338    def validate(self):
1339        """
1340        Validate a complete GRIB2 message.
1341
1342        The g2c library does its own internal validation when g2_gribend() is called, but
1343        we will check in grib2io also. The validation checks if the first 4 bytes in
1344        self._msg is 'GRIB' and '7777' as the last 4 bytes and that the message length in
1345        section 0 equals the length of the packed message.
1346
1347        Returns
1348        -------
1349        `True` if the packed GRIB2 message is complete and well-formed, `False` otherwise.
1350        """
1351        valid = False
1352        if hasattr(self,'_msg'):
1353            if self._msg[0:4]+self._msg[-4:] == b'GRIB7777':
1354                if self.section0[-1] == len(self._msg):
1355                    valid = True
1356        return valid

Validate a complete GRIB2 message.

The g2c library does its own internal validation when g2_gribend() is called, but we will check in grib2io also. The validation checks if the first 4 bytes in self._msg is 'GRIB' and '7777' as the last 4 bytes and that the message length in section 0 equals the length of the packed message.

Returns
  • True if the packed GRIB2 message is complete and well-formed, False otherwise.
@dataclass
class Grib2GridDef:
1808@dataclass
1809class Grib2GridDef:
1810    """
1811    Class for Grid Definition Template Number and Template as attributes.
1812
1813    This allows for cleaner looking code when passing these metadata around.
1814    For example, the `grib2io._Grib2Message.interpolate` method and
1815    `grib2io.interpolate` function accepts these objects.
1816    """
1817    gdtn: int
1818    gdt: NDArray
1819
1820    @classmethod
1821    def from_section3(cls, section3):
1822        return cls(section3[4],section3[5:])
1823
1824    @property
1825    def nx(self):
1826        return self.gdt[7]
1827
1828    @property
1829    def ny(self):
1830        return self.gdt[8]
1831
1832    @property
1833    def npoints(self):
1834        return self.gdt[7] * self.gdt[8]
1835
1836    @property
1837    def shape(self):
1838        return (self.ny, self.nx)

Class for Grid Definition Template Number and Template as attributes.

This allows for cleaner looking code when passing these metadata around. For example, the grib2io._Grib2Message.interpolate method and grib2io.interpolate function accepts these objects.

Grib2GridDef( gdtn: int, gdt: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]])
gdtn: int
gdt: numpy.ndarray[typing.Any, numpy.dtype[+_ScalarType_co]]
@classmethod
def from_section3(cls, section3):
1820    @classmethod
1821    def from_section3(cls, section3):
1822        return cls(section3[4],section3[5:])
nx
1824    @property
1825    def nx(self):
1826        return self.gdt[7]
ny
1828    @property
1829    def ny(self):
1830        return self.gdt[8]
npoints
1832    @property
1833    def npoints(self):
1834        return self.gdt[7] * self.gdt[8]
shape
1836    @property
1837    def shape(self):
1838        return (self.ny, self.nx)