Protocolo automático para el tratamiento de imágenes Landsat

Jornadas SIG Libre

Geotech/Spatial Data Science

17-18 de septiembre de 2025 | Girona

Importante!!

 

  • GEE Account

  • Anaconda

  • SACO/NextCloud:

https://saco.csic.es/s/PNpRW7Bz9366mFJ

 

  • Introducción teórica sobre el Protocolo

  • Práctica con la librería Ndvi2Gif para obtención de datos auxiliares de las escenas

  • Procesado de una escena

  • Obtención del hidroperiodo del ciclo 2024/2025

Esquema previsto (dios mediante) del Taller:

Laboratorio de SIG Y Teledetección

Estación Biológica de Doñana

http://www.ebd.csic.es/web/last/inicio

http://last-ebd.blogspot.com/

  • Geoportal y Servicios Web Mapping
  • Asistencia científica en SIG y Teledetección
  • Asistencia equipamiento tecnológico (GPS, radiometría, ...)
  • Cursos (GEE, GeoPython, SIG, Teledetección, drones)
  • Programa de seguimiento Doñana e ICTS
  • Vuelos de drones (Multi/Hiperespectrales, térmicas, LiDAR)
  • Impresión 3D

Javier Bustamante Díaz

Responsable científico

(Investigador Científico)

Ricardo Díaz-Delgado Hernández

Dr. Biología

Isabel Afán Asencio

Responsable Técnico

Dra. Biología

David Aragonés Borrego

Ingeniero Forestal

Diego García Díaz

Geógrafo

Pedro Gómez Giráldez

Dr. Ingeniería de Montes

Gabriela Patricia Romero Olivos

Geógrafa

Protocolo Automático v2

Máquina Virtual

  • CentOS 7
  • 64 bits
  • 64 Gigas RAM
  • 32 núcleos
  • Conectado a un NAS (32 TB)

Características

Software

  • MongoDB
  • GDAL
  • Git
  • Anaconda (python 3.8)
    • Numpy
    • Rasterio
    • Fiona
    • Geopandas
    • Scipy
https://github.com/Digdgeo/ProtocoloV2

auto_download.py

protocolov2.py

productos.py

coast.py

hidroperiodo.py

utils.py

crontab

GitHub Repository

Código

Carpetas

/ori

/geo

/rad

/nor

/pro

/hyd

/data

/mongo

auto_download.py

download_landsat_scenes(username,
password,
latitude,
longitude,
days_back=15,
end_date=None, 
process=True,
max_cloud_cover=100, 
output_dir='path/to/save/rars')

Es donde sucede "la magia". Se trata del script que se encarga de buscar escenas nuevas para el rango de fechas seleccionado, cotejarlas con la base de datos y en caso de no estar, realizar la descarga y llamar a las clases Landsat() y Product() para procesar la escena. Si todo va bien, termina enviando un mail llamando a las funciones de utils.py

  • Desarrollado en el LAST-EBD como respuesta a un problema complejo, como es el estudio de la inundación en la marisma de Doñana.
  • Las condiciones cambiantes de la marisma hacen muy difícil definir una metodología para estimar la cubierta de agua (en una serie temporal) por medio de índices de agua u otro tipo de algoritmos.
  • Por ello, se decidió llevar a cabo una normalización (basada en Áreas Pseudo Invariantes (PIAs)) de las imágenes en base a la reflectividad de unas áreas en las que a priori debería ser homogenéa y constante a lo largo del tiempo (arenas, mar y embalses, cascos históricos, etc...)
  • Versión 2 del Protocolo automático para el tratamiento de imágenes Landsat.
- Uso de imágenes ori en L2 (Reflectividad en superficie)

- Nuevos productos ofrecidos

- Integración con Geoportal

- Nueva metodología para la máscara de agua (NoData)

- Integramente en Python

 

  • Basado en las campañas de campo realizadas por personal del LAST y de la Estación Biológica de Doñana.
  • Reajustado dentro del Proyecto Ecopotential con unos 1440 puntos de campo para Landsat 8.
  • Descarga Automática (SR-Colección 2 L2)
  • Geo (recorte 202/34)
  • Rad (aplicar coeficientes escalado)
  • Normalización
  • Productos 
  • Actualización de MongoDB y PostgreSQL

 

Busqueda y descarga 

Hillshade

MDT

Geo & Rad

Productos

  • NDVI
  • NDWI
  • MNDWI
  • FLOOD
  • DEPTH
  • Turbidity
  • Coastline
  • Hydroperiod

Normalización

Escena

Landsat 

Pias 

PIAs Ref

FMASK 

Data

No

GeoServer

ICTS

LAST

├── 20250615l9oli202_34_censo_aereo_l3.csv
├── 20250615l9oli202_34_flood.png
├── 20250615l9oli202_34_lagunas_carola.csv
├── 20250615l9oli202_34_lagunas_principales.csv
├── 20250615l9oli202_34_resumen_lagunas_carola.csv
├── 20250615l9oli202_34_resumen_lagunas.csv
├── 20250615l9oli202_34_rgb.png
└── 20250615l9oli202_34_superficie_inundada.csv

 

20250615l9oli202_34

{"_id": "20240815l9oli202_34",  
"usgs_id": "LC92020342024228LGN01",  
"tier_id": "LC09_L1TP_202034_20240815_20240815_02_T1",  
"lpgs": "LPGS_16.4.0",  
"category": "T1",  
"Clouds": {"cloud_scene": 0.06,    
			"land cloud cover": 0.06,    
            "cloud_PN": 0  },  
"Info": {    
	"Tecnico": "LAST-EBD Auto",    
	"Iniciada": {"$date": {"$numberLong": "1724690285816"}},    
	"Pasos": {"rad": "",      
	"nor": {"Normalize": "True",        
	"Nor-Values": {
    		"blue": 
				{"Parametros": 
                	{"slope": 1.020172035444719,              
                	"intercept": -0.00711237939135306,              
                	"std": 0.016095327213406563,              
                	"r": 0.98906989518987,              
                	"N": 54886,              
                	"iter": 1},            
                "Tipo_Area": 
                	{"Mar": 49347,              
                	"Embalses": 793,              
                    "Pinar": 1564,              
                    "Urbano-1": 575,              
                    "Urbano-2": 944,              
                    "Aeropuertos": 547,              
                    "Arena": 664,              
                    "Pastizales": 266,              
                    "Mineria": 186}},  
    "Productos": 
         ["NDVI",    "NDWI",    "MNDWI",    
         {"Flood": {        
                 "El Rincon del Pescador": 72.72,        
                 "Marismillas": 0.99,        
                 "Caracoles": 0,        
                 "FAO": 14.94,        
                 "Marisma Occidental": 1.26,        
                 "Marisma Oriental": 6.57,        
                 "Entremuros": 115.29}    },    
          "Turbidity",   
          "Depth"  ]}
{"_id": "hidroperiodo_2023-2024",  
"escenas": [    
	{"escena_id": "20240409l9oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 19896.300000000003},    
    {"escena_id": "20240503l8oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 11161.53},    
    {"escena_id": "20240511l9oli202_34",      
    "nubes_marismas": 8.31,      
    "ha_inundacion": 8354.789999999999},    
    {"escena_id": "20240527l9oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 3754.53},    
    {"escena_id": "20240604l8oli202_34",      
    "nubes_marismas": 0.06,      
    "ha_inundacion": 2661.84},    
    {"escena_id": "20240714l9oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 249.93},    
    {"escena_id": "20240722l8oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 215.1},    
    {"escena_id": "20240807l8oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 209.43},    
    {"escena_id": "20240815l9oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 211.76999999999998},    
    {"escena_id": "20231008l8oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 305.91},    
    {"escena_id": "20231024l8oli202_34",      
    "nubes_marismas": 7.69,      
    "ha_inundacion": 1131.8400000000001},    
    {"escena_id": "20231117l9oli202_34",      
    "nubes_marismas": 0,      
    "ha_inundacion": 408.51},    
    {"escena_id": "20231125l8oli202_34",      
    "nubes_marismas": 0.03,      
    "ha_inundacion": 421.83},    
    {"escena_id": "20231219l9oli202_34",      
    "nubes_marismas": 4.76,      
    "ha_inundacion": 589.5899999999999},    
    {"escena_id": "20231227l8oli202_34",      
    "nubes_marismas": 0.63,      
    "ha_inundacion": 512.8199999999999},    
    {"escena_id": "20240112l8oli202_34",      
    "nubes_marismas": 0.09,      
    "ha_inundacion": 696.69},    
    {"escena_id": "20240120l9oli202_34",      
    "nubes_marismas": 0.01,      
    "ha_inundacion": 1138.59},    
    {"escena_id": "20240221l9oli202_34",      
    "nubes_marismas": 7.73,      
    "ha_inundacion": 4288.32},    
    {"escena_id": "20240229l8oli202_34",      
    "nubes_marismas": 0.01,      
    "ha_inundacion": 3901.3200000000006},    
    {"escena_id": "20240316l8oli202_34",      
    "nubes_marismas": 6.05,      
    "ha_inundacion": 4150.98}]}

Landsat

Hidroperiodo

Productos

  • Inundación
  • Profundidad
  • Hidroperiodo
  • Línea de Costa

  • dtm

  • fmask_escena

  • hillshade_escena

  • ndvi_escena

  • ndwi_escena

  • mndwi_escena

  • cob_veg

  • ndvi_p10

  • ndvi_mean

  • fmask

  • mndwi

  • ndwi

  • swir1

Criterios:

Swir1 <= 0.12

slope > 8 (**embalses NDWI_scene & MNDWI_scene**)

hillshade < p30

ndvip10 > 0.3 & ndvimean > 0.5

cobveg > 75

ndvi_scene > 0.6 & dtm > 2.5

fmask_scene clouds & shadows

sum(reclass((fmask_scene, ndwi_scene, mndwi_scene)) >= 2

 

 

Flood

hidroperiodo.py

funciones:

get_escenas_values()

get_hydroperiod()

get_products()

Gini() & IRT()**

get_escenas_values() se encarga de calcular los días "que vale" cada escena del ciclo hidrológico (1 de septiembre - 31 de agosto).

get_hydroperiod() calcula el número de días que un pixel ha estado inundado o seco o ha sido NoData.

get_products() genera los rasters con los hidroperiodos y el número de días válidos del ciclo para cada pixel. 

Valid days

Hydroperiod

 Primer día inundado     Último día inundado        

 

  • Integración en el Protocolo del LAST
  • Facilidad para obtener la línea de costa a partir de la máscara de agua del Protocolo v2
  • Rapidez del procesado (toda la serie Landsat en ~20 minutos)
  • Proceso inmediatamente transferible a Sentinel 2 y otros satélites/sensores
  • Incapacidad de detectar la duna embrionaria por la resolución espacial de Landsat.
  • Las líneas de costa anuales obtenidas no permiten generar un dtm, ya que se superponen los valores.

¡Toda esta metodología es transferible a otras escenas!

Solo se necesitan una imagen de referencia y unas PIAs sobre la misma. Junto con una serie de capas auxiliares que son fácilmente obtenibles con Ndvi2Gif.

NDVI2GIF -> GitHubtoColab

Creación del entorno de Conda:

  1. conda env create -f environment.yml
  2. conda activate landsat
  3. python -m ipykernel install --user --name landsat --display-name "Landsat (jslg)"
  4. jupyter lab
  5. opcional para ndvi2gif si falla ipyleaflet: pip install -U anywidget jupyterlab_widgets ipywidgets jupyter-leaflet