miércoles, 7 de diciembre de 2011

Análisis Geomorfológico e Hidrológico de Modelos Digitales de Elevaciones con MATLAB

Análisis Geomorfológico e Hidrológico de Modelos Digitales de Elevaciones con MATLAB
Geomorphological and Hydrological Analysis of Digital Elevation Models with MATLAB
Por Eber Risco Sence


El análisis físico del terreno a partir de los valores de elevación contenidos en las celdas del Modelo Digital del Terreno, quizás por la inherente buena disposición de los mismos para el estudio mediante procedimientos y algoritmos diversos, así como por el mayor nivel de detalle y precisión respecto a las representaciones clásicas, permite obtener resultados de mucha mayor amplitud, tanto en cantidad y calidad como en el propio significado de los mismos (Olaya 2004). 

Caracterización Matemática
Evans (1972) menciona que el terreno puede ser expresado mediante una superficie cuadrática de la forma:
 
Al hacer que esta función se ajuste en la mayor medida posible a la superficie real del terreno tal y como ésta se encuentra recogida en el MDT con una menor o mayor precisión, resulta poco lógico, intuyéndose de antemano que la fidelidad para con el relieve real de tal función será poco menos que nula. No obstante, si en lugar de trabajar a nivel global con la totalidad de la malla lo hacemos a nivel local tomando el entorno inmediato de un punto, esta simplificación cobra mayor sentido y los parámetros que puedan calcularse a partir de la ecuación que se deduzca en dicha porción de la malla serán más representativos de las características reales del punto tomado (Olaya 2004).
De esa manera podemos determinar el nivel local a una malla de 9 celdas centrada alrededor de la celda denominada z5.
 
A partir de lo anterior la definición de los parámetros que configuran la superficie cuadrática a la que se pretende aproximar el entorno local de z5 de acuerdo con las siguiente expresiones, todas ellas en función de las 9 celdas consideradas en dicho entorno y de la distancia entre ellas, que coincide lógicamente con el tamaño de celda y que representamos de forma simbólica como delta de s.


De acuerdo a la formula polinómica, los coeficientes p, q, r, s ,t representan las siguientes derivadas parciales:
Figura 1: Modelo digital de elevaciones


Figura 2: Modelo digital de elevaciones mostrando curvas de nivel

Figura 3: Modelo digital de elevaciones en 3D


Parámetros geomorfológicos

Pendiente
Partiendo de la función matemática descrita anteriormente, y a la cual asimilábamos la forma del MDE en un entorno local de un punto dado, la pendiente puede calcularse a través de las primeras derivadas de dicha función (Olaya 2004).
Así la pendiente se define de la siguiente manera:
La pendiente deseada representa el ángulo entre dicho vector gradiente y la vertical, el cual, haciendo uso de conceptos básicos de cálculo de ángulos entre vectores, se obtiene según la expresión
Figura 4: Modelo digital de pendientes

Orientación
La orientación en un punto puede definirse como el ángulo existente entre el vector que señala el norte y la proyección sobre el plano horizontal del vector gradiente (Felicísimo 1994).
La forma habitual de utilizar este parámetro es expresado en grados sexagesimales, de acuerdo a la siguiente expresión (Olaya 2004).

 
Figura 5: Modelo digital de orientaciones
Sombreado
El sombreado en una superficie se obtiene de una iluminación hipotética, mediante la determinación de los valores de iluminación para cada celda. Esto se logra mediante el establecimiento de una fuente de luz hipotética y el cálculo de los valores de iluminación de cada celda en relación a las celdas vecinas (ESRI 2008).
Figura 6: Modelo digital de sombreado

Rugosidad
Upton y Fingleton, citados por Felicísimo (1994) mencionan que las coordenadas rectangulares de un vector unitario perpendicular a la superficie en un punto i vienen dadas por las expresiones:

El módulo del vector suma de un conjunto de vectores es un indicador de agrupación y, por tanto, inversamente proporcional a la rugosidad. El módulo del vector suma se calcula, para un conjunto de n datos vecinos al punto problema (los 8 más próximos, por ejemplo) mediante la expresión:
Band, citado por Felicísimo (1994) menciona que resulta conveniente estandarizar el valor de R dividiéndolo por el tamaño muestral obteniendo así el módulo medio. El resultado puede variar entre los valores extremos de 0 (dispersión máxima) y 1 (alineamiento completo). El módulo medio es complementario del parámetro estadístico denominado varianza esférica.
 Figura 7: Modelo digital de rugosidad
 Índice Elevación-Relieve
Olaya (2010) en la ayuda del programa SEXTANTE para gvSIG describe este índice según la siguiente expresión:
 Siendo Zm, Zmin, Zmax los valores medio, mínimo y máximo de elevación en la ventana de análisis 3x3 alrededor de cada celda.
 Figura 8: Índice de elevación-relieve
Índice Varianza Pendiente
Felicísimo (1994) describe este índice como el análisis mediante una ventana móvil sobre el modelo digital de pendientes. El cálculo de la varianza se realiza sobre los valores incluidos en ella y el resultado puede utilizarse como estimación de la rugosidad.
 Figura 9: Índice de varianza pendiente


Parámetros geomorfológicos (Concavidad/Convexidad)
Curvatura
Felicísimo (1994) describe la curvatura en un punto, h, puede definirse como la tasa de cambio en la pendiente y depende, por tanto, de las derivadas de segundo grado de la altitud (es decir, de los cambios de pendiente en el entorno del punto).
 Figura 10: Curvatura


Figura 11: Fuente: Hengl  y Reuter (2008)

Curvatura vertical
Krcho y Young, citados por Hengl  y Reuter (2008) definen la curvatura vertical de acuerdo a la siguiente expresión:

Figura 12: Curvatura vertical

Curvatura tangencial (horizontal)
Krcho y Hofierka, citados por Hengl  y Reuter (2008) definen la curvatura tangencial de acuerdo a la siguiente expresión:
Figura 13: Curvatura tangencial

Curvatura planta
Esta curvatura es perpendicular a la dirección de la máxima pendiente. Evans, citado por Hengl  y Reuter (2008) definen esta curvatura con la siguiente expresión:
Figura 14: Plan curvature

Curvatura inesférica (Unsphericity curvature)
Shary, citado por Hengl  y Reuter (2008) introduce un sistema de 12 curvaturas, donde las tres curvaturas básicas son: mean, unsphericity y difference curvaturas, son usadas para derivar las restantes nueve curvaturas. La Unsphericity curvatura (UNSPHC) puede definirse como:
Estos parámetros describen como la superficie difiere de una esfera y consecuentemente tiene valores de cero solo si la superficie es una esfera. 

Figura 15: Unsphericity curvature

Curvatura diferencial (Difference curvature)
Es la segunda curva independiente de las doce curvas, se define como la mitad de la diferencia entre la curvatura vertical y horizontal (Shary, citado por Hengl  y Reuter 2008).
DIFFC indica cuál de las dos curvas es formalmente la más fuerte.
  Figura 16: Difference curvature

Las tres curvas independientes (MEANC, DIFFC y UNSPHC), las restantes curvas pueden ser derivadas: la curvatura mínima y máxima (MINC, MAXC), curvatura horizontal y vertical (HEXC, VEXC), la curvatura Gaussiana total (TOTGC), curvatura acumulación total (TOTAC) y curvatura anillo total (TOTRC) (Shary, citado por Hengl  y Reuter 2008).

Curvatura Gaussiana total (Total Gaussian curvature)
Puede determinarse usando la siguiente expresión (Gauss, citado por Hengl  y Reuter 2008).
Figura 17: Total Gaussian curvature

Curvatura ROTOR (ROTOR curvature)
Florinsky, citado por Hengl  y Reuter (2008) usa también ROTOR (una curvatura de las líneas de flujo que son perpendiculares a los contornos).

ROTOR es una medida de la torsión de las líneas de flujo. Las líneas de flujo giran hacia la derecha (sentido de las agujas del reloj) cuando ROTOR>0, mientras las líneas de flujo giran en sentido anti-horario cuando ROTOR<0.
 Figura 18: ROTOR curvature

Parámetros hidrológicos
Acumulación de flujo
Hengl  y Reuter (2008) describen el principio de la acumulación de flujo: cuando la proporción d de drenaje de una celda a sus vecinos (debe sumar 1) son conocidos, también la proporción de recepción r que drena hacia una celda. Las proporciones de recepción determinan que fracción de cada vecino son recibidas. La cantidad de masa (o volumen, área u otra propiedad) A que se acumula en la celda i está dada por la suma de A en cada celda vecina multiplicada por la respectiva fracción de recepción r más la masa u otra cualidad ingresante I en la celda misma.
 Figura 19: Flujo acumulado (Matriz simple)

Figura 20: Flujo acumulado (Matriz múltiple)

Cuencas hidrográficas
La dirección de flujo es el elemento clave para la delimitación de las cuencas hidrográficas, es con el cual se da forma con toda precisión a las cuencas vertientes.
El conocimiento de esta dirección de flujo de cada celda nos permitirá deducir para cada una de ellas, y siguiendo el recorrido desde las mismas hacia aguas abajo, si pertenecen o no a la cuenca que pretendemos definir (Olaya 2004).
 Figura 21: Cuencas hidrográficas

Índice de humedad
El indice de humedad para cada celda está representado mediante la siguiente expresión (Trauth 2010)
 
 Figura 22: Índice de humedad

Índice de potencia de cauce
El índice de potencia de cauce es frecuentemente usado en geomorfología, ciencias del suelo y disciplinas relacionadas. Como una media de potencia de cauce es un indicador del potencial transporte de sedimentos y de erosión del agua. Se define mediante la siguiente expresión ((Trauth 2010)
 Figura 23: Índice de potencia de cauce

Parámetros geomorfológicos (Factores de pendiente)
Factor de profundidad S
Kinnell (s.f.) menciona que el factor de profundidad de pendiente puede ser estimado mediante la siguiente formula:
 Donde β es el ángulo de la pendiente.
Nearing, citado por Sanchez (s.f.) propone la siguiente fórmula para determinar el factor de  la profundidad de la pendiente.
 Donde θ es el ángulo de la pendiente.
Figura 24: Factor profundidad de pendiente
 Figura 25: Factor profundidad de pendiente (Nearing)
Factor Longitud L
Kinnell (s.f.) menciona que el factor de longitud de pendiente puede ser estimado mediante la siguiente formula:
Donde  Ai,j-in es el área contribuyente de escorrentía para la celda i,j y D es el tamaño  de celda, x es un factor que determina las variaciones en el ancho del flujo producto de la orientación de la celda con respecto a las curvas de nivel, m es el exponente de la longitud de la pendiente, definido por: m=0.1342*Ln(θ)+0.192
 Figura 26: Factor longitud de pendiente

Factor longitud y profundidad de Pendiente LS 
Se representa mediante el producto de los dos factores anteriores.
 
Figura 27: Factor profundidad- longitud de pendiente

Parámetros hidro-geomorfológicos (Transporte de sedimentos)
Capacidad de transporte de sedimentos
Es un índice utilizado para estimar el potencial topográfico para la erosión o deposición por medio de una expresión que representa el cambio en la capacidad de transporte de sedimentos en la dirección del flujo (Martínez 1999).
 Donde
LS=índice de capacidad de transporte de sedimentos.
AS=área de drenaje específica.
β = ángulo de la pendiente.
 Figura 28: Capacidad de transporte de sedimentos

Referencias:
1. Chavez, Richard. Sf. Modelo de Riesgo de Erosión de Suelo en cuatro micro-cuencas de la Cordillera los Maribios, Nicaragua.
2.    ESRI (Environmental Systems Research Institute, US). 2008. Manual de Ayuda de ArcGIS.
3.  Evans, I. S., General geomorphometry, derivatives of altitude, and descriptive statistics. En Chorley, R. J. 1972. (ed.) Spatial Analysis in Geomorphology, Methuen, London. pp.17-90.
4.   Felicísimo, Ángel. 1994. Modelos Digitales del Terreno: Introducción y aplicaciones en las ciencias ambientales. Pentalfa. 122 p.
5. Hengl, T., Reuter, H.I. (eds) 2008. Geomorphometry: Concepts, Software, Applications. Developments in Soil Science, vol. 33, Elsevier, 772 p.
6.    Kinnell, P. sf. The miscalculation of the USLE topographic factors in GIS. University of Camberra, Australia.
7.     Martínez-Casasnovas, J.A., 1999. Modelos digitales de terreno: Estructuras de datos y aplicaciones en el análisis de formas del terreno y en Edafología. QUADERNS DMACS Núm. 25, Departament de Medi Ambient i Ciències del Sòl, Universitat de Lleida, Lleida.
8.  Olaya Ferrero, V.2004. Hidrología Computacional y Modelos Digitales del Terreno. Madrid, ES, Víctor Olaya.365 p.
9.   Olaya Ferrero. 2010. Manual de Ayuda de SEXTANTE para gvSIG. 

Video Aplicación en MATLAB: