# install.packages("raster")
# install.packages("sp")
# install.packages("rgdal")
# install.packages("dplyr")
# install.packages("sf")
library(raster)
## Loading required package: sp
library(sp)
library(rgdal)
## Please note that rgdal will be retired by the end of 2023,
## plan transition to sf/stars/terra functions using GDAL and PROJ
## at your earliest convenience.
##
## rgdal: version: 1.5-27, (SVN revision 1148)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.2.1, released 2020/12/29
## Path to GDAL shared files: C:/Users/USUARIO/Documents/R/win-library/4.1/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 7.2.1, January 1st, 2021, [PJ_VERSION: 721]
## Path to PROJ shared files: C:/Users/USUARIO/Documents/R/win-library/4.1/rgdal/proj
## PROJ CDN enabled: FALSE
## Linking to sp version:1.4-5
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal.
## Overwritten PROJ_LIB was C:/Users/USUARIO/Documents/R/win-library/4.1/rgdal/proj
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:raster':
##
## intersect, select, union
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(sf)
## Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1
setwd("C:/Users/USUARIO/Documents/Geomatica/CLASE 10")
b4 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b04.tif")
b5 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b05.tif")
b6 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b06.tif")
b7 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b07.tif")
b8 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b08.tif")
b9 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b09.tif")
b10 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b10.tif")
b11 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b11.tif")
b12 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b12.tif")
b13 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b13.tif")
b14 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b14.tif")
b15 <- raster("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/facade_b15.tif")
# Informacion de las imagenes
b4
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b04.tif
## names : facade_b04
## values : 0, 66 (min, max)
b5
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b05.tif
## names : facade_b05
## values : 0, 66 (min, max)
b6
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b06.tif
## names : facade_b06
## values : 0, 66 (min, max)
b7
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b07.tif
## names : facade_b07
## values : 0, 66 (min, max)
b8
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b08.tif
## names : facade_b08
## values : 0, 76 (min, max)
b9
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b09.tif
## names : facade_b09
## values : 0, 78 (min, max)
b10
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b10.tif
## names : facade_b10
## values : 1, 79 (min, max)
b11
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b11.tif
## names : facade_b11
## values : 1, 63 (min, max)
b12
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b12.tif
## names : facade_b12
## values : 0, 43 (min, max)
b13
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b13.tif
## names : facade_b13
## values : 0, 66 (min, max)
b14
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b14.tif
## names : facade_b14
## values : 0, 66 (min, max)
b15
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b15.tif
## names : facade_b15
## values : 0, 76 (min, max)
par(mfrow = c(3,4))
plot(b4, main = "facade_b04", col = gray(0:100 / 100))
plot(b5, main = "facade_b05", col = gray(0:100 / 100))
plot(b6, main = "facade_b06", col = gray(0:100 / 100))
plot(b7, main = "facade_b05", col = gray(0:100 / 100))
plot(b8, main = "facade_b08", col = gray(0:100 / 100))
plot(b9, main = "facade_b09", col = gray(0:100 / 100))
plot(b10, main = "facade_b10", col = gray(0:100 / 100))
plot(b11, main = "facade_b11", col = gray(0:100 / 100))
plot(b12, main = "facade_b12", col = gray(0:100 / 100))
plot(b13, main = "facade_b13", col = gray(0:100 / 100))
plot(b14, main = "facade_b14", col = gray(0:100 / 100))
plot(b15, main = "facade_b15", col = gray(0:100 / 100))
landsat <- stack(b4,b5,b6,b7,b8,b9,b10,b11,b12,b13,b14,b15)
landsat
## class : RasterStack
## dimensions : 800, 1000, 8e+05, 12 (nrow, ncol, ncell, nlayers)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## names : facade_b04, facade_b05, facade_b06, facade_b07, facade_b08, facade_b09, facade_b10, facade_b11, facade_b12, facade_b13, facade_b14, facade_b15
## min values : 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0
## max values : 66, 66, 66, 66, 76, 78, 79, 63, 43, 66, 66, 76
puntos_de_control <- read.csv("C:/Users/USUARIO/Documents/Geomatica/CLASE 10/training.csv")
puntos_de_control <- puntos_de_control %>%
st_as_sf(coords = c("x", "y"))
plot(puntos_de_control$geometry)
df <- extract(landsat, puntos_de_control); df
## facade_b04 facade_b05 facade_b06 facade_b07 facade_b08 facade_b09
## [1,] 0 0 0 0 0 0
## [2,] 0 0 0 0 0 0
## [3,] 0 0 0 0 0 0
## [4,] 0 0 0 0 0 0
## [5,] 0 0 0 0 0 0
## [6,] 0 0 0 0 0 0
## [7,] 0 0 0 0 0 0
## [8,] 0 0 0 0 0 0
## [9,] 0 0 0 0 0 0
## [10,] 0 0 0 0 0 1
## [11,] 0 0 0 0 0 1
## [12,] 0 0 0 0 0 1
## [13,] 0 0 0 0 0 1
## [14,] 0 0 0 0 0 1
## [15,] 0 0 0 0 0 1
## [16,] 0 0 0 0 0 1
## [17,] 0 0 0 0 0 1
## [18,] 0 0 0 0 0 2
## [19,] 0 0 0 0 0 2
## [20,] 0 0 0 0 0 2
## [21,] 0 0 0 0 0 2
## [22,] 0 0 0 0 0 2
## [23,] 0 0 0 0 0 2
## [24,] 0 0 0 0 0 2
## [25,] 0 0 0 0 0 2
## [26,] 65 65 65 65 75 77
## [27,] 62 62 62 62 72 74
## [28,] 65 65 65 65 75 77
## [29,] 65 65 65 65 75 77
## [30,] 65 65 65 65 75 77
## [31,] 65 65 65 65 75 77
## [32,] 21 21 21 21 31 33
## [33,] 65 65 65 65 75 77
## [34,] 22 22 22 22 32 34
## [35,] 65 65 65 65 75 77
## [36,] 65 65 65 65 75 77
## [37,] 22 22 22 22 32 34
## [38,] 22 22 22 22 32 34
## [39,] 65 65 65 65 75 77
## [40,] 65 65 65 65 75 77
## [41,] 23 23 23 23 33 35
## [42,] 24 24 24 24 34 36
## [43,] 65 65 65 65 75 77
## [44,] 25 25 25 25 35 37
## [45,] 65 65 65 65 75 77
## [46,] 65 65 65 65 75 77
## [47,] 65 65 65 65 75 77
## [48,] 65 65 65 65 75 77
## [49,] 65 65 65 65 75 77
## [50,] 62 62 62 62 72 74
## [51,] 65 65 65 65 75 77
## [52,] 65 65 65 65 75 77
## [53,] 62 62 62 62 72 74
## [54,] 65 65 65 65 75 77
## [55,] 24 24 24 24 34 36
## [56,] 65 65 65 65 75 77
## [57,] 65 65 65 65 75 77
## [58,] 65 65 65 65 75 77
## [59,] 62 62 62 62 72 74
## [60,] 62 62 62 62 72 74
## [61,] 62 62 62 62 72 74
## [62,] 65 65 65 65 75 77
## [63,] 62 62 62 62 72 74
## [64,] 63 63 63 63 73 75
## [65,] 63 63 63 63 73 75
## [66,] 63 63 63 63 73 75
## [67,] 65 65 65 65 75 77
## [68,] 64 64 64 64 74 76
## [69,] 65 65 65 65 75 77
## [70,] 65 65 65 65 75 77
## [71,] 64 64 64 64 74 76
## [72,] 65 65 65 65 75 77
## [73,] 65 65 65 65 75 77
## [74,] 64 64 64 64 74 76
## [75,] 64 64 64 64 74 76
## [76,] 64 64 64 64 74 76
## [77,] 65 65 65 65 75 77
## [78,] 64 64 64 64 74 76
## [79,] 65 65 65 65 75 77
## [80,] 65 65 65 65 75 77
## [81,] 65 65 65 65 75 77
## [82,] 65 65 65 65 75 77
## [83,] 65 65 65 65 75 77
## [84,] 65 65 65 65 75 77
## [85,] 65 65 65 65 75 77
## [86,] 65 65 65 65 75 77
## [87,] 65 65 65 65 75 77
## [88,] 65 65 65 65 75 77
## [89,] 65 65 65 65 75 77
## [90,] 0 0 0 0 0 2
## [91,] 65 65 65 65 75 77
## [92,] 65 65 65 65 75 77
## [93,] 65 65 65 65 75 77
## [94,] 65 65 65 65 75 77
## [95,] 65 65 65 65 75 77
## [96,] 65 65 65 65 75 77
## [97,] 65 65 65 65 75 77
## [98,] 66 66 66 66 76 78
## [99,] 65 65 65 65 75 77
## [100,] 66 66 66 66 76 78
## facade_b10 facade_b11 facade_b12 facade_b13 facade_b14 facade_b15
## [1,] 1 2 0 0 0 0
## [2,] 1 1 0 0 0 0
## [3,] 1 1 0 0 0 0
## [4,] 1 1 0 0 0 0
## [5,] 1 1 0 0 0 0
## [6,] 1 2 0 0 0 0
## [7,] 1 2 0 0 0 0
## [8,] 1 2 0 0 0 0
## [9,] 1 2 0 0 0 0
## [10,] 2 2 0 0 0 0
## [11,] 2 2 0 0 0 0
## [12,] 2 1 0 0 0 0
## [13,] 2 1 0 0 0 0
## [14,] 2 1 0 0 0 0
## [15,] 2 2 0 0 0 0
## [16,] 2 2 0 0 0 0
## [17,] 2 2 0 0 0 0
## [18,] 3 2 0 0 0 0
## [19,] 3 1 0 0 0 0
## [20,] 3 2 0 0 0 0
## [21,] 3 2 0 0 0 0
## [22,] 3 2 0 0 0 0
## [23,] 3 2 0 0 0 0
## [24,] 3 2 0 0 0 0
## [25,] 3 2 0 0 0 0
## [26,] 78 62 42 65 65 75
## [27,] 75 59 39 62 62 72
## [28,] 78 62 42 65 65 75
## [29,] 78 62 42 65 65 75
## [30,] 78 62 42 65 65 75
## [31,] 78 62 42 65 65 75
## [32,] 34 18 0 21 21 31
## [33,] 78 62 42 65 65 75
## [34,] 35 19 0 22 22 32
## [35,] 78 62 42 65 65 75
## [36,] 78 62 42 65 65 75
## [37,] 35 19 0 22 22 32
## [38,] 35 19 0 22 22 32
## [39,] 78 62 42 65 65 75
## [40,] 78 62 42 65 65 75
## [41,] 36 20 0 23 23 33
## [42,] 37 21 1 24 24 34
## [43,] 78 62 42 65 65 75
## [44,] 38 22 2 25 25 35
## [45,] 78 62 42 65 65 75
## [46,] 78 62 42 65 65 75
## [47,] 78 62 42 65 65 75
## [48,] 78 62 42 65 65 75
## [49,] 78 62 42 65 65 75
## [50,] 75 59 39 62 62 72
## [51,] 78 62 42 65 65 75
## [52,] 78 62 42 65 65 75
## [53,] 75 59 39 62 62 72
## [54,] 78 62 42 65 65 75
## [55,] 37 21 1 24 24 34
## [56,] 78 62 42 65 65 75
## [57,] 78 62 42 65 65 75
## [58,] 78 62 42 65 65 75
## [59,] 75 59 39 62 62 72
## [60,] 75 59 39 62 62 72
## [61,] 75 59 39 62 62 72
## [62,] 78 62 42 65 65 75
## [63,] 75 59 39 62 62 72
## [64,] 76 60 40 63 63 73
## [65,] 76 60 40 63 63 73
## [66,] 76 60 40 63 63 73
## [67,] 78 62 42 65 65 75
## [68,] 77 61 41 64 64 74
## [69,] 78 62 42 65 65 75
## [70,] 78 62 42 65 65 75
## [71,] 77 61 41 64 64 74
## [72,] 78 62 42 65 65 75
## [73,] 78 62 42 65 65 75
## [74,] 77 61 41 64 64 74
## [75,] 77 61 41 64 64 74
## [76,] 77 61 41 64 64 74
## [77,] 78 62 42 65 65 75
## [78,] 77 61 41 64 64 74
## [79,] 78 62 42 65 65 75
## [80,] 78 62 42 65 65 75
## [81,] 78 62 42 65 65 75
## [82,] 78 62 42 65 65 75
## [83,] 78 62 42 65 65 75
## [84,] 78 62 42 65 65 75
## [85,] 78 62 42 65 65 75
## [86,] 78 62 42 65 65 75
## [87,] 78 62 42 65 65 75
## [88,] 78 62 42 65 65 75
## [89,] 78 62 42 65 65 75
## [90,] 3 2 0 0 0 0
## [91,] 78 62 42 65 65 75
## [92,] 78 62 42 65 65 75
## [93,] 78 62 42 65 65 75
## [94,] 78 62 42 65 65 75
## [95,] 78 62 42 65 65 75
## [96,] 78 62 42 65 65 75
## [97,] 78 62 42 65 65 75
## [98,] 79 63 43 66 66 76
## [99,] 78 62 42 65 65 75
## [100,] 79 63 43 66 66 76
ms <- aggregate(df, list(puntos_de_control$class), mean)
ms
## Group.1 facade_b04 facade_b05 facade_b06 facade_b07 facade_b08 facade_b09
## 1 fissure 45.266667 45.266667 45.266667 45.266667 55.26667 57.26667
## 2 gap 8.862069 8.862069 8.862069 8.862069 10.24138 11.34483
## 3 wall 62.607143 62.607143 62.607143 62.607143 72.42857 74.42857
## facade_b10 facade_b11 facade_b12 facade_b13 facade_b14 facade_b15
## 1 58.26667 42.266667 22.600000 45.266667 45.266667 55.26667
## 2 12.34483 9.896552 5.689655 8.862069 8.862069 10.24138
## 3 75.42857 59.696429 40.017857 62.607143 62.607143 72.42857
# install.packages("writexl")
library(writexl)
# write_xlsx(ms, "C:/Users/USUARIO/Documents/Geomatica/CLASE 10/Ultimo excel.xlsx")
#Preguntas
# ¿cual es la resolucion espacial de las bandas?, la resolucion seria: 0.1, 0.1 (x, y), como se muestra en los siguientes ejemplos
b4
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b04.tif
## names : facade_b04
## values : 0, 66 (min, max)
b5
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b05.tif
## names : facade_b05
## values : 0, 66 (min, max)
b6
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b06.tif
## names : facade_b06
## values : 0, 66 (min, max)
b7
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b07.tif
## names : facade_b07
## values : 0, 66 (min, max)
b8
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b08.tif
## names : facade_b08
## values : 0, 76 (min, max)
b9
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b09.tif
## names : facade_b09
## values : 0, 78 (min, max)
b15
## class : RasterLayer
## dimensions : 800, 1000, 8e+05 (nrow, ncol, ncell)
## resolution : 0.1, 0.1 (x, y)
## extent : 40000, 40100, 2e+05, 200080 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : facade_b15.tif
## names : facade_b15
## values : 0, 76 (min, max)
# ¿cual banda da mejor respuesta a cada clase? segun la grafica la banda 10 tiene mejor reflectancia que las demas
# Conclusiones
# 1. obtener informacion sobre gap es dificil debido a la reflectancia de las 12 bandas.
# 2. algunas bandas tienen reflectancia cero