inconsistencies between readGDAL and other software (GrADS) when reading GRIB1 files
On Thu, 17 Nov 2011, Ariel Ortiz-Bobea wrote:
Hi, I've been importing GRIB1 files into R through the readGDAL command in the GDAL package.
The readGDAL() function in the rgdal package uses GDAL and its drivers, so see: http://www.gdal.org/frmt_grib.html for details. Note also: https://stat.ethz.ch/pipermail/r-sig-geo/2011-January/010760.html http://grass.osgeo.org/wiki/GRIB for a bug in the GDAL driver and a workaround, and other GRIB threads on this list.
However, the GRIB1 files are multi-band and during the import process I "lose" the variable names and get generic names: band1, band2, etc.
As far as I can see, band names are not supported in GDAL in: http://www.gdal.org/classGDALDataset.html or http://www.gdal.org/classGDALRasterBand.html but may be in http://www.gdal.org/classGDALRasterAttributeTable.html. Could you please try running: GDALinfo("mygrib.grb", returnRAT=TRUE) to see if you get the band names? Roger
I have tried to circumvent this by using other software (GrADS) to obtain summary statistics for each variable (of which I know the names) and then compare these to summary statistics of each variable/band in R. Matching summary statistics would allow me to get back the real names. However, there are discrepancies between summary statistics. More precisely, there are differences between minimums, means and maximums between both approaches AND these differences are not of the same magnitude (ruling out a homogenous shift in values). Moreover, GrADS gives temperatures in Kelvin and readGDAL seems to be transforming them to Celsius BUT the difference between minimum, maximum and means is not consistently 273.15 as one would expect. For instance, for soil temperature the difference is 271.441 (for the mean), 270.514 (for the max) and 273.149 (for the min), but for Dew point it is 269.436 (for the mean), 272.58 (for the max) and 273.349 (for the min). So it is not even consistent across different temperature variables. I wonder if this has to do with how GrADS or readGDAL read these GRIB1 files (or the way I'm importing this)... would there be any hidden processing going on that would explain this? Would using netCDF format files rather than GRIB1 help in any way? Any suggestions or comments on how to import this data in a manner in a proper way? Thanks a lot in advance for any guidance on this, Ariel Here is my example file: http://dl.dropbox.com/u/45311184/US48_merged_AWIP32.1985120100.RS.flx (in GRIB1 format) R code: ------- # import flx<- readGDAL("1985/US48_merged_AWIP32.1985120100.RS.flx") # fix missing values (imported as 9999) for (n in 1:10 ){ is.na(flx[[n]]) <- flx[[n]] > 9000 } # print summary statistics for each variable/band for (n in 1:10 ){ print(summary(flx[[n]])[c(4,6,1)]) } GrADS: -------- code: http://dl.dropbox.com/u/45311184/gradscode.txt files to be included in the same folder with the .flx file http://dl.dropbox.com/u/45311184/US48_merged_AWIP32.1985120100.RS.flx.ctl http://dl.dropbox.com/u/45311184/US48_merged_AWIP32.1985120100.RS.flx.idx ----- Ariel Ortiz-Bobea PhD Candidate in Agricultural & Resource Economics University of Maryland - College Park -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/inconsistencies-between-readGDAL-and-other-software-GrADS-when-reading-GRIB1-files-tp7006160p7006160.html Sent from the R-sig-geo mailing list archive at Nabble.com.
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand Department of Economics, NHH Norwegian School of Economics, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: Roger.Bivand at nhh.no