Creación de vídeos al vuelo en el notebook

Vamos a ver cómo poder ver un vídeo en el notebook sin necesidad de crear una copia en disco del mismo. En realidad es una excusa para, en el camino, hablar de web map services (WMS), de imágenes satelitales, sistemas de coordenadas CRS, Cartopy, Matplotlib, meteorología, predicción, servidores de datos THREDDS,… Como todo lo anterior no cabe en un título lo dejaremos como ‘Creación de vídeos al vuelo en el notebook’.

Antes de nada, vamos a necesitar las siguientes bibliotecas, además de cosas de la stdlib:

  • Numpy (para manejar datos)
  • Matplotlib (para pintar en gráficos)
  • Cartopy (para pintar en mapas)
  • IPython + Jupyter notebook (para representar los datos en el notebook)
  • owslib (para acceder a Web Map Services)
  • XArray (para manejar ficheros NetCDF de forma sencilla)
  • NetCDF4 (driver para leer y escribir ficheros NetCDF)
  • Siphon (para acceder a datos en servidores THREDDS)
  • ffmpeg (para utilizar como respaldo de Matplotlib en la generación del vídeo).

Todo lo anterior es fácilmente instalable vía conda usando el canal de conda-forge tanto en Linux como en Windows. Los siguientes comandos en una consola o en el Anaconda Prompt (Windows) harán el trabajo:

Luego activamos el entorno virtual. Elige una sola de las líneas de debajo, la que corresponda con tu conda y sistema operativo:

Y, por último, dentro de nuestra terminal nos vamos a la carpeta en la que se encuentre este notebook y arrancamos el notebook:

Como siempre, vamos a importar todo lo necesario:

WMS

Un servicio de mapas web es un estándar definido por el Open Geospatial Consortium que ofrece, grosso modo, poder hacer determinadas peticiones al servicio y este nos devolverá información geográfica que podremos representar en un mapa.

Existen otra serie de servicios complementarios en los que no vamos a entrar puesto que lo que vamos a hacer hoy es bastante simple.

¿Para qué vamos a usar un WMS? Primero de todo, vamos a descargarnos imágenes de satélite del MSG. En este caso del MSG que está colocado a 0º de longitud actualmente. Si no ando muy equivocado este es el Meteosat 10 o Meteosat Second Generation 3 (el décimo de todos los Meteosat o el tercero de los Meteosat de segunda generación).

[Curiosidad, hace unos meses se cumplieron 40 años desde el lanzamiento del Meteosat 1].

Para el que no lo sepa, el MSG3 es un satélite de órbita geoestacionaria. Esto quiere decir que orbita desde el plano que define el ecuador terrestre mirando siempre al mismo punto de la tierra.

Bueno, al lío, vamos a descargar imágenes del MSG3 del canal infrarrojo 10.8 \(\mu m\).

[Breve inciso: ese canal nos permite diferenciar características térmicas y así ver diferencias entre cuerpos más calientes o más fríos como la superficie terrestre o la parte alta de las nubes].

Accedemos a un WMS de EUMETSAT.

Vemos los contenidos del WMS y destacamos del resto la capa que usaremos. También guardamos el nombre de la capa de interés en la variable layer.

Veamos las coordenadas disponibles para nuestra capa de interés:

Y los estilos disponibles para la capa:

Y los posibles sistemas de coordenadas a usar:

[Breve inciso: resumidísimo y sin entrar en mucho detalle, estamos analizando objetos tridimensionales que queremos representar en un espacio bidimensional como la pantalla de nuestro PC, teléfono, en un papel,… Necesitamos poder transformar posiciones de un espacio a otro. Para ello se han inventado los sistemas de referencia espaciales que nos permiten hacer esas transformaciones. Estas transformaciones son imperfectas ya que hemos de sacrificar algunas de las características (ángulo, área, distancia) del original ya que solo podemos mantener alguna de ellas sin distorsiones].

Vamos a usar WGS84, EPSG 4326 o Plate Carree un poco de forma indistinta aunque NO son lo mismo pero la extensión que tenemos no permite entrar en muchos detalles. Esto es lo que se suele usar para los mapas de OpenStreetMap.

Bueno, centrémonos, vamos a descargar las imágenes del canal IR 10.8 cada 6 horas desde hace 5 días hasta hoy usando un WMS de Eumetsat:

Definimos las fechas para las que descargaremos información:

Le decimos al WMS cómo queremos las cosas (qué capa, en qué área,…):

Hacemos la petición al WMS para cada fecha con nuestros requerimientos:

Obviamente, dependiendo del momento en que ejecutes el notebook o lo leas en nuestro blog las fechas no tienen porqué ser las mismas ni estar actualizadas…

Ahora mismo ya podríamos crear un vídeo, pero de momento solo vamos a usar matplotlib para visualizar las tres primeras imágenes:

Meteorología y predicción

Ahora vamos a usar los datos de un modelo global de predicción atmosférica. El modelo en cuestión es el GFS (Global Forecast System) y lo gestiona el National Center for Environmental Prediction (NCEP). Estos son los modelos que nos dicen si lloverá mañana y son muy sofisticados pero también hacen muchas concesiones para poder correrlos en un supercomputador en tiempo y forma. No voy a entrar mucho más en ello, si alguien tiene alguna pregunta concreta que se dirija a @Dr_Meteo, mi alter ego tuitero o use los comentarios en el blog.

[INCISO: esa cuenta de twitter está en modo mantenimiento por eso de la vigilancia masiva, si no contesto usad @pybonacci].

Bien, vamos a descargar datos del GFS cuya ejecución empezase en el momento de la primera imagen disponible del Meteosat (hace cinco días) y ver cómo el modelo ha predicho la nubosidad y precipitación hasta hoy viendo conjuntamente la imagen del Meteosat y los campos del modelo.

El modelo produce muchísimas variables, todas ellas interesantísimas, en muchos niveles verticales y para un horizonte de predicción de hasta 384 horas. Como he comentado anteriormente me voy a centrar en estas cinco:

NumberLevel/LayerParameterForecast ValidDescription
292surfaceAPCP0-6 hour accTotal Precipitation [kg/m^2]
317low cloud layerTCDC0-6 hour aveTotal Cloud Cover [%]
318middle cloud layerTCDC0-6 hour aveTotal Cloud Cover [%]
319high cloud layerTCDC0-6 hour aveTotal Cloud Cover [%]
320entire atmosphereTCDC0-6 hour aveTotal Cloud Cover [%]

Las variables anteriores son un promedio en un periodo de 6 horas por lo que no es exactamente lo mismo que nos ofrece una imagen de satélite, que es un instante, pero para lo que queremos hacer nos vale.

Vamos a descargar lo que necesitamos de un servidor THREDDS, que es un tipo de servidor que ofrece la posibilidad de descargar datos geofísicos siguiendo una serie de estándares, siguiendo la misma filosofía a lo que hemos visto con WMS.

Generamos la URL para el día en que comenzará la predicción, que será el mismo momento en que empiezan las imágenes del Meteosat:

Creamos una instancia de un objeto que nos permitirá descargar los datos del servidor THREDDS.

Veamos todos los datos disponibles:

Ahora vamos a usar un servicio especial de THREDDS que nos permite descargar solo un subconjunto de datos y además nos permite hacerlo en formato netCDF, que es más amigable para los que estén en Windows. El servicio en cuestión se llama NetCDFSubset (NCS). Como solo voy a descargar los datos en común con las imágenes de Meteosat voy, además, a seleccionar solo un subconjunto de todos los ficheros de GFS disponibles en THREDDS:

Generamos las url’s para cada momento de la predicción del modelo.

Empezamos a pedir los datos del modelo y los guardamos en memoria. Los datos a las 00.00 del primer dia no existen puesto que no puede haber un promedio en las 6 horas antes a empezar. Por tanto, esos datos del primer momento no se descargan.

AVISO: Esto puede tardar un rato, es mucha información y está metido en un while que seguirá intentando descargarlo todo hasta que no queden ficheros para descargar. Si, por la razón que sea, veis que tarda demasiado y los ficheros que está intentando descargar se repiten mucho cortadlo porque algo está pasando (el servidor está caído, está en mantenimiento,…).

Ya hemos descargado los ficheros y los tenemos en objetos NetCDF en memoria. Vamos a usar xarray para ver el primero (el de las 06:00 UTC desde la hora de inicio).

Cartopy

Vamos a usar cartopy para crear cada frame ya que nos permite trabajar de forma simple con mapas.

Un ejemplo de cómo crear un frame:

En los mapas con datos del modelo veis una raya en longitud 0º. Es debido a que los datos terminan antes de longitud 360º. Esto es fácilmente corregible pero os lo dejo como ejercicio.

Vídeo

Ya tenemos todos los ingredientes para crear el vídeo que queríamos hacer y poder visualizar directamente en el notebook. Tendremos que crear primero los frames que vamos a usar en el vídeo y después crear la animación.

Para lo primero podemos hacer:

Generamos la animación:

Y, por último, lo visualizamos en el notebook:

En caso de que queráis una copia en disco podéis usar:

Si el vídeo no se ve en la web siempre lo podéis descargar de aquí.

Saludos.

Deja un comentario

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

three + one =