[R-sig-Geo] reading HDF5 files
Paul Hiemstra
p.hiemstra at geo.uu.nl
Tue Aug 14 11:38:24 CEST 2007
Dear Roger/Mike,
Thanks for all the feedback! This is the full gdal_translate command I used:
gdal_translate -of GTiff
HDF5:"RAD_NL21_PCP_NA_200707202345.h5"://image1/image_data
../regenradar_tif/RAD_NL21_PCP_NA_20070720234.tif
And the gdalinfo results
HDF5:
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)
GTiff
E:\LSO data\Regenbui 20-7\regenradar_rasters\regenradar_tif>gdalinfo
RAD_NL21_PCP_NA_200707200000.tif
Driver: GTiff/GeoTIFF
Files: RAD_NL21_PCP_NA_200707200000.tif
RAD_NL21_PCP_NA_200707200000.aux
Size is 256, 256
Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.2572235630016,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Metadata:
AREA_OR_POINT=Area
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_product_corners=0.000000 49.769001 0.000000 55.296001
9.743000 54.818001 8.337000 49.373001
calibration=calibration_out_of_image=255
statistics=stat_max_value=45.000000
image1=image_bytes_per_pixel=1d
overview=number_vector_groups=0d
radar1=radar_location=5.179000 52.103001
radar2=radar_location=4.790000 52.955002
Corner Coordinates:
Upper Left ( 0.0, 0.0)
Lower Left ( 0.0, 256.0)
Upper Right ( 256.0, 0.0)
Lower Right ( 256.0, 256.0)
Center ( 128.0, 128.0)
Band 1 Block=256x32 Type=Byte, ColorInterp=Gray
Metadata:
image_data=DISPLAY_ORIGIN=UL
thanks,
Paul
Roger Bivand schreef:
> On Tue, 14 Aug 2007, Paul Hiemstra 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
>
> Please state the gdalinfo of the input HDF5, and the output gTiff.
> Your use of gdal_translate seems to be mangling the (IIRC) projected
> coordinate reference system of the HDF5 to a geographical coordinate
> reference system, with the result that the coordinate values are "out
> of bounds". I think you need to show the full gdal_translate command,
> and to look at the help page:
>
> http://www.gdal.org/gdal_translate.html
>
> and try out the -a_srs option. Maybe also look at the -sds option if
> you are not sure that the subset is being chosen correctly.
>
> Once gdalinfo on the output gTiff shows a geographical coordinate
> reference system, rgdal will work (unless this is near-polar and you
> have a raster cell boundary > |90|, which does happen sometimes - then
> you need to adjust the input data).
>
> Roger
>
>
>>> 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
>>
>> 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
More information about the R-sig-Geo
mailing list