2

I need to correct the reflectance values of LANDSAT-8 data for steep terrain and find it too tedious and clunky in ArcGIS. The data I have is: LANSDAT-8 DEM (ASTER) Slope and aspect (Derived using DEM in ArcGIS) Sun Elevation and Azimuth (Obtained from L8 MTL file)

All of which are clipped to a polygon of my study area.

Two methods I wanted to use were the cosine and minneart corrections. I found a package for R called "landsat" - http://cran.r-project.org/web/packages/landsat/index.html thatcontains a function named "topocorr" - http://www.inside-r.org/packages/cran/landsat/docs/topocorr which takes the uncorrected image, slope, aspect, azimuth and elevation then produces a corrected image.

I'm struggling to use it however as the input data type for the input image is "matrix, data frame, or SpatialGridDataFrame" formats. How can I convert L8 data into this format?

Thanks.

  • 1
    Read it it with the readGDAL function in the rgdal package, which returns a SpatialGridDataFrame. – mdsumner Oct 19 '14 at 12:52
1

As pointed out by @mdsumner, You can read a SpatialGridDataFrame directly using readGDAL in the rgdal package. You can also easily coerce raster objects to a SpatialGridDataFrame or matrix. If your imagery is in separate files I would recommend reading them in as a stack and then coercing to a SpatialGridDataFrame.

require(raster)
require(rgdal)

# Using readGDAL to read raster (if in multi-band stack format) 
y = readGDAL(system.file("pictures/Rlogo.jpg", package = "rgdal")[1])
  head(y@data)

# Using raster to create stack from individual bands and coerce to SpatialGridDataFrame
y <- stack( raster(system.file("pictures/Rlogo.jpg", package = "rgdal")[1], band=1),
            raster(system.file("pictures/Rlogo.jpg", package = "rgdal")[1], band=2),
            raster(system.file("pictures/Rlogo.jpg", package = "rgdal")[1], band=3))
class(y)

y <- as(y, "SpatialGridDataFrame")
  class(y)

# do something to the data and write raster to disk  
y@data <- y@data * 0.01  
writeGDAL(y, "corrected.tif", drivername="GTiff", type="Float32") 
0

install raster package, then load band by band. Then you can use landsat package. Ex:

library(raster)
setwd(L8 folder)
b1=raster("L8b1.tif")
b2=raster("L8b2.tif")
image=stack(b1,b2)
  • 1
    Reading a raster or stack object does not result in a SpatialGridDataFrame object. You can coerce to a sp or matrix class object but your current solution does not provide this information. – Jeffrey Evans Feb 9 '15 at 20:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.