Converting .hdf file to raster layer

Hello All,

Can any one can help me to convert this Lidar data latitude and longitude
in Degree Decimal format using R.

  latitude longitude elevation intensity

1 204140.7   4733835    334.14       6.4
2 204150.1   4733849    334.86      11.4

Appreciate your help.


Vinit



On Thu, Aug 14, 2014 at 3: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

  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]]

--
Michael Sumner
Software and Database Engineer
Australian Antarctic Division
Hobart, Australia
e-mail: mdsumner at gmail.com

        [[alternative HTML version deleted]]