library(raster)
## Loading required package: sp
#Blue
b2 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B2.tif')
#Green
b3 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B3.tif')
#Red
b4 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B4.tif')
#Near Infrared (NIR)
b5 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B5.tif')
# BANDA 6
b6 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B6.tif')
# BANDA 7
b7 <- raster('C:/PRACTOCA #3/rs/LC08_044034_20170614_B7.tif')
b2
## class : RasterLayer
## dimensions : 1245, 1497, 1863765 (nrow, ncol, ncell)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## source : C:/PRACTOCA #3/rs/LC08_044034_20170614_B2.tif
## names : LC08_044034_20170614_B2
## values : 0.0748399, 0.7177562 (min, max)
## class : RasterLayer
## dimensions : 1245, 1497, 1863765 (nrow, ncol, ncell)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## source : C:/PRACTOCA #3/rs/LC08_044034_20170614_B2.tif
## names : LC08_044034_20170614_B2
## values : 0.0748399, 0.7177562 (min, max)
# coordinate reference system (crs)
crs(b2)
## CRS arguments:
## +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## CRS arguments:
## +proj=utm +zone=10 +datum=WGS84 +units=m
## +no_defs
# Number of cells, row, columns
ncell(b2)
## [1] 1863765
## [1] 1863765
dim(b2)
## [1] 1245 1497 1
## [1] 1245 1497 1
# spatial resolution
res(b2)
## [1] 30 30
## [1] 30 30
# Number of bands
nlayers(b2)
## [1] 1
## [1] 1
# Do the bands have the same extent, number of rows and columns, projection, resolution, and origin
compareRaster(b2,b3)
## [1] TRUE
## [1] TRUE
s <- stack(b5, b4, b3)
# Check the properties of the RasterStack
s
## class : RasterStack
## dimensions : 1245, 1497, 1863765, 3 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## names : LC08_044034_20170614_B5, LC08_044034_20170614_B4, LC08_044034_20170614_B3
## min values : 0.0008457669, 0.0208406653, 0.0425921641
## max values : 1.0124315, 0.7861769, 0.6924697
## class : RasterStack
## dimensions : 1245, 1497, 1863765, 3 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## names : LC08_044034_20170614_B5, LC08_044034_20170614_B4, LC08_044034_20170614_B3
## min values : 0.0008457669, 0.0208406653, 0.0425921641
## max values : 1.0124315, 0.7861769, 0.6924697
par(mfrow = c(2,2))
plot(b2, main = "Blue", col = gray(0:100/100))
plot(b3, main = "Green", col = gray(0:100/100))
plot(b4, main = "Red", col = gray(0:100/100))
plot(b5, main = "NIR", col = gray(0:100/100))
landsatRGB <- stack(b4, b3, b2)
plotRGB(landsatRGB,axes = TRUE, stretch = "lin", main = "Landsat True Color Composite")
par(mfrow = c(1,2))
plotRGB(landsatRGB, axes=TRUE, stretch="lin", main="Landsat True Color Composite")
landsatFCC <- stack(b5, b4, b3)
plotRGB(landsatFCC, axes=TRUE, stretch="lin", main="Landsat False Color Composite")
par(mfrow = c(1,2))
plotRGB(landsatRGB, axes=TRUE, stretch="lin", main="Landsat True Color Composite")
landsatURB <- stack(b7, b6, b4)
plotRGB(landsatURB, axes=TRUE, stretch="lin", main="Landsat False Color Composite")
# El area urbana se muestra de un color en un color magenta mientras que la vegetacion en distintos tonos de verde y amarillo
par(mfrow = c(1,2))
plotRGB(landsatRGB, axes=TRUE, stretch="lin", main="Landsat True Color Composite")
landsatAGRI <- stack(b6, b5, b2)
plotRGB(landsatAGRI, axes=TRUE, stretch="lin", main="Landsat False Color Composite")
# Con la combinacion de las bandas b6, b5 y b2, se puede detectar las zonas de uso agricola. donde las areas destinadas para esta actividad se representan de un verde brillante
par(mfrow = c(1,2))
plotRGB(landsatRGB, axes=TRUE, stretch="lin", main="Landsat True Color Composite")
landsatVVG <- stack(b5, b6, b2)
plotRGB(landsatVVG, axes=TRUE, stretch="lin", main="Landsat False Color Composite")
# La vegetacion en estado saludable se puede observar de color naranja, verificando aquellas areas con vegetacion primitiva
par(mfrow = c(1,2))
plotRGB(landsatRGB, axes=TRUE, stretch="lin", main="Landsat True Color Composite")
landsatATM <- stack(b7, b6, b5)
plotRGB(landsatATM, axes=TRUE, stretch="lin", main="Landsat False Color Composite")
filenames <- paste0('C:/PRACTOCA #3/rs/LC08_044034_20170614_B', 1:11, ".tif")
filenames
## [1] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B1.tif"
## [2] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B2.tif"
## [3] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B3.tif"
## [4] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B4.tif"
## [5] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B5.tif"
## [6] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B6.tif"
## [7] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B7.tif"
## [8] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B8.tif"
## [9] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B9.tif"
## [10] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B10.tif"
## [11] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B11.tif"
## [1] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B1.tif"
## [2] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B2.tif"
## [3] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B3.tif"
## [4] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B4.tif"
## [5] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B5.tif"
## [6] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B6.tif"
## [7] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B7.tif"
## [8] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B8.tif"
## [9] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B9.tif"
## [10] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B10.tif"
## [11] "C:/PRACTOCA #3/rs/LC08_044034_20170614_B11.tif"
landsat <- stack(filenames)
landsat
## class : RasterStack
## dimensions : 1245, 1497, 1863765, 11 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## names : LC08_044034_20170614_B1, LC08_044034_20170614_B2, LC08_044034_20170614_B3, LC08_044034_20170614_B4, LC08_044034_20170614_B5, LC08_044034_20170614_B6, LC08_044034_20170614_B7, LC08_044034_20170614_B8, LC08_044034_20170614_B9, LC08_044034_20170614_B10, LC08_044034_20170614_B11
## min values : 9.641791e-02, 7.483990e-02, 4.259216e-02, 2.084067e-02, 8.457669e-04, -7.872183e-03, -5.052945e-03, 3.931751e-02, -4.337332e-04, 2.897978e+02, 2.885000e+02
## max values : 0.73462820, 0.71775615, 0.69246972, 0.78617686, 1.01243150, 1.04320455, 1.11793602, 0.82673049, 0.03547901, 322.43139648, 317.99530029
## class : RasterStack
## dimensions : 1245, 1497, 1863765, 11 (nrow, ncol, ncell, nlayers)
## resolution : 30, 30 (x, y)
## extent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs
## names : LC08_044034_20170614_B1, LC08_044034_20170614_B2, LC08_044034_20170614_B3, LC08_044034_20170614_B4, LC08_044034_20170614_B5, LC08_044034_20170614_B6, LC08_044034_20170614_B7, LC08_044034_20170614_B8, LC08_044034_20170614_B9, LC08_044034_20170614_B10, LC08_044034_20170614_B11
## min values : 9.641791e-02, 7.483990e-02, 4.259216e-02, 2.084067e-02, 8.457669e-04, -7.872183e-03, -5.052945e-03, 3.931751e-02, -4.337332e-04, 2.897978e+02, 2.885000e+02
## max values : 0.73462820, 0.71775615, 0.69246972, 0.78617686, 1.01243150, 1.04320455, 1.11793602, 0.82673049, 0.03547901, 322.43139648, 317.99530029
landsat <- subset(landsat, 1:7)
names(landsat)
## [1] "LC08_044034_20170614_B1" "LC08_044034_20170614_B2"
## [3] "LC08_044034_20170614_B3" "LC08_044034_20170614_B4"
## [5] "LC08_044034_20170614_B5" "LC08_044034_20170614_B6"
## [7] "LC08_044034_20170614_B7"
## [1] "LC08_044034_20170614_B1"
## [2] "LC08_044034_20170614_B2"
## [3] "LC08_044034_20170614_B3"
## [4] "LC08_044034_20170614_B4"
## [5] "LC08_044034_20170614_B5"
## [6] "LC08_044034_20170614_B6"
## [7] "LC08_044034_20170614_B7"
names(landsat) <- c('ultra-blue', 'blue', 'green', 'red', 'NIR','SWIR1', 'SWIR2')
names(landsat)
## [1] "ultra.blue" "blue" "green" "red" "NIR"
## [6] "SWIR1" "SWIR2"
## [1] "ultra.blue" "blue"
## [3] "green" "red"
## [5] "NIR" "SWIR1"
## [7] "SWIR2"
# Using extent
extent(landsat)
## class : Extent
## xmin : 594090
## xmax : 639000
## ymin : 4190190
## ymax : 4227540
## class : Extent
## xmin : 594090
## xmax : 639000
## ymin : 4190190
## ymax : 4227540
e <- extent(624387, 635752, 4200047, 4210939)
# crop landsat by the extent
landsatcrop <- crop(landsat, e)
x <- writeRaster(landsatcrop, filename = "cropped-landsat.tif", overwrite=TRUE)
pairs(landsatcrop[[1:2]], main = "Ultra-blue versus Blue")
pairs(landsatcrop[[4:5]], main = "Red versus NIR")
pairs(landsatcrop[[2:3]], main = "Blue versus Green")
pairs(landsatcrop[[6:7]], main = "SWIR1 versus SWIR2")
