[R-sig-Geo] Converting .hdf file to raster layer
Loïc
loic.dutrieux at wur.nl
Thu Aug 21 12:56:07 CEST 2014
Hi Arnold,
GDAL provided with the rgdal binaries for windows isn't built with the
hdf4 driver. Theoretically it is possible to build rgdal from source
against GDAL contained in OSGeo4W or QGIS but from what I've heard it is
not trivial at all.
One alternative for reading hdf layer into R when working in windows is
by converting it to another format (e.g.: GeoTiff) using gdal_translate
(either from command line or via gdalUtils) and reading that new layer
into R.
The lines below illustrate briefly the two step process.
library(raster)
library(gdalUtils)
sds <- get_subdatasets('/full/path/to/file.hdf')
# translate first subdataset of hdf file to tiff
gdal_translate(sds[1], dst_dataset = 'out.tif')
r <- raster('out.tif')
Cheers,
Loïc
On 21/08/2014 11:31, Ahmadou Dicko wrote:
> I suspect that your GDAL was not build with the HDF4 driver.
>
> What is the output of :
> grep("hdf4", gdalDrivers()$name, ignore.case = TRUE, value = TRUE)
>
>
>
> On Thu, Aug 21, 2014 at 7:52 AM, Arnold Salvacion <
> arnold_salvacion at yahoo.com> wrote:
>
>> Hi Ahmadou,
>>
>> Have tried your code but still no success. However, I have actually tried
>> to download the .hdf file directly from net. What I did is import the data
>> from local folder. Below is the error that I got:
>>
>> Error in .local(.Object, ...) :
>> `C:\Users\Owner\Documents\VHP.G16.C07.NC.P1982006.VH.hdf' not recognised as a supported file format.
>>
>>
>> Best regards,
>>
>> Arnold
>>
>>
>>
>>
>> On Friday, August 15, 2014 2:33 AM, Ahmadou Dicko <
>> dicko.ahmadou at gmail.com> wrote:
>>
>>
>> You can use the gdalUtils package to get information on datasets inside
>> your hdf file and then stack them using the raster package
>>
>> library(raster)
>> library(rgdal) ## built with HDF4 support
>> library(gdalUtils)
>>
>> fsrc <- "
>> ftp://ftp.star.nesdis.noaa.gov/pub/corp/scsb/wguo/data/VHP_16km/VH/VHP.G16.C07.NC.P1981035.ND.hdf
>> "
>> f <- file.path("/tmp", basename(fsrc))
>> if (!file.exists(f)) download.file(fsrc, f, mode = "wb")
>>
>> info <- gdalinfo(f)
>> data <- grep("HDF4_SDS", info, value = TRUE)
>> data <- gsub("SUBDATASET_\\d+_NAME=|\\s+", "", data)
>>
>>
>> sds <- stack(data)
>> sds
>> ## class : RasterStack
>> ## dimensions : 904, 2500, 2260000, 3 (nrow, ncol, ncell, nlayers)
>> ## resolution : 1, 1 (x, y)
>> ## extent : 0, 2500, 0, 904 (xmin, xmax, ymin, ymax)
>> ## coord. ref. : NA
>> ## names : VHP.G16.C07.NC.P1981035.ND.1,
>> VHP.G16.C07.NC.P1981035.ND.2, VHP.G16.C07.NC.P1981035.ND.3
>> ## min values : -32768,
>> -32768, 0
>> ## max values : 32767,
>> 32767, 255
>>
>> ##You can use RasterBrick too, by doing `brick(sds)`
>>
>> sessionInfo()
>> ## R version 3.1.1 Patched (2014-08-12 r66349)
>> ## Platform: x86_64-unknown-linux-gnu (64-bit)
>>
>> ## locale:
>> ## [1] LC_CTYPE=en_US.utf8 LC_NUMERIC=C
>> ## [3] LC_TIME=en_US.utf8 LC_COLLATE=en_US.utf8
>> ## [5] LC_MONETARY=en_US.utf8 LC_MESSAGES=en_US.utf8
>> ## [7] LC_PAPER=en_US.utf8 LC_NAME=C
>> ## [9] LC_ADDRESS=C LC_TELEPHONE=C
>> ## [11] LC_MEASUREMENT=en_US.utf8 LC_IDENTIFICATION=C
>>
>> ## attached base packages:
>> ## [1] stats graphics grDevices utils datasets methods
>> ## [7] base
>>
>> ## other attached packages:
>> ## [1] gdalUtils_0.3.2 rgdal_0.8-16 raster_2.2-32
>> ## [4] sp_1.0-15
>>
>> ## loaded via a namespace (and not attached):
>> ## [1] codetools_0.2-8 compiler_3.1.1 foreach_1.4.2
>> ## [4] grid_3.1.1 iterators_1.0.7 lattice_0.20-29
>> ## [7] R.methodsS3_1.6.1 R.oo_1.18.0 R.utils_1.32.4
>> ## [10] tools_3.1.1
>>
>> grep("hdf4", gdalDrivers()$name, ignore.case = TRUE, value = TRUE)
>> ## [1] "HDF4" "HDF4Image"
>>
>>
>>
>> On Thu, Aug 14, 2014 at 8:52 AM, Michael Sumner <mdsumner at gmail.com>
>> wrote:
>>
>> Just had a look, the file has subdatasets and seemingly raster can't
>> handle those. (I haven't explored this much but might have a chance
>> to).
>>
>> library(raster)
>> library(rgdal) ## built with HDF4, HDF5 and NetCDF4
>>
>> fsrc <- "
>> ftp://ftp.star.nesdis.noaa.gov/pub/corp/scsb/wguo/data/VHP_16km/VH/VHP.G16.C07.NC.P1981035.ND.hdf
>> "
>> f <- basename(fsrc)
>> if (!file.exists(f)) download.file(fsrc, f, mode = "wb")
>>
>> ## found with system GDAL, see below
>> sds <- c('HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":0',
>> 'HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":1',
>> 'HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":2'
>> )
>>
>> ##r <- brick(sds)
>> ##Error in GDALinfo(filename, silent = silent, returnRAT = RAT,
>> returnCategoryNames = RAT) :
>> ## object 'RATlist' not found
>> ##Error in .rasterObjectFromFile(x, band = band, objecttype =
>> "RasterLayer", :
>> ## Cannot create a RasterLayer object
>> from this file.
>> ## In addition: Warning message:
>> ## In dim(x) : no bands in dataset
>>
>>
>> ## try with rgdal directly, so far so good but
>> ## no spatial-reference (again see gdalinfo output below)
>> raster(readGDAL(sds[1]) )
>>
>> raster(readGDAL(sds[1]) )
>> HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":0 has GDAL driver
>> HDF4Image
>> and has 904 rows and 2500 columns
>> class : RasterLayer
>> dimensions : 904, 2500, 2260000 (nrow, ncol, ncell)
>> resolution : 1, 1 (x, y)
>> extent : 0, 2500, 0, 904 (xmin, xmax, ymin, ymax)
>> coord. ref. : NA
>> data source : in memory
>> names : band1
>> values : -99.9, 322 (min, max)
>>
>> Warning message:
>> In readGDAL(sds[1]) : GeoTransform values not available
>>
>>
>> ## proceed without caution, again values from gdalinfo
>> b <- brick(stack(sapply(sds, function(x) raster(readGDAL(x)))))
>> proj4string(b) <- CRS(" +proj=longlat +ellps=WGS84 +datum=WGS84
>> +no_defs +towgs84=0,0,0")
>> extent(b) <- c(-180, 180, -55.15200043, 75.02400208)
>>
>> library(maptools
>> data(wrld_simpl)
>> plot(b, addfun = function() plot(wrld_simpl, add = TRUE))
>>
>> That looks ok to me, I don't know if it's WGS84 or the sphere or
>> something else.
>>
>> Cheers, Mike.
>>
>>
>>
>>
>>
>> system(sprintf("gdalinfo %s", f))
>>
>> Driver: HDF4/Hierarchical Data Format Release 4
>> Files: VHP.G16.C07.NC.P1981035.ND.hdf
>> Size is 512, 512
>> Coordinate System is `'
>> Metadata:
>> ANCILLARY_FILES=FILE_CONFIGURE:vh.config
>> FILE_PRELAUNCH_CALIBRATION:../ancillary/AVHRR_calibration_prelaunch.txt
>> FILE_POSTLAUNCH_CALIBRATION:../ancillary/AVHRR_calibration_postlaunch.txt
>> FILE_IGBP_LANDTYPE:../ancillary/igbp_landtype_GVIx.hdf
>> FILE_METADATA_REGIONS:../ancillary/regions_for_metadata.txt
>> FILE_EDF_NDVI:../ancillary/NVI_counts_ByLine_G04.hdf
>> FILE_EDF_BT4:../ancillary/BT4_counts_ByLine_G04.hdf
>>
>> CITATION_TO_DOCUMENTS=User Guide of Vegetation Health(VH) system
>> (version 1.3, March 21 2012)
>> CONFIGURE_FILE_CONTENT=[Options for vh.exe]
>> DIR_Ancillary= ../ancillary
>> DIR_GVI= data/VH_unitTest1/weekly
>> DIR_VH= data/VH_unitTest1/G04
>> DIR_VH_META= data/VH_unitTest1/G04/meta
>> DIR_CLIMAT= data/VH_EDF_v2/climate_G04
>> DIR_CLIMAT_META= data/VH_EDF_v2/climate_G04/meta
>> FILE_PREFIX= VHP
>> Input_Data_Type= VHP
>> ResolutionString= G04
>> Days_Per_Period= 7
>> FilterSize= 15
>> applyEDFonNDVI= 1 # 0: none, 1: linebyline for NVI;
>> applyEDFonBT= 1
>> Instrument= AVHRR # so far AVHRR is the only option
>> FORMAT_GVI= NETCDF #or HDF4
>> FORMAT_CLIMAT= HDF4
>> FORMAT_ND= NETCDF
>> FORMAT_SM= NETCDF
>> FORMAT_VH= NETCDF
>>
>> [Periods of GVI data used for VH]
>> # this section control which satellite will be used for calculating
>> ND, SM and VH
>> #satID satNumber yearWeek1 yearWeek2
>> NC 07 198135 198449
>> NF 09 198509 198844
>> NH 11 198846 199436
>> NJ 14 199504 200052
>> NL 16 200101 200401
>> NL 16 200405 200410
>> NL 16 200425 200428
>> NL 16 200430 200523
>> NN 18 200524 201052
>> NP 19 201101 399999
>>
>> [Periods of AVHRR data used for GVI climatology]
>> #this section controls which satellite will be used for creating VH
>> climatology
>> #satID satNumber yearWeek1 yearWeek2
>> NC 07 198142 198450
>> NF 09 198515 198752
>> NH 11 198920 199252
>> NJ 14 199520 199952
>> NL 16 200120 200252
>>
>> [END]
>>
>>
>> CONTACT=NOAA/NESDIS/STAR/EMB
>> DATE_BEGIN=239
>> DATE_END=245
>> DAYS_PER_PERIOD=7
>> END_LATITUDE_RANGE=-55.15200043
>> END_LONGITUDE_RANGE=180
>> INPUT_FILENAMES=data/VH_unitTest1/weekly/
>> VHP.G04.C07.NC.P1981035.S239.E245.nc
>> <http://vhp.g04.c07.nc.p1981035.s239.e245.nc/>
>>
>> INPUT_FILES=1
>> INSTRUMENT=AVHRR
>> PERIOD_OF_YEAR=35
>> PRODUCT_NAME=Vegetation Health
>> PROJECTION=Plate_Carree
>> SATELLITE=NC
>> START_LATITUDE_RANGE=75.02400208
>> START_LONGITUDE_RANGE=-180
>> TIME_BEGIN=00:00 UTC (use day time data only)
>> TIME_END=23:59 UTC (use day time data only)
>> VERSION=VH (vh.exe,version 1.3, March 21 2012)
>> YEAR=1981
>> Subdatasets:
>> SUBDATASET_1_NAME=HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":0
>> SUBDATASET_1_DESC=[904x2500] BT4 (16-bit integer)
>> SUBDATASET_2_NAME=HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":1
>> SUBDATASET_2_DESC=[904x2500] NDVI (16-bit integer)
>> SUBDATASET_3_NAME=HDF4_SDS:UNKNOWN:"VHP.G16.C07.NC.P1981035.ND.hdf":2
>> SUBDATASET_3_DESC=[904x2500] QA (8-bit unsigned integer)
>> 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)
>>
>> On Thu, Aug 14, 2014 at 10:27 AM, Arnold Salvacion
>> <arnold_salvacion at yahoo.com> wrote:
>>> Dear Colleagues,
>>>
>>> Does any here have already experience converting NOAA AVHRR VHP Product
>> in .hdf (
>> ftp://ftp.star.nesdis.noaa.gov/pub/corp/scsb/wguo/data/VHP_16km/VH/ )
>> format to a raster layer in R?
>>>
>>>
>>> Best regards,
>>>
>>> Arnold
>>> [[alternative HTML version deleted]]
>>>
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>
>>
>>
>> --
>> Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> Hobart, Australia
>> e-mail: mdsumner at gmail.com
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>>
>>
>>
>> --
>> Ahmadou H. DICKO
>> statistician economist (Ingénieur Statisticien Économiste)
>> PhD candidate in Climate change economics
>> Faculty of economics and managment - Cheikh Anta Diop University
>> West African Science Service Center on Climate Change and Adaptated Land
>> Use (WASCAL)
>> Center for Development Research (ZEF) - University of Bonn
>>
>> email :
>> ahmadou.dicko at ucad.edu.sn
>> twitter : @dickoah
>> github : github/dickoa <https://github.com/dickoa>
>> tel : +221 33 827 55 16
>> portable: +221 77 123 81 69
>>
>>
>>
>>
>
>
More information about the R-sig-Geo
mailing list