Cómo unir varios rasters en RStudio para ver paisajes extensos

Afuera llueve, para variar. No es esa lluvia suave que te deja caminar tranquila, sino ese aguacero valdiviano que golpea los vidrios y te hace sentir que el mundo se va a inundar en cualquier momento. Estoy acá sentada, con mi mate ya un poco lavado, mirando la pantalla de mi laptop que, pobrecita, hace lo que puede. Hace unos meses, durante una semana de lluvias intensas en agosto, me puse en la cabeza que quería ver todo el cordón de bosque nativo que rodea la cuenca, no solo el pedacito que alcanzaba a ver desde mi ventana, sino todo el paisaje extendido. Y ahí fue cuando me choqué con la realidad: los satélites no sacan fotos a medida de lo que uno quiere ver.

Resulta que yo estaba súper entusiasmada con mis imágenes de Sentinel-2, que tienen esa resolución espacial de 10 metros en las bandas visibles que te permite ver casi cada árbol, pero cuando quise cargar la zona que me interesaba, me di cuenta de que mi bosque quedaba justo en la costura. Sí, justo en el borde donde termina una imagen y empieza la otra. Es como tener un rompecabezas de mil piezas y darte cuenta de que compraste dos cajas distintas que no sabes cómo encajar.

El drama de los bordes y las cajas invisibles

A finales de la primavera pasada, me pasé varias tardes frustrada porque tenía dos archivos .tif que, en mi cabeza, deberían haberse pegado solos al abrirlos en RStudio. Pero no. En la pantalla se veían como dos cuadrados separados y, lo peor de todo, es que cuando intentaba cargarlos juntos, mi computador se ponía a silbar. Literalmente. El sonido del ventilador de mi laptop acelerándose al máximo mientras la consola de RStudio procesa el mosaico de tres gigabytes es algo que ya me pone nerviosa, como si el aparato fuera a salir volando o a derretirse en mis manos.

Pantalla de laptop mostrando dos imágenes satelitales separadas en el software RStudio.

Lo que pasa es que cada imagen de Sentinel-2 es un 'tile' o mosaico estándar de 100x100 km. Para el satélite, el mundo está dividido en esas cuadrículas perfectas, pero para los que vivimos abajo, la naturaleza no sabe de cuadrículas. Mi bosque favorito empezaba en el kilómetro 99 de una imagen y seguía en el kilómetro 1 de la siguiente. Y ahí estaba yo, sin saber qué era una proyección ni un raster, tratando de entender por qué si yo veía las fotos 'pegadas' en el navegador, en R se veían como dos entes que no se hablaban.

Al principio, mi miedo era romper algo. Pensaba que si usaba la función equivocada iba a desordenar los píxeles o, peor, que los datos iban a quedar desplazados y el mapa me iba a mentir. Como te contaba en aquel post sobre por qué elegí RStudio para estudiar el bosque nativo chileno, no soy programadora, así que cada paso que doy es como caminar por un humedal: con cuidado de no hundirme en el barro de los errores de código.

Descubriendo el paquete 'terra' (mi salvavidas)

Después de mucho buscar en foros y cerrar pestañas que me daban dolor de cabeza, entendí que mi mejor aliado iba a ser el paquete terra. Antes la gente usaba uno que se llamaba raster, pero parece que terra es más rápido y maneja mejor la memoria, algo vital cuando tu compu tiene ya sus añitos y el procesador sufre con cada mapa.

Aprendí que los rasters no son fotos mágicas, sino matrices gigantes de números. Y unir dos matrices no es como pegar dos fotos con masking tape; hay que asegurarse de que hablen el mismo idioma. Ese idioma es el Sistema de Referencia de Coordenadas (CRS). Una tarde gris de hace un par de meses, me di cuenta de que mis imágenes no se unían porque una estaba en una zona UTM y la otra en una distinta, o algo así entendí. Si no tienen el mismo CRS, R simplemente no sabe cómo poner una al lado de la otra.

Para simplificarme la vida, decidí que todo debía estar en el estándar internacional, el código EPSG para WGS 84, que es el 4326. Es como el 'inglés' de los mapas. Una vez que logré que ambos archivos tuvieran la misma proyección usando la función project(), sentí que por fin estaba avanzando, aunque todavía me faltaba el paso final: el mosaico.

El error que me hizo cerrar la laptop

Aquí es donde viene la parte donde casi me rindo. El fin de semana pasado, después de haber preparado todo, escribí el comando para unir las imágenes. Le di 'Enter' y me fui a hacerme un té. Volví, y ahí estaba: ese texto rojo en la consola diciendo 'extents do not overlap' después de haber esperado diez minutos a que terminara el proceso. ¡Casi me pongo a llorar! ¿Cómo que no se solapan si yo las veo una al lado de la otra?

Error en la consola de RStudio indicando que las extensiones de los rasters no coinciden.

Resulta que a veces, por unos pocos metros de diferencia en los bordes, R se confunde. O quizás yo había descargado mal una de las imágenes. Si quieres evitar este tipo de dolores de cabeza desde el principio, te recomiendo mirar cómo lo hago yo en mi guía sobre cómo descargar imágenes Sentinel para procesar en RStudio fácilmente, porque ahí explico cómo elegir bien el área para que no te pasen estas cosas.

El problema era que estaba usando la función merge(), que es un poco más estricta. Al final, descubrí que mosaic() es mucho más amable. La diferencia, según lo que entendí de un manual que leí a medias, es que mosaic() te permite decidir qué hacer en las zonas donde las imágenes se solapan (donde se enciman una con otra). Puedes pedirle que saque el promedio de los valores de los píxeles, lo que hace que la costura sea casi invisible. Con merge(), simplemente pone una encima de la otra y a veces se nota un salto de color medio feo.

El truco de los 'VRT' o cómo no matar tu RAM

Aquí te va el secreto que de verdad me cambió la vida y que es un poco distinto a lo que te dicen los tutoriales básicos. La mayoría de la gente te dice: 'une todos tus rasters y guárdalos en un archivo .tif gigante'. ¡No lo hagas! Bueno, a menos que tengas un computador de la NASA. Si unes tres o cuatro imágenes de Sentinel-2, vas a terminar con un archivo de varios gigabytes que va a hacer que RStudio se cuelgue cada vez que intentes hacer un zoom.

Lo que yo aprendí a hacer es crear un 'Mosaico Virtual' o VRT. Es como una instrucción que le dice a R: 'Mira, estos cinco archivos son en realidad uno solo, pero no los pegues todavía, solo muéstramelos juntos'. Es mucho más eficiente porque procesa los datos bajo demanda. No saturas la memoria RAM porque el archivo físico gigante nunca se llega a crear en tu disco duro, solo existe en la 'imaginación' de R mientras trabajas. Para alguien como yo, que aprende de a poco los fines de semana, esto fue como descubrir el fuego.

Al final, después de pelear con el código y escuchar el ventilador sufrir, logré ver la cuenca completa. Ahí estaba el río Valdivia serpenteando, el verde intenso del bosque nativo que se volvía más oscuro hacia la cordillera, y todo sin esa raya molesta en el medio que me dividía el paisaje. Me quedé un buen rato mirando la pantalla, con el mate ya frío, sintiendo esa satisfacción de cuando algo por fin cuadra. Afuera seguía lloviendo, pero en mi pantalla, por fin, el sol de los datos satelitales iluminaba todo mi mapa sin costuras.

Artículos relacionados