Introducción

este cuaderno es para ilustrar el cálculo de la pendiente y el 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 avión. 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, es la dirección de la brújula a la que se dirige una pendiente.

Creación de Raster

en esta sección creamos un raster de ficción para emular un DEM usando la función de raster de la biblioteca homónima

rm(list=ls())
library(raster)
## Loading required package: sp
dem <- raster(ncol=3, nrow=3, xmn=100, xmx=115, ymn=100, ymx=115)
ncell(dem)
## [1] 9
res(dem)
## [1] 5 5

asignar valores al DEM:

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

tracemos el DEM a lo largo de los valores de elevación:

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

asignemos un sistema de referencia de coordenadas a DEM:

crs(dem) <- CRS("+init=epsg:3115")
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on GRS80 ellipsoid in CRS definition,
##  but +towgs84= values preserved

La función del SRC incluye un término llamado EPSG, está relacionado con un código que proviene del EPSG Geodetic Parameter Dataset (también registro EPSG) que es un registro público de sistemas de referencia espacial, elipsoides terrestres, transformaciones de coordenadas y unidades de medida relacionadas. A cada entidad se le asigna un código EPSG entre 1024-32767, junto con una representación estándar de texto conocido legible por máquina (WKT). El conjunto de datos es mantenido activamente por el Comité de Geomática del IOGP.

La mayoría de los sistemas de información geográfica (SIG) y las bibliotecas de SIG utilizan los códigos del EPSG como identificadores de sistemas de referencia espacial (SRID) y datos de definición del EPSG para identificar las proyecciones y realizar transformaciones entre estos sistemas, mientras que algunos también apoyan los SRID emitidos por otras organizaciones (como Esri). (Tomado de https://en.wikipedia.org/wiki/EPSG_Geodetic_Parameter_Dataset)

Estos códigos se pueden encontrar en la página web https://epsg.io, el código para el MAGNA de Bogotá es EPSG:6247 para la red urbana y EPSG:3116 para la zona de Bogotá

Cálculo de la pendiente y el aspecto

Usaremos la función terreno para el cálculo de la pendiente y el aspecto y los graficaremos, podemos ver que la pendiente para es de 75 grados y el aspecto es de 181 grados, lo que significa que está orientado hacia el sur

(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        : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-77.0775079166667 +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= "Slope")
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        : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-77.0775079166667 +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 = "Aspect")
text(aspecto)

Asignación

Ahora calculamos la pendiente y el aspecto para un DEM con dimensión de celda de 16 y resolución de 5, usaremos la función CRS para asignar el código para el código ESPG de la red urbana de Bogotá MAGNA, que es 6247

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

valuesA <- c(50, 45, 50, 48, 30, 29, 30, 29, 10, 9, 9, 10, 25, 23, 19, 21)

(values(demA) <- valuesA) 
##  [1] 50 45 50 48 30 29 30 29 10  9  9 10 25 23 19 21
plot(demA, main = "DEM")
text(demA)

crs(demA) <- CRS('+init=epsg:3115')
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded datum Unknown based on GRS80 ellipsoid in CRS definition,
##  but +towgs84= values preserved

Ahora calculamos la pendiente y el aspecto y los trazamos

(slope1= terrain(demA, '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        : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-77.0775079166667 +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(slope1, main = "Slope")
text(slope1)

(aspectA = terrain(demA, '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        : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-77.0775079166667 +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)
(aspectA = terrain(demA, '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        : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-77.0775079166667 +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(aspectA, main = "Aspect")
text(aspectA)

Aquí podemos ver que la pendiente de las celdas superiores es de 75 y 76 grados, para las celdas inferiores es de 36 y 42 grados, lo que demuestra que es un terreno muy irregular con pendientes muy fuertes; por el aspecto podemos ver que tres de las celdas están orientadas hacia el sur (180-181 grados) y una está orientada hacia el sureste (164 grados)