NOAA-MDL / grib2io

Python interface to the NCEP G2C Library for reading and writing GRIB2 messages.
https://noaa-mdl.github.io/grib2io
MIT License
33 stars 12 forks source link

unexpected "GRIB message already complete" #100

Closed rafatoshio closed 1 year ago

rafatoshio commented 1 year ago

Dear, I'm testing the grib2io module and I've encountered a curious problem. In the example below, I can run it for several files, but there is one that yields an error. If I alter one grid point value (see the commented line), the program runs without any problems.

Here you can find my files.

import grib2io
import numpy as np

with open("newdata", "rb") as f:
    newdata = np.load(f)

grb = grib2io.open("apcp_template.grib2")
olddata = grb.select(shortName="APCP")[0]

newmsg = grib2io.Grib2Message(
    gdtn=olddata.gdtn, pdtn=olddata.pdtn, drtn=olddata.drtn
)

newmsg.discipline = 0
newmsg.section1 = olddata.section1
newmsg.section3 = olddata.section3
newmsg.section4 = olddata.section4
newmsg.section5 = olddata.section5

# newdata[-1, -1] = 1   # <-- REPLACING THE LAST ELEMENT WORKS
newmsg.data = newdata
newmsg.pack()

the error message:

GRIB message already complete.  Cannot add new section.
Traceback (most recent call last):
  File "test_grib2io.py", line 23, in <module>
    newmsg.pack()
  File ".../lib/python3.11/site-packages/grib2io/_grib2io.py", line 856, in pack
    self._msg, self._pos = g2clib.grib2_end(self._msg)
                           ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "g2clib.pyx", line 545, in grib2io.g2clib.grib2_end
RuntimeError: error in grib2_end, error code = -2

I'm using manjaro. Python 3.11.5 grib2io version: 2.0.0

EricEngle-NOAA commented 1 year ago

Hello. I spent some time working through this issue and came to an interesting result. This is a pretty rare occurrence. You are using simple packing scheme with no decimal or binary scaling which means you are packing the integer part of the floating point values...no decimal places are used. And because it is simple packing, the grid of integer values are packed into the GRIB2 message as is.

It just so happens that the last 4 values of the grid round to the value of 55. The value of 55 when encoded into ASCII resolves to the character "7". The end section of GRIB2, section 8, is just a string of 4, sevens, "7777". So the g2c library thinks the message has already been completed.

I will reach out to the team that manages the library and probably submit an issue to them. Here is your original code with some lines added to show the last 4 values packed.

import grib2io
import numpy as np

with open("newdata", "rb") as f:
    newdata = np.load(f)

grb = grib2io.open("apcp_template.grib2")
olddata = grb.select(shortName="APCP")[0]

newmsg = grib2io.Grib2Message(
    gdtn=olddata.gdtn, pdtn=olddata.pdtn, drtn=olddata.drtn
)

newmsg.discipline = 0
newmsg.section1 = olddata.section1
newmsg.section3 = olddata.section3
newmsg.section4 = olddata.section4
newmsg.section5 = olddata.section5

#newdata[-1, -1] = 1   # <-- REPLACING THE LAST ELEMENT WORKS
last4 = newdata[-1:,-4:].reshape((4))
print(f'{last4=}')
last4int = np.round(last4).astype(np.int8)
print(f'{last4int=}')
last4int_asbytes = (np.round(last4int).tobytes())
print(f'{last4int_asbytes=}')
newmsg.data = newdata
newmsg.decScaleFactor = 0
newmsg.binScaleFactor = 0
newmsg.nBitsPacking = 0
newmsg.pack()

Output:

last4=array([55.08  , 55.325 , 54.9575, 54.59  ])
last4int=array([55, 55, 55, 55], dtype=int8)
last4int_asbytes=b'7777'
GRIB message already complete.  Cannot add new section.
Traceback (most recent call last):
  File "/Users/ericengle/Downloads/grib2io_issue/./test_grib2io.py", line 33, in <module>
    newmsg.pack()
  File "/Users/ericengle/.pyenv/versions/3.10.13/lib/python3.10/site-packages/grib2io/_grib2io.py", line 857, in pack
    self._msg, self._pos = g2clib.grib2_end(self._msg)
  File "g2clib.pyx", line 545, in grib2io.g2clib.grib2_end
RuntimeError: error in grib2_end, error code = -2

Solution

To "fix" this, set the decimal scale factor to something greater than 0 before packing.

newmsg.decScaleFactor = 1
rafatoshio commented 1 year ago

Thank you very much; I suppose I was very lucky, though.