cuaderno 5 elevacion
Juan Guillermo Narvaez Diaz
17/4/2021
Introduccion
En este cuaderno 5 es ta creado usando RMarkdown para general un modelo de elevacion del departamento de Guainia, los datos seran utilizados para ampliar la informacion que se tiene del departamento, tales como la hidrografia, la visualizacion del terreno. exploraremoslos datos de elevacion raster, del paquete elevatr.
library(rasterVis)## Warning: package 'rasterVis' was built under R version 4.0.5## Loading required package: raster## Warning: package 'raster' was built under R version 4.0.5## Loading required package: sp## Loading required package: terra## Warning: package 'terra' was built under R version 4.0.5## terra version 1.1.4## Loading required package: lattice## Loading required package: latticeExtra## Warning: package 'latticeExtra' was built under R version 4.0.5library(raster)
library(rgl)## Warning: package 'rgl' was built under R version 4.0.5library(rgdal)## Warning: package 'rgdal' was built under R version 4.0.5## rgdal: version: 1.5-23, (SVN revision 1121)
## 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/Juan Narvaez/Documents/R/win-library/4.0/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/Juan Narvaez/Documents/R/win-library/4.0/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 rgdal.
## Overwritten PROJ_LIB was C:/Users/Juan Narvaez/Documents/R/win-library/4.0/rgdal/proj##
## Attaching package: 'rgdal'## The following object is masked from 'package:terra':
##
## projectlibrary(elevatr)## Warning: package 'elevatr' was built under R version 4.0.5library(sf)## Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1library(tidyverse)## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --## v ggplot2 3.3.3 v purrr 0.3.4
## v tibble 3.1.0 v dplyr 1.0.5
## v tidyr 1.1.3 v stringr 1.4.0
## v readr 1.4.0 v forcats 0.5.1## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::collapse() masks terra::collapse()
## x dplyr::desc() masks terra::desc()
## x tidyr::expand() masks terra::expand()
## x tidyr::extract() masks terra::extract(), raster::extract()
## x tidyr::fill() masks terra::fill()
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
## x ggplot2::layer() masks latticeExtra::layer()
## x dplyr::near() masks terra::near()
## x tidyr::pack() masks terra::pack()
## x dplyr::select() masks terra::select(), raster::select()
## x tidyr::separate() masks terra::separate()
## x purrr::transpose() masks terra::transpose()library(mapview)## Warning: package 'mapview' was built under R version 4.0.5las diferentes fuentes de datos de elevacion como Shuttle Radar Topography Mission (SRTM), USGS National Elevation Dataset (NED), Global DEM (GDEM), la forma de utilizarlos es atravez de get_elev_raster () es un data.frame con ubicaciones x (longitud) e y (latitud). para utilizar el marco de datos requiere un CRS (es decir, un archivo de proyecc. ión). El nivel de zoom (z) predeterminado es 9. Pude obtener datos de elevación con zoom = 10, hasta un zoom de z=14 si son zonas muy pequeñas. primer se cargara el shp del geoportal del DANE `
(munici <- st_read("94_GUAINIA/ADMINISTRATIVO/MGN_MPIO_POLITICO.shp"))## Reading layer `MGN_MPIO_POLITICO' from data source `C:\Users\Juan Narvaez\Documents\geomatica\94_GUAINIA\ADMINISTRATIVO\MGN_MPIO_POLITICO.shp' using driver `ESRI Shapefile'
## Simple feature collection with 9 features and 9 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -70.94249 ymin: 1.165633 xmax: -66.84722 ymax: 4.045026
## geographic CRS: WGS 84## Simple feature collection with 9 features and 9 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -70.94249 ymin: 1.165633 xmax: -66.84722 ymax: 4.045026
## geographic CRS: WGS 84
## DPTO_CCDGO MPIO_CCDGO MPIO_CNMBR
## 1 94 94001 INÍRIDA
## 2 94 94343 BARRANCO MINA
## 3 94 94663 MAPIRIPANA
## 4 94 94883 SAN FELIPE
## 5 94 94884 PUERTO COLOMBIA
## 6 94 94885 LA GUADALUPE
## 7 94 94886 CACAHUAL
## 8 94 94887 PANÁ-PANÁ (Campo Alegre)
## 9 94 94888 MORICHAL (Morichal Nuevo)
## MPIO_CRSLC MPIO_NAREA
## 1 Decreto 1593 de Agosto 5 de 1974 15970.268
## 2 Resolución 83 del 1 de Febrero de 1988 que aprueba con modi 9467.804
## 3 ACUERDO COMISARIAL 018 DEL 25 DE AGOSTO DE 1990 4927.971
## 4 Resolución 83 2926.365
## 5 Resolución 83 15700.545
## 6 Resolución 83 1178.868
## 7 Resolución 83 2335.323
## 8 Resolución 83 de 1988 10227.139
## 9 1988 8554.996
## MPIO_NANO DPTO_CNMBR Shape_Leng Shape_Area geometry
## 1 2017 GUAINÍA 9.256685 1.2869895 POLYGON ((-67.67638 3.91228...
## 2 2017 GUAINÍA 6.825914 0.7652557 POLYGON ((-68.91332 3.68215...
## 3 2017 GUAINÍA 5.010789 0.3990354 POLYGON ((-70.10453 3.38436...
## 4 2017 GUAINÍA 4.029719 0.2346116 POLYGON ((-67.34976 2.50451...
## 5 2017 GUAINÍA 9.166487 1.2632740 POLYGON ((-67.5385 3.177568...
## 6 2017 GUAINÍA 1.750810 0.0943385 POLYGON ((-66.96243 1.66858...
## 7 2017 GUAINÍA 2.364174 0.1876089 POLYGON ((-67.49529 3.75681...
## 8 2017 GUAINÍA 7.805305 0.8250362 POLYGON ((-68.77132 2.42182...
## 9 2017 GUAINÍA 6.916544 0.6917881 POLYGON ((-69.46493 2.88960...son los atributos de los municipios.
head(munici)## Simple feature collection with 6 features and 9 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -70.94249 ymin: 1.165633 xmax: -66.84722 ymax: 4.045026
## geographic CRS: WGS 84
## DPTO_CCDGO MPIO_CCDGO MPIO_CNMBR
## 1 94 94001 INÍRIDA
## 2 94 94343 BARRANCO MINA
## 3 94 94663 MAPIRIPANA
## 4 94 94883 SAN FELIPE
## 5 94 94884 PUERTO COLOMBIA
## 6 94 94885 LA GUADALUPE
## MPIO_CRSLC MPIO_NAREA
## 1 Decreto 1593 de Agosto 5 de 1974 15970.268
## 2 Resolución 83 del 1 de Febrero de 1988 que aprueba con modi 9467.804
## 3 ACUERDO COMISARIAL 018 DEL 25 DE AGOSTO DE 1990 4927.971
## 4 Resolución 83 2926.365
## 5 Resolución 83 15700.545
## 6 Resolución 83 1178.868
## MPIO_NANO DPTO_CNMBR Shape_Leng Shape_Area geometry
## 1 2017 GUAINÍA 9.256685 1.2869895 POLYGON ((-67.67638 3.91228...
## 2 2017 GUAINÍA 6.825914 0.7652557 POLYGON ((-68.91332 3.68215...
## 3 2017 GUAINÍA 5.010789 0.3990354 POLYGON ((-70.10453 3.38436...
## 4 2017 GUAINÍA 4.029719 0.2346116 POLYGON ((-67.34976 2.50451...
## 5 2017 GUAINÍA 9.166487 1.2632740 POLYGON ((-67.5385 3.177568...
## 6 2017 GUAINÍA 1.750810 0.0943385 POLYGON ((-66.96243 1.66858...tail(munici)## Simple feature collection with 6 features and 9 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -70.63844 ymin: 1.165633 xmax: -66.84722 ymax: 3.775236
## geographic CRS: WGS 84
## DPTO_CCDGO MPIO_CCDGO MPIO_CNMBR MPIO_CRSLC
## 4 94 94883 SAN FELIPE Resolución 83
## 5 94 94884 PUERTO COLOMBIA Resolución 83
## 6 94 94885 LA GUADALUPE Resolución 83
## 7 94 94886 CACAHUAL Resolución 83
## 8 94 94887 PANÁ-PANÁ (Campo Alegre) Resolución 83 de 1988
## 9 94 94888 MORICHAL (Morichal Nuevo) 1988
## MPIO_NAREA MPIO_NANO DPTO_CNMBR Shape_Leng Shape_Area
## 4 2926.365 2017 GUAINÍA 4.029719 0.2346116
## 5 15700.545 2017 GUAINÍA 9.166487 1.2632740
## 6 1178.868 2017 GUAINÍA 1.750810 0.0943385
## 7 2335.323 2017 GUAINÍA 2.364174 0.1876089
## 8 10227.139 2017 GUAINÍA 7.805305 0.8250362
## 9 8554.996 2017 GUAINÍA 6.916544 0.6917881
## geometry
## 4 POLYGON ((-67.34976 2.50451...
## 5 POLYGON ((-67.5385 3.177568...
## 6 POLYGON ((-66.96243 1.66858...
## 7 POLYGON ((-67.49529 3.75681...
## 8 POLYGON ((-68.77132 2.42182...
## 9 POLYGON ((-69.46493 2.88960...descargaremos los datos de elevacion de Guainia con un zoom =10 ya que el departamento es muy extenso.
elevacion<- get_elev_raster(munici,z=10)## Note: Your request will download approximately 124.2Mb.## Mosaicing & Projecting## Note: Elevation units are in meters.
## Note: The coordinate reference system is:
## GEOGCRS["WGS 84 (with axis order normalized for visualization)",
## DATUM["World Geodetic System 1984",
## ELLIPSOID["WGS 84",6378137,298.257223563,
## LENGTHUNIT["metre",1]]],
## PRIMEM["Greenwich",0,
## ANGLEUNIT["degree",0.0174532925199433]],
## CS[ellipsoidal,2],
## AXIS["geodetic longitude (Lon)",east,
## ORDER[1],
## ANGLEUNIT["degree",0.0174532925199433,
## ID["EPSG",9122]]],
## AXIS["geodetic latitude (Lat)",north,
## ORDER[2],
## ANGLEUNIT["degree",0.0174532925199433,
## ID["EPSG",9122]]]]La función get_elev_raster descarga datos de elevación como un ráster (es decir, como un DEM ráster), el raster al ser temporal, con la siguiente la volvemos algo mas maliable y atemporal.
if (require(rgdal)){
rf <- writeRaster(elevacion, filename=file.path("C:/Users/Juan Narvaez/Documents/geomatica/Guainia_Dem_z10.tif"), format="GTiff", overwrite=TRUE)
}verificamos los datos del raster
elevacion## class : RasterLayer
## dimensions : 5121, 6665, 34131465 (nrow, ncol, ncell)
## resolution : 0.0006857925, 0.0006857925 (x, y)
## extent : -71.01562, -66.44482, 0.7029995, 4.214943 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : C:/Users/Juan Narvaez/AppData/Local/Temp/Rtmp616Mkg/file1d286ac31678.tif
## names : file1d286ac31678
## values : -32768, 32767 (min, max)El funcionamiento de as()cualquiera de las formas depende de la definición de métodos de coerción. Los métodos se definen automáticamente cuando las dos clases están relacionadas por herencia; es decir, cuando una de las clases es una subclase de la otra.
munici_sp <- as(munici, "Spatial")hacemos el mapa con el raster de elevacion del departamento de Guainia.
plot(elevacion, main="DEscarga del DEM [metros]")
plot(munici_sp,col="NA",border="black", add=TRUE)
text(coordinates(munici_sp), labels=as.character(munici$MPIO_CNMBR),
col="black", cex=0.20)
### recorte de datos de elvacion para que coincidan con la zona de estudio El DEM cubre una extension mayor a la zona de estudio, en este caso el departamento de guainia. con este codigo recortamos los datos del raster para el municipio y corregimientos.
elev_crop = crop(elevacion, munici_sp)
plot(elev_crop, main="Modelo digital de elevacion recortado ")
plot(munici_sp, add=TRUE)
text(coordinates(munici_sp), labels=as.character(munici_sp$MPIO_CNMBR), cex=0.2)
visualizamos el nuevo elementos
elev_crop## class : RasterLayer
## dimensions : 4198, 5971, 25066258 (nrow, ncol, ncell)
## resolution : 0.0006857925, 0.0006857925 (x, y)
## extent : -70.94225, -66.84738, 1.165909, 4.044867 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : C:/Users/Juan Narvaez/AppData/Local/Temp/Rtmp616Mkg/raster/r_tmp_2021-05-02_164328_7464_14716.grd
## names : file1d286ac31678
## values : -461, 972 (min, max)se proyectan nuevamente los datos de elevacion
es una buena idea utilizar coordenadas del mapa en lugar de coordenadas geográficas. proyectamos los datos de elevacion.
crs(elev_crop)## CRS arguments: +proj=longlat +datum=WGS84 +no_defsspatialref <- CRS(SRS_string="EPSG:32618")proyectamos los datos de elevación de las coordenadas geográficas WGS84 en el UTM
#creamos un objeto raster, para darnos el control
#projectExtent (no se transfieren valores)
pr3 <- projectExtent(elev_crop, crs=spatialref)
#tamaño de la celda
res(pr3) <- 180
rep_elev <- projectRaster(elev_crop, pr3)rep_elev## class : RasterLayer
## dimensions : 1792, 2552, 4573184 (nrow, ncol, ncell)
## resolution : 180, 180 (x, y)
## extent : 950785.1, 1410145, 129107, 451667 (xmin, xmax, ymin, ymax)
## crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
## source : memory
## names : file1d286ac31678
## values : -268.236, 960.6074 (min, max)olvemos a representar el SpatialPolygonsDataFrame en nuestros departamento
(rep_munic = spTransform(munici_sp,spatialref))## class : SpatialPolygonsDataFrame
## features : 9
## extent : 951413.4, 1410054, 130083.7, 450762.9 (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 : 94, 94001, BARRANCO MINA, 1988, 1178.86832903, 2017, GUAINÍA, 1.75081007547, 0.0943384971244
## max values : 94, 94888, SAN FELIPE, Resolución 83 del 1 de Febrero de 1988 que aprueba con modi, 15970.2680471, 2017, GUAINÍA, 9.25668481051, 1.28698950156plot(rep_elev, main="reproyeccion del modelo digital de elevacion ")
plot(rep_munic, lwd=0.5, add=TRUE)
text(coordinates(rep_munic), labels=as.character(rep_munic$MPIO_CNMBR), cex=0.2)
guardamos el nuevo dem, como una reproyeccion
writeRaster(rep_elev, filename="C:/Users/Juan Narvaez/Documents/geomatica/rep_Guainia_Dem_z10.tif", datatype='INT4S', overwrite=TRUE)sessionInfo()## R version 4.0.4 (2021-02-15)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 19042)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=Spanish_Colombia.1252 LC_CTYPE=Spanish_Colombia.1252
## [3] LC_MONETARY=Spanish_Colombia.1252 LC_NUMERIC=C
## [5] LC_TIME=Spanish_Colombia.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] mapview_2.9.0 forcats_0.5.1 stringr_1.4.0
## [4] dplyr_1.0.5 purrr_0.3.4 readr_1.4.0
## [7] tidyr_1.1.3 tibble_3.1.0 ggplot2_3.3.3
## [10] tidyverse_1.3.0 sf_0.9-7 elevatr_0.3.4
## [13] rgdal_1.5-23 rgl_0.105.22 rasterVis_0.50.1
## [16] latticeExtra_0.6-29 lattice_0.20-41 terra_1.1-4
## [19] raster_3.4-5 sp_1.4-5
##
## loaded via a namespace (and not attached):
## [1] fs_1.5.0 satellite_1.0.2 lubridate_1.7.10
## [4] progress_1.2.2 webshot_0.5.2 RColorBrewer_1.1-2
## [7] httr_1.4.2 tools_4.0.4 backports_1.2.1
## [10] bslib_0.2.4 utf8_1.1.4 R6_2.5.0
## [13] KernSmooth_2.23-18 DBI_1.1.1 colorspace_2.0-0
## [16] manipulateWidget_0.10.1 withr_2.4.1 prettyunits_1.1.1
## [19] tidyselect_1.1.0 leaflet_2.0.4.1 curl_4.3
## [22] compiler_4.0.4 leafem_0.1.3 cli_2.3.1
## [25] rvest_1.0.0 xml2_1.3.2 sass_0.3.1
## [28] scales_1.1.1 classInt_0.4-3 hexbin_1.28.2
## [31] proxy_0.4-25 digest_0.6.27 rmarkdown_2.7
## [34] base64enc_0.1-3 jpeg_0.1-8.1 pkgconfig_2.0.3
## [37] htmltools_0.5.1.1 highr_0.8 dbplyr_2.1.0
## [40] fastmap_1.1.0 htmlwidgets_1.5.3 rlang_0.4.10
## [43] readxl_1.3.1 rstudioapi_0.13 shiny_1.6.0
## [46] jquerylib_0.1.3 generics_0.1.0 zoo_1.8-9
## [49] jsonlite_1.7.2 crosstalk_1.1.1 magrittr_2.0.1
## [52] Rcpp_1.0.6 munsell_0.5.0 fansi_0.4.2
## [55] lifecycle_1.0.0 stringi_1.5.3 yaml_2.2.1
## [58] grid_4.0.4 parallel_4.0.4 promises_1.2.0.1
## [61] crayon_1.4.1 miniUI_0.1.1.1 haven_2.3.1
## [64] hms_1.0.0 knitr_1.31 pillar_1.5.1
## [67] stats4_4.0.4 codetools_0.2-18 reprex_1.0.0
## [70] glue_1.4.2 evaluate_0.14 modelr_0.1.8
## [73] png_0.1-7 vctrs_0.3.6 httpuv_1.5.5
## [76] cellranger_1.1.0 gtable_0.3.0 assertthat_0.2.1
## [79] xfun_0.22 mime_0.10 xtable_1.8-4
## [82] broom_0.7.5 e1071_1.7-5 later_1.1.0.1
## [85] class_7.3-18 viridisLite_0.3.0 units_0.7-0
## [88] ellipsis_0.3.1