[R-sig-Geo] reading HDF5 files

Tim Keitt tkeitt at gmail.com
Tue Aug 14 19:05:54 CEST 2007 

On 8/14/07, Paul Hiemstra <p.hiemstra at geo.uu.nl> wrote:
> Mike,
>
> I tried, and failed, to directly load the HDF5 files into R using rgdal.
> The next thing I tried was to convert the .h5 file to a GTiff using
> FWTools. This seemed to work, but when I tried to load the file into R
> through rgdal I got the following error message:
>
> test = readGDAL("regenradar_hdf5/RAD_NL21_PCP_NA_200707200000.tif")
>
> > test = readGDAL("RAD_NL21_PCP_NA_200707200000.tif")
> RAD_NL21_PCP_NA_200707200000.tif has GDAL driver GTiff
> and has 256 rows and 256 columns
> Closing GDAL dataset handle 0x0226a220...  destroyed ... done.
> Error in `proj4string<-`(`*tmp*`, value = <S4 object of class "CRS">) :
>         Geographical CRS given to non-conformant data
> > traceback()
> 8: stop("Geographical CRS given to non-conformant data")
> 7: `proj4string<-`(`*tmp*`, value = <S4 object of class "CRS">)
> 6: SpatialGrid(grid, proj4string)
> 5: initialize(value, ...)
> 4: initialize(value, ...)
> 3: new("SpatialGridDataFrame", SpatialGrid(grid, proj4string), data = data)
> 2: SpatialGridDataFrame(grid = grid, data = data.frame(df), proj4string = CRS(p4s))
> 1: readGDAL("RAD_NL21_PCP_NA_200707200000.tif")
> > sessionInfo()
> R version 2.5.1 (2007-06-27)
> i386-pc-mingw32

This seems a case where failing gracefully would save a lot of
trouble. Perhaps a warning and a NULL CRS? Then things can be patched
up in R. I am not a fan of the defensive programming style prevalent
in R. (Although I fully recognize that it is a matter of taste and
many people would disagree with me.)

THK


> locale:
> LC_COLLATE=Dutch_Netherlands.1252;LC_CTYPE=Dutch_Netherlands.1252;LC_MONETARY=Dutch_Netherlands.1252;LC_NUMERIC=C;LC_TIME=Dutch_Netherlands.1252
>
> attached base packages:
> [1] "stats"     "graphics"  "grDevices" "utils"     "datasets"  "methods"   "base"
>
> other attached packages:
>    rgdal       sp
> "0.5-14" "0.9-14"
>
> Do you have any idea what caused this error and how I can circumvent it?
>
> Many thanks in advance,
>
> Paul
>
> Michael Sumner schreef:
> > Note that readGDAL can read the subdataset directly (escaping the quotes
> > appropriately):
> >
> > d <- readGDAL("HDF5:\"RAD_NL21_PCP_NA_200707200000.h5\"://image1/image_data")
> >
> > It also seems that your gdal_translate to TIFF output has the extension
> > .img, but then you attempt to read from a tif.
> >
> > Just a thought.
> >
> > Cheers, Mike.
> >
> > ==============Original message text===============
> > On Fri, 10 Aug 2007 19:09:00 +1000 Roger Bivand wrote:
> >
> > Paul,
> >
> > I have added a CC to the R-sig-geo list, since I myself have no experiece
> > of using this format in practice. If you are not a member, please join at:
> >
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
> > so that you will see any replies.
> >
> > Did you try gdalinfo on the output file?
> >
> > The error message suggests that the output file has a geographical
> > coordinate reference system, with coordinate values outside the feasible
> > range. From the HDF5 description, however, it ought to be projected,
> > right? Do you need to make more use of the arguments taken by gdalinfo to
> > give the output file the correct projected CRS?
> >
> > Roger
> >
> > On Fri, 10 Aug 2007, Paul Hiemstra wrote:
> >
> >
> >> Dear Roger,
> >>
> >> I tried the route through a GTiff before I sent you the e-mail but
> >>
> > somehow it
> >
> >> did not work. When I use gdalinfo on one of the .h5 files I get this:
> >>
> >> PS E:\LSO data\Regenbui 20-7\regenradar_rasters\regenradar_hdf5> gdalinfo
> >> RAD_NL21_PCP_NA_200707200000.h5
> >> Driver: HDF5/Hierarchical Data Format Release 5
> >> Files: RAD_NL21_PCP_NA_200707200000.h5
> >> Size is 512, 512
> >> Coordinate System is `'
> >> Metadata:
> >>  map_projection:projection_indication=Y
> >>  map_projection:projection_name=STEREOGRAPHIC
> >>  map_projection:projection_proj4_params=+proj=stere +lat_0=90
> >> +lon_0=0.0 +lat_ts=60.0 +a=6378.388 +b=6356.912 +x_0=0 +y_0=0
> >> geographic: geo_number_rows=256
> >> geographic: geo_number_columns=256
> >> geographic: geo_pixel_size_x=2.500000
> >> geographic: geo_pixel_size_y=-2.500000
> >> geographic: geo_par_pixel=X,Y
> >> geographic: geo_dim_pixel=KM,KM
> >> geographic: geo_column_offset=0.000000
> >> geographic: geo_row_offset=1490.906006
> >> geographic: geo_pixel_def=LU
> >> geographic: geo_product_corners=0.000000 49.769001 0.000000 55.296001
> >> 9.743000 54.818001 8.337000 49.373001
> >> calibration: calibration_flag=Y
> >> calibration: calibration_formulas=GEO = 0.5 * PV + -31.5
> >> calibration: calibration_missing_data=255
> >> calibration: calibration_out_of_image=255
> >> statistics: stat_min_value=-1.500000
> >> statistics: stat_max_value=45.000000
> >> image1: image_product_name=RAD_NL21_PCP_H0.8_NA
> >> image1: image_geo_parameter=REFLECTIVITY_[DBZ]
> >> image1: image_size=65536d
> >> image1: image_bytes_per_pixel=1d
> >> overview: product_group_name=RAD_NL21_PCP_NA
> >> overview: products_missing=N/A
> >> overview: hdftag_version_number=3.5
> >> overview: number_image_groups=1
> >> overview: number_satellite_groups=0d
> >> overview: number_radar_groups=2d
> >> overview: number_station_groups=0
> >> overview: product_datetime_start=20-JUL-2007;00:00:11.000
> >> overview: product_datetime_end=20-JUL-2007;00:00:11.000
> >> overview: number_lightning_groups=0d
> >> overview: number_classification_groups=0d
> >> overview: number_grid_groups=0d
> >> overview: number_point_groups=0d
> >> overview: number_vector_groups=0d
> >> radar1: radar_name=De_Bilt
> >> radar1: radar_location=5.179000 52.103001
> >> radar2: radar_name=Den_Helder
> >> radar2: radar_location=4.790000 52.955002
> >> Subdatasets:
> >> SUBDATASET_0_NAME=HDF5:"RAD_NL21_PCP_NA_200707200000.h5"://image1/image_data
> >> SUBDATASET_0_DESC=[256x256] //image1/image_data (8-bit unsigned character)
> >> Corner Coordinates:
> >> Upper Left  (    0.0,    0.0)
> >> Lower Left  (    0.0,  512.0)
> >> Upper Right (  512.0,    0.0)
> >> Lower Right (  512.0,  512.0)
> >> Center      (  256.0,  256.0)
> >>
> >> Then I use gdal_translate:
> >> PS E:\LSO data\Regenbui 20-7\regenradar_rasters\regenradar_hdf5>
> >> gdal_translate -of GTiff
> >> HDF5:"RAD_NL21_PCP_NA_200707200000.h5"://image1/image_data
> >> RAD_NL21_PCP_NA_200707200000.img
> >>
> >> Then I try to read the file using rgdal produces the following error:
> >>
> >>
> >>>  test = readGDAL("regenradar_hdf5/RAD_NL21_PCP_NA_200707200000.tif")
> >>>
> >> regenradar_hdf5/RAD_NL21_PCP_NA_200707200000.tif has GDAL driver GTiff
> >> and has 256 rows and 256 columns
> >> Closing GDAL dataset handle 0x0130b898...  destroyed ... done.
> >> Error in `proj4string<-`(`*tmp*`, value = <S4 object of class "CRS">) :
> >>        Geographical CRS given to non-conformant data
> >>
> >> Do you have an idea how to solve this?
> >>
> >> kind regards,
> >>
> >> Paul
> >>
> >>
> >> Roger Bivand schreef:
> >>
> >>>  Dear Paul,
> >>>
> >>>  There is information in the README.windows file in the package, this
> >>>  should display it:
> >>>
> >>>  file.show(system.file("README.windows", package="rgdal"))
> >>>
> >>>  If you do not want to build against the FWTools binary, then just install
> >>>  FWTools, and use the gdal_translate command line tool to convert your data
> >>>  to a format that the standard rgdal windows binary can read, such as
> >>>  gTiff. You'll need to play around a bit to find out what the data are
> >>>  called.
> >>>
> >>>  On Linux, this is easier, because you can build GDAL against HDF5 before
> >>>  installing rgdal, and there it "just works".
> >>>
> >>>  Hope this helps,
> >>>
> >>>  Roger
> >>>
> >>>  On Thu, 9 Aug 2007, Paul Hiemstra wrote:
> >>>
> >>>
> >>>>  Dear mr. Bivand
> >>>>
> >>>>  My name is Paul Hiemstra, a PhD student of Edzer. I am trying to read
> >>>>  HDF5 files through rgdal. If I use GDALinfo() on one of the files it
> >>>>  produces the following error:
> >>>>
> >>>>  Error in .local(.Object, ...) :
> >>>>        GDAL Error 4: `RAD_NL21_PCP_NA_200707200000.h5' not recognised as
> >>>>  a supported file format.
> >>>>
> >>>>  On the mailing list
> >>>>  (http://finzi.psych.upenn.edu/R/Rhelp02a/archive/92270.html) you >> >  metioned using the dll files from FWTools to get support for HDF5 into
> >>>>  rgdal. How can this be done?
> >>>>
> >>>>  I also tried the hdf5 package on CRAN, but that does not seem to work,
> >>>>  returning an empty object.
> >>>>
> >>>>  many thanks,
> >>>>
> >>>>  Paul
> >>>>
> >>>>
> >>>>
> >>
> >>
> >
> >
>
>
> --
> Drs. Paul Hiemstra
> Department of Physical Geography
> Faculty of Geosciences
> University of Utrecht
> Heidelberglaan 2
> P.O. Box 80.115
> 3508 TC Utrecht
> Phone:  +31302535773
> Fax:    +31302531145
> http://intamap.geo.uu.nl/~paul
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>


-- 
Timothy H. Keitt, University of Texas at Austin
Contact info and schedule at http://www.keittlab.org/tkeitt/
Reprints at http://www.keittlab.org/tkeitt/papers/
ODF attachment? See http://www.openoffice.org/

More information about the R-sig-Geo mailing list