Cómo calcular la temperatura superficial del suelo con Landsat y R

Afuera llueve, obvio, porque es Valdivia y estamos en pleno invierno, y el sonido de las gotas contra el techo de zinc es lo único que me acompaña mientras peleo con una línea de código que no quiere funcionar. Les cuento esto porque hace unos meses, a finales del verano pasado, me puse en la cabeza que quería entender por qué mi casa se sentía tan fresca comparada con el centro de la ciudad, y así fue como terminé metida en el mundo de la temperatura superficial del suelo (o LST, por sus siglas en inglés), tratando de descifrar qué dicen realmente esos píxeles que bajan desde el espacio.

Una tarde de lluvia y píxeles rojos

Todo empezó una tarde de lluvia en junio, de esas en las que el café se enfría antes de que alcances a darle el segundo sorbo y el ventilador de mi laptop empieza a zumbar como si estuviera a punto de despegar (pobre, procesar toda la región de Los Ríos en un raster no es tarea fácil). Yo estaba toda emocionada porque ya sabía crear composiciones de color en RStudio para ver el bosque nativo, pero pasar de ver colores a entender el calor es otra historia; y entonces, cuando por fin logré que R me mostrara un mapa, me salió un cuadro completamente rojo, con valores que decían que el suelo estaba a más de mil grados, ¡una locura!

Me sentí derrotada, la verdad. Estuve a punto de cerrar la tapa y ponerme a ver una serie, pero después de mirar fijo la pantalla un buen rato, me di cuenta de que había cometido el error más típico de principiante: olvidé aplicar el factor de escala multiplicativo en la fórmula de radiancia. Los datos de Landsat no vienen listos para leer en grados Celsius; vienen como Niveles Digitales (DN) que son solo números enteros, y si no usas las constantes que vienen en ese archivo de texto gigante llamado MTL, los resultados no tienen ningún sentido. Para pasar de esos números raros a algo real, primero hay que llegar a Kelvin y luego restar 273.15, que es la constante estándar para convertir a Celsius. Sin ese pequeño detalle, mi mapa parecía el centro de un volcán en vez de una tarde nublada.

Detalle de archivo de metadatos MTL de Landsat en una pantalla de computadora

El laberinto de EarthExplorer y el archivo MTL

Para hacer esto, lo primero fue ir a EarthExplorer (que siempre me marea un poco con tantos filtros) y buscar una imagen de Landsat 8 que no tuviera nubes, algo que aquí en el sur es casi un milagro de fin de semana. Cuando descargas la carpeta, te encuentras con un montón de archivos .TIF. Lo importante es saber que para la temperatura vamos a usar la Banda 10, que es la térmica. Aunque Landsat 8 y 9 traen la Banda 11 también, en los foros de gente que sabe mucho dicen que es mejor quedarse con la 10 por unos problemas de calibración que tiene la otra (algo de una luz dispersa, no me preguntes el detalle técnico porque todavía me pierdo).

Lo que me salvó la vida fue encontrar el archivo .MTL. Es un archivo de texto que la mayoría ignora, pero ahí están los coeficientes de reescalado y las constantes K1 y K2 que son específicas para cada imagen. Sin eso, estás adivinando. Una vez que cargué la Banda 10 en R y leí los metadatos, el proceso empezó a tener forma. Es un alivio cuando ves que los números en la consola empiezan a bajar de 40.000 a valores humanos de dos cifras.

De números raros a temperatura de brillo

A finales de marzo, cuando todavía quedaba algo de sol, empecé a notar que no bastaba con la temperatura de brillo. Verán, la temperatura de brillo es lo que el satélite 've' asumiendo que la Tierra es un espejo perfecto (un cuerpo negro, le dicen en los libros), pero la realidad es que el suelo no emite calor de forma perfecta. Aquí es donde entra la longitud de onda central de la Banda 10 de Landsat 8, que anda por los 10.89 micrómetros. Es una medida minúscula pero vital para los cálculos físicos.

En R, el código se ve más o menos así (después de muchas pruebas y errores): primero transformas los DN a radiancia espectral usando los multiplicadores del MTL, y luego usas esa radiancia para calcular la temperatura de brillo. Pero ojo, que esa temperatura todavía es 'mentirosa' porque no considera si lo que está mirando el satélite es agua, cemento o mi querido bosque valdiviano. El sensor TIRS de Landsat tiene una resolución espacial nativa de 100 metros (aunque cuando lo descargamos ya viene remuestreado a 30 metros para que calce con las otras bandas), y en esos 100 metros pueden pasar muchas cosas.

El truco de la emisividad (y por qué no es tan perfecto)

Para arreglar esa temperatura 'mentirosa', necesitamos calcular la emisividad. Lo que casi todos los tutoriales te dicen es que calcules el NDVI y que, a partir de ahí, saques la proporción de vegetación para estimar la emisividad. Yo lo hice así por varias semanas, sintiéndome una experta, hasta que un domingo por la mañana me puse a leer con más calma y me di cuenta de una cosa: corregir la emisividad solo con el NDVI es una simplificación excesiva. El problema es que el NDVI asume que si no hay verde, todo se comporta igual, y eso ignora totalmente lo complejo que es el suelo, especialmente en las ciudades.

En áreas urbanas densas, como el centro de Valdivia o las zonas industriales, este método provoca errores significativos porque el cemento, el asfalto y los techos de zinc emiten el calor de formas muy distintas a un suelo desnudo con tierra. Al final, mi mapa me decía que ciertas zonas estaban más frías de lo que yo sentía al caminar por ahí. Aprendí que, aunque para empezar está bien usar el NDVI, si uno quiere hacer algo serio para la ciudad, hay que tener cuidado. Si te interesa limpiar un poco tus datos antes de estos pasos, a veces ayuda saber pasos para reclasificar un raster en R y simplificar mapas de vegetación para separar mejor lo que es bosque de lo que es construcción.

Cuaderno con fórmulas matemáticas junto a un mapa térmico satelital en RStudio

Ver el calor del cemento desde mi ventana

Después de varios meses de errores en la consola y de mirar el zumbido de mi pobre notebook, por fin logré renderizar el mapa final. Fue una satisfacción increíble. Ver el mapa en la pantalla y notar cómo el bosque cerca de mi casa, ese verde profundo que veo por la ventana, resaltaba en tonos azules y celestes (frío), mientras que el cemento de la costanera y los techos del centro se veían en naranjos intensos, me hizo sentir que todo el tiempo invertido valió la pena.

No soy geógrafa, y a veces me pregunto por qué estudiar teledetección satelital con R para ver cambios ambientales si a veces es tan frustrante, pero luego veo esos mapas y entiendo que es la única forma de ver lo que a simple vista no notamos. Entender que la temperatura del suelo no es solo un número, sino el resultado de cómo tratamos la tierra y cuántos árboles dejamos en pie, es algo que me cambió la forma de mirar mi propio barrio. Ahora, cada vez que salgo y siento el calor que rebota del asfalto, me acuerdo de mi código en R y de esos 10.89 micrómetros que, desde allá arriba, nos están contando la historia de nuestro clima.

Artículos relacionados