GRIB File – Getting Time Dimension Using Python or Java

grib-2javanoaapython

I have a GRIB file from NOAA that contains several variables for a specific subregion in Norht-East Atlantic.

I'm using ECMWF grip-api Python library to parse the file. While I'm able to get the values for each {latitude, longitude} pair for each GRIB message with this code:

def iterate_latitude_longitude():
f = open(INPUT)

while 1:
    gid = grib_new_from_file(f)
    if gid is None: break

    iterid = grib_iterator_new(gid,0)

    #Get the missingValue for this GRIB message (e.g. NaN)
    missingValue = grib_get_double(gid,"missingValue")

    i=0
    while 1:
        #Get the value from the geo iterator obtained previously, value is a tuple with [latitude, longitude, value]
        result = grib_iterator_next(iterid)
        if not result: break
        [lat,lon,value] = result

        sys.stdout.write("-Point # %d - lat=%.6f lon=%.6f value=" % (i,lat,lon))
        #Check if missing value is present or not
        if value == missingValue:
            print "missing"
        else:
            print "%.6f" % value
        #increase iterator
        i += 1

    grib_iterator_delete(iterid)
    grib_release(gid)

f.close()

While this is great, this file contains forecast data, it was extracted from this place, more precisely, this file:

multi_1.nww3.t00z.grib2

Which contains world data. There are other files for different times: 00-60-12-18-00, thus, every 6 hours.

What I'm missing is the relation between the value obtained for each latitude,longitude and the time dimension. Does anybody know how to extract this info from a grib file, whether using Python or Java?

I know I can use things like wgrib2 or pygrib to convert to CSV and do the other way around…but it just looks so inefficient

Best Answer

You probably didn't want to just iterate over the whole file, but rather iterate over each of the (interesting) values, for each of the time steps you want.

Here is some code adapted from the sample code

import traceback
import sys

from gribapi import *

INPUT='multi_1.nww3.t00z.grib2'
VERBOSE=1 # verbose error reporting

def example():
    f = open(INPUT)

    keys = [
        'stepRange',
        'shortName',
        ]

    while 1:
        gid = grib_new_from_file(f)
        if gid is None: break

        for key in keys:
            if not grib_is_defined(gid, key): raise Exception("Key was not defined")
            print '%s=%s' % (key, grib_get(gid, key))

        missingValue = grib_get_double(gid,"missingValue")

        iterid = grib_iterator_new(gid,0)
        i=0
        while 1:
            result = grib_iterator_next(iterid)
            if not result: break

            [lat,lon,value] = result

            sys.stdout.write("- %d - lat=%.6f lon=%.6f value=" % (i,lat,lon))

            if value == missingValue:
                print "missing"
            else:
                print "%.6f" % value

            i += 1

        grib_iterator_delete(iterid) 
        grib_release(gid)

    f.close()

def main():
    try:
        example()
    except GribInternalError,err:
        if VERBOSE:
            traceback.print_exc(file=sys.stderr)
        else:
            print >>sys.stderr,err.msg

        return 1

if __name__ == "__main__":
    sys.exit(main())

And some sample output:

stepRange=0
shortName=unknown
....
- 288 - lat=77.000000 lon=0.000000 value=4.890000
- 289 - lat=77.000000 lon=1.250000 value=5.720000
- 290 - lat=77.000000 lon=2.500000 value=6.690000
- 291 - lat=77.000000 lon=3.750000 value=7.620000
- 292 - lat=77.000000 lon=5.000000 value=8.350000
- 293 - lat=77.000000 lon=6.250000 value=8.810000
- 294 - lat=77.000000 lon=7.500000 value=9.090000
- 295 - lat=77.000000 lon=8.750000 value=9.250000
....
stepRange=3
shortName=ws
....
- 288 - lat=77.000000 lon=0.000000 value=34.030000
- 289 - lat=77.000000 lon=1.250000 value=20.160000
- 290 - lat=77.000000 lon=2.500000 value=4.290000
- 291 - lat=77.000000 lon=3.750000 value=351.300000
- 292 - lat=77.000000 lon=5.000000 value=342.490000
- 293 - lat=77.000000 lon=6.250000 value=337.310000
- 294 - lat=77.000000 lon=7.500000 value=334.990000
- 295 - lat=77.000000 lon=8.750000 value=334.460000
- 296 - lat=77.000000 lon=10.000000 value=334.800000
- 297 - lat=77.000000 lon=11.250001 value=335.990000
- 298 - lat=77.000000 lon=12.500001 value=337.970000
- 299 - lat=77.000000 lon=13.750001 value=338.940000

Note that per your comments on the other answer, the results are in 3 hour blocks (not 6 hour as noted in the question).