Diego García Díaz PRO
GIS & Remote Sensing. Python & Statistics. Surf & MTB :)
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
- 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
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:
- conda env create -f environment.yml
- conda activate landsat
- python -m ipykernel install --user --name landsat --display-name "Landsat (jslg)"
- jupyter lab
- opcional para ndvi2gif si falla ipyleaflet: pip install -U anywidget jupyterlab_widgets ipywidgets jupyter-leaflet
JSLG
By Diego García Díaz
