(This article was first published on geo-affine » R, and kindly contributed to R-bloggers)
In my last post I was explaining the usage of QGis to do a layerstack of a Landsat-scene. Due to the fact that further research and trying out resulted in frustration I decided to stick with a software I know well: R.
So download the needed layers here and open up your flavoured version of R (in my case RStudio).
What do we need is a package called raster
install.packages("raster")
library("raster")
Now define the working directory where your ETM-bands are stored which is “ETM” in my case:
setwd("~/wd_r/ETM")
Let’s do the last step and create the stack using one line and store this raster object using a second line:
A=stack(c("p134r027_7dt20020722.SR.b01.tif","p134r027_7dt20020722.SR.b02.tif","p134r027_7dt20020722.SR.b03.tif","p134r027_7dt20020722.SR.b04.tif","p134r027_7dt20020722.SR.b05.tif", "p134r027_7dt20020722.SR.b07.tif"))
writeRaster(A, filename="merge.tif")
plotRGB(A, 3, 2, 1, stretch="hist")
That’s it. Easy, wasn’t it?
To leave a comment for the author, please follow the link and comment on his blog: geo-affine » R.
R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series,ecdf, trading) and more...


Zero Inflated Models and Generalized Linear Mixed Models with R.
Zuur, Saveliev, Ieno (2012).