[gdal-dev] grib driver

Thu Jan 17 12:12:32 EST 2008

Hello,

I am happy that you were able to build GDAL with the GRIB driver and try 
it. I also thank you for the feedback, which I will use to improve the 
driver.

You are right, the choice for "metric" is accidentally fixed in the 
driver (in case you are interested, this is in file gribdataset.cpp, 
function GRIBRasterBand::ReadGribData, you will find a line with "sChar 
f_unit = 2; /* None = 0, English = 1, Metric = 2 */").
Consider that it was quite a headache to integrate degrib with GDAL, and 
I was so happy to get it to work after a huge effort, that I forgot that 
this was left hardcoded.
I guess this could be offered as an option to the user (e.g. with a 
command-line parameter for gdal_translate and gdalwarp).

To check what happens with the image's description and the 
georeferencing, I downloaded the image gfs.2008011612/gfs.t12.pgrb2f00.

Description:
Band 2 of this image contains the following metadata that could be used 
in a description:
- msgNum = 2
- subgNum = 0
- refTime = 1200484800.0000000
- validTime = 1200484800.0000000
- element = "TMP"
- comment = "Temperature [K]"
- unitName = "[K]"
- foreSec = 0.00000000000000000
- shortFstLevel = "1000-ISBL"
- longFstLevel = "1000[Pa] ISBL="Isobaric surface""
When I was implementing the GDAL-GRIB driver (a year ago), I had checked 
those fields in all GRIB images available to me at that time, and 
concluded that longFstLevel was sufficient. With your images it makes 
sense to compose a combination of these fields (note that GDAL only has 
place for one description string).
This also goes for bands 8 and 9: For both bands, msgNum = 8, but 
subgNum = 0,1. This can be incorporated in the description, but for GDAL 
they will remain separate bands.

Georeferencing:
This image contains (among others) the following info in the metadata:
- Scan direction = GRIB2BIT_2 (meaning from south to north)
- (lat1,lon1) = (90,0) .. considering the scan direction, this is the 
lower-left point
- (lat2,lon2) = (-90,359.5) .. this is the upper-right point of the image
Well, I was not happy with this, as it differs from all images I tried 
til now. At first sight I would conclude that the image is upside-down. 
Only by opening a band I can see that it is straight-up, and that the 
lower-left corner must be assigned (lat,lon)=(-90,0). The image contains 
the whole world, with the greenwich meridian on the left.
To solve this, I must add some intelligence to the code that passes the 
georeferencing parameters to GDAL (file gribdataset.cpp, function 
GRIBDataset::SetMetaData). However, I wish I had some documentation of 
the GRIB format, as I am not sure if I can simply "sort" the corner 
coordinates of the image (why is the scan direction given?).
For now, you can force gdal_translate to assign the correct georeference 
to the image during import, by adding a parameter: -a_ullr 0 90 360 -90 
. Unless your GIS expects longitudes to be between -180 and 180, this 
will solve the problem.

When I find some time, I will implement above improvements (user's 
choice of unit, better description, and more robust georeferencing). But 
if you have some experience in C++, feel free to assist.

Cheers,

Bas.

-- 
Ir. V. (Bas) Retsios
Software Developer
Geo-information Processing Department
International Institute for Geo-information Science and Earth Observation (ITC)
P.O. Box 6,  7500 AA Enschede, The Netherlands
Phone +31 (0)53 4874 573, telefax +31 (0)53 4874 335
E-mail retsios at itc.nl, Internet http://www.itc.nl

winkey wrote:

>After trying a second time I finally got the grib driver to work for me.
>
>the grib file i am testing with is the gfs forecast model .5 degree
>global grid available from
>http://www.ftp.ncep.noaa.gov/data/nccf/com/gfs/prod/gfs.2008011412/gfs.t12z.pgrb2f00
>
>replace 2008011412 with current date, last 2 digits are hour 00, 06, 12,
>18
>
>this is a list of my experiences.
>
>first with gdalinfo I noticed that not enough info was reported to
>choose the band(s) you want
>
>degrib reports:
> 2.0, 158408, 2, TMP="Temperature [K]", 1000-ISBL, 11/05/2007 00:00,
>11/05/2007 00:00, 0.00
>
>gdalinfo reports:
>Band 2 Block=720x1 Type=Float64, ColorInterp=Undefined
>  Description = 1000[Pa] ISBL="Isobaric surface"
>
>
>There seems to be no way to specify the units. degrib has a feature to
>convert the values in the grid to english or metric units. For instance
>in the band above default would have the values in kelvin, specifying
>english would convert those values to Fahrenheit, and metric would
>convert to celsius.
>
>
>Two part messages seem to be represented by two different bands. This is
>not necessarily a problem but without the variable name being reported
>by gdalinfo it is hard to find what your looking for.
>
>degrib reports:
>8.0, 901490, 2, UGRD="u-component of wind [m/s]", 2000-ISBL, 11/05/2007
>00:00, 11/05/2007 00:00, 0.00
>8.1, 901490, 2, VGRD="v-component of wind [m/s]", 2000-ISBL, 11/05/2007
>00:00, 11/05/2007 00:00, 0.00
>
>gdalinfo reports:
>Band 8 Block=720x1 Type=Float64, ColorInterp=Undefined
>  Description = 2000[Pa] ISBL="Isobaric surface"
>Band 9 Block=720x1 Type=Float64, ColorInterp=Undefined
>  Description = 2000[Pa] ISBL="Isobaric surface"
>
>
>I tried contouring band 169 (aka grib msg 147) 925 mb temps it seemed to
>contour only the northern hemisphere and nothing seemed in the right
>place. once again I may have the wrong band, but that would not effect
>hemisphere only the contours
>
>
>thanks
>
>brian
>
>
>
>_______________________________________________
>gdal-dev mailing list
>gdal-dev at lists.osgeo.org
>http://lists.osgeo.org/mailman/listinfo/gdal-dev
>  
>

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/gdal-dev/attachments/20080117/96c65d28/attachment-0001.html