title: "EXPLORACION"

LIBRERIA

library(terra)
## terra version 1.0.10
library(raster)
## Loading required package: sp
library(sp)

ENTORNO DE TRABAJO

ruta <- getwd()
carpeta <- 'Datos'
file <- paste0(ruta,'/',carpeta,'/','LC08_044034_20170614_B', 1:11, ".tif")
landsat <- rast(file)
landsat <- subset(landsat,1:7)#eliminar banda del 8 al 11
names(landsat) <- c('ultra-blue','blue','green','red','NIR','SWIR1','SWIR2')

UNA SOLA BANDA

par(mfrow = c(2,2))
plot(landsat$blue, main = "Blue", col = rgb((0:100)/100, red = 0, green = 0))
plot(landsat$green, main = "Green", col = rgb((0:100 / 100), red = 0, blue = 0))
plot(landsat$red, main = "Red", col = rgb((0:100 / 100),green = 0, blue = 0))
plot(landsat$NIR, main = "NIR", col = rgb((0:100 / 100), red = 1, green=0,blue=0))

Tres Bandas RGB color natural

landsatRGB <- c(landsat$red, landsat$green, landsat$blue)
# plot 3-band image
plotRGB(landsatRGB, axes = TRUE, stretch = "lin", main = "Landsat Composicion de color natural")

Falso Color

landsatNirRG <- c(landsat$NIR , landsat$red, landsat$green)
# plot 3-band image
plotRGB(landsatNirRG, axes = TRUE, stretch = "lin", main = "Landsat Composicion de Falso Color")

Recorte

xmx <- 625752
xmn <- 614387
ymx <- 4210939
ymn <- 4200047
ext(landsat)
## SpatExtent : 594090, 639000, 4190190, 4227540 (xmin, xmax, ymin, ymax)
e <- ext(xmn, xmx, ymn, ymx)
# crop landsat by the extent
landsatcrop <- crop(landsat, e)
par(mfrow = c(1,2))
landsatcropRGB <- c(landsatcrop$red, landsatcrop$green, landsatcrop$blue)
landsatcropNirRG <- c(landsatcrop$NIR , landsatcrop$red, landsatcrop$green)

plotRGB(landsatcropRGB, axes = TRUE, stretch = "lin", main = "Landsat Color natural")
plotRGB(landsatcropNirRG, axes = TRUE, stretch = "lin", main = "Landsat Falso Color")

Guardar

writeRaster(landsatcropRGB, filename="Recorte.tif", overwrite=TRUE)

Comparar Bandas

pairs(landsatcrop[[1:2]], main = "Ultra-blue versus Blue")

pairs(landsatcrop[[4:5]], main = "Red versus NIR")

Extraer Valores de Pixeles

Lon <- c(xmx, xmx, xmn, xmn)
Lat <- c(ymx, ymn,  ymn, ymx)
LonLat <- cbind(Lon, Lat)
df <- extract(landsat, LonLat)
head(df)

PUNTOS DE PERFIL ESPECTRAL

######puntos#######
agua <- extract(landsat, cbind(622749, 4210129))
edificio <- extract(landsat, cbind(615411,4205858))
cultivo <- extract(landsat, cbind(624173, 4200370))
descubierto <-extract(landsat, cbind(624408, 4207517))
######matriz#######
pt <- rbind(agua,edificio,cultivo,descubierto)
pt$head <- c('agua','edificio','cultivo','descubierto')
rownames(pt) <- pt[,8]
pt <- pt[,-8]
pt

PERFIL ESPECTRAL

Colores <- c('Blue', 'gray', 'green', 'darkred')
pt <- as.matrix(pt)
plot(0, ylim=c(0,0.6), xlim = c(1,7), type='n', xlab="Bandas", ylab = "Reflectancia")
for (i in 1:nrow(pt)){
  lines(pt[i,], type = "l", lwd = 3, lty = 1, col = Colores[i])
}
# Title
title(main="Perfil Espectral", font.main = 2)
# Legend
legend("topleft", rownames(pt),
       cex=0.8, col=Colores, lty = 1, lwd =3, bty = "n")