Saltar al contenido

Animación del iceberg A68-A usando imágenes MODIS de AQUA en Python

Etiquetas:

Hace poco hemos estado dibujando la trayectoria del iceberg A68 usando cartopy y folium. En el artículo de hoy vamos a ver una parte de esa trayectoria desde el espacio usando información del instrumento MODIS a bordo del satélite AQUA, xarray, matplotlib, cartopy y Python.

La parte de trayectoria que vamos a procesar corresponderá al fragmento A68-A en su aproximación a las islas Georgias del sur después de haber vagado más de mil kilómetros durante meses por las gélidas aguas más meridionales del océano atlántico.

MODIS y AQUA

El AQUA es un satélite de observación de la tierra operado por la NASA y lleva, entre otros, un instrumento llamado MODIS. MODIS viene de moderate resolution imaging spectroradiometer.

El MODIS recoge información en 36 bandas distintas con resoluciones horizontales que varían entre los 250 metros y los 1000 metros dependiendo de la banda o canal que estemos usando.

Para el ejercicio de hoy voy a usar las bandas 1 (rojo a 250 metros), la banda 3 (azul a 500 metros) y la banda 4 (verde a 500 metros). Voy a usar un producto con algo de preproceso que me dan las imágenes agregadas a una resolución de 1 kilómetro para todas las bandas (producto L1B).

Los tres canales o bandas que voy a usar están dentro del visible (rojo, azul y verde) y con ellas crearé una imagen visible.

La información del producto L1B la podéis descargar de aquí recortadas para una zona de interés, que es mucho más práctico que descargar las imágenes completas de cada pasada del satélite.

Si quieres reproducir lo que voy a hacer aquí, este es el listado de imágenes usado:

All files subset to the following area of interest:
GEOLOCATION:
             NorthBound = -52.7
             SouthBound = -60.7
             WestBound  = -45.5
             EastBound  = -35.2

All output will be reprojected as follows:
    Projection = GEO,    ResampleType = nearest,
The following files were specified for post-processing:
FILENAME: MYD021KM.A2020356.1845.061.2020357152436.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020356.1705.061.2020357152223.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020356.1710.061.2020357151613.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020356.1525.061.2020357151624.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020356.1530.061.2020357150952.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020355.1940.061.2020356152932.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020355.1800.061.2020356152931.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020355.1620.061.2020356152042.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020355.1805.061.2020356151936.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020355.1625.061.2020356151413.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020354.1855.061.2020355152125.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020354.1715.061.2020355152104.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020354.1540.061.2020355152104.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020354.1720.061.2020355151627.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020354.1545.061.2020355151604.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020353.1950.061.2020354152209.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020353.1815.061.2020354151637.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020353.1635.061.2020354151415.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020353.1640.061.2020354150823.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020351.1650.061.2020352152406.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020351.1830.061.2020352152338.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020349.1700.061.2020350154442.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020348.1755.061.2020349152953.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020348.1620.061.2020349151636.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020347.1715.061.2020348161921.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020346.1630.061.2020347165526.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020346.1810.061.2020347165416.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020345.1725.061.2020346152510.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020345.1550.061.2020346151738.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020344.1645.061.2020345152147.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020343.1600.061.2020344185609.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020343.1740.061.2020344185543.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020342.1655.061.2020343153112.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020341.1750.061.2020342152916.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020341.1615.061.2020342151135.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020340.1710.061.2020341151109.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020339.1805.061.2020340152426.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

FILENAME: MYD021KM.A2020339.1625.061.2020340151611.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020338.1720.061.2020339155700.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020337.1635.061.2020338161833.hdf
     SDS: MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020337.1815.061.2020338161638.hdf
     SDS: MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 01 (645 nm)

FILENAME: MYD021KM.A2020336.1730.061.2020337221013.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 04 (555 nm), MODIS Spectral Band 03 (469 nm)

FILENAME: MYD021KM.A2020336.1555.061.2020337220910.hdf
     SDS: MODIS Spectral Band 01 (645 nm), MODIS Spectral Band 03 (469 nm), MODIS Spectral Band 04 (555 nm)

Todas las imágenes tienen formato HDF. He hablado anteriormente sobre este formato en varias ocasiones.

Bibliotecas a usar

Como siempre, importamos una serie de cosas de las que haremos uso:

from pathlib import Path
from io import BytesIO

import numpy as np
import matplotlib.pyplot as plt
import rasterio
import xarray as xr
import imageio
import cartopy.crs as ccrs

rasterio nos ayuda para poder leer el formato HDF. Lo llamamos desde xarray. cartopy lo uso para crear el mapa mientras que imageio me ayudará a crear la animación.

Creación de cada fotograma

La siguiente función incluye toda la funcionalidad que necesito para generar cada uno de los fotogramas de la animación final:

def create_frame(r_file, g_file, b_file):
    # read data
    r = xr.open_rasterio(r_file)
    g = xr.open_rasterio(g_file)
    b = xr.open_rasterio(b_file)
    
    lon = r.coords['x'].values
    lat = r.coords['y'].values
    
    # play with histogram
    R = np.clip(r.data[0,...] * 1 / 65535, 0, 1)
    G = np.clip(g.data[0,...] * 1 / 65535, 0, 1)
    B = np.clip(b.data[0,...] * 1 / 65535, 0, 1)
    
    # raw
    RGB = np.dstack([R, G, B])
    
    # more histogram stuff
    R_ = np.power(R, 1 / 2)
    G_ = np.power(G, 1 / 2)
    B_ = np.power(B, 1 / 2)
    
    RGB_ = np.dstack([R_, G_, B_])
    
    xmin = -49
    xmax = -32
    ymin = -60
    ymax = -50
    
    fig, (ax1, ax2) = plt.subplots(
        ncols=2, 
        nrows=1, 
        figsize=(20, 6),
        subplot_kw={'projection': ccrs.PlateCarree()}
    )
    
    ax1.stock_img()
    ax1.coastlines(resolution='50m', color='yellow', linewidth=2)
    ax1.set_xlim(xmin, xmax)
    ax1.set_ylim(ymin, ymax)
    ax1.imshow(
        RGB, 
        origin='upper',
        extent=(lon.min(), lon.max(), lat.min(), lat.max())
    )
    
    ax2.stock_img()
    ax2.coastlines(resolution='50m', color='yellow', linewidth=2)
    ax2.set_xlim(xmin, xmax)
    ax2.set_ylim(ymin, ymax)
    ax2.imshow(
        RGB_, 
        origin='upper',
        extent=(lon.min(), lon.max(), lat.min(), lat.max())
    )
    
    fig.suptitle('.'.join(str(r_file).split('.')[:-2]), fontsize=15)
    fig.tight_layout()
    
    buffer = BytesIO()
    fig.savefig(buffer, format='png', dpi=300)
    
    buffer.seek(0)
    
    return buffer

La función recibe la ruta de cada fichero de cada banda (rojo, verde y azul) para un momento dado de pasada del satélite sobre la zona de interés y hace varias cosas:

  • Lee cada uno de los canales usando xarray y rasterio.
  • Extraemos las posiciones de los datos.
  • Toqueteamos un poco los valores para tener un histograma un poco modificado.
  • Guardamos la imagen con los tres canales en RGB y en RGB_
  • Finalmente, usamos cartopy y matplotlib para crear la imagen final.

El tema de los histogramas, no quiero entrar ahora en mucho detalle. Si alguien tiene un poco de interés le puede echar un ojo a este enlace.

En la función se crean dos imágenes, una para RGB y otra para RGB_. En la segunda se ve ‘más’ color. Lo dicho, no voy a entrar mucho en ello ahora e igual en un próximo artículo me meto más en detalle.

Una vez que tenemos la función solo nos queda gestionar todo el tema de recolectar las imágenes, guardar los fotogramas y crear la animación. Eso lo puedes ver en la siguiente porción de código:

# Asumo que las imágenes están en el mismo sitio
# desde el que se ejecuta el código
r_files = sorted(list(Path('.').glob('*EV_250_Aggr1km_RefSB_1.hdf')))
g_files = sorted(list(Path('.').glob('*EV_500_Aggr1km_RefSB_2.hdf')))
b_files = sorted(list(Path('.').glob('*EV_500_Aggr1km_RefSB_1.hdf')))

# creamos los frames
frames = []
for r_file, g_file, b_file in zip(r_files, g_files, b_files):
    buffer = create_frame(r_file, g_file, b_file)
    frames.append(imageio.imread(buffer))

# creo el vídeo/animación
imageio.mimsave('A68_MODIS_SouthGeorgia.mp4', frames, fps=1)

El anterior código no tiene mucha historia. Si tienes alguna duda usa los comentarios (de más abajo) para preguntar. Sobre imageio ya hablé aquí.

El resultado final sería (ponlo a pantalla completa y verás que algo, el iceberg, se va acercando desde el oeste a la isla para luego girar hacia el sur):

Un saludo.

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *

− four = two

Pybonacci