On Fri, 16 Jan 2009, Sebastian P. Luque wrote:
Hi,
I'm trying to use HDF5 files with this structure:
Sorry for the delay. Having now set up a reprodicable route (OSGeo4W GDAL
and drivers, rgdal built against OSGeo4W as in a posting a couple of days
ago), I can see the problem, but have no resolution.
Doing gdalinfo on the data itself, I see geographical coordinates, not
projected. GDAL does not retrieve the grid metadata correctly, and assigns
default 1x1 grid cells starting at c(0.5, 0.5). Using gdal_translate to
convert it to a GTiff, I see that I have to set p4s to NA to read at all,
and the problems remain. I would be interested in knowing who or what
wrote the HDF5 file, because it doesn't seem to be what it says it is (no
grid metadata, different CRS from your assumption), and so the data may
actually be in the wrong order too.
Roger
PS - I can't see that you've asked on the GDAL list, hence the reply here.
---<--------------------cut here---------------start------------------->---
$ gdalinfo conc_200901011200.hdf
Driver: HDF5/Hierarchical Data Format Release 5
Files: conc_200901011200.hdf
Size is 512, 512
Coordinate System is `'
Subdatasets:
SUBDATASET_0_NAME=HDF5:"conc_200901011200.hdf"://Data/data[00]
SUBDATASET_0_DESC=[760x1120] //Data/data[00] (32-bit floating-point)
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)
---<--------------------cut here---------------end--------------------->---
The file (467 kb placed in
http://members.shaw.ca/sluque/conc_200901011200.hdf) is a grid of a
single variable (ice concentration) in a well defined polar
stereographic projection, that loads ok in R, via rgdal's readGDAL():
---<--------------------cut here---------------start------------------->---
R> p4s <- "+proj=stere +a=6378273 +b=6356889.44891 +lat_0=90 +lat_ts=70
+lon_0=-45"
R> ice <- readGDAL("HDF5:\"conc_200901011200.hdf\"://Data/data[00]", p4s=p4s)
HDF5:"conc_200901011200.hdf"://Data/data[00] has GDAL driver HDF5Image
and has 760 rows and 1120 columns
R> summary(ice)
Object of class SpatialGridDataFrame
Coordinates:
min max
x 0 1120
y 0 760
Is projected: TRUE
proj4string :
[+proj=stere +a=6378273 +b=6356889.44891 +lat_0=90 +lat_ts=70 +lon_0=-45]
Number of points: 2
Grid attributes:
cellcentre.offset cellsize cells.dim
x 0.5 1 1120
y 0.5 1 760
Data attributes:
Min. 1st Qu. Median Mean 3rd Qu. Max.
-32800 -99 -99 -117 0 100
R> ## We are concerned only with the 0-100 range; other values are codes with
R> ## various meanings, so we remove them
R> i...@data[ice@data[[1]] < 0 | i...@data[[1]] > 100, ] <- NA
R> summary(ice)
Object of class SpatialGridDataFrame
Coordinates:
min max
x 0 1120
y 0 760
Is projected: TRUE
proj4string :
[+proj=stere +a=6378273 +b=6356889.44891 +lat_0=90 +lat_ts=70 +lon_0=-45]
Number of points: 2
Grid attributes:
cellcentre.offset cellsize cells.dim
x 0.5 1 1120
y 0.5 1 760
Data attributes:
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
0.0 0.0 0.0 23.3 36.7 100.0 478521.0
## but something is wrong with visualizing these data
pdf("~/tmp/ice.pdf")
spplot(ice)
dev.off()
R> sessionInfo()
R version 2.8.1 (2008-12-22)
x86_64-pc-linux-gnu
locale:
LC_CTYPE=en_CA.UTF-8;LC_NUMERIC=C;LC_TIME=en_CA.UTF-8;LC_COLLATE=en_CA.UTF-8;LC_MONETARY=C;LC_MESSAGES=en_CA.UTF-8;LC_PAPER=en_CA.UTF-8;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHONE=C;LC_MEASUREMENT=en_CA.UTF-8;LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] rgdal_0.5-35 sp_0.9-29 slmisc_0.7.0 lattice_0.17-20
loaded via a namespace (and not attached):
[1] grid_2.8.1
---<--------------------cut here---------------end--------------------->---
I'm not attaching the resulting file because it's too large. However,
the plot should look like http://members.shaw.ca/sluque/full.ps I'd
appreciate any feedback on what the problem might be. Thanks.
Cheers,
--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: roger.biv...@nhh.no
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
