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 extension module. The decoding/encoding of GRIB2 metadata is translated into more descriptive, plain language metadata by providing a mapping of the integer coded metadata to the appropriate GRIB2 code tables. The NCEP GRIB2 code tables are a part of the grib2io package.

Interpolation

As of grib2io v2.4.0, spatial interpolation via NCEPLIPS-ip Fortran library is now a part of the grib2io package. The separate component package, grib2io-interp, has been deprecated. grib2io-interp provided interpolation via F2PY interface to NCEPLIBS-ip, which has become difficult since the removal of distutils from Python 3.12+.

NCEPLIBS-ip interpolation Fortran subroutines contain the BIND(C) attribute which provides an equivalent C-interface. grib2io now provides a Cython-based interface, iplib, to these Fortran subroutines via their C-interface. If NCEPLIBS-ip was built with OpenMP support, iplib will provide functions for getting and setting the number of OpenMP threads.

Xarray Backend

grib2io provides a Xarray backend engine so that many GRIB2 messages can be represented as N-dimensional DataArray objects and collected along common coordinates as Datasets or DataTrees. The Xarray backend engine API is experimental and is subject to change without backward compatibility.

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', 'msgs_from_index']
 8
 9try:
10    from . import __config__
11    __version__ = __config__.grib2io_version
12    has_interpolation = __config__.has_interpolation
13    has_openmp_support = __config__.has_openmp_support
14    g2c_static = __config__.g2c_static
15    ip_static = __config__.ip_static
16    extra_objects = __config__.extra_objects
17except(ImportError):
18    pass
19
20from .g2clib import __version__ as __g2clib_version__
21from .g2clib import _has_jpeg
22from .g2clib import _has_png
23from .g2clib import _has_aec
24
25has_jpeg_support = bool(_has_jpeg)
26has_png_support  = bool(_has_png)
27has_aec_support = bool(_has_aec)
28
29from .tables.originating_centers import _ncep_grib2_table_version
30ncep_grib2_table_version = _ncep_grib2_table_version
31g2c_version = __g2clib_version__
32
33def show_config():
34    """Print grib2io build configuration information."""
35    print(f'grib2io version {__version__} Configuration:')
36    print(f'')
37    print(f'NCEPLIBS-g2c library version: {__g2clib_version__}')
38    print(f'\tStatic library: {g2c_static}')
39    print(f'\tJPEG compression support: {has_jpeg_support}')
40    print(f'\tPNG compression support: {has_png_support}')
41    print(f'\tAEC compression support: {has_aec_support}')
42    print(f'')
43    print(f'NCEPLIPS-ip support: {has_interpolation}')
44    print(f'\tStatic library: {ip_static}')
45    print(f'\tOpenMP support: {has_openmp_support}')
46    print(f'')
47    print(f'Static libs:')
48    for lib in extra_objects:
49        print(f'\t{lib}')
50    print(f'')
51    print(f'NCEP GRIB2 Table Version: {_ncep_grib2_table_version}')
class open:
 95class open():
 96    """
 97    GRIB2 File Object.
 98
 99    This class can accommodate a physical file with one or more GRIB2
100    messages or a bytes object containing a GRIB2 messages.
101
102    A physical file can contain one or more GRIB2 messages.  When instantiated,
103    class `grib2io.open`, the file named `filename` is opened for reading (`mode
104    = 'r'`) and is automatically indexed.  The indexing procedure reads some of
105    the GRIB2 metadata for all GRIB2 Messages.  A GRIB2 Message may contain
106    submessages whereby Section 2-7 can be repeated.  grib2io accommodates for
107    this by flattening any GRIB2 submessages into multiple individual messages.
108
109    It is important to note that GRIB2 files from some Meteorological agencies
110    contain other data than GRIB2 messages.  GRIB2 files from ECMWF can contain
111    GRIB1 and GRIB2 messages.  grib2io checks for these and safely ignores them.
112
113    Attributes
114    ----------
115    closed : bool
116        `True` is file handle is close; `False` otherwise.
117    current_message : int
118        Current position of the file in units of GRIB2 Messages. (read only)
119    levels : tuple
120        Tuple containing a unique list of wgrib2-formatted level/layer strings.
121    messages : int
122        Count of GRIB2 Messages contained in the file.
123    mode : str
124        File IO mode of opening the file.
125    name : str
126        Full path name of the GRIB2 file.
127    save_index : bool
128        Whether to save a pickle-based index file for the GRIB2 file. Default is `True`.
129    size : int
130        Size of the file in units of bytes.
131    use_index
132        Whether to use an existing pickle-based index file for the GRIB2 file. Default is `True`.
133    variables : tuple
134        Tuple containing a unique list of variable short names (i.e. GRIB2
135        abbreviation names).
136    """
137
138    __slots__ = ('_fileid', '_filehandle', '_hasindex', '_index',
139                 '_pos', 'closed', 'current_message', 'messages', 'mode',
140                 'name', 'size', 'save_index', 'use_index', '_msgs')
141
142    def __init__(
143        self,
144        filename: Union[bytes, str, Path],
145        mode: Literal["r", "w", "x"] = "r",
146        *,
147        save_index = True,
148        use_index = True,
149        _xarray_backend = False,
150        **kwargs,
151        ):
152        """
153        Initialize GRIB2 File object instance.
154
155        Parameters
156        ----------
157        filename
158            File path containing GRIB2 messages OR bytes.
159        mode
160            File access mode where "r" opens the files for reading only; "w"
161            opens the file for overwriting and "x" for writing to a new file.
162        save_index
163            Whether to save a pickle-based index file for the GRIB2 file. Default is True.
164        use_index
165            Whether to use an existing pickle-based index file for the GRIB2 file. Default is True.
166        """
167
168        # All write modes are read/write.
169        # All modes are binary.
170        if mode in ("a", "x", "w"):
171            mode += "+"
172        mode = mode + "b"
173
174        self.mode = mode
175        self.save_index = save_index
176        self.use_index = use_index
177
178        if isinstance(filename, bytes):
179            if 'r' not in self.mode:
180                raise ValueError("bytes grib2 can only be opened as read")
181            self.current_message = 0
182            if filename[:2] == _GZIP_HEADER:
183                filename = gzip.decompress(filename)
184            if filename[:4]+filename[-4:] != b'GRIB7777':
185                raise ValueError("Invalid GRIB bytes")
186            self._filehandle = BytesIO(filename)
187            self.name = "<in-memory-file>"
188            self.size = len(filename)
189            self._fileid = hashlib.sha1((self.name+str(self.size)).encode('ASCII')).hexdigest()
190            self._index = build_index(self._filehandle)
191            self._msgs = msgs_from_index(self._index, filehandle=self._filehandle)
192            self.messages = len(self._msgs)
193
194        else:
195            self.current_message = 0
196            if 'r' in mode:
197                self._filehandle = builtins.open(filename, mode=mode)
198                # Some GRIB2 files are gzipped, so check for that here, but
199                # raise error when using xarray backend.
200                # Gzip files contain a 2-byte header b'\x1f\x8b'.
201                if self._filehandle.read(2) == _GZIP_HEADER:
202                    self._filehandle.close()
203                    if _xarray_backend:
204                        raise RuntimeError('Gzip GRIB2 files are not supported by the Xarray backend.')
205                    self._filehandle = gzip.open(filename, mode=mode)
206                else:
207                    self._filehandle = builtins.open(filename, mode=mode)
208            else:
209                self._filehandle = builtins.open(filename, mode=mode)
210                self.name = os.path.abspath(filename)
211            self.name = os.path.abspath(filename)
212            fstat = os.stat(self.name)
213            self.size = fstat.st_size
214            self._fileid = hashlib.sha1((self.name+str(fstat.st_ino)+
215                                     str(self.size)).encode('ASCII')).hexdigest()
216
217            if 'r' in self.mode:
218
219                indexfile = f"{self.name}_{self._fileid}.grib2ioidx"
220                if self.use_index and os.path.exists(indexfile):
221                    try:
222                        with builtins.open(indexfile, 'rb') as file:
223                            index = pickle.load(file)
224                        self._index = index
225                    except Exception as e:
226                        print(f"found indexfile: {indexfile}, but unable to load it: {e}")
227                        print(f"re-forming index from grib2file, but not writing indexfile")
228                        self._index = build_index(self._filehandle)
229                else:
230                    self._index = build_index(self._filehandle)
231                    if self.save_index:
232                        try:
233                            serialize_index(self._index, indexfile)
234                        except Exception as e:
235                            print(f"index was not serialized for future use: {e}")
236
237                self._msgs = msgs_from_index(self._index, filehandle=self._filehandle)
238
239                self.messages = len(self._msgs)
240            elif 'w' or 'x' in self.mode:
241                self.messages = 0
242                self.current_message = None
243
244        self.closed = self._filehandle.closed
245
246
247    def __delete__(self, instance):
248        self.close()
249        del self._index
250
251
252    def __enter__(self):
253        return self
254
255
256    def __exit__(self, atype, value, traceback):
257        self.close()
258
259
260    def __iter__(self):
261        yield from self._msgs
262
263
264    def __len__(self):
265        return self.messages
266
267
268    def __repr__(self):
269        strings = []
270        for k in self.__slots__:
271            if k.startswith('_'): continue
272            strings.append('%s = %s\n'%(k,eval('self.'+k)))
273        return ''.join(strings)
274
275
276    def __getitem__(self, key):
277        if isinstance(key,int):
278            if abs(key) >= len(self._msgs):
279                raise IndexError("index out of range")
280            else:
281                return self._msgs[key]
282        elif isinstance(key,str):
283            return self.select(shortName=key)
284        elif isinstance(key,slice):
285            return self._msgs[key]
286        else:
287            raise KeyError('Key must be an integer, slice, or GRIB2 variable shortName.')
288
289
290
291    @property
292    def levels(self):
293        """Provides a unique tuple of level strings."""
294        if self._hasindex:
295            return tuple(sorted(set([msg.level for msg in self._msgs])))
296        else:
297            return None
298
299
300    @property
301    def variables(self):
302        """Provides a unique tuple of variable shortName strings."""
303        if self._hasindex:
304            return tuple(sorted(set([msg.shortName for msg in self._msgs])))
305        else:
306            return None
307
308
309    def close(self):
310        """Close the file handle."""
311        if not self._filehandle.closed:
312            self.messages = 0
313            self.current_message = 0
314            self._filehandle.close()
315            self.closed = self._filehandle.closed
316
317
318    def read(self, size: Optional[int]=None):
319        """
320        Read size amount of GRIB2 messages from the current position.
321
322        If no argument is given, then size is None and all messages are returned
323        from the current position in the file. This read method follows the
324        behavior of Python's builtin open() function, but whereas that operates
325        on units of bytes, we operate on units of GRIB2 messages.
326
327        Parameters
328        ----------
329        size: default=None
330            The number of GRIB2 messages to read from the current position. If
331            no argument is give, the default value is None and remainder of
332            the file is read.
333
334        Returns
335        -------
336        read
337            ``Grib2Message`` object when size = 1 or a list of Grib2Messages
338            when size > 1.
339        """
340        if size is not None and size < 0:
341            size = None
342        if size is None or size > 1:
343            start = self.tell()
344            stop = self.messages if size is None else start+size
345            if size is None:
346                self.current_message = self.messages-1
347            else:
348                self.current_message += size
349            return self._msgs[slice(start,stop,1)]
350        elif size == 1:
351            self.current_message += 1
352            return self._msgs[self.current_message]
353        else:
354            None
355
356
357    def seek(self, pos: int):
358        """
359        Set the position within the file in units of GRIB2 messages.
360
361        Parameters
362        ----------
363        pos
364            The GRIB2 Message number to set the file pointer to.
365        """
366        if self._hasindex:
367            self._filehandle.seek(self._index['sectionOffset'][0][pos])
368            self.current_message = pos
369
370
371    def tell(self):
372        """Returns the position of the file in units of GRIB2 Messages."""
373        return self.current_message
374
375
376    def select(self, **kwargs):
377        """Select GRIB2 messages by `Grib2Message` attributes."""
378        # TODO: Added ability to process multiple values for each keyword (attribute)
379        idxs = []
380        nkeys = len(kwargs.keys())
381        for k,v in kwargs.items():
382            for m in self._msgs:
383                if hasattr(m,k) and getattr(m,k) == v: idxs.append(m._msgnum)
384        idxs = np.array(idxs,dtype=np.int32)
385        return [self._msgs[i] for i in [ii[0] for ii in collections.Counter(idxs).most_common() if ii[1] == nkeys]]
386
387
388    def write(self, msg):
389        """
390        Writes GRIB2 message object to file.
391
392        Parameters
393        ----------
394        msg
395            GRIB2 message objects to write to file.
396        """
397        if isinstance(msg,list):
398            for m in msg:
399                self.write(m)
400            return
401
402        if issubclass(msg.__class__,_Grib2Message):
403            # TODO: We can consider letting pack return packed bytes instead of associating with message object
404            if hasattr(msg,'_msg'):
405                # write already packed bytes
406                self._filehandle.write(msg._msg)
407            else:
408                if msg._signature == msg._generate_signature() and \
409                   msg._data is None and \
410                   hasattr(msg._ondiskarray,'filehandle'):
411                   # write unchanged message from input
412                   offset = msg._ondiskarray.filehandle.tell()
413                   msg._ondiskarray.filehandle.seek(msg._ondiskarray.offset)
414                   self._filehandle.write(msg._ondiskarray.filehandle.read(msg.section0[-1]))
415                   msg._ondiskarray.filehandle.seek(offset)
416                else:
417                    msg.pack()
418                    self._filehandle.write(msg._msg)
419            self.size = os.path.getsize(self.name)
420            self.messages += 1
421        else:
422            raise TypeError("msg must be a Grib2Message object.")
423        return
424
425
426    def flush(self):
427        """Flush the file object buffer."""
428        self._filehandle.flush()
429
430
431    def levels_by_var(self, name: str):
432        """
433        Return a list of level strings given a variable shortName.
434
435        Parameters
436        ----------
437        name
438            Grib2Message variable shortName
439
440        Returns
441        -------
442        levels_by_var
443            A list of unique level strings.
444        """
445        return list(sorted(set([msg.level for msg in self.select(shortName=name)])))
446
447
448    def vars_by_level(self, level: str):
449        """
450        Return a list of variable shortName strings given a level.
451
452        Parameters
453        ----------
454        level
455            Grib2Message variable level
456
457        Returns
458        -------
459        vars_by_level
460            A list of unique variable shortName strings.
461        """
462        return list(sorted(set([msg.shortName for msg in self.select(level=level)])))

GRIB2 File Object.

This class can accommodate a physical file with one or more GRIB2 messages or a bytes object containing a GRIB2 messages.

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. (read only)
  • 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.
  • save_index (bool): Whether to save a pickle-based index file for the GRIB2 file. Default is True.
  • size (int): Size of the file in units of bytes.
  • use_index: Whether to use an existing pickle-based index file for the GRIB2 file. Default is True.
  • variables (tuple): Tuple containing a unique list of variable short names (i.e. GRIB2 abbreviation names).
open( filename: Union[bytes, str, pathlib.Path], mode: Literal['r', 'w', 'x'] = 'r', *, save_index=True, use_index=True, _xarray_backend=False, **kwargs)
142    def __init__(
143        self,
144        filename: Union[bytes, str, Path],
145        mode: Literal["r", "w", "x"] = "r",
146        *,
147        save_index = True,
148        use_index = True,
149        _xarray_backend = False,
150        **kwargs,
151        ):
152        """
153        Initialize GRIB2 File object instance.
154
155        Parameters
156        ----------
157        filename
158            File path containing GRIB2 messages OR bytes.
159        mode
160            File access mode where "r" opens the files for reading only; "w"
161            opens the file for overwriting and "x" for writing to a new file.
162        save_index
163            Whether to save a pickle-based index file for the GRIB2 file. Default is True.
164        use_index
165            Whether to use an existing pickle-based index file for the GRIB2 file. Default is True.
166        """
167
168        # All write modes are read/write.
169        # All modes are binary.
170        if mode in ("a", "x", "w"):
171            mode += "+"
172        mode = mode + "b"
173
174        self.mode = mode
175        self.save_index = save_index
176        self.use_index = use_index
177
178        if isinstance(filename, bytes):
179            if 'r' not in self.mode:
180                raise ValueError("bytes grib2 can only be opened as read")
181            self.current_message = 0
182            if filename[:2] == _GZIP_HEADER:
183                filename = gzip.decompress(filename)
184            if filename[:4]+filename[-4:] != b'GRIB7777':
185                raise ValueError("Invalid GRIB bytes")
186            self._filehandle = BytesIO(filename)
187            self.name = "<in-memory-file>"
188            self.size = len(filename)
189            self._fileid = hashlib.sha1((self.name+str(self.size)).encode('ASCII')).hexdigest()
190            self._index = build_index(self._filehandle)
191            self._msgs = msgs_from_index(self._index, filehandle=self._filehandle)
192            self.messages = len(self._msgs)
193
194        else:
195            self.current_message = 0
196            if 'r' in mode:
197                self._filehandle = builtins.open(filename, mode=mode)
198                # Some GRIB2 files are gzipped, so check for that here, but
199                # raise error when using xarray backend.
200                # Gzip files contain a 2-byte header b'\x1f\x8b'.
201                if self._filehandle.read(2) == _GZIP_HEADER:
202                    self._filehandle.close()
203                    if _xarray_backend:
204                        raise RuntimeError('Gzip GRIB2 files are not supported by the Xarray backend.')
205                    self._filehandle = gzip.open(filename, mode=mode)
206                else:
207                    self._filehandle = builtins.open(filename, mode=mode)
208            else:
209                self._filehandle = builtins.open(filename, mode=mode)
210                self.name = os.path.abspath(filename)
211            self.name = os.path.abspath(filename)
212            fstat = os.stat(self.name)
213            self.size = fstat.st_size
214            self._fileid = hashlib.sha1((self.name+str(fstat.st_ino)+
215                                     str(self.size)).encode('ASCII')).hexdigest()
216
217            if 'r' in self.mode:
218
219                indexfile = f"{self.name}_{self._fileid}.grib2ioidx"
220                if self.use_index and os.path.exists(indexfile):
221                    try:
222                        with builtins.open(indexfile, 'rb') as file:
223                            index = pickle.load(file)
224                        self._index = index
225                    except Exception as e:
226                        print(f"found indexfile: {indexfile}, but unable to load it: {e}")
227                        print(f"re-forming index from grib2file, but not writing indexfile")
228                        self._index = build_index(self._filehandle)
229                else:
230                    self._index = build_index(self._filehandle)
231                    if self.save_index:
232                        try:
233                            serialize_index(self._index, indexfile)
234                        except Exception as e:
235                            print(f"index was not serialized for future use: {e}")
236
237                self._msgs = msgs_from_index(self._index, filehandle=self._filehandle)
238
239                self.messages = len(self._msgs)
240            elif 'w' or 'x' in self.mode:
241                self.messages = 0
242                self.current_message = None
243
244        self.closed = self._filehandle.closed

Initialize GRIB2 File object instance.

Parameters
  • filename: File path containing GRIB2 messages OR bytes.
  • mode: 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.
  • save_index: Whether to save a pickle-based index file for the GRIB2 file. Default is True.
  • use_index: Whether to use an existing pickle-based index file for the GRIB2 file. Default is True.
mode
save_index
use_index
closed
levels
291    @property
292    def levels(self):
293        """Provides a unique tuple of level strings."""
294        if self._hasindex:
295            return tuple(sorted(set([msg.level for msg in self._msgs])))
296        else:
297            return None

Provides a unique tuple of level strings.

variables
300    @property
301    def variables(self):
302        """Provides a unique tuple of variable shortName strings."""
303        if self._hasindex:
304            return tuple(sorted(set([msg.shortName for msg in self._msgs])))
305        else:
306            return None

Provides a unique tuple of variable shortName strings.

def close(self):
309    def close(self):
310        """Close the file handle."""
311        if not self._filehandle.closed:
312            self.messages = 0
313            self.current_message = 0
314            self._filehandle.close()
315            self.closed = self._filehandle.closed

Close the file handle.

def read(self, size: Optional[int] = None):
318    def read(self, size: Optional[int]=None):
319        """
320        Read size amount of GRIB2 messages from the current position.
321
322        If no argument is given, then size is None and all messages are returned
323        from the current position in the file. This read method follows the
324        behavior of Python's builtin open() function, but whereas that operates
325        on units of bytes, we operate on units of GRIB2 messages.
326
327        Parameters
328        ----------
329        size: default=None
330            The number of GRIB2 messages to read from the current position. If
331            no argument is give, the default value is None and remainder of
332            the file is read.
333
334        Returns
335        -------
336        read
337            ``Grib2Message`` object when size = 1 or a list of Grib2Messages
338            when size > 1.
339        """
340        if size is not None and size < 0:
341            size = None
342        if size is None or size > 1:
343            start = self.tell()
344            stop = self.messages if size is None else start+size
345            if size is None:
346                self.current_message = self.messages-1
347            else:
348                self.current_message += size
349            return self._msgs[slice(start,stop,1)]
350        elif size == 1:
351            self.current_message += 1
352            return self._msgs[self.current_message]
353        else:
354            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):
357    def seek(self, pos: int):
358        """
359        Set the position within the file in units of GRIB2 messages.
360
361        Parameters
362        ----------
363        pos
364            The GRIB2 Message number to set the file pointer to.
365        """
366        if self._hasindex:
367            self._filehandle.seek(self._index['sectionOffset'][0][pos])
368            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):
371    def tell(self):
372        """Returns the position of the file in units of GRIB2 Messages."""
373        return self.current_message

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

def select(self, **kwargs):
376    def select(self, **kwargs):
377        """Select GRIB2 messages by `Grib2Message` attributes."""
378        # TODO: Added ability to process multiple values for each keyword (attribute)
379        idxs = []
380        nkeys = len(kwargs.keys())
381        for k,v in kwargs.items():
382            for m in self._msgs:
383                if hasattr(m,k) and getattr(m,k) == v: idxs.append(m._msgnum)
384        idxs = np.array(idxs,dtype=np.int32)
385        return [self._msgs[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):
388    def write(self, msg):
389        """
390        Writes GRIB2 message object to file.
391
392        Parameters
393        ----------
394        msg
395            GRIB2 message objects to write to file.
396        """
397        if isinstance(msg,list):
398            for m in msg:
399                self.write(m)
400            return
401
402        if issubclass(msg.__class__,_Grib2Message):
403            # TODO: We can consider letting pack return packed bytes instead of associating with message object
404            if hasattr(msg,'_msg'):
405                # write already packed bytes
406                self._filehandle.write(msg._msg)
407            else:
408                if msg._signature == msg._generate_signature() and \
409                   msg._data is None and \
410                   hasattr(msg._ondiskarray,'filehandle'):
411                   # write unchanged message from input
412                   offset = msg._ondiskarray.filehandle.tell()
413                   msg._ondiskarray.filehandle.seek(msg._ondiskarray.offset)
414                   self._filehandle.write(msg._ondiskarray.filehandle.read(msg.section0[-1]))
415                   msg._ondiskarray.filehandle.seek(offset)
416                else:
417                    msg.pack()
418                    self._filehandle.write(msg._msg)
419            self.size = os.path.getsize(self.name)
420            self.messages += 1
421        else:
422            raise TypeError("msg must be a Grib2Message object.")
423        return

Writes GRIB2 message object to file.

Parameters
  • msg: GRIB2 message objects to write to file.
def flush(self):
426    def flush(self):
427        """Flush the file object buffer."""
428        self._filehandle.flush()

Flush the file object buffer.

def levels_by_var(self, name: str):
431    def levels_by_var(self, name: str):
432        """
433        Return a list of level strings given a variable shortName.
434
435        Parameters
436        ----------
437        name
438            Grib2Message variable shortName
439
440        Returns
441        -------
442        levels_by_var
443            A list of unique level strings.
444        """
445        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):
448    def vars_by_level(self, level: str):
449        """
450        Return a list of variable shortName strings given a level.
451
452        Parameters
453        ----------
454        level
455            Grib2Message variable level
456
457        Returns
458        -------
459        vars_by_level
460            A list of unique variable shortName strings.
461        """
462        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.
current_message
messages
name
size
def show_config():
34def show_config():
35    """Print grib2io build configuration information."""
36    print(f'grib2io version {__version__} Configuration:')
37    print(f'')
38    print(f'NCEPLIBS-g2c library version: {__g2clib_version__}')
39    print(f'\tStatic library: {g2c_static}')
40    print(f'\tJPEG compression support: {has_jpeg_support}')
41    print(f'\tPNG compression support: {has_png_support}')
42    print(f'\tAEC compression support: {has_aec_support}')
43    print(f'')
44    print(f'NCEPLIPS-ip support: {has_interpolation}')
45    print(f'\tStatic library: {ip_static}')
46    print(f'\tOpenMP support: {has_openmp_support}')
47    print(f'')
48    print(f'Static libs:')
49    for lib in extra_objects:
50        print(f'\t{lib}')
51    print(f'')
52    print(f'NCEP 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):
1828def interpolate(a,
1829    method: Union[int, str],
1830    grid_def_in,
1831    grid_def_out,
1832    method_options=None,
1833    num_threads=1,
1834):
1835    """
1836    This is the module-level interpolation function.
1837
1838    This interfaces with the 4-byte (32-bit float) [NCEPLIBS-ip library](https://github.com/NOAA-EMC/NCEPLIBS-ip)
1839    through grib2io's internal iplib Cython extension module.
1840
1841    Parameters
1842    ----------
1843    a : numpy.ndarray or tuple
1844        Input data.  If `a` is a `numpy.ndarray`, scalar interpolation will be
1845        performed.  If `a` is a `tuple`, then vector interpolation will be
1846        performed with the assumption that u = a[0] and v = a[1] and are both
1847        `numpy.ndarray`.
1848
1849        These data are expected to be in 2-dimensional form with shape (ny, nx)
1850        or 3-dimensional (:, ny, nx) where the 1st dimension represents another
1851        spatial, temporal, or classification (i.e. ensemble members) dimension.
1852        The function will properly flatten the (ny,nx) dimensions into (nx * ny)
1853        acceptable for input into the interpolation subroutines. If needed, these
1854        data will be converted to `np.float32`.
1855    method
1856        Interpolate method to use. This can either be an integer or string using
1857        the following mapping:
1858
1859        | Interpolate Scheme | Integer Value |
1860        | :---:              | :---:         |
1861        | 'bilinear'         | 0             |
1862        | 'bicubic'          | 1             |
1863        | 'neighbor'         | 2             |
1864        | 'budget'           | 3             |
1865        | 'spectral'         | 4             |
1866        | 'neighbor-budget'  | 6             |
1867
1868    grid_def_in : grib2io.Grib2GridDef
1869        Grib2GridDef object for the input grid.
1870    grid_def_out : grib2io.Grib2GridDef
1871        Grib2GridDef object for the output grid or station points.
1872    method_options : list of ints, optional
1873        Interpolation options. See the NCEPLIBS-ip documentation for
1874        more information on how these are used.
1875    num_threads : int, optional
1876        Number of OpenMP threads to use for interpolation. The default
1877        value is 1. If NCEPLIBS-ip and grib2io's iplib extension module
1878        was not built with OpenMP, then this keyword argument and value
1879        will have no impact.
1880
1881    Returns
1882    -------
1883    interpolate
1884        Returns a `numpy.ndarray` of dtype `np.float32` when scalar interpolation
1885        is performed or a `tuple` of `numpy.ndarray`s when vector interpolation is
1886        performed with the assumptions that 0-index is the interpolated u and
1887        1-index is the interpolated v.
1888    """
1889
1890    from . import iplib
1891
1892    prev_num_threads = 1
1893    try:
1894        prev_num_threads = iplib.openmp_get_num_threads()
1895        iplib.openmp_set_num_threads(num_threads)
1896    except(AttributeError):
1897        pass
1898
1899    print(f"grib2io.interpolate thread report: OpenMP num threads = {iplib.openmp_get_num_threads()}")
1900
1901    if isinstance(method,int) and method not in _interp_schemes.values():
1902        raise ValueError('Invalid interpolation method.')
1903    elif isinstance(method,str):
1904        if method in _interp_schemes.keys():
1905            method = _interp_schemes[method]
1906        else:
1907            raise ValueError('Invalid interpolation method.')
1908
1909    if method_options is None:
1910        method_options = np.zeros((20),dtype=np.int32)
1911        if method in {3,6}:
1912            method_options[0:2] = -1
1913
1914    mi = grid_def_in.npoints
1915    mo = grid_def_out.npoints
1916
1917    # Adjust shape of input array(s)
1918    a, newshp = _adjust_array_shape_for_interp(a, grid_def_in, grid_def_out)
1919
1920    # Call interpolation subroutines according to type of a.
1921    if isinstance(a,np.ndarray):
1922        # Scalar
1923        km = a.shape[0]
1924        if np.any(np.isnan(a)):
1925            ibi = np.ones((km), dtype=np.int32)
1926            li = np.where(np.isnan(a),0,1).astype(np.uint8)
1927        else:
1928            ibi = np.zeros((km), dtype=np.int32)
1929            li = np.zeros(a.shape,dtype=np.uint8)
1930        no, rlat, rlon, ibo, lo, go, iret = iplib.interpolate_scalar(
1931            method,
1932            method_options,
1933            grid_def_in.gdtn,
1934            np.array(grid_def_in.gdt, dtype=np.int32),
1935            grid_def_out.gdtn,
1936            np.array(grid_def_out.gdt, dtype=np.int32),
1937            mi,
1938            mo,
1939            km,
1940            ibi,
1941            li,
1942            a,
1943        )
1944        out = np.where(lo==0,np.nan,go).reshape(newshp)
1945    elif isinstance(a,tuple):
1946        # Vector
1947        km = a[0].shape[0]
1948        if np.any(np.isnan(a)):
1949            ibi = np.ones((km), dtype=np.int32)
1950            li = np.where(np.isnan(a),0,1).astype(np.uint8)
1951        else:
1952            ibi = np.zeros((km), dtype=np.int32)
1953            li = np.zeros(a[0].shape,dtype=np.uint8)
1954        no, rlat, rlon, crot, srot, ibo, lo, uo, vo, iret = iplib.interpolate_vector(
1955            method,
1956            method_options,
1957            grid_def_in.gdtn,
1958            np.array(grid_def_in.gdt, dtype=np.int32),
1959            grid_def_out.gdtn,
1960            np.array(grid_def_out.gdt, dtype=np.int32),
1961            mi,
1962            mo,
1963            km,
1964            ibi,
1965            li,
1966            a[0],
1967            a[1],
1968        )
1969        uo = np.where(lo==0,np.nan,uo).reshape(newshp)
1970        vo = np.where(lo==0,np.nan,vo).reshape(newshp)
1971        out = (uo,vo)
1972
1973    try:
1974        iplib.openmp_set_num_threads(prev_num_threads)
1975    except(AttributeError):
1976        pass
1977
1978    return out

This is the module-level interpolation function.

This interfaces with the 4-byte (32-bit float) NCEPLIBS-ip library through grib2io's internal iplib Cython extension module.

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. If needed, these data will be converted to np.float32.

  • 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 NCEPLIBS-ip and grib2io's iplib extension module was not built with OpenMP, then this keyword argument and value will have no impact.
Returns
  • interpolate: Returns a numpy.ndarray of dtype np.float32 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: Union[int, str], grid_def_in, lats, lons, method_options=None, num_threads=1):
1981def interpolate_to_stations(
1982    a,
1983    method: Union[int, str],
1984    grid_def_in,
1985    lats,
1986    lons,
1987    method_options=None,
1988    num_threads=1
1989):
1990    """
1991    Module-level interpolation function for interpolation to stations.
1992
1993    Interfaces with the 4-byte (32-bit float) [NCEPLIBS-ip library](https://github.com/NOAA-EMC/NCEPLIBS-ip)
1994    via grib2io's iplib Cython exntension module. It supports scalar and
1995    vector interpolation according to the type of object a.
1996
1997    Parameters
1998    ----------
1999    a : numpy.ndarray or tuple
2000        Input data.  If `a` is a `numpy.ndarray`, scalar interpolation will be
2001        performed.  If `a` is a `tuple`, then vector interpolation will be
2002        performed with the assumption that u = a[0] and v = a[1] and are both
2003        `numpy.ndarray`.
2004
2005        These data are expected to be in 2-dimensional form with shape (ny, nx)
2006        or 3-dimensional (:, ny, nx) where the 1st dimension represents another
2007        spatial, temporal, or classification (i.e. ensemble members) dimension.
2008        The function will properly flatten the (ny,nx) dimensions into (nx * ny)
2009        acceptable for input into the interpolation subroutines. If needed, these
2010        data will be converted to `np.float32`.
2011    method
2012        Interpolate method to use. This can either be an integer or string using
2013        the following mapping:
2014
2015        | Interpolate Scheme | Integer Value |
2016        | :---:              | :---:         |
2017        | 'bilinear'         | 0             |
2018        | 'bicubic'          | 1             |
2019        | 'neighbor'         | 2             |
2020        | 'budget'           | 3             |
2021        | 'spectral'         | 4             |
2022        | 'neighbor-budget'  | 6             |
2023
2024    grid_def_in : grib2io.Grib2GridDef
2025        Grib2GridDef object for the input grid.
2026    lats : numpy.ndarray or list
2027        Latitudes for station points
2028    lons : numpy.ndarray or list
2029        Longitudes for station points
2030    method_options : list of ints, optional
2031        Interpolation options. See the NCEPLIBS-ip documentation for
2032        more information on how these are used.
2033    num_threads : int, optional
2034        Number of OpenMP threads to use for interpolation. The default
2035        value is 1. If NCEPLIBS-ip and grib2io's iplib extension module
2036        was not built with OpenMP, then this keyword argument and value
2037        will have no impact.
2038
2039    Returns
2040    -------
2041    interpolate_to_stations
2042        Returns a `numpy.ndarray` of dtype `np.float32` when scalar
2043        interpolation is performed or a `tuple` of `numpy.ndarray`s
2044        when vector interpolation is performed with the assumptions
2045        that 0-index is the interpolated u and 1-index is the
2046        interpolated v.
2047    """
2048    from . import iplib
2049
2050    prev_num_threads = 1
2051    try:
2052        prev_num_threads = iplib.openmp_get_num_threads()
2053        iplib.openmp_set_num_threads(num_threads)
2054    except(AttributeError):
2055        pass
2056
2057    if isinstance(method,int) and method not in _interp_schemes.values():
2058        raise ValueError('Invalid interpolation method.')
2059    elif isinstance(method,str):
2060        if method in _interp_schemes.keys():
2061            method = _interp_schemes[method]
2062        else:
2063            raise ValueError('Invalid interpolation method.')
2064
2065    if method_options is None:
2066        method_options = np.zeros((20),dtype=np.int32)
2067        if method in {3,6}:
2068            method_options[0:2] = -1
2069
2070    # Check lats and lons
2071    if isinstance(lats,list):
2072        nlats = len(lats)
2073    elif isinstance(lats,np.ndarray) and len(lats.shape) == 1:
2074        nlats = lats.shape[0]
2075    else:
2076        raise ValueError('Station latitudes must be a list or 1-D NumPy array.')
2077    if isinstance(lons,list):
2078        nlons = len(lons)
2079    elif isinstance(lons,np.ndarray) and len(lons.shape) == 1:
2080        nlons = lons.shape[0]
2081    else:
2082        raise ValueError('Station longitudes must be a list or 1-D NumPy array.')
2083    if nlats != nlons:
2084        raise ValueError('Station lats and lons must be same size.')
2085
2086    mi = grid_def_in.npoints
2087    mo = nlats
2088
2089    # Adjust shape of input array(s)
2090    a, newshp = _adjust_array_shape_for_interp_stations(a, grid_def_in, mo)
2091
2092    # Use gdtn = -1 for stations and an empty template array
2093    gdtn = -1
2094    gdt = np.zeros((200),dtype=np.int32)
2095
2096    # Call interpolation subroutines according to type of a.
2097    if isinstance(a,np.ndarray):
2098        # Scalar
2099        km = a.shape[0]
2100        ibi = np.zeros((km), dtype=np.int32)
2101        li = np.zeros(a.shape,dtype=np.uint8)
2102        no, rlat, rlon, ibo, lo, go, iret = iplib.interpolate_scalar(
2103            method,
2104            method_options,
2105            grid_def_in.gdtn,
2106            np.array(grid_def_in.gdt, dtype=np.int32),
2107            gdtn,
2108            gdt,
2109            mi,
2110            mo,
2111            km,
2112            ibi,
2113            li,
2114            a,
2115            lats=np.array(lats, dtype=np.float32),
2116            lons=np.array(lons, dtype=np.float32),
2117        )
2118        out = go.reshape(newshp)
2119
2120    elif isinstance(a,tuple):
2121        # Vector
2122        km = a[0].shape[0]
2123        ibi = np.zeros((km), dtype=np.int32)
2124        li = np.zeros(a[0].shape,dtype=np.uint8)
2125        no, rlat, rlon, crot, srot, ibo, lo, uo, vo, iret = iplib.interpolate_vector(
2126            method,
2127            method_options,
2128            grid_def_in.gdtn,
2129            np.array(grid_def_in.gdt, dtype=np.int32),
2130            gdtn,
2131            gdt,
2132            mi,
2133            mo,
2134            km,
2135            ibi,
2136            li,
2137            a[0],
2138            a[1],
2139            lats=np.array(lats, dtype=np.float32),
2140            lons=np.array(lons, dtype=np.float32),
2141        )
2142        out = (uo.reshape(newshp),vo.reshape(newshp))
2143
2144    try:
2145        iplib.openmp_set_num_threads(prev_num_threads)
2146    except(AttributeError):
2147        pass
2148
2149    return out

Module-level interpolation function for interpolation to stations.

Interfaces with the 4-byte (32-bit float) NCEPLIBS-ip library via grib2io's iplib Cython exntension module. 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. If needed, these data will be converted to np.float32.

  • 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 NCEPLIBS-ip and grib2io's iplib extension module was not built with OpenMP, then this keyword argument and value will have no impact.
Returns
  • interpolate_to_stations: Returns a numpy.ndarray of dtype np.float32 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:
683class Grib2Message:
684    """
685    Creation class for a GRIB2 message.
686
687    This class returns a dynamically-created Grib2Message object that
688    inherits from `_Grib2Message` and grid, product, data representation
689    template classes according to the template numbers for the respective
690    sections. If `section3`, `section4`, or `section5` are omitted, then
691    the appropriate keyword arguments for the template number `gdtn=`,
692    `pdtn=`, or `drtn=` must be provided.
693
694    Parameters
695    ----------
696    section0
697        GRIB2 section 0 array.
698    section1
699        GRIB2 section 1 array.
700    section2
701        Local Use section data.
702    section3
703        GRIB2 section 3 array.
704    section4
705        GRIB2 section 4 array.
706    section5
707        GRIB2 section 5 array.
708
709    Returns
710    -------
711    Msg
712        A dynamically-create Grib2Message object that inherits from
713        _Grib2Message, a grid definition template class, product
714        definition template class, and a data representation template
715        class.
716    """
717    def __new__(self, section0: NDArray = np.array([struct.unpack('>I',b'GRIB')[0],0,0,2,0]),
718                      section1: NDArray = np.zeros((13),dtype=np.int64),
719                      section2: Optional[bytes] = None,
720                      section3: Optional[NDArray] = None,
721                      section4: Optional[NDArray] = None,
722                      section5: Optional[NDArray] = None, *args, **kwargs):
723
724        if np.all(section1==0):
725            try:
726                # Python >= 3.10
727                section1[5:11] = datetime.datetime.fromtimestamp(0, datetime.UTC).timetuple()[:6]
728            except(AttributeError):
729                # Python < 3.10
730                section1[5:11] = datetime.datetime.utcfromtimestamp(0).timetuple()[:6]
731
732        bases = list()
733        if section3 is None:
734            if 'gdtn' in kwargs.keys():
735                gdtn = kwargs['gdtn']
736                Gdt = templates.gdt_class_by_gdtn(gdtn)
737                bases.append(Gdt)
738                section3 = np.zeros((Gdt._len+5),dtype=np.int64)
739                section3[4] = gdtn
740            else:
741                raise ValueError("Must provide GRIB2 Grid Definition Template Number or section 3 array")
742        else:
743            gdtn = section3[4]
744            Gdt = templates.gdt_class_by_gdtn(gdtn)
745            bases.append(Gdt)
746
747        if section4 is None:
748            if 'pdtn' in kwargs.keys():
749                pdtn = kwargs['pdtn']
750                Pdt = templates.pdt_class_by_pdtn(pdtn)
751                bases.append(Pdt)
752                section4 = np.zeros((Pdt._len+2),dtype=np.int64)
753                section4[1] = pdtn
754            else:
755                raise ValueError("Must provide GRIB2 Production Definition Template Number or section 4 array")
756        else:
757            pdtn = section4[1]
758            Pdt = templates.pdt_class_by_pdtn(pdtn)
759            bases.append(Pdt)
760
761        if section5 is None:
762            if 'drtn' in kwargs.keys():
763                drtn = kwargs['drtn']
764                Drt = templates.drt_class_by_drtn(drtn)
765                bases.append(Drt)
766                section5 = np.zeros((Drt._len+2),dtype=np.int64)
767                section5[1] = drtn
768            else:
769                raise ValueError("Must provide GRIB2 Data Representation Template Number or section 5 array")
770        else:
771            drtn = section5[1]
772            Drt = templates.drt_class_by_drtn(drtn)
773            bases.append(Drt)
774
775        # attempt to use existing Msg class if it has already been made with gdtn,pdtn,drtn combo
776        try:
777            Msg = _msg_class_store[f"{gdtn}:{pdtn}:{drtn}"]
778        except KeyError:
779            @dataclass(init=False, repr=False)
780            class Msg(_Grib2Message, *bases):
781                pass
782            _msg_class_store[f"{gdtn}:{pdtn}:{drtn}"] = Msg
783
784        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:
 787@dataclass
 788class _Grib2Message:
 789    """
 790    GRIB2 Message base class.
 791    """
 792    # GRIB2 Sections
 793    section0: NDArray = field(init=True,repr=False)
 794    section1: NDArray = field(init=True,repr=False)
 795    section2: bytes = field(init=True,repr=False)
 796    section3: NDArray = field(init=True,repr=False)
 797    section4: NDArray = field(init=True,repr=False)
 798    section5: NDArray = field(init=True,repr=False)
 799    bitMapFlag: templates.Grib2Metadata = field(init=True,repr=False,default=255)
 800
 801    # Section 0 looked up attributes
 802    indicatorSection: NDArray = field(init=False,repr=False,default=templates.IndicatorSection())
 803    discipline: templates.Grib2Metadata = field(init=False,repr=False,default=templates.Discipline())
 804
 805    # Section 1 looked up attributes
 806    identificationSection: NDArray = field(init=False,repr=False,default=templates.IdentificationSection())
 807    originatingCenter: templates.Grib2Metadata = field(init=False,repr=False,default=templates.OriginatingCenter())
 808    originatingSubCenter: templates.Grib2Metadata = field(init=False,repr=False,default=templates.OriginatingSubCenter())
 809    masterTableInfo: templates.Grib2Metadata = field(init=False,repr=False,default=templates.MasterTableInfo())
 810    localTableInfo: templates.Grib2Metadata = field(init=False,repr=False,default=templates.LocalTableInfo())
 811    significanceOfReferenceTime: templates.Grib2Metadata = field(init=False,repr=False,default=templates.SignificanceOfReferenceTime())
 812    year: int = field(init=False,repr=False,default=templates.Year())
 813    month: int = field(init=False,repr=False,default=templates.Month())
 814    day: int = field(init=False,repr=False,default=templates.Day())
 815    hour: int = field(init=False,repr=False,default=templates.Hour())
 816    minute: int = field(init=False,repr=False,default=templates.Minute())
 817    second: int = field(init=False,repr=False,default=templates.Second())
 818    refDate: datetime.datetime = field(init=False,repr=False,default=templates.RefDate())
 819    productionStatus: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ProductionStatus())
 820    typeOfData: templates.Grib2Metadata = field(init=False,repr=False,default=templates.TypeOfData())
 821
 822    # Section 3 looked up common attributes.  Other looked up attributes are available according
 823    # to the Grid Definition Template.
 824    gridDefinitionSection: NDArray = field(init=False,repr=False,default=templates.GridDefinitionSection())
 825    sourceOfGridDefinition: int = field(init=False,repr=False,default=templates.SourceOfGridDefinition())
 826    numberOfDataPoints: int = field(init=False,repr=False,default=templates.NumberOfDataPoints())
 827    interpretationOfListOfNumbers: templates.Grib2Metadata = field(init=False,repr=False,default=templates.InterpretationOfListOfNumbers())
 828    gridDefinitionTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.GridDefinitionTemplateNumber())
 829    gridDefinitionTemplate: list = field(init=False,repr=False,default=templates.GridDefinitionTemplate())
 830    _earthparams: dict = field(init=False,repr=False,default=templates.EarthParams())
 831    _dxsign: float = field(init=False,repr=False,default=templates.DxSign())
 832    _dysign: float = field(init=False,repr=False,default=templates.DySign())
 833    _llscalefactor: float = field(init=False,repr=False,default=templates.LLScaleFactor())
 834    _lldivisor: float = field(init=False,repr=False,default=templates.LLDivisor())
 835    _xydivisor: float = field(init=False,repr=False,default=templates.XYDivisor())
 836    shapeOfEarth: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ShapeOfEarth())
 837    earthShape: str = field(init=False,repr=False,default=templates.EarthShape())
 838    earthRadius: float = field(init=False,repr=False,default=templates.EarthRadius())
 839    earthMajorAxis: float = field(init=False,repr=False,default=templates.EarthMajorAxis())
 840    earthMinorAxis: float = field(init=False,repr=False,default=templates.EarthMinorAxis())
 841    resolutionAndComponentFlags: list = field(init=False,repr=False,default=templates.ResolutionAndComponentFlags())
 842    ny: int = field(init=False,repr=False,default=templates.Ny())
 843    nx: int = field(init=False,repr=False,default=templates.Nx())
 844    scanModeFlags: list = field(init=False,repr=False,default=templates.ScanModeFlags())
 845    projParameters: dict = field(init=False,repr=False,default=templates.ProjParameters())
 846
 847    # Section 4
 848    productDefinitionTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.ProductDefinitionTemplateNumber())
 849    productDefinitionTemplate: NDArray = field(init=False,repr=False,default=templates.ProductDefinitionTemplate())
 850
 851    # Section 5 looked up common attributes.  Other looked up attributes are
 852    # available according to the Data Representation Template.
 853    numberOfPackedValues: int = field(init=False,repr=False,default=templates.NumberOfPackedValues())
 854    dataRepresentationTemplateNumber: templates.Grib2Metadata = field(init=False,repr=False,default=templates.DataRepresentationTemplateNumber())
 855    dataRepresentationTemplate: list = field(init=False,repr=False,default=templates.DataRepresentationTemplate())
 856    typeOfValues: templates.Grib2Metadata = field(init=False,repr=False,default=templates.TypeOfValues())
 857
 858    def __post_init__(self):
 859        """Set some attributes after init."""
 860        self._auto_nans = _AUTO_NANS
 861        self._coordlist = np.zeros((0), dtype=np.float32)
 862        self._data = None
 863        self._deflist = np.zeros((0), dtype=np.int64)
 864        self._msgnum = -1
 865        self._ondiskarray = None
 866        self._orig_section5 = np.copy(self.section5)
 867        self._signature = self._generate_signature()
 868        try:
 869            self._sha1_section3 = hashlib.sha1(self.section3).hexdigest()
 870        except(TypeError):
 871            pass
 872        self.bitMapFlag = templates.Grib2Metadata(self.bitMapFlag,table='6.0')
 873        self.bitmap = None
 874
 875    @property
 876    def _isNDFD(self):
 877        """Check if GRIB2 message is from NWS NDFD"""
 878        return np.all(self.section1[0:2]==[8,65535])
 879
 880    @property
 881    def _isAerosol(self):
 882        """Check if GRIB2 message contains aerosol data"""
 883        is_aero_template = self.productDefinitionTemplateNumber.value in tables.AEROSOL_PDTNS
 884        is_aero_param = ((str(self.parameterCategory) == '13') |
 885                     (str(self.parameterCategory) == '20')) and str(self.parameterNumber) in tables.AEROSOL_PARAMS
 886        # Check table 4.205 aerosol presence
 887        is_aero_type = (str(self.parameterCategory) == '205' and
 888                        str(self.parameterNumber) == '1')
 889        return is_aero_template or is_aero_param or is_aero_type
 890
 891    @property
 892    def gdtn(self):
 893        """Return Grid Definition Template Number"""
 894        return self.section3[4]
 895
 896
 897    @property
 898    def gdt(self):
 899        """Return Grid Definition Template."""
 900        return self.gridDefinitionTemplate
 901
 902
 903    @property
 904    def pdtn(self):
 905        """Return Product Definition Template Number."""
 906        return self.section4[1]
 907
 908
 909    @property
 910    def pdt(self):
 911        """Return Product Definition Template."""
 912        return self.productDefinitionTemplate
 913
 914
 915    @property
 916    def drtn(self):
 917        """Return Data Representation Template Number."""
 918        return self.section5[1]
 919
 920
 921    @property
 922    def drt(self):
 923        """Return Data Representation Template."""
 924        return self.dataRepresentationTemplate
 925
 926
 927    @property
 928    def pdy(self):
 929        """Return the PDY ('YYYYMMDD')."""
 930        return ''.join([str(i) for i in self.section1[5:8]])
 931
 932
 933    @property
 934    def griddef(self):
 935        """Return a Grib2GridDef instance for a GRIB2 message."""
 936        return Grib2GridDef.from_section3(self.section3)
 937
 938
 939    @property
 940    def lats(self):
 941        """Return grid latitudes."""
 942        return self.latlons()[0]
 943
 944
 945    @property
 946    def lons(self):
 947        """Return grid longitudes."""
 948        return self.latlons()[1]
 949
 950    @property
 951    def min(self):
 952        """Return minimum value of data."""
 953        return np.nanmin(self.data)
 954
 955
 956    @property
 957    def max(self):
 958        """Return maximum value of data."""
 959        return np.nanmax(self.data)
 960
 961
 962    @property
 963    def mean(self):
 964        """Return mean value of data."""
 965        return np.nanmean(self.data)
 966
 967
 968    @property
 969    def median(self):
 970        """Return median value of data."""
 971        return np.nanmedian(self.data)
 972
 973
 974    @property
 975    def shape(self):
 976        """Return shape of data."""
 977        return self.griddef.shape
 978
 979
 980    def __repr__(self):
 981        """
 982        Return an unambiguous string representation of the object.
 983
 984        Returns
 985        -------
 986        repr
 987            A string representation of the object, including information from
 988            sections 0, 1, 3, 4, 5, and 6.
 989        """
 990        info = ''
 991        for sect in [0,1,3,4,5,6]:
 992            for k,v in self.attrs_by_section(sect,values=True).items():
 993                info += f'Section {sect}: {k} = {v}\n'
 994        return info
 995
 996    def __str__(self):
 997        """
 998        Return a readable string representation of the object.
 999
1000        Returns
1001        -------
1002        str
1003            A formatted string representation of the object, including
1004            selected attributes.
1005        """
1006        return (f'{self._msgnum}:d={self.refDate}:{self.shortName}:'
1007                f'{self.fullName} ({self.units}):{self.level}:'
1008                f'{self.leadTime}')
1009
1010
1011    def _generate_signature(self):
1012        """Generature SHA-1 hash string from GRIB2 integer sections."""
1013        return hashlib.sha1(np.concatenate((self.section0,self.section1,
1014                                            self.section3,self.section4,
1015                                            self.section5))).hexdigest()
1016
1017
1018    def attrs_by_section(self, sect: int, values: bool=False):
1019        """
1020        Provide a tuple of attribute names for the given GRIB2 section.
1021
1022        Parameters
1023        ----------
1024        sect
1025            The GRIB2 section number.
1026        values
1027            Optional (default is `False`) argument to return attributes values.
1028
1029        Returns
1030        -------
1031        attrs_by_section
1032            A list of attribute names or dict of name:value pairs if `values =
1033            True`.
1034        """
1035        if sect in {0,1,6}:
1036            attrs = templates._section_attrs[sect]
1037        elif sect in {3,4,5}:
1038            def _find_class_index(n):
1039                _key = {3:'Grid', 4:'Product', 5:'Data'}
1040                for i,c in enumerate(self.__class__.__mro__):
1041                    if _key[n] in c.__name__:
1042                        return i
1043                else:
1044                    return []
1045            if sys.version_info.minor <= 8:
1046                attrs = templates._section_attrs[sect]+\
1047                        [a for a in dir(self.__class__.__mro__[_find_class_index(sect)]) if not a.startswith('_')]
1048            else:
1049                attrs = templates._section_attrs[sect]+\
1050                        self.__class__.__mro__[_find_class_index(sect)]._attrs()
1051        else:
1052            attrs = []
1053        if values:
1054            return {k:getattr(self,k) for k in attrs}
1055        else:
1056            return attrs
1057
1058
1059    def pack(self):
1060        """
1061        Pack GRIB2 section data into a binary message.
1062
1063        It is the user's responsibility to populate the GRIB2 section
1064        information with appropriate metadata.
1065        """
1066        # Create beginning of packed binary message with section 0 and 1 data.
1067        self._sections = []
1068        self._msg,self._pos = g2clib.grib2_create(self.indicatorSection[2:4],self.identificationSection)
1069        self._sections += [0,1]
1070
1071        # Add section 2 if present.
1072        if isinstance(self.section2,bytes) and len(self.section2) > 0:
1073            self._msg,self._pos = g2clib.grib2_addlocal(self._msg,self.section2)
1074            self._sections.append(2)
1075
1076        # Add section 3.
1077        self.section3[1] = self.nx * self.ny
1078        self._msg,self._pos = g2clib.grib2_addgrid(self._msg,self.gridDefinitionSection,
1079                                                   self.gridDefinitionTemplate,
1080                                                   self._deflist)
1081        self._sections.append(3)
1082
1083        # Prepare data.
1084        if self._data is None:
1085            if self._ondiskarray is None:
1086                raise ValueError("Grib2Message object has no data, thus it cannot be packed.")
1087        field = np.copy(self.data)
1088        if self.scanModeFlags is not None:
1089            if self.scanModeFlags[3]:
1090                fieldsave = field.astype('f') # Casting makes a copy
1091                field[1::2,:] = fieldsave[1::2,::-1]
1092        fld = field.astype('f')
1093
1094        # Prepare bitmap, if necessary
1095        bitmapflag = self.bitMapFlag.value
1096        if bitmapflag == 0:
1097            if self.bitmap is not None:
1098                bmap = np.ravel(self.bitmap).astype(DEFAULT_NUMPY_INT)
1099            else:
1100                bmap = np.ravel(np.where(np.isnan(fld),0,1)).astype(DEFAULT_NUMPY_INT)
1101        else:
1102            bmap = None
1103
1104        # Prepare data for packing if nans are present
1105        fld = np.ravel(fld)
1106        if bitmapflag in {0,254}:
1107            fld = np.where(np.isnan(fld),0,fld)
1108        else:
1109            if np.isnan(fld).any():
1110                if hasattr(self,'priMissingValue'):
1111                    fld = np.where(np.isnan(fld),self.priMissingValue,fld)
1112            if hasattr(self,'_missvalmap'):
1113                if hasattr(self,'priMissingValue'):
1114                    fld = np.where(self._missvalmap==1,self.priMissingValue,fld)
1115                if hasattr(self,'secMissingValue'):
1116                    fld = np.where(self._missvalmap==2,self.secMissingValue,fld)
1117
1118        # Add sections 4, 5, 6, and 7.
1119        self._msg,self._pos = g2clib.grib2_addfield(self._msg,self.pdtn,
1120                                                    self.productDefinitionTemplate,
1121                                                    self._coordlist,
1122                                                    self.drtn,
1123                                                    self.dataRepresentationTemplate,
1124                                                    fld,
1125                                                    bitmapflag,
1126                                                    bmap)
1127        self._sections.append(4)
1128        self._sections.append(5)
1129        self._sections.append(6)
1130        self._sections.append(7)
1131
1132        # Finalize GRIB2 message with section 8.
1133        self._msg, self._pos = g2clib.grib2_end(self._msg)
1134        self._sections.append(8)
1135        self.section0[-1] = len(self._msg)
1136
1137
1138    @property
1139    def data(self) -> np.array:
1140        """Access the unpacked data values."""
1141        if self._data is None:
1142            if self._auto_nans != _AUTO_NANS:
1143                self._data = self._ondiskarray
1144            self._data = np.asarray(self._ondiskarray)
1145        return self._data
1146
1147
1148    @data.setter
1149    def data(self, arr):
1150        """
1151        Set the internal data array, enforcing shape (ny, nx) and dtype float32.
1152
1153        If the Grid Definition Section (section 3) of Grib2Message object is
1154        not fully formed (i.e. nx, ny = 0, 0), then the shape of the data array
1155        will be used to set nx and ny of the Grib2Message object. It will be the
1156        responsibility of the user to populate the rest of the Grid Definition
1157        Section attributes.
1158
1159        Parameters
1160        ----------
1161        arr : array_like
1162            A 2D array whose shape must match ``(self.ny, self.nx)``.
1163            It will be converted to ``float32`` and C-contiguous if needed.
1164
1165        Raises
1166        ------
1167        ValueError
1168            If the shape of `arr` does not match the expected dimensions.
1169        """
1170        if not isinstance(arr, np.ndarray):
1171            raise ValueError(f"Grib2Message data only supports numpy arrays")
1172        if self.nx == 0 and self.ny == 0:
1173            self.ny = arr.shape[0]
1174            self.nx = arr.shape[1]
1175        if arr.shape != (self.ny, self.nx):
1176            raise ValueError(
1177                f"Data shape mismatch: expected ({self.ny}, {self.nx}), "
1178                f"got {arr.shape}"
1179            )
1180        # Ensure contiguous memory layout (important for C interoperability)
1181        if not arr.flags["C_CONTIGUOUS"]:
1182            arr = np.ascontiguousarray(arr, dtype=np.float32)
1183        self._data = arr
1184
1185
1186    def flush_data(self):
1187        """
1188        Flush the unpacked data values from the Grib2Message object.
1189
1190        Notes
1191        -----
1192        If the Grib2Message object was constructed from "scratch" (i.e.
1193        not read from file), this method will remove the data array from
1194        the object and it cannot be recovered.
1195        """
1196        self._data = None
1197        self.bitmap = None
1198
1199
1200    def __getitem__(self, item):
1201        return self.data[item]
1202
1203
1204    def __setitem__(self, item):
1205        raise NotImplementedError('assignment of data not supported via setitem')
1206
1207
1208    def latlons(self, *args, **kwrgs):
1209        """Alias for `grib2io.Grib2Message.grid` method."""
1210        return self.grid(*args, **kwrgs)
1211
1212
1213    def grid(self, unrotate: bool=True):
1214        """
1215        Return lats,lons (in degrees) of grid.
1216
1217        Currently can handle reg. lat/lon,cglobal Gaussian, mercator,
1218        stereographic, lambert conformal, albers equal-area, space-view and
1219        azimuthal equidistant grids.
1220
1221        Parameters
1222        ----------
1223        unrotate
1224            If `True` [DEFAULT], and grid is rotated lat/lon, then unrotate the
1225            grid, otherwise `False`, do not.
1226
1227        Returns
1228        -------
1229        lats, lons : numpy.ndarray
1230            Returns two numpy.ndarrays with dtype=numpy.float32 of grid
1231            latitudes and longitudes in units of degrees.
1232        """
1233        if self._sha1_section3 in _latlon_datastore.keys():
1234            return (_latlon_datastore[self._sha1_section3]['latitude'],
1235                    _latlon_datastore[self._sha1_section3]['longitude'])
1236        gdtn = self.gridDefinitionTemplateNumber.value
1237        gdtmpl = self.gridDefinitionTemplate
1238        reggrid = self.gridDefinitionSection[2] == 0 # This means regular 2-d grid
1239        if gdtn == 0:
1240            # Regular lat/lon grid
1241            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1242            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1243            dlon = self.gridlengthXDirection
1244            dlat = self.gridlengthYDirection
1245            if lon2 < lon1 and dlon < 0: lon1 = -lon1
1246            lats = np.linspace(lat1,lat2,self.ny)
1247            if reggrid:
1248                lons = np.linspace(lon1,lon2,self.nx)
1249            else:
1250                lons = np.linspace(lon1,lon2,self.ny*2)
1251            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1252        elif gdtn == 1: # Rotated Lat/Lon grid
1253            pj = pyproj.Proj(self.projParameters)
1254            lat1,lon1 = self.latitudeFirstGridpoint,self.longitudeFirstGridpoint
1255            lat2,lon2 = self.latitudeLastGridpoint,self.longitudeLastGridpoint
1256            if lon1 > 180.0: lon1 -= 360.0
1257            if lon2 > 180.0: lon2 -= 360.0
1258            lats = np.linspace(lat1,lat2,self.ny)
1259            lons = np.linspace(lon1,lon2,self.nx)
1260            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1261            if unrotate:
1262                from grib2io.utils import rotated_grid
1263                lats,lons = rotated_grid.unrotate(lats,lons,self.anglePoleRotation,
1264                                                  self.latitudeSouthernPole,
1265                                                  self.longitudeSouthernPole)
1266        elif gdtn == 40: # Gaussian grid (only works for global!)
1267            from grib2io.utils.gauss_grid import gaussian_latitudes
1268            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1269            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1270            nlats = self.ny
1271            if not reggrid: # Reduced Gaussian grid.
1272                nlons = 2*nlats
1273                dlon = 360./nlons
1274            else:
1275                nlons = self.nx
1276                dlon = self.gridlengthXDirection
1277            lons = np.linspace(lon1,lon2,nlons)
1278            # Compute Gaussian lats (north to south)
1279            lats = gaussian_latitudes(nlats)
1280            if lat1 < lat2:  # reverse them if necessary
1281                lats = lats[::-1]
1282            lons,lats = np.meshgrid(lons,lats)
1283        elif gdtn in {10,20,30,31,110}:
1284            # Mercator, Lambert Conformal, Stereographic, Albers Equal Area,
1285            # Azimuthal Equidistant
1286            dx,dy = self.gridlengthXDirection, self.gridlengthYDirection
1287            lon1,lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1288            pj = pyproj.Proj(self.projParameters)
1289            llcrnrx, llcrnry = pj(lon1,lat1)
1290            x = llcrnrx+dx*np.arange(self.nx)
1291            y = llcrnry+dy*np.arange(self.ny)
1292            x,y = np.meshgrid(x, y)
1293            lons,lats = pj(x, y, inverse=True)
1294        elif gdtn == 90:
1295            # Satellite Projection
1296            dx = self.gridlengthXDirection
1297            dy = self.gridlengthYDirection
1298            pj = pyproj.Proj(self.projParameters)
1299            x = dx*np.indices((self.ny,self.nx),'f')[1,:,:]
1300            x -= 0.5*x.max()
1301            y = dy*np.indices((self.ny,self.nx),'f')[0,:,:]
1302            y -= 0.5*y.max()
1303            lons,lats = pj(x,y,inverse=True)
1304            # Set lons,lats to 1.e30 where undefined
1305            abslons = np.fabs(lons)
1306            abslats = np.fabs(lats)
1307            lons = np.where(abslons < 1.e20, lons, 1.e30)
1308            lats = np.where(abslats < 1.e20, lats, 1.e30)
1309        elif gdtn == 32769:
1310            # Special NCEP Grid, Rotated Lat/Lon, Arakawa E-Grid (Non-Staggered)
1311            from grib2io.utils import arakawa_rotated_grid
1312            from grib2io.utils.rotated_grid import DEG2RAD
1313            di, dj = 0.0, 0.0
1314            do_180 = False
1315            idir = 1 if self.scanModeFlags[0] == 0 else -1
1316            jdir = -1 if self.scanModeFlags[1] == 0 else 1
1317            grid_oriented = 0 if self.resolutionAndComponentFlags[4] == 0 else 1
1318            do_rot = 1 if grid_oriented == 1 else 0
1319            la1 = self.latitudeFirstGridpoint
1320            lo1 = self.longitudeFirstGridpoint
1321            clon = self.longitudeCenterGridpoint
1322            clat = self.latitudeCenterGridpoint
1323            lasp = clat - 90.0
1324            losp = clon
1325            llat, llon = arakawa_rotated_grid.ll2rot(la1,lo1,lasp,losp)
1326            la2, lo2 = arakawa_rotated_grid.rot2ll(-llat,-llon,lasp,losp)
1327            rlat = -llat
1328            rlon = -llon
1329            if self.nx == 1:
1330                di = 0.0
1331            elif idir == 1:
1332                ti = rlon
1333                while ti < llon:
1334                    ti += 360.0
1335                di = (ti - llon)/float(self.nx-1)
1336            else:
1337                ti = llon
1338                while ti < rlon:
1339                    ti += 360.0
1340                di = (ti - rlon)/float(self.nx-1)
1341            if self.ny == 1:
1342               dj = 0.0
1343            else:
1344                dj = (rlat - llat)/float(self.ny-1)
1345                if dj < 0.0:
1346                    dj = -dj
1347            if idir == 1:
1348                if llon > rlon:
1349                    llon -= 360.0
1350                if llon < 0 and rlon > 0:
1351                    do_180 = True
1352            else:
1353                if rlon > llon:
1354                    rlon -= 360.0
1355                if rlon < 0 and llon > 0:
1356                    do_180 = True
1357            xlat1d = llat + (np.arange(self.ny)*jdir*dj)
1358            xlon1d = llon + (np.arange(self.nx)*idir*di)
1359            xlons, xlats = np.meshgrid(xlon1d,xlat1d)
1360            rot2ll_vectorized = np.vectorize(arakawa_rotated_grid.rot2ll)
1361            lats, lons = rot2ll_vectorized(xlats,xlons,lasp,losp)
1362            if do_180:
1363                lons = np.where(lons>180.0,lons-360.0,lons)
1364            vector_rotation_angles_vectorized = np.vectorize(arakawa_rotated_grid.vector_rotation_angles)
1365            rots = vector_rotation_angles_vectorized(lats, lons, clat, losp, xlats)
1366            del xlat1d, xlon1d, xlats, xlons
1367        else:
1368            raise ValueError('Unsupported grid')
1369
1370        _latlon_datastore[self._sha1_section3] = dict(latitude=lats,longitude=lons)
1371        try:
1372            _latlon_datastore[self._sha1_section3]['vector_rotation_angles'] = rots
1373        except(NameError):
1374            pass
1375
1376        return lats, lons
1377
1378
1379    def map_keys(self):
1380        """
1381        Unpack data grid replacing integer values with strings.
1382
1383        These types of fields are categorical or classifications where data
1384        values do not represent an observable or predictable physical quantity.
1385        An example of such a field would be [Dominant Precipitation Type -
1386        DPTYPE](https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-201.shtml)
1387
1388        Returns
1389        -------
1390        map_keys
1391            numpy.ndarray of string values per element.
1392        """
1393        hold_auto_nans = _AUTO_NANS
1394        set_auto_nans(False)
1395
1396        if (np.all(self.section1[0:2] == [7, 14]) and self.shortName == "PWTHER") or (
1397            self._isNDFD and self.shortName in {"WX", "WWA"}
1398        ):
1399            keys = utils.decode_wx_strings(self.section2)
1400            if hasattr(self,'priMissingValue') and self.priMissingValue not in [None,0]:
1401                keys[int(self.priMissingValue)] = 'Missing'
1402            if hasattr(self,'secMissingValue') and self.secMissingValue not in [None,0]:
1403                keys[int(self.secMissingValue)] = 'Missing'
1404            u,inv = np.unique(self.data,return_inverse=True)
1405            fld = np.array([keys[x] for x in u])[inv].reshape(self.data.shape)
1406        else:
1407            # For data whose units are defined in a code table (i.e. classification or mask)
1408            tblname = re.findall(r'\d\.\d+',self.units,re.IGNORECASE)[0]
1409            fld = self.data.astype(np.int32).astype(str)
1410            tbl = tables.get_table(tblname,expand=True)
1411            for val in np.unique(fld):
1412                fld = np.where(fld==val,tbl[val],fld)
1413        set_auto_nans(hold_auto_nans)
1414        return fld
1415
1416
1417    def to_bytes(self, validate: bool=True):
1418        """
1419        Return packed GRIB2 message in bytes format.
1420
1421        This will be useful for exporting data in non-file formats. For example,
1422        can be used to output grib data directly to S3 using the boto3 client
1423        without the need to write a temporary file to upload first.
1424
1425        Parameters
1426        ----------
1427        validate: default=True
1428            If `True` (DEFAULT), validates first/last four bytes for proper
1429            formatting, else returns None. If `False`, message is output as is.
1430
1431        Returns
1432        -------
1433        to_bytes
1434            Returns GRIB2 formatted message as bytes.
1435        """
1436        if hasattr(self,'_msg'):
1437            if validate:
1438                if self.validate():
1439                    return self._msg
1440                else:
1441                    return None
1442            else:
1443                return self._msg
1444        else:
1445            return None
1446
1447
1448    def interpolate(self, method, grid_def_out, method_options=None, drtn=None,
1449                    num_threads=1):
1450        """
1451        Grib2Message Interpolator
1452
1453        Performs spatial interpolation via the [NCEPLIBS-ip
1454        library](https://github.com/NOAA-EMC/NCEPLIBS-ip). This interpolate
1455        method only supports scalar interpolation. If you need to perform
1456        vector interpolation, use the module-level `grib2io.interpolate`
1457        function.
1458
1459        Parameters
1460        ----------
1461        method
1462            Interpolate method to use. This can either be an integer or string
1463            using the following mapping:
1464
1465            | Interpolate Scheme | Integer Value |
1466            | :---:              | :---:         |
1467            | 'bilinear'         | 0             |
1468            | 'bicubic'          | 1             |
1469            | 'neighbor'         | 2             |
1470            | 'budget'           | 3             |
1471            | 'spectral'         | 4             |
1472            | 'neighbor-budget'  | 6             |
1473
1474        grid_def_out : grib2io.Grib2GridDef
1475            Grib2GridDef object of the output grid.
1476        method_options : list of ints, optional
1477            Interpolation options. See the NCEPLIBS-ip documentation for
1478            more information on how these are used.
1479        drtn
1480            Data Representation Template to be used for the returned
1481            interpolated GRIB2 message. When `None`, the data representation
1482            template of the source GRIB2 message is used. Once again, it is the
1483            user's responsibility to properly set the Data Representation
1484            Template attributes.
1485        num_threads : int, optional
1486            Number of OpenMP threads to use for interpolation. The default
1487            value is 1. If NCEPLIBS-ip and grib2io's iplib extension module
1488            was not built with OpenMP, then this keyword argument and value
1489            will have no impact.
1490
1491        Returns
1492        -------
1493        interpolate
1494            If interpolating to a grid, a new Grib2Message object is returned.
1495            The GRIB2 metadata of the new Grib2Message object is identical to
1496            the input except where required to be different because of the new
1497            grid specs and possibly a new data representation template.
1498
1499            If interpolating to station points, the interpolated data values are
1500            returned as a numpy.ndarray.
1501        """
1502        section0 = self.section0
1503        section0[-1] = 0
1504        gds = [0, grid_def_out.npoints, 0, 255, grid_def_out.gdtn]
1505        section3 = np.concatenate((gds,grid_def_out.gdt))
1506        drtn = self.drtn if drtn is None else drtn
1507
1508        msg = Grib2Message(section0,self.section1,self.section2,section3,
1509                           self.section4,None,self.bitMapFlag.value,drtn=drtn)
1510
1511        msg._msgnum = -1
1512        msg._deflist = self._deflist
1513        msg._coordlist = self._coordlist
1514        shape = (msg.ny,msg.nx)
1515        ndim = 2
1516        if msg.typeOfValues == 0:
1517            dtype = 'float32'
1518        elif msg.typeOfValues == 1:
1519            dtype = 'int32'
1520        msg._data = interpolate(self.data,method,Grib2GridDef.from_section3(self.section3),grid_def_out,
1521                                method_options=method_options,num_threads=num_threads).reshape(msg.ny,msg.nx)
1522        msg.section5[0] = grid_def_out.npoints
1523        return msg
1524
1525
1526    def subset(self, lats, lons):
1527        """
1528        Return a spatial subset.
1529
1530        Currently only supports regular grids of the following types:
1531
1532        | Grid Type                                                    | gdtn  |
1533        | :---:                                                        | :---: |
1534        | Latitude/Longitude, Equidistant Cylindrical, or Plate Carree | 0     |
1535        | Rotated Latitude/Longitude                                   | 1     |
1536        | Mercator                                                     | 10    |
1537        | Polar Stereographic                                          | 20    |
1538        | Lambert Conformal                                            | 30    |
1539        | Albers Equal-Area                                            | 31    |
1540        | Gaussian Latitude/Longitude                                  | 40    |
1541        | Equatorial Azimuthal Equidistant Projection                  | 110   |
1542
1543        Parameters
1544        ----------
1545        lats
1546            List or tuple of latitudes.  The minimum and maximum latitudes will
1547            be used to define the southern and northern boundaries.
1548
1549            The order of the latitudes is not important.  The function will
1550            determine which is the minimum and maximum.
1551
1552            The latitudes should be in decimal degrees with 0.0 at the equator,
1553            positive values in the northern hemisphere increasing to 90, and
1554            negative values in the southern hemisphere decreasing to -90.
1555        lons
1556            List or tuple of longitudes.  The minimum and maximum longitudes
1557            will be used to define the western and eastern boundaries.
1558
1559            The order of the longitudes is not important.  The function will
1560            determine which is the minimum and maximum.
1561
1562            GRIB2 longitudes should be in decimal degrees with 0.0 at the prime
1563            meridian, positive values increasing eastward to 360.  There are no
1564            negative GRIB2 longitudes.
1565
1566            The typical west longitudes that start at 0.0 at the prime meridian
1567            and decrease to -180 westward, are converted to GRIB2 longitudes by
1568            '360 - (absolute value of the west longitude)' where typical
1569            eastern longitudes are unchanged as GRIB2 longitudes.
1570
1571        Returns
1572        -------
1573        subset
1574            A spatial subset of a GRIB2 message.
1575        """
1576        if self.gdtn not in [0, 1, 10, 20, 30, 31, 40, 110]:
1577            raise ValueError(
1578                """
1579
1580Subset only works for
1581    Latitude/Longitude, Equidistant Cylindrical, or Plate Carree (gdtn=0)
1582    Rotated Latitude/Longitude (gdtn=1)
1583    Mercator (gdtn=10)
1584    Polar Stereographic (gdtn=20)
1585    Lambert Conformal (gdtn=30)
1586    Albers Equal-Area (gdtn=31)
1587    Gaussian Latitude/Longitude (gdtn=40)
1588    Equatorial Azimuthal Equidistant Projection (gdtn=110)
1589
1590"""
1591            )
1592
1593        if self.nx == 0 or self.ny == 0:
1594            raise ValueError(
1595                """
1596
1597Subset only works for regular grids.
1598
1599"""
1600            )
1601
1602        newmsg = Grib2Message(
1603            np.copy(self.section0),
1604            np.copy(self.section1),
1605            np.copy(self.section2),
1606            np.copy(self.section3),
1607            np.copy(self.section4),
1608            np.copy(self.section5),
1609        )
1610
1611        msglats, msglons = self.grid()
1612
1613        la1 = np.max(lats)
1614        lo1 = np.min(lons)
1615        la2 = np.min(lats)
1616        lo2 = np.max(lons)
1617
1618        # Find the indices of the first and last grid points to the nearest
1619        # lat/lon values in the grid.
1620        first_lat = np.abs(msglats - la1)
1621        first_lon = np.abs(msglons - lo1)
1622        max_idx = np.maximum(first_lat, first_lon)
1623        first_j, first_i = np.where(max_idx == np.min(max_idx))
1624
1625        last_lat = np.abs(msglats - la2)
1626        last_lon = np.abs(msglons - lo2)
1627        max_idx = np.maximum(last_lat, last_lon)
1628        last_j, last_i = np.where(max_idx == np.min(max_idx))
1629
1630        setattr(newmsg, "latitudeFirstGridpoint", msglats[first_j[0], first_i[0]])
1631        setattr(newmsg, "longitudeFirstGridpoint", msglons[first_j[0], first_i[0]])
1632        setattr(newmsg, "nx", np.abs(first_i[0] - last_i[0]))
1633        setattr(newmsg, "ny", np.abs(first_j[0] - last_j[0]))
1634
1635        # Set *LastGridpoint attributes even if only used for gdtn=[0, 1, 40].
1636        # This information is used to subset xarray datasets and even though
1637        # unnecessary for some supported grid types, it won't affect a grib2io
1638        # message to set them.
1639        setattr(newmsg, "latitudeLastGridpoint", msglats[last_j[0], last_i[0]])
1640        setattr(newmsg, "longitudeLastGridpoint", msglons[last_j[0], last_i[0]])
1641
1642        setattr(
1643            newmsg,
1644            "data",
1645            self.data[
1646                min(first_j[0], last_j[0]) : max(first_j[0], last_j[0]),
1647                min(first_i[0], last_i[0]) : max(first_i[0], last_i[0]),
1648            ].copy(),
1649        )
1650
1651        # Need to set the newmsg._sha1_section3 to a blank string so the grid
1652        # method ignores the cached lat/lon values.  This will force the grid
1653        # method to recompute the lat/lon values for the subsetted grid.
1654        newmsg._sha1_section3 = ""
1655        newmsg.grid()
1656
1657        return newmsg
1658
1659    def validate(self):
1660        """
1661        Validate a complete GRIB2 message.
1662
1663        The g2c library does its own internal validation when g2_gribend() is called, but
1664        we will check in grib2io also. The validation checks if the first 4 bytes in
1665        self._msg is 'GRIB' and '7777' as the last 4 bytes and that the message length in
1666        section 0 equals the length of the packed message.
1667
1668        Returns
1669        -------
1670        `True` if the packed GRIB2 message is complete and well-formed, `False` otherwise.
1671        """
1672        valid = False
1673        if hasattr(self,'_msg'):
1674            if self._msg[0:4]+self._msg[-4:] == b'GRIB7777':
1675                if self.section0[-1] == len(self._msg):
1676                    valid = True
1677        return valid

GRIB2 Message base class.

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

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[tuple[typing.Any, ...], numpy.dtype[~_ScalarT]]

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[tuple[typing.Any, ...], numpy.dtype[~_ScalarT]]

Product Definition Template

numberOfPackedValues: int

Number of Packed Values

dataRepresentationTemplate: list

Data Representation Template

gdtn
891    @property
892    def gdtn(self):
893        """Return Grid Definition Template Number"""
894        return self.section3[4]

Return Grid Definition Template Number

gdt
897    @property
898    def gdt(self):
899        """Return Grid Definition Template."""
900        return self.gridDefinitionTemplate

Return Grid Definition Template.

pdtn
903    @property
904    def pdtn(self):
905        """Return Product Definition Template Number."""
906        return self.section4[1]

Return Product Definition Template Number.

pdt
909    @property
910    def pdt(self):
911        """Return Product Definition Template."""
912        return self.productDefinitionTemplate

Return Product Definition Template.

drtn
915    @property
916    def drtn(self):
917        """Return Data Representation Template Number."""
918        return self.section5[1]

Return Data Representation Template Number.

drt
921    @property
922    def drt(self):
923        """Return Data Representation Template."""
924        return self.dataRepresentationTemplate

Return Data Representation Template.

pdy
927    @property
928    def pdy(self):
929        """Return the PDY ('YYYYMMDD')."""
930        return ''.join([str(i) for i in self.section1[5:8]])

Return the PDY ('YYYYMMDD').

griddef
933    @property
934    def griddef(self):
935        """Return a Grib2GridDef instance for a GRIB2 message."""
936        return Grib2GridDef.from_section3(self.section3)

Return a Grib2GridDef instance for a GRIB2 message.

lats
939    @property
940    def lats(self):
941        """Return grid latitudes."""
942        return self.latlons()[0]

Return grid latitudes.

lons
945    @property
946    def lons(self):
947        """Return grid longitudes."""
948        return self.latlons()[1]

Return grid longitudes.

min
950    @property
951    def min(self):
952        """Return minimum value of data."""
953        return np.nanmin(self.data)

Return minimum value of data.

max
956    @property
957    def max(self):
958        """Return maximum value of data."""
959        return np.nanmax(self.data)

Return maximum value of data.

mean
962    @property
963    def mean(self):
964        """Return mean value of data."""
965        return np.nanmean(self.data)

Return mean value of data.

median
968    @property
969    def median(self):
970        """Return median value of data."""
971        return np.nanmedian(self.data)

Return median value of data.

shape
974    @property
975    def shape(self):
976        """Return shape of data."""
977        return self.griddef.shape

Return shape of data.

def attrs_by_section(self, sect: int, values: bool = False):
1018    def attrs_by_section(self, sect: int, values: bool=False):
1019        """
1020        Provide a tuple of attribute names for the given GRIB2 section.
1021
1022        Parameters
1023        ----------
1024        sect
1025            The GRIB2 section number.
1026        values
1027            Optional (default is `False`) argument to return attributes values.
1028
1029        Returns
1030        -------
1031        attrs_by_section
1032            A list of attribute names or dict of name:value pairs if `values =
1033            True`.
1034        """
1035        if sect in {0,1,6}:
1036            attrs = templates._section_attrs[sect]
1037        elif sect in {3,4,5}:
1038            def _find_class_index(n):
1039                _key = {3:'Grid', 4:'Product', 5:'Data'}
1040                for i,c in enumerate(self.__class__.__mro__):
1041                    if _key[n] in c.__name__:
1042                        return i
1043                else:
1044                    return []
1045            if sys.version_info.minor <= 8:
1046                attrs = templates._section_attrs[sect]+\
1047                        [a for a in dir(self.__class__.__mro__[_find_class_index(sect)]) if not a.startswith('_')]
1048            else:
1049                attrs = templates._section_attrs[sect]+\
1050                        self.__class__.__mro__[_find_class_index(sect)]._attrs()
1051        else:
1052            attrs = []
1053        if values:
1054            return {k:getattr(self,k) for k in attrs}
1055        else:
1056            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):
1059    def pack(self):
1060        """
1061        Pack GRIB2 section data into a binary message.
1062
1063        It is the user's responsibility to populate the GRIB2 section
1064        information with appropriate metadata.
1065        """
1066        # Create beginning of packed binary message with section 0 and 1 data.
1067        self._sections = []
1068        self._msg,self._pos = g2clib.grib2_create(self.indicatorSection[2:4],self.identificationSection)
1069        self._sections += [0,1]
1070
1071        # Add section 2 if present.
1072        if isinstance(self.section2,bytes) and len(self.section2) > 0:
1073            self._msg,self._pos = g2clib.grib2_addlocal(self._msg,self.section2)
1074            self._sections.append(2)
1075
1076        # Add section 3.
1077        self.section3[1] = self.nx * self.ny
1078        self._msg,self._pos = g2clib.grib2_addgrid(self._msg,self.gridDefinitionSection,
1079                                                   self.gridDefinitionTemplate,
1080                                                   self._deflist)
1081        self._sections.append(3)
1082
1083        # Prepare data.
1084        if self._data is None:
1085            if self._ondiskarray is None:
1086                raise ValueError("Grib2Message object has no data, thus it cannot be packed.")
1087        field = np.copy(self.data)
1088        if self.scanModeFlags is not None:
1089            if self.scanModeFlags[3]:
1090                fieldsave = field.astype('f') # Casting makes a copy
1091                field[1::2,:] = fieldsave[1::2,::-1]
1092        fld = field.astype('f')
1093
1094        # Prepare bitmap, if necessary
1095        bitmapflag = self.bitMapFlag.value
1096        if bitmapflag == 0:
1097            if self.bitmap is not None:
1098                bmap = np.ravel(self.bitmap).astype(DEFAULT_NUMPY_INT)
1099            else:
1100                bmap = np.ravel(np.where(np.isnan(fld),0,1)).astype(DEFAULT_NUMPY_INT)
1101        else:
1102            bmap = None
1103
1104        # Prepare data for packing if nans are present
1105        fld = np.ravel(fld)
1106        if bitmapflag in {0,254}:
1107            fld = np.where(np.isnan(fld),0,fld)
1108        else:
1109            if np.isnan(fld).any():
1110                if hasattr(self,'priMissingValue'):
1111                    fld = np.where(np.isnan(fld),self.priMissingValue,fld)
1112            if hasattr(self,'_missvalmap'):
1113                if hasattr(self,'priMissingValue'):
1114                    fld = np.where(self._missvalmap==1,self.priMissingValue,fld)
1115                if hasattr(self,'secMissingValue'):
1116                    fld = np.where(self._missvalmap==2,self.secMissingValue,fld)
1117
1118        # Add sections 4, 5, 6, and 7.
1119        self._msg,self._pos = g2clib.grib2_addfield(self._msg,self.pdtn,
1120                                                    self.productDefinitionTemplate,
1121                                                    self._coordlist,
1122                                                    self.drtn,
1123                                                    self.dataRepresentationTemplate,
1124                                                    fld,
1125                                                    bitmapflag,
1126                                                    bmap)
1127        self._sections.append(4)
1128        self._sections.append(5)
1129        self._sections.append(6)
1130        self._sections.append(7)
1131
1132        # Finalize GRIB2 message with section 8.
1133        self._msg, self._pos = g2clib.grib2_end(self._msg)
1134        self._sections.append(8)
1135        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>
1138    @property
1139    def data(self) -> np.array:
1140        """Access the unpacked data values."""
1141        if self._data is None:
1142            if self._auto_nans != _AUTO_NANS:
1143                self._data = self._ondiskarray
1144            self._data = np.asarray(self._ondiskarray)
1145        return self._data

Access the unpacked data values.

def flush_data(self):
1186    def flush_data(self):
1187        """
1188        Flush the unpacked data values from the Grib2Message object.
1189
1190        Notes
1191        -----
1192        If the Grib2Message object was constructed from "scratch" (i.e.
1193        not read from file), this method will remove the data array from
1194        the object and it cannot be recovered.
1195        """
1196        self._data = None
1197        self.bitmap = None

Flush the unpacked data values from the Grib2Message object.

Notes

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):
1208    def latlons(self, *args, **kwrgs):
1209        """Alias for `grib2io.Grib2Message.grid` method."""
1210        return self.grid(*args, **kwrgs)

Alias for grib2io.Grib2Message.grid method.

def grid(self, unrotate: bool = True):
1213    def grid(self, unrotate: bool=True):
1214        """
1215        Return lats,lons (in degrees) of grid.
1216
1217        Currently can handle reg. lat/lon,cglobal Gaussian, mercator,
1218        stereographic, lambert conformal, albers equal-area, space-view and
1219        azimuthal equidistant grids.
1220
1221        Parameters
1222        ----------
1223        unrotate
1224            If `True` [DEFAULT], and grid is rotated lat/lon, then unrotate the
1225            grid, otherwise `False`, do not.
1226
1227        Returns
1228        -------
1229        lats, lons : numpy.ndarray
1230            Returns two numpy.ndarrays with dtype=numpy.float32 of grid
1231            latitudes and longitudes in units of degrees.
1232        """
1233        if self._sha1_section3 in _latlon_datastore.keys():
1234            return (_latlon_datastore[self._sha1_section3]['latitude'],
1235                    _latlon_datastore[self._sha1_section3]['longitude'])
1236        gdtn = self.gridDefinitionTemplateNumber.value
1237        gdtmpl = self.gridDefinitionTemplate
1238        reggrid = self.gridDefinitionSection[2] == 0 # This means regular 2-d grid
1239        if gdtn == 0:
1240            # Regular lat/lon grid
1241            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1242            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1243            dlon = self.gridlengthXDirection
1244            dlat = self.gridlengthYDirection
1245            if lon2 < lon1 and dlon < 0: lon1 = -lon1
1246            lats = np.linspace(lat1,lat2,self.ny)
1247            if reggrid:
1248                lons = np.linspace(lon1,lon2,self.nx)
1249            else:
1250                lons = np.linspace(lon1,lon2,self.ny*2)
1251            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1252        elif gdtn == 1: # Rotated Lat/Lon grid
1253            pj = pyproj.Proj(self.projParameters)
1254            lat1,lon1 = self.latitudeFirstGridpoint,self.longitudeFirstGridpoint
1255            lat2,lon2 = self.latitudeLastGridpoint,self.longitudeLastGridpoint
1256            if lon1 > 180.0: lon1 -= 360.0
1257            if lon2 > 180.0: lon2 -= 360.0
1258            lats = np.linspace(lat1,lat2,self.ny)
1259            lons = np.linspace(lon1,lon2,self.nx)
1260            lons,lats = np.meshgrid(lons,lats) # Make 2-d arrays.
1261            if unrotate:
1262                from grib2io.utils import rotated_grid
1263                lats,lons = rotated_grid.unrotate(lats,lons,self.anglePoleRotation,
1264                                                  self.latitudeSouthernPole,
1265                                                  self.longitudeSouthernPole)
1266        elif gdtn == 40: # Gaussian grid (only works for global!)
1267            from grib2io.utils.gauss_grid import gaussian_latitudes
1268            lon1, lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1269            lon2, lat2 = self.longitudeLastGridpoint, self.latitudeLastGridpoint
1270            nlats = self.ny
1271            if not reggrid: # Reduced Gaussian grid.
1272                nlons = 2*nlats
1273                dlon = 360./nlons
1274            else:
1275                nlons = self.nx
1276                dlon = self.gridlengthXDirection
1277            lons = np.linspace(lon1,lon2,nlons)
1278            # Compute Gaussian lats (north to south)
1279            lats = gaussian_latitudes(nlats)
1280            if lat1 < lat2:  # reverse them if necessary
1281                lats = lats[::-1]
1282            lons,lats = np.meshgrid(lons,lats)
1283        elif gdtn in {10,20,30,31,110}:
1284            # Mercator, Lambert Conformal, Stereographic, Albers Equal Area,
1285            # Azimuthal Equidistant
1286            dx,dy = self.gridlengthXDirection, self.gridlengthYDirection
1287            lon1,lat1 = self.longitudeFirstGridpoint, self.latitudeFirstGridpoint
1288            pj = pyproj.Proj(self.projParameters)
1289            llcrnrx, llcrnry = pj(lon1,lat1)
1290            x = llcrnrx+dx*np.arange(self.nx)
1291            y = llcrnry+dy*np.arange(self.ny)
1292            x,y = np.meshgrid(x, y)
1293            lons,lats = pj(x, y, inverse=True)
1294        elif gdtn == 90:
1295            # Satellite Projection
1296            dx = self.gridlengthXDirection
1297            dy = self.gridlengthYDirection
1298            pj = pyproj.Proj(self.projParameters)
1299            x = dx*np.indices((self.ny,self.nx),'f')[1,:,:]
1300            x -= 0.5*x.max()
1301            y = dy*np.indices((self.ny,self.nx),'f')[0,:,:]
1302            y -= 0.5*y.max()
1303            lons,lats = pj(x,y,inverse=True)
1304            # Set lons,lats to 1.e30 where undefined
1305            abslons = np.fabs(lons)
1306            abslats = np.fabs(lats)
1307            lons = np.where(abslons < 1.e20, lons, 1.e30)
1308            lats = np.where(abslats < 1.e20, lats, 1.e30)
1309        elif gdtn == 32769:
1310            # Special NCEP Grid, Rotated Lat/Lon, Arakawa E-Grid (Non-Staggered)
1311            from grib2io.utils import arakawa_rotated_grid
1312            from grib2io.utils.rotated_grid import DEG2RAD
1313            di, dj = 0.0, 0.0
1314            do_180 = False
1315            idir = 1 if self.scanModeFlags[0] == 0 else -1
1316            jdir = -1 if self.scanModeFlags[1] == 0 else 1
1317            grid_oriented = 0 if self.resolutionAndComponentFlags[4] == 0 else 1
1318            do_rot = 1 if grid_oriented == 1 else 0
1319            la1 = self.latitudeFirstGridpoint
1320            lo1 = self.longitudeFirstGridpoint
1321            clon = self.longitudeCenterGridpoint
1322            clat = self.latitudeCenterGridpoint
1323            lasp = clat - 90.0
1324            losp = clon
1325            llat, llon = arakawa_rotated_grid.ll2rot(la1,lo1,lasp,losp)
1326            la2, lo2 = arakawa_rotated_grid.rot2ll(-llat,-llon,lasp,losp)
1327            rlat = -llat
1328            rlon = -llon
1329            if self.nx == 1:
1330                di = 0.0
1331            elif idir == 1:
1332                ti = rlon
1333                while ti < llon:
1334                    ti += 360.0
1335                di = (ti - llon)/float(self.nx-1)
1336            else:
1337                ti = llon
1338                while ti < rlon:
1339                    ti += 360.0
1340                di = (ti - rlon)/float(self.nx-1)
1341            if self.ny == 1:
1342               dj = 0.0
1343            else:
1344                dj = (rlat - llat)/float(self.ny-1)
1345                if dj < 0.0:
1346                    dj = -dj
1347            if idir == 1:
1348                if llon > rlon:
1349                    llon -= 360.0
1350                if llon < 0 and rlon > 0:
1351                    do_180 = True
1352            else:
1353                if rlon > llon:
1354                    rlon -= 360.0
1355                if rlon < 0 and llon > 0:
1356                    do_180 = True
1357            xlat1d = llat + (np.arange(self.ny)*jdir*dj)
1358            xlon1d = llon + (np.arange(self.nx)*idir*di)
1359            xlons, xlats = np.meshgrid(xlon1d,xlat1d)
1360            rot2ll_vectorized = np.vectorize(arakawa_rotated_grid.rot2ll)
1361            lats, lons = rot2ll_vectorized(xlats,xlons,lasp,losp)
1362            if do_180:
1363                lons = np.where(lons>180.0,lons-360.0,lons)
1364            vector_rotation_angles_vectorized = np.vectorize(arakawa_rotated_grid.vector_rotation_angles)
1365            rots = vector_rotation_angles_vectorized(lats, lons, clat, losp, xlats)
1366            del xlat1d, xlon1d, xlats, xlons
1367        else:
1368            raise ValueError('Unsupported grid')
1369
1370        _latlon_datastore[self._sha1_section3] = dict(latitude=lats,longitude=lons)
1371        try:
1372            _latlon_datastore[self._sha1_section3]['vector_rotation_angles'] = rots
1373        except(NameError):
1374            pass
1375
1376        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):
1379    def map_keys(self):
1380        """
1381        Unpack data grid replacing integer values with strings.
1382
1383        These types of fields are categorical or classifications where data
1384        values do not represent an observable or predictable physical quantity.
1385        An example of such a field would be [Dominant Precipitation Type -
1386        DPTYPE](https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table4-201.shtml)
1387
1388        Returns
1389        -------
1390        map_keys
1391            numpy.ndarray of string values per element.
1392        """
1393        hold_auto_nans = _AUTO_NANS
1394        set_auto_nans(False)
1395
1396        if (np.all(self.section1[0:2] == [7, 14]) and self.shortName == "PWTHER") or (
1397            self._isNDFD and self.shortName in {"WX", "WWA"}
1398        ):
1399            keys = utils.decode_wx_strings(self.section2)
1400            if hasattr(self,'priMissingValue') and self.priMissingValue not in [None,0]:
1401                keys[int(self.priMissingValue)] = 'Missing'
1402            if hasattr(self,'secMissingValue') and self.secMissingValue not in [None,0]:
1403                keys[int(self.secMissingValue)] = 'Missing'
1404            u,inv = np.unique(self.data,return_inverse=True)
1405            fld = np.array([keys[x] for x in u])[inv].reshape(self.data.shape)
1406        else:
1407            # For data whose units are defined in a code table (i.e. classification or mask)
1408            tblname = re.findall(r'\d\.\d+',self.units,re.IGNORECASE)[0]
1409            fld = self.data.astype(np.int32).astype(str)
1410            tbl = tables.get_table(tblname,expand=True)
1411            for val in np.unique(fld):
1412                fld = np.where(fld==val,tbl[val],fld)
1413        set_auto_nans(hold_auto_nans)
1414        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):
1417    def to_bytes(self, validate: bool=True):
1418        """
1419        Return packed GRIB2 message in bytes format.
1420
1421        This will be useful for exporting data in non-file formats. For example,
1422        can be used to output grib data directly to S3 using the boto3 client
1423        without the need to write a temporary file to upload first.
1424
1425        Parameters
1426        ----------
1427        validate: default=True
1428            If `True` (DEFAULT), validates first/last four bytes for proper
1429            formatting, else returns None. If `False`, message is output as is.
1430
1431        Returns
1432        -------
1433        to_bytes
1434            Returns GRIB2 formatted message as bytes.
1435        """
1436        if hasattr(self,'_msg'):
1437            if validate:
1438                if self.validate():
1439                    return self._msg
1440                else:
1441                    return None
1442            else:
1443                return self._msg
1444        else:
1445            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):
1448    def interpolate(self, method, grid_def_out, method_options=None, drtn=None,
1449                    num_threads=1):
1450        """
1451        Grib2Message Interpolator
1452
1453        Performs spatial interpolation via the [NCEPLIBS-ip
1454        library](https://github.com/NOAA-EMC/NCEPLIBS-ip). This interpolate
1455        method only supports scalar interpolation. If you need to perform
1456        vector interpolation, use the module-level `grib2io.interpolate`
1457        function.
1458
1459        Parameters
1460        ----------
1461        method
1462            Interpolate method to use. This can either be an integer or string
1463            using the following mapping:
1464
1465            | Interpolate Scheme | Integer Value |
1466            | :---:              | :---:         |
1467            | 'bilinear'         | 0             |
1468            | 'bicubic'          | 1             |
1469            | 'neighbor'         | 2             |
1470            | 'budget'           | 3             |
1471            | 'spectral'         | 4             |
1472            | 'neighbor-budget'  | 6             |
1473
1474        grid_def_out : grib2io.Grib2GridDef
1475            Grib2GridDef object of the output grid.
1476        method_options : list of ints, optional
1477            Interpolation options. See the NCEPLIBS-ip documentation for
1478            more information on how these are used.
1479        drtn
1480            Data Representation Template to be used for the returned
1481            interpolated GRIB2 message. When `None`, the data representation
1482            template of the source GRIB2 message is used. Once again, it is the
1483            user's responsibility to properly set the Data Representation
1484            Template attributes.
1485        num_threads : int, optional
1486            Number of OpenMP threads to use for interpolation. The default
1487            value is 1. If NCEPLIBS-ip and grib2io's iplib extension module
1488            was not built with OpenMP, then this keyword argument and value
1489            will have no impact.
1490
1491        Returns
1492        -------
1493        interpolate
1494            If interpolating to a grid, a new Grib2Message object is returned.
1495            The GRIB2 metadata of the new Grib2Message object is identical to
1496            the input except where required to be different because of the new
1497            grid specs and possibly a new data representation template.
1498
1499            If interpolating to station points, the interpolated data values are
1500            returned as a numpy.ndarray.
1501        """
1502        section0 = self.section0
1503        section0[-1] = 0
1504        gds = [0, grid_def_out.npoints, 0, 255, grid_def_out.gdtn]
1505        section3 = np.concatenate((gds,grid_def_out.gdt))
1506        drtn = self.drtn if drtn is None else drtn
1507
1508        msg = Grib2Message(section0,self.section1,self.section2,section3,
1509                           self.section4,None,self.bitMapFlag.value,drtn=drtn)
1510
1511        msg._msgnum = -1
1512        msg._deflist = self._deflist
1513        msg._coordlist = self._coordlist
1514        shape = (msg.ny,msg.nx)
1515        ndim = 2
1516        if msg.typeOfValues == 0:
1517            dtype = 'float32'
1518        elif msg.typeOfValues == 1:
1519            dtype = 'int32'
1520        msg._data = interpolate(self.data,method,Grib2GridDef.from_section3(self.section3),grid_def_out,
1521                                method_options=method_options,num_threads=num_threads).reshape(msg.ny,msg.nx)
1522        msg.section5[0] = grid_def_out.npoints
1523        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 NCEPLIBS-ip and grib2io's iplib extension module 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 subset(self, lats, lons):
1526    def subset(self, lats, lons):
1527        """
1528        Return a spatial subset.
1529
1530        Currently only supports regular grids of the following types:
1531
1532        | Grid Type                                                    | gdtn  |
1533        | :---:                                                        | :---: |
1534        | Latitude/Longitude, Equidistant Cylindrical, or Plate Carree | 0     |
1535        | Rotated Latitude/Longitude                                   | 1     |
1536        | Mercator                                                     | 10    |
1537        | Polar Stereographic                                          | 20    |
1538        | Lambert Conformal                                            | 30    |
1539        | Albers Equal-Area                                            | 31    |
1540        | Gaussian Latitude/Longitude                                  | 40    |
1541        | Equatorial Azimuthal Equidistant Projection                  | 110   |
1542
1543        Parameters
1544        ----------
1545        lats
1546            List or tuple of latitudes.  The minimum and maximum latitudes will
1547            be used to define the southern and northern boundaries.
1548
1549            The order of the latitudes is not important.  The function will
1550            determine which is the minimum and maximum.
1551
1552            The latitudes should be in decimal degrees with 0.0 at the equator,
1553            positive values in the northern hemisphere increasing to 90, and
1554            negative values in the southern hemisphere decreasing to -90.
1555        lons
1556            List or tuple of longitudes.  The minimum and maximum longitudes
1557            will be used to define the western and eastern boundaries.
1558
1559            The order of the longitudes is not important.  The function will
1560            determine which is the minimum and maximum.
1561
1562            GRIB2 longitudes should be in decimal degrees with 0.0 at the prime
1563            meridian, positive values increasing eastward to 360.  There are no
1564            negative GRIB2 longitudes.
1565
1566            The typical west longitudes that start at 0.0 at the prime meridian
1567            and decrease to -180 westward, are converted to GRIB2 longitudes by
1568            '360 - (absolute value of the west longitude)' where typical
1569            eastern longitudes are unchanged as GRIB2 longitudes.
1570
1571        Returns
1572        -------
1573        subset
1574            A spatial subset of a GRIB2 message.
1575        """
1576        if self.gdtn not in [0, 1, 10, 20, 30, 31, 40, 110]:
1577            raise ValueError(
1578                """
1579
1580Subset only works for
1581    Latitude/Longitude, Equidistant Cylindrical, or Plate Carree (gdtn=0)
1582    Rotated Latitude/Longitude (gdtn=1)
1583    Mercator (gdtn=10)
1584    Polar Stereographic (gdtn=20)
1585    Lambert Conformal (gdtn=30)
1586    Albers Equal-Area (gdtn=31)
1587    Gaussian Latitude/Longitude (gdtn=40)
1588    Equatorial Azimuthal Equidistant Projection (gdtn=110)
1589
1590"""
1591            )
1592
1593        if self.nx == 0 or self.ny == 0:
1594            raise ValueError(
1595                """
1596
1597Subset only works for regular grids.
1598
1599"""
1600            )
1601
1602        newmsg = Grib2Message(
1603            np.copy(self.section0),
1604            np.copy(self.section1),
1605            np.copy(self.section2),
1606            np.copy(self.section3),
1607            np.copy(self.section4),
1608            np.copy(self.section5),
1609        )
1610
1611        msglats, msglons = self.grid()
1612
1613        la1 = np.max(lats)
1614        lo1 = np.min(lons)
1615        la2 = np.min(lats)
1616        lo2 = np.max(lons)
1617
1618        # Find the indices of the first and last grid points to the nearest
1619        # lat/lon values in the grid.
1620        first_lat = np.abs(msglats - la1)
1621        first_lon = np.abs(msglons - lo1)
1622        max_idx = np.maximum(first_lat, first_lon)
1623        first_j, first_i = np.where(max_idx == np.min(max_idx))
1624
1625        last_lat = np.abs(msglats - la2)
1626        last_lon = np.abs(msglons - lo2)
1627        max_idx = np.maximum(last_lat, last_lon)
1628        last_j, last_i = np.where(max_idx == np.min(max_idx))
1629
1630        setattr(newmsg, "latitudeFirstGridpoint", msglats[first_j[0], first_i[0]])
1631        setattr(newmsg, "longitudeFirstGridpoint", msglons[first_j[0], first_i[0]])
1632        setattr(newmsg, "nx", np.abs(first_i[0] - last_i[0]))
1633        setattr(newmsg, "ny", np.abs(first_j[0] - last_j[0]))
1634
1635        # Set *LastGridpoint attributes even if only used for gdtn=[0, 1, 40].
1636        # This information is used to subset xarray datasets and even though
1637        # unnecessary for some supported grid types, it won't affect a grib2io
1638        # message to set them.
1639        setattr(newmsg, "latitudeLastGridpoint", msglats[last_j[0], last_i[0]])
1640        setattr(newmsg, "longitudeLastGridpoint", msglons[last_j[0], last_i[0]])
1641
1642        setattr(
1643            newmsg,
1644            "data",
1645            self.data[
1646                min(first_j[0], last_j[0]) : max(first_j[0], last_j[0]),
1647                min(first_i[0], last_i[0]) : max(first_i[0], last_i[0]),
1648            ].copy(),
1649        )
1650
1651        # Need to set the newmsg._sha1_section3 to a blank string so the grid
1652        # method ignores the cached lat/lon values.  This will force the grid
1653        # method to recompute the lat/lon values for the subsetted grid.
1654        newmsg._sha1_section3 = ""
1655        newmsg.grid()
1656
1657        return newmsg

Return a spatial subset.

Currently only supports regular grids of the following types:

Grid Type gdtn
Latitude/Longitude, Equidistant Cylindrical, or Plate Carree 0
Rotated Latitude/Longitude 1
Mercator 10
Polar Stereographic 20
Lambert Conformal 30
Albers Equal-Area 31
Gaussian Latitude/Longitude 40
Equatorial Azimuthal Equidistant Projection 110
Parameters
  • lats: List or tuple of latitudes. The minimum and maximum latitudes will be used to define the southern and northern boundaries.

The order of the latitudes is not important. The function will determine which is the minimum and maximum.

The latitudes should be in decimal degrees with 0.0 at the equator, positive values in the northern hemisphere increasing to 90, and negative values in the southern hemisphere decreasing to -90.

  • lons: List or tuple of longitudes. The minimum and maximum longitudes will be used to define the western and eastern boundaries.

The order of the longitudes is not important. The function will determine which is the minimum and maximum.

GRIB2 longitudes should be in decimal degrees with 0.0 at the prime meridian, positive values increasing eastward to 360. There are no negative GRIB2 longitudes.

The typical west longitudes that start at 0.0 at the prime meridian and decrease to -180 westward, are converted to GRIB2 longitudes by '360 - (absolute value of the west longitude)' where typical eastern longitudes are unchanged as GRIB2 longitudes.

Returns
  • subset: A spatial subset of a GRIB2 message.
def validate(self):
1659    def validate(self):
1660        """
1661        Validate a complete GRIB2 message.
1662
1663        The g2c library does its own internal validation when g2_gribend() is called, but
1664        we will check in grib2io also. The validation checks if the first 4 bytes in
1665        self._msg is 'GRIB' and '7777' as the last 4 bytes and that the message length in
1666        section 0 equals the length of the packed message.
1667
1668        Returns
1669        -------
1670        `True` if the packed GRIB2 message is complete and well-formed, `False` otherwise.
1671        """
1672        valid = False
1673        if hasattr(self,'_msg'):
1674            if self._msg[0:4]+self._msg[-4:] == b'GRIB7777':
1675                if self.section0[-1] == len(self._msg):
1676                    valid = True
1677        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:
2152@dataclass
2153class Grib2GridDef:
2154    """
2155    Class for Grid Definition Template Number and Template as attributes.
2156
2157    This allows for cleaner looking code when passing these metadata around.
2158    For example, the `grib2io._Grib2Message.interpolate` method and
2159    `grib2io.interpolate` function accepts these objects.
2160    """
2161    gdtn: int
2162    gdt: NDArray
2163
2164    @classmethod
2165    def from_section3(cls, section3):
2166        return cls(section3[4],section3[5:])
2167
2168    @property
2169    def nx(self):
2170        """Number of grid points in x-direction."""
2171        return int(self.gdt[7])
2172
2173    @property
2174    def ny(self):
2175        """Number of grid points in y-direction."""
2176        return int(self.gdt[8])
2177
2178    @property
2179    def npoints(self):
2180        """Total number of grid points."""
2181        return int(self.gdt[7] * self.gdt[8])
2182
2183    @property
2184    def shape(self):
2185        """Shape of the grid."""
2186        return (int(self.ny), int(self.nx))
2187
2188    def to_section3(self):
2189        """Return a full GRIB2 section3 array."""
2190        return np.array(
2191            [0, self.npoints, 0, 0, self.gdtn] + list(self.gdt)
2192        ).astype(np.int64)

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[tuple[typing.Any, ...], numpy.dtype[~_ScalarT]])
gdtn: int
gdt: numpy.ndarray[tuple[typing.Any, ...], numpy.dtype[~_ScalarT]]
@classmethod
def from_section3(cls, section3):
2164    @classmethod
2165    def from_section3(cls, section3):
2166        return cls(section3[4],section3[5:])
nx
2168    @property
2169    def nx(self):
2170        """Number of grid points in x-direction."""
2171        return int(self.gdt[7])

Number of grid points in x-direction.

ny
2173    @property
2174    def ny(self):
2175        """Number of grid points in y-direction."""
2176        return int(self.gdt[8])

Number of grid points in y-direction.

npoints
2178    @property
2179    def npoints(self):
2180        """Total number of grid points."""
2181        return int(self.gdt[7] * self.gdt[8])

Total number of grid points.

shape
2183    @property
2184    def shape(self):
2185        """Shape of the grid."""
2186        return (int(self.ny), int(self.nx))

Shape of the grid.

def to_section3(self):
2188    def to_section3(self):
2189        """Return a full GRIB2 section3 array."""
2190        return np.array(
2191            [0, self.npoints, 0, 0, self.gdtn] + list(self.gdt)
2192        ).astype(np.int64)

Return a full GRIB2 section3 array.

def msgs_from_index(index, filehandle=None):
653def msgs_from_index(index, filehandle=None):
654    """
655    return list of Grib2Messages from index
656    msgs can get the data so long as an open filehandle is provided
657    """
658
659    zipped = zip(
660        index["section0"],
661        index["section1"],
662        index["section2"],
663        index["section3"],
664        index["section4"],
665        index["section5"],
666        index["bmapflag"]
667    )
668    msgs = [Grib2Message(*sections) for sections in zipped]
669
670    if filehandle is not None:
671        for n, (msg, offset, secpos) in enumerate(zip(msgs, index["offset"], index["sectionOffset"])):
672            msg._ondiskarray = Grib2MessageOnDiskArray(shape=(msg.ny,msg.nx),
673                                                       ndim=2,
674                                                       dtype=TYPE_OF_VALUES_DTYPE[msg.typeOfValues],
675                                                       filehandle=filehandle,
676                                                       msg=msg,
677                                                       offset=offset,
678                                                       bitmap_offset=secpos[6],
679                                                       data_offset=secpos[7])
680            msg._msgnum = n
681    return msgs

return list of Grib2Messages from index msgs can get the data so long as an open filehandle is provided