Omitir los comandos de cinta
Saltar al contenido principal
Inicio de sesión

Saber más

Administrar permisosAdministrar permisos
Caso 3: Imágenes Sentinel2 - Seguimiento de zonas inundadas (zonas húmedas y balsas) - (post 2/2)
18/04/2020 - Producto, Formación, Noticia

Segundo Post en el que, una vez evaluado el estado hídrico de zonas húmedas y balsas, según una imagen en Falso Color Natural de 3 de febrero de 2020 (que se realiza trimestralmente), prepararemos a los servicios contraincendios un archivo compatible con el sistema GPS y los dispositivos móviles, con sus coordenadas y una etiqueta con la disponibilidad de agua para repostar, que permita elegir al mando, bomberos forestales y pilotos dónde acudir entre los lugares más cercanos al siniestro.

Proceso vectorial

2. Ahora descargamos la cartografía de detalle con embalses, balsas, piscinas, etc. Lo haremos desde el mapa Topográfico de Navarra a escala 1.5.000, modelo Base Topográfica Armonizada (MTNa5-BTA), fecha de referencia 2017 o más actualizada (cada 3 años) en formato SIG. Interesa el shapefile "Hidrografia_pol_2017" (o 2020, 2023…). Solo el Embalse de Yesa (casi todo en territorio aragonés), está parcialmente contenido.

Interesa la información de los atributos id_fenomen, ESTADO y Descrip, especialmente éste último y el área de las entidades por debajo del tamaño del píxel (10x10 m más el previsible "desencuadre" geográfico entre ambas, con lo que nos daremos un margen de seguridad de 3 x 3 píxeles y eliminaremos aquellas menores a 900 m2).

Desde la tabla de atributos, en selección por expresión, recogemos aquellos en los que "id_fenomen" no contenga los tipos que interesan: 'Embalse', 'Estanque' y 'Laguna'.

"id_fenomen" LIKE 'Corriente%' OR "id_fenomen" = 'Isla' OR  "id_fenomen"  =  'Piscina' OR $area<900

 

Activamos la sesión de edición y borramos estas features. Cerramos y guardamos. Ya tenemos las zonas inundadas que nos interesa analizar.

 

Igual que hicimos en el ejercicio del caso 2, utilizaremos: Procesos – Caja de herramientas – Geoalgoritmos de QGIS – Análisis ráster – Estadísticas de zona. Como en la figura:

 

Calculando los estadísticos Número, Suma y Mayoría (moda) para cada entidad geográfica vectorial, intentaremos, no solo conocer si tienen o no agua (utilizando la moda sería muy simplista) sino si, respecto al polígono (que son las huellas de las zonas inundadas en "máximo nivel de ocupación"), podemos ofrecer una pequeña clasificación de su estado según sea el régimen en el momento de la imagen. Para ello nos ayudaremos del total de píxeles analizados (Número) con valor 1 (Suma).

Con la calculadora de campos, creamos un nuevo atributo "PRES_AGUA" de tipo Número entero, que proporcione el % de píxeles con agua según el ráster "B4yB8-129".

 

round("_sum" * 100 / "_count" ,0)

 3. También sería un dato relevante el volumen de agua disponible, aunque el cálculo no puede hacerse de forma ni siquiera aproximada, ya que, aunque el sensor utilizado en el LiDAR de 2017 tenía cierta capacidad batimétrica, ésta no superaba los 3 m. Podremos conocer con cierta exactitud la altitud de la superficie, en la mayoría de los casos será inexacta la del fondo.

El objetivo de este dato no es otro que guiar al equipo de emergencias para distinguir, de las fuentes de agua más próximas al incendio, cuales permitirán recargar, incluso varias veces, el bambi (cesta – bomba para recoger y arrojar el agua).

Así que nos permitiremos simplificar mucho y tomarlo como una aproximación (suponiendo que tratamos con un cilindro, cuando, si acaso, es un cono truncado), como el área / 2 x altura * PRES_AGUA / 100. No vaya a ser que el helicóptero trate de cargar la cesta de agua y no haya volumen suficiente; así que seremos conservadores con las distintas tipologías de "vasos" que pueden tener embalses, lagunas, estanques y balsas, aplicando un coeficiente corrector del -20% de inclinación de los bordes.

De todos estos datos, tenemos el área (función $area) y el % de presencia de agua; solo nos falta la altura de la lámina de agua. Para ello nos ayudaremos de los Modelos Digitales de Elevación (MDE) que sirve IDENA.

Damos de alta el servicio WCS de IDENA

 

 https://idena.navarra.es/ogc/wcs

E incorporamos al proyecto: MDT_2M_VE2017 y MDS_2M_VE2017.

El MDT es el Modelo Digital del Terreno (la superficie de la tierra sin edificios, vegetación, etc. y el MDS es el Modelo Digital de Superficies, tal como lo capturó el LiDAR.

La diferencia entre ambos es el valor de "altura"*, que necesitamos.

* Entre los servicios WCS está Alturas_2M_VE2017 pero su objetivo es mostrar la de edificaciones y vegetación, con lo que no sirve para este propósito.

La belleza de ambas coberturas se puede apreciar, modificando la simbología a Mapa de Sombras (Hillshade), pero es lo que tenemos en los mapas de relieve en blanco y negro por WMS. La diferencia aquí, radica en que el ráster tiene valores reales de altitud para una mala de 2 x 2 m/píxel y no es una imagen como el caso del WMS. Por esto tarda más en el proceso de carga.

Para la operación posterior, en cada capa: Propiedades – Fuente, simplificaremos el nombre (en el proyecto) a MDT y MDS respectivamente.

Como ya hicimos anteriormente, con el Geoalgoritmo: Análisis ráster – Estadísticas de zona. Como en la figura:

 

Calculamos en este caso la Media de altitud del MDT (suelo) y MDS (vuelo) desde ambos MDE en cada zona inundada, utilizando ese nombre para el prefijo.

Atención al tiempo de procesado, ¡que estamos trabajando con ráster por WCS con una resolución de 2 x 2 m!

Ya con toda la información necesaria, mejor guardamos la sesión de edición para que se consoliden los nuevos atributos del MDE. Luego, con la Calculadora de campos, creamos el nuevo atributo "VOLUMENm3", de tipo Número entero (64 bit), con la expresión:

ROUND($area / 2 * ("MDS_mean"-"MDT_mean") * ("PRES_AGUA" / 100) * ("PRES_AGUA" * 0.2),0)

 Recordamos la fórmula y vemos la expresión: área / 2 x altura * PRES_AGUA / 100 con un coeficiente corrector del -20% por la inclinación de los bordes (20% de pendiente media)

  • a) área: función $area (superficie en m2 de la entidad geográfica poligonal);
  • b) altura: ("MDS_mean"-"MDT_mean"), que son medias de los MDE para todo el polígono, en este caso imprescindible entre () para que tenga prioridad la resta sobre el resto de la expresión;
  • c) porcentaje del polígono ocupado por agua según el análisis ráster: ("PRES_AGUA" / 100);
  • d) coeficiente corrector por pendiente del vaso: ("PRES_AGUA" * 0.2)
  • e) el atributo se creó como de tipo entero (grande, de 64 bit), para evitar incómodos decimales, para un dato que ya declaramos poco preciso: round(valor,nº decimales), realiza el redondeo al entero más próximo.
  • f) recordemos que el único polígono incompleto es el Embalse de Yesa (<1/5 parte), pero la existencia de agua está prácticamente asegurada al sr la zona próxima a la presa.

Cerramos la sesión de edición y guardamos los cambios.

4. Con estos valores, ajustamos una simbología graduada por este atributo, 4 o 5 clases es suficiente, una que recoja =0 (atención a que la primera clase es >= y <=, resto > y <=, con enlazar contornos de clase activado), que proporcione una valoración – guía durante las tareas de extinción.

Como el archivo debe ser compatible con aplicaciones que utilizan información GPS, mantener la simbología que hemos preparado y aceptar geometrías de polígonos, optamos por el formato KML.

Una vez conformes con la clasificación según la simbología, aprovechamos los valores de corte para crear con la calculadora de campos un nuevo atributo "ETIQUETA", de tipo Texto, tamaño 20. Cada uno ajusta los valores a sus resultados, con la expresión:

CASE

WHEN  "VOLUMENm3" =0 THEN 'Vacío'

WHEN "VOLUMENm3" >0 AND "VOLUMENm3" <= 1000 THEN 'Nivel bajo'

WHEN "VOLUMENm3" > 1000 AND "VOLUMENm3" <= 10000 THEN 'Nivel medio'

WHEN "VOLUMENm3" >10000 AND "VOLUMENm3" <= 50000THEN 'Nivel alto'

ELSE 'Completo'

END

Cerramos la sesión de edición y guardamos los cambios.

5. La mejor creación de un KML nos la brinda el propio Google Earth Pro. Todos lo tendréis ya instalado, pero si no es así, descargadlo desde:

https://www.google.es/earth/download/gep/agree.html. El código de activación es GEPFREE

Como la capa "Hidrografia_pol_2017" tiene codificación de caracteres ISO-8859-15 (latin-9) y GE trabaja en UTF-8, prepararemos la capa: Exportar - Guardar objetos como ´formato shape, nombre de shapefile "Repostaje_AGUA_CONTRA_INCENDIOS", SRC EPSG:4326 (WGS84), codificación UTF-8. Añadir el archivo al mapa.

Con la tabla de atributos de esta nueva capa a la vista, activamos la sesión de edición, para eliminar aquellos campos que no interesan. Dejamos NOMBRE, Descrip, VOLUMENm3 y ETIQUETA. Cerramos y guardamos los cambios.

 

6. Ya en Google Earth Pro – Archivo – Importar – Tipo Shape – Archivo:

Repostaje_AGUA_CONTRA_INCENDIOS.shp

¿Desea aplicar una plantilla de estilo a los elementos introducidos? SI – Crear nueva plantilla (Aceptar si es la primera vez, si no, podemos aplicar la del trimestre anterior).

En el cuadro de diálogo "Configuración de plantilla de estilo"

- Pestaña Nombre: campo: NOMBRE

- Pestaña Color: Establecer color desde campo: ETIQUETA. Rampa de rojo a azul. Modificar la asignación de colores entre los niveles medio y bajo. Check de "Crear subcarpetas para cada segmento" activado. Rellenar los nombres que harán de carpetas agrupadoras.

- Pestaña Altura: Sujetar elementos al suelo

 

Aceptamos, guardamos la plantilla de estilo y se crea el archivo KML, que puede llamarse igual que el shapefile, porque, aunque lo importemos en QGIS, el archivo prj que se crea recoge el mismo SRC. Se sitúa en el panel de contenidos en "Sitios temporales" y, por defecto, con la visibilidad desactivada y el árbol de contenidos compactado.

Activamos la visibilidad y desplegamos el árbol:

 

Como las carpetas aparecen desordenadas y bajo la que replica el nombre del archivo, con todas compactadas, pueden reordenarse, solo arrastrándolas, eliminando finalmente la que agrupaba a todas, de forma muy similar al trabajo en el Panel de capas de QGIS con los grupos:

 

Ahora, finalmente, sobre "Repostaje_AGUA_CONTRA_INCENDIOS.shp", botón derecho – Propiedades, eliminamos ".shp" y aceptamos. Y de nuevo con botón derecho – Guardar sitio como, sobrescribimos el KML o creamos un KMZ, según el software de los dispositivos dónde vaya a ser utilizado.

Si lo incorporamos en QGIS, ofrecerá la carga de las carpetas creadas en GEP, incorporando cada una de ellas como una capa, pero con la plantilla por defecto de los archivos KML, con múltiples campos en valor NULL y perdiendo la simbología.

También se puede importar al Visor de IDENA, dónde mantendrá simbología y modelo de datos:

 

Lo valioso es que ahora tenemos a nuestra disposición un procedimiento y un archivo final con información actualizable cada trimestre, iniciativa del Departamento de Desarrollo Rural y Medio Ambiente, que poner a disposición de los servicios de emergencias.

Ah. Y recordad que se puede guardar definitivamente el proyecto en QGIS.

 

Esperamos que os haya gustado este ejercicio.

En teoría, cerramos esta serie. Ya tenéis a vuestra disposición:

Caso 1: Imágenes Sentinel2 - Detección de cambios en zonas rurales (parte 1 y parte 2)

Caso 2: Seguimiento del ciclo de cultivos en una serie de parcelas (parte 1 y parte 2)

Caso 3: Seguimiento de zonas inundadas (zonas húmedas y balsas) - (parte 1 y parte 2)

Más aquellas peticiones que nos hagáis llegar a sitna@navarra.es

SITNA

Datos adjuntos