Ana María Montaño Hernández

29 de marzo de 2020

1. Introducción:

En este cuaderno de R se va a ilustrar el cálculo de pendiente y aspecto. En matemáticas, la pendiente o gradiente de un plano es un número que describe tanto la dirección, como la inclinación del plano. Por otro lado, el aspecto se refiere a la orientación de la pendiente del plano en relación con los puntos cardinales. En resumen, el aspecto es la dirección de la brújula que enfrenta una pendiente.

2. Creación del Raster

rm(list = ls())
library(raster)
## Loading required package: sp

Se va a crear un DEM. Hay que tener en cuenta, que es necesario pasar las coordenadas del cuadro delimitador (mínimo y máximo a lo largo de los dos ejes, x e y)

dem <- raster(ncol= 3, nrow= 3, xmn= 100, xmx= 115, ymn= 100, ymx= 115)

¿Cuántas celdas componen el DEM?

ncell(dem)
## [1] 9

¿Cuál es la resolución espacial (es decir,el tamaño de la celda ) del DEM?

res(dem)
## [1] 5 5

Asignar valores de elevación al DEM:

valores <- c(50, 45, 50, 30, 30, 30, 8, 10, 10)
(values(dem) <- valores)
## [1] 50 45 50 30 30 30  8 10 10

Ahora, se va a trazar el DEM a lo largo de los valores de elevación:

plot(dem, main = "DEM")
text(dem)

Después, se le asigna un sistema de referencia de coordenadas al DEM:

crs(dem) <- CRS("+init=epsg:3115")

¿Qué significa el código EPSG?

EPSG es el acrónimo de European Petroleum Survey Group y es una organización que mantiene una base de datos de parámetros geodésicos con códigos estándar, los códigos EPSG, para sistemas de coordenadas, datos, esferoides, proyecciones cartográficas y elipsoides de todo el mundo. Esta base de datos contiene los parámetros para estos objetos geográficos, y a cada uno de estos, se le asigna un número único.

Fuente: Support Portal. En: http://support.virtual-surveyor.com/es/support/solutions/articles/1000261353-que-es-un-c%C3%B3digo-epsg-

¿Dónde se pueden encontrar los códigos EPSG?

Los códigos EPSG se pueben buscar en epsg.io

¿Cuál es el código EPSG para MAGNA Ciudad de Bogotá?

Según, epsg.io el código EPSG para MAGNA-SIRGAS/ Colombia Bogota zone es EPSG: 3116

3. Cálculo de pendiente y aspecto.

(slope = terrain(dem, "slope", unit= "degrees", neighbors= 8))
## class      : RasterLayer 
## dimensions : 3, 3, 9  (nrow, ncol, ncell)
## resolution : 5, 5  (x, y)
## extent     : 100, 115, 100, 115  (xmin, xmax, ymin, ymax)
## crs        : +init=epsg:3115 +proj=tmerc +lat_0=4.596200416666666 +lon_0=-77.07750791666666 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs 
## source     : memory
## names      : slope 
## values     : 75.25766, 75.25766  (min, max)
plot(slope, main= "Pendiente")
text(slope)

(aspecto = terrain(dem, "aspect", unit= "degrees", neighbors= 8 ))
## class      : RasterLayer 
## dimensions : 3, 3, 9  (nrow, ncol, ncell)
## resolution : 5, 5  (x, y)
## extent     : 100, 115, 100, 115  (xmin, xmax, ymin, ymax)
## crs        : +init=epsg:3115 +proj=tmerc +lat_0=4.596200416666666 +lon_0=-77.07750791666666 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs 
## source     : memory
## names      : aspect 
## values     : 180.7538, 180.7538  (min, max)
plot(aspecto, main = "Aspecto")
text(aspecto)

4. Asignación:

Se debe escribir el código R para calcular la pendiente y el aspecto para un DEM que tiene 16 celdas y una resolución espacial de 5 5. - Para crear este DEM, hay que tener en cuenta que es necesario pasar las coordenadas del cuadro delimitador (mínimo y máximo a lo largo de los dos ejes, x e y)

Ndem <- raster(ncol= 4, nrow=4, xmn=100, xmx=120, ymn=100, ymx= 120 )

Número de celdas:

ncell(Ndem)
## [1] 16

Resolución espacial:

res(Ndem)
## [1] 5 5

Asignación de valores de elevación del DEM:

valores <- c(50, 45, 50, 48, 30, 29, 30, 29, 10, 9, 9, 10, 25, 23, 19, 21)
(values(Ndem) <- valores)
##  [1] 50 45 50 48 30 29 30 29 10  9  9 10 25 23 19 21

Ahora, se plotea:

plot(Ndem, main = "Nuevo DEM")
text(Ndem)

Ahora, se le asigna un sistema de coordenadas de referencia al DEM

crs(Ndem) <- CRS("+init=epsg:3115")
(slope = terrain(Ndem, "slope", unit= "degrees", neighbors = 8))
## class      : RasterLayer 
## dimensions : 4, 4, 16  (nrow, ncol, ncell)
## resolution : 5, 5  (x, y)
## extent     : 100, 120, 100, 120  (xmin, xmax, ymin, ymax)
## crs        : +init=epsg:3115 +proj=tmerc +lat_0=4.596200416666666 +lon_0=-77.07750791666666 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs 
## source     : memory
## names      : slope 
## values     : 36.05503, 75.62313  (min, max)
plot(slope, main= "Pendiente del Nuevo DEM")
text(slope)

(aspecto = terrain(Ndem, "aspect", unit = "degrees", neighbors = 8))
## class      : RasterLayer 
## dimensions : 4, 4, 16  (nrow, ncol, ncell)
## resolution : 5, 5  (x, y)
## extent     : 100, 120, 100, 120  (xmin, xmax, ymin, ymax)
## crs        : +init=epsg:3115 +proj=tmerc +lat_0=4.596200416666666 +lon_0=-77.07750791666666 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs 
## source     : memory
## names      : aspect 
## values     : 164.0546, 181.4688  (min, max)
plot(aspecto, main= "Aspecto del Nuevo DEM")
text(aspecto)