[gdal-dev] question about reading grib2 file

[Cassanova, Bill]

Hi Matt,

Can you be a little more explicit about what you are asking?  Generally speaking GRIB2 files do not contain explicit lat-lon information.
GRIB2's contain a projection and within that projection the bounding box including the number of points along each axis is contained.

As an experiment run gdalinfo on your grib2 file and you will see output similar to the below:

gdalinfo foo.grib2
Driver: GRIB/GRIdded Binary (.grb)
Files: foo.grib2
Size is 720, 361
Coordinate System is:
GEOGCS["Coordinate System imported from GRIB file",
    DATUM["unknown",
        SPHEROID["Sphere",6371229,0]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433]]
Origin = (0.000000000000000,270.000000000000000)
Pixel Size = (0.500000000000000,-0.500000000000000)
Corner Coordinates:
Upper Left  (       0.000,     270.000) (  0d 0'0.01"E,270d 0'0.00"N)
Lower Left  (   0.0000000,  89.5000000) (  0d 0'0.01"E, 89d30'0.00"N)
Upper Right (     360.000,     270.000) (360d 0'0.00"E,270d 0'0.00"N)
Lower Right (     360.000,      89.500) (360d 0'0.00"E, 89d30'0.00"N)
Center      (     180.000,     179.750) (180d 0'0.00"E,179d45'0.00"N)
Band 1 Block=720x1 Type=Float64, ColorInterp=Undefined
  Description = 0[-] SFC="Ground or water surface"
  Metadata:
    GRIB_UNIT=[K]
    GRIB_COMMENT=Temperature [K]
    GRIB_ELEMENT=TMP
    GRIB_SHORT_NAME=0-SFC
    GRIB_REF_TIME=  1295481600 sec UTC
    GRIB_VALID_TIME=  1296151200 sec UTC
    GRIB_FORECAST_SECONDS=669600 sec
    GRIB_PDS_PDTN=0
    GRIB_PDS_TEMPLATE_NUMBERS=0 0 2 0 96 0 0 0 1 0 0 0 186 1 0 0 0 0 0 255 0 0 0 0 0

If you are trying to read the values contained in Upper Left, Lower Left, etc?

If this is the case you need to get the geo_transform using code similar to the C-code below:

GDALAllRegister();


   GDALDatasetH my_data_set;

   int my_x_size, my_y_size;
   int my_raster_bands;

   GDALRasterBandH my_raster_band;

   double my_geo_transform[6];

   std::string my_projection_ref;

   my_data_set = NULL;
   my_raster_band = NULL;
   my_x_size = my_y_size = my_raster_bands = 0;

   ACE_LOG_MSG->log(LM_INFO, "%D %n Opening %s \n", input_grib2_file.c_str());

   my_data_set = GDALOpen(input_grib2_file.c_str(), GA_ReadOnly);


   if (NULL == my_data_set)
   {
      ACE_LOG_MSG->log(LM_INFO, "%D %n NULL GDALDatasetH pointer for file %s, "
                       "file invalid, quitting\n", input_grib2_file.c_str());
      return (1);
   }

   GDALGetGeoTransform(my_data_set, my_geo_transform);

   my_x_size = GDALGetRasterXSize(my_data_set);

   my_y_size = GDALGetRasterYSize(my_data_set);

   my_raster_bands = GDALGetRasterCount(my_data_set);

   ACE_LOG_MSG->log(LM_INFO,
                    "%D %n Data Set: x_size = %d y_size = %d raster_bands = %d origin = %f, %f, pixel_size=%f, %f\n",
                    my_x_size, my_y_size, my_raster_bands, my_geo_transform[0],
                    my_geo_transform[3], my_geo_transform[1],
                    my_geo_transform[5]);

>From this information you can calculate the lat-lon of each point within the file.

To answer your other question about wgrib2 -spread....wgrib2 calculates the lat-lon coordinates on the fly and does not store these explicitly in the file...Unless they are put there as a wgrib2 data segment which most people don't do.


Hope this helps,
Bill

-----Original Message-----
From: gdal-dev-bounces at lists.osgeo.org [mailto:gdal-dev-bounces at lists.osgeo.org] On Behalf Of Matt Funk
Sent: Thursday, January 20, 2011 4:38 PM
To: gdal-dev at lists.osgeo.org
Subject: [gdal-dev] question about reading grib2 file

Hi,

i am somewhat new to gdal and the grib format.
I need to read in a grib2 file which contains latitude longitude and
elevation.
So i do:
ds = gdal.Open(filename)
and
print "ds.GetRasterCount(): %s" % ds.GetRasterCount()
gives: 1

I can proceed to read in this band, but all the data that i get is the
elevation.
Not sure if this is an incredibly stupid question, but how do i access
the latitude/longitude data?

I know that there is latitude/longitude since when outputting things to
a csv file via the wgrib2 utility i get lat/lon and elevation?

thanks for any advice
matt
_______________________________________________
gdal-dev mailing list
gdal-dev at lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev