Datos de Elevación Antioquia

library(rasterVis)
library(raster)
library(rgl)
library(rgdal)
library(elevatr)
library(sf)
library(tidyverse)
library(mapview)

munic <- sf::st_read("C:/Users/hinca/Documents/(A) Hinca/UNAL/Semestre 5/Geomática/05_ANTIOQUIA/ADMINISTRATIVO/MGN_MPIO_POLITICO.shp")

Reading layer `MGN_MPIO_POLITICO' from data source 
  `C:\Users\hinca\Documents\(A) Hinca\UNAL\Semestre 5\Geomática\05_ANTIOQUIA\ADMINISTRATIVO\MGN_MPIO_POLITICO.shp' 
  using driver `ESRI Shapefile'
Simple feature collection with 125 features and 9 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -77.12783 ymin: 5.418558 xmax: -73.88128 ymax: 8.873974
Geodetic CRS:  WGS 84

elevation <- get_elev_raster(munic, z = 10)

 Accessing raster elevation [-------------------------]   0%
 Accessing raster elevation [-------------------------]   1%
 Accessing raster elevation [-------------------------]   2%
 Accessing raster elevation [>------------------------]   2%
 Accessing raster elevation [>------------------------]   3%
 Accessing raster elevation [>------------------------]   4%
 Accessing raster elevation [>------------------------]   5%

Error in curl::curl_fetch_disk(url, x$path, handle = handle) : 
  Operation was aborted by an application callback

NA

if (require(rgdal)) {
 rf <- writeRaster(elevation, filename=file.path("C:/Users/hinca/Documents/(A) Hinca/UNAL/Semestre 5/Geomática/antioquia_dem_z10.tif"), format="GTiff", overwrite=TRUE)
}

elevation <-  raster("C:/Users/hinca/Documents/(A) Hinca/UNAL/Semestre 5/Geomática/antioquia_dem_z10.tif")

elevation

class      : RasterLayer 
dimensions : 6134, 5667, 34761378  (nrow, ncol, ncell)
resolution : 0.0006824976, 0.0006824976  (x, y)
extent     : -77.34375, -73.47604, 4.915657, 9.102097  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : antioquia_dem_z10.tif 
names      : antioquia_dem_z10 
values     : -1934, 4758  (min, max)

munic_sp <- as(munic, "Spatial")

plot(elevation, main="Datos de DEM descargados sin tratamiento [meters]")
plot(munic_sp,col="NA",border="black", add=TRUE)

text(coordinates(munic_sp), labels=as.character(munic$MPIO_CNMBR), 
     col="black", cex=0.20)

elev_crop = crop(elevation, munic_sp)
plot(elev_crop, main="DEM Recortado para Antioquia")
plot(munic_sp, add=TRUE)

text(coordinates(munic_sp), labels=as.character(munic_sp$MPIO_CNMBR), cex=0.2)

elev_crop

class      : RasterLayer 
dimensions : 5063, 4757, 24084691  (nrow, ncol, ncell)
resolution : 0.0006824976, 0.0006824976  (x, y)
extent     : -77.12808, -73.88144, 5.418657, 8.874143  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : memory
names      : antioquia_dem_z10 
values     : -1609, 4758  (min, max)

library(sp)
spatialref <- CRS(SRS_string="EPSG:32618")

pr3 <- projectExtent(elev_crop, crs=spatialref)
res(pr3) <- 180
rep_elev <- projectRaster(elev_crop, pr3)

rep_elev

class      : RasterLayer 
dimensions : 2126, 1998, 4247748  (nrow, ncol, ncell)
resolution : 180, 180  (x, y)
extent     : 264195.4, 623835.4, 598928.5, 981608.5  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs 
source     : memory
names      : antioquia_dem_z10 
values     : -692.154, 4490.655  (min, max)

(rep_munic = spTransform(munic_sp,spatialref))

class       : SpatialPolygonsDataFrame 
features    : 125 
extent      : 265363.8, 623567.1, 598938.4, 981215  (xmin, xmax, ymin, ymax)
crs         : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs 
variables   : 9
names       : DPTO_CCDGO, MPIO_CCDGO, MPIO_CNMBR,                          MPIO_CRSLC,    MPIO_NAREA, MPIO_NANO, DPTO_CNMBR,     Shape_Leng,      Shape_Area 
min values  :         05,      05001,  ABEJORRAL,                                1541,   15.83597275,      2017,  ANTIOQUIA, 0.172962481259, 0.0012928296672 
max values  :         05,      05895,   ZARAGOZA, Ordenanza 9 de Diciembre 15 de 1983, 2959.36315089,      2017,  ANTIOQUIA,  6.61177822771,   0.23894871119 

plot(rep_elev, main="Reprojected digital elevation model")
plot(rep_munic, lwd=0.5, add=TRUE)

text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)

writeRaster(rep_elev, filename="C:/Users/hinca/Documents/(A) Hinca/UNAL/Semestre 5/Geomática/rep_boyaca_dem.tif", datatype='INT4S', overwrite=TRUE)

hist(rep_elev)

promedio <- cellStats(rep_elev, 'mean')
minimo <- cellStats(rep_elev, 'min')
maximo <- cellStats(rep_elev, 'max')
desviacion  <- cellStats(rep_elev, 'sd')

metricas <- c('mean', 'min', 'max', 'std')
valores <- c(promedio, minimo, maximo, desviacion)

(df_estadisticas <- data.frame(metricas, valores))

slope = terrain(rep_elev,opt='slope', unit='degrees')
aspect = terrain(rep_elev,opt='aspect',unit='degrees')
hill = hillShade(slope,aspect,40,315)

plot(rep_elev,main="DEM for Antioquia municipalities [meters]", col=terrain.colors(25,alpha=0.8))
plot(rep_munic, add=TRUE, lwd=0.5)

text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)

plot(slope,main="Slope for Antioquia municipalities [degrees]", col=topo.colors(6,alpha=0.6))
plot(rep_munic, add=TRUE, lwd=0.5)

text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)

plot(aspect,main="Aspect for several municipalities [degrees]", col=rainbow(10,alpha=0.7))
plot(rep_munic, add=TRUE, lwd=0.5)

text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)

plot(hill,
        col=grey(1:100/100),  
        legend=FALSE,         
        main="DEM for Antioquia municipalities [m]",
        axes=FALSE)           

plot(rep_elev, 
        axes=FALSE,
        col=terrain.colors(12, alpha=0.6), add=TRUE)

plot(rep_munic, add=TRUE, lwd=0.5)
text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)

elmat %>%
  sphere_shade(texture = "desert") %>%
  add_water(detect_water(elmat), color = "desert") %>%
  add_shadow(ray_shade(elmat), 0.5) %>%
  add_shadow(ambient_shade(elmat), 0) %>%
  plot_map()

Error in hillshade[, , 1] * scales::rescale(shadowmap, c(max_darken, 1)) : 
  non-conformable arrays

lista <- list("MUTATÁ" , "CHIGORODÓ" , "CAREPA" , "APARTADÓ" , "TURBO" , "NECOCLÍ" , "ARBOLETES" , "SAN JUAN DE URABÁ")
(algunosmunicipios <- munic %>% filter(MPIO_CNMBR %in% unlist(lista)))

Simple feature collection with 8 features and 9 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -77.12783 ymin: 7.132796 xmax: -76.21803 ymax: 8.873974
Geodetic CRS:  WGS 84
  DPTO_CCDGO MPIO_CCDGO        MPIO_CNMBR
1         05      05051         ARBOLETES
2         05      05045          APARTADÓ
3         05      05837             TURBO
4         05      05172         CHIGORODÓ
5         05      05147            CAREPA
6         05      05659 SAN JUAN DE URABÁ
7         05      05490           NECOCLÍ
8         05      05480            MUTATÁ
                                           MPIO_CRSLC
1 Decreto departamental 340 del 30 de Agosto  de 1958
2                 Ordenanza 7 de Noviembre 30 de 1967
3                                                1840
4                Ordenanza 52 del 27 de Abril de 1915
5                 Ordenanza 7 de Diciembre 15 de 1983
6                     Ordenanza 58 de Marzo 3 de 1986
7                Ordenanza 23 de Noviembre 28 de 1977
8                                                1951
  MPIO_NAREA MPIO_NANO DPTO_CNMBR Shape_Leng
1   754.2525      2017  ANTIOQUIA  1.8586821
2   535.1801      2017  ANTIOQUIA  1.2876523
3  2959.3632      2017  ANTIOQUIA  6.6117782
4   721.6826      2017  ANTIOQUIA  1.6103856
5   387.2706      2017  ANTIOQUIA  1.3099406
6   251.9969      2017  ANTIOQUIA  0.8078772
7  1255.5439      2017  ANTIOQUIA  2.4147100
8  1198.2352      2017  ANTIOQUIA  1.6610843
  Shape_Area                       geometry
1 0.06185481 POLYGON ((-76.41355 8.87382...
2 0.04380111 POLYGON ((-76.42487 8.07299...
3 0.23894871 POLYGON ((-76.40525 8.47647...
4 0.05902320 POLYGON ((-76.80956 7.77069...
5 0.03168377 POLYGON ((-76.76941 7.87896...
6 0.02066743 POLYGON ((-76.47999 8.81150...
7 0.10289352 POLYGON ((-76.62422 8.70221...
8 0.08789726 POLYGON ((-76.49619 7.55139...

algunosmunicipios <- as(algunosmunicipios, "Spatial")

elev_cropA = crop(elevation, algunosmunicipios)
elev_cropA = mask(elev_cropA, algunosmunicipios)
plot(elev_cropA, main=" DEM Urabá de Antioquia ")
plot(algunosmunicipios, add=TRUE)

text(coordinates(algunosmunicipios), labels=as.character(algunosmunicipios$MPIO_CNMBR), cex=0.3)

slopeA = terrain(elev_cropA, opt='slope', unit='degrees')
aspectA = terrain(elev_cropA, opt='aspect',unit='degrees')
hillM = hillShade(slopeA, aspectA, 40, 315)

plot(slopeA,main="Pendiente Urabá de Antioquia [grados]", col=topo.colors(6,alpha=0.6))
plot(algunosmunicipios, add=TRUE, lwd=0.5)

text(coordinates(algunosmunicipios), labels=as.character(algunosmunicipios$MPIO_CNMBR), cex=0.3)

plot(aspectA,main="Aspecto Urabá de Antioquia", col=rainbow(10,alpha=0.7))
plot(algunosmunicipios, add=TRUE, lwd=0.5)

text(coordinates(algunosmunicipios), labels=as.character(algunosmunicipios$MPIO_CNMBR), cex=0.3)

lista <- list("RIONEGRO" , "PEÑOL" , "GUARNE" , "SANTUARIO" , "MARINILLA" , "SAN VICENTE" , "GRANADA" , "CONCEPCIÓN" , "BARBOSA")
(algunosmunicipiosm <- munic %>% filter(MPIO_CNMBR %in% unlist(lista)))

Simple feature collection with 9 features and 9 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -75.50126 ymin: 6.018431 xmax: -75.05946 ymax: 6.513189
Geodetic CRS:  WGS 84
  DPTO_CCDGO MPIO_CCDGO  MPIO_CNMBR MPIO_CRSLC
1         05      05079     BARBOSA       1812
2         05      05697   SANTUARIO       1838
3         05      05206  CONCEPCIÓN       1773
4         05      05674 SAN VICENTE       1814
5         05      05615    RIONEGRO       1783
6         05      05541       PEÑOL       1774
7         05      05440   MARINILLA       1787
8         05      05318      GUARNE       1817
9         05      05313     GRANADA       1817
  MPIO_NAREA MPIO_NANO DPTO_CNMBR Shape_Leng
1  205.66617      2017  ANTIOQUIA  0.7700180
2   83.45109      2017  ANTIOQUIA  0.6097303
3  202.05449      2017  ANTIOQUIA  0.7150559
4  230.49117      2017  ANTIOQUIA  0.7710366
5  195.94020      2017  ANTIOQUIA  0.8907217
6  118.31435      2017  ANTIOQUIA  0.6187559
7  114.87120      2017  ANTIOQUIA  0.6129668
8  151.81382      2017  ANTIOQUIA  0.7063743
9  190.22512      2017  ANTIOQUIA  0.7639613
   Shape_Area                       geometry
1 0.016803898 POLYGON ((-75.32148 6.51265...
2 0.006814663 POLYGON ((-75.24933 6.17184...
3 0.016508022 POLYGON ((-75.21393 6.45258...
4 0.018827491 POLYGON ((-75.33483 6.39938...
5 0.015999436 POLYGON ((-75.37699 6.23349...
6 0.009663738 POLYGON ((-75.23488 6.32415...
7 0.009381138 POLYGON ((-75.29652 6.26197...
8 0.012398701 POLYGON ((-75.40888 6.34896...
9 0.015535174 POLYGON ((-75.122 6.22429, ...

algunosmunicipiosm <- as(algunosmunicipiosm, "Spatial")

elev_cropM = crop(elevation, algunosmunicipiosm)
elev_cropM = mask(elev_cropM, algunosmunicipiosm)
plot(elev_cropM, main=" DEM Zona Montañosa de Antioquia ")
plot(algunosmunicipiosm, add=TRUE)

text(coordinates(algunosmunicipiosm), labels=as.character(algunosmunicipiosm$MPIO_CNMBR), cex=0.3)

slopeM = terrain(elev_cropM, opt='slope', unit='degrees')
aspectM = terrain(elev_cropM, opt='aspect',unit='degrees')
hillM = hillShade(slopeM, aspectM, 40, 315)

plot(slopeM,main="Pendiente Zona Montañosa Antioquia [grados]", col=topo.colors(6,alpha=0.6))
plot(algunosmunicipiosm, add=TRUE, lwd=0.5)

text(coordinates(algunosmunicipiosm), labels=as.character(algunosmunicipiosm$MPIO_CNMBR), cex=0.3)

plot(aspectM,main="Aspecto Zona Montañosa de Antioquia", col=rainbow(10,alpha=0.7))
plot(algunosmunicipiosm, add=TRUE, lwd=0.5)

text(coordinates(algunosmunicipiosm), labels=as.character(algunosmunicipiosm$MPIO_CNMBR), cex=0.3)