💡 OBJETIVOS DE LA ACTIVIDAD:

► Familiarizarse con los procedimientos necesarios para realizar mosaicos de imágenes de satélite cuando una única escena o gránulo es insuficiente para el análisis de cualquier fenómeno sobre el territorio.

💡 MATERIALES PARA LA ACTIVIDAD:

Los materiales para el desarrollo de esta actividad son los mismos de la actividad pasada.

► Imágenes para mosaico. Contiene los siguientres ficheros:

► script_laboratorio_3_mosaico_imágenes

INTRODUCCIÓN

Un mosaico es la combinación de dos o más imágenes raster en una única imagen. En el caso de las imágenes de satélite, una única escena no cubre completamente el área de estudio deseada, por lo que puede requerirse combinar dos o más escenas. Las órbitas de los satélites están diseñadas para permitir el solape o superposición parcial de las imágenes, compartiendo franjas comunes que son utilizadas como guías a tal efecto.

CARGA DE LOS PAQUETES NECESARIOS PARA LA ACTIVIDAD Y DEFINICIÓN DE LA ESTRUCTURA DE CARPETAS

Como en casos anteriores, se trabajará con terra.

library(terra)
## Warning: package 'terra' was built under R version 4.3.3
## terra 1.8.21

En esta actividad es conveniente crear una carpeta denominada lab_3 que incluya una subcarpeta denominada datos/mosaico_imagenes ey ésta, a su vez, dos carpetas adicionales:

A continuación, debe establecerse la carpeta dónde se encuentran los ficheros que serán utilizados en esta sesión.

Una vez creado dicha carpeta, se descargarán los ficheros necesarios para esta actividad directamente desde el servidor web.

Y se descomprimirán

Un procedimiento alternativo es realizar todo el proceso manualmente. Para ello, los ficheros necesarios para esta práctica se descargan aquí

IMPORTACIÓN DE FICHEROS RÁSTER

Se crea el objeto ráster raster1 a partir de un fichero en formato .tif

raster1 <- rast("D:/G174_2025/LABORATORIO_3_Preprocesamiento_imagenes_satelite/datos/mosaico_imagenes/raster1.tif")

A continuación se revisa el contenido de este objeto.

raster1 
## class       : SpatRaster 
## dimensions  : 1488, 1142, 3  (nrow, ncol, nlyr)
## resolution  : 30, 30  (x, y)
## extent      : 414855, 449115, 4438995, 4483635  (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 30N (EPSG:32630) 
## source      : raster1.tif 
## names       : raster1_1, raster1_2, raster1_3 
## min values  :      8881,      7609,      6917 
## max values  :     65535,     65535,     65535

Es muy habitual que los mosaicos se realicen a partir de dos o más escenas de Landsat contiguas. Es habitual que las escenas incluyan valores NA correspondientes al marco espacial sin datos de cada imagen. En consecuencia, es conveniente asignar el valor 0 a los píxeles con un valor NA de ambas escenas. Este valor sin datos se ignora al asignar las funciones. Esto se consigue con la función NAflag()[https://rspatial.github.io/terra/reference/NAflag.html]

NAflag(raster1) <- 0  # Asigna 0 a los valores NA

A continuación se cartografía la primera imagen en falso (FCC), que luega será utilizada en la composición.

plotRGB(raster1, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Primera imagen")

A continuación se cargan las restantes imágenes.

raster2 <- rast("D:/G174_2025/LABORATORIO_3_Preprocesamiento_imagenes_satelite/datos/mosaico_imagenes/raster2.tif")
raster3 <- rast("D:/G174_2025/LABORATORIO_3_Preprocesamiento_imagenes_satelite/datos/mosaico_imagenes/raster3.tif")
raster4 <- rast("D:/G174_2025/LABORATORIO_3_Preprocesamiento_imagenes_satelite/datos/mosaico_imagenes/raster4.tif")

Y se cartografía todas ellas en un único gráfico

par(mfrow = c(2, 2))
plotRGB(raster1, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Primera imagen")
plotRGB(raster2, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Segunda imagen")
plotRGB(raster3, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Tercera imagen")
plotRGB(raster4, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Cuarta imagen")

Volvemos al estado original

par(mfrow = c(1, 1))

CREACIÓN DE UN MOSAICO DE IMÁGENES

Cuando el paquete terra crea un mosaico con varias imágenes, aquella parte que es común (normalmente, en uno o varios bordes) se superponen entre sí. Además, todos los objetos deben tener un origen, una resolución y un sistema de coordenadas de referencia iguales.

Para crear un nuevo objeto raster con unas dimensiones mayores que la de los objetos individuales debe usarse la función terra::mosaic()[https://rspatial.github.io/terra/reference/mosaic.html]. Esta función aplica a su vez una serie de rutinas para computar el valor de los píxeles en las zonas de superposición. En R, se pueden aplicar tres funciones a la región superpuesta: la media, el mínimo y el máximo.

Mosaico aplicando la función promedio a la región superpuesta

Para aplicar la función promedio a la región superpuesta se asigna el argumento fun = mean. Esta función calcula la media de los valores de los píxel de la región superpuesta y asigna ese valor medio al ráster de salida.

mosaico_promedio <- mosaic(raster1, raster2, raster3, raster4, 
                          fun = mean)
mosaico_promedio 
## class       : SpatRaster 
## dimensions  : 2595, 1994, 3  (nrow, ncol, nlyr)
## resolution  : 30, 30  (x, y)
## extent      : 414855, 474675, 4438845, 4516695  (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 30N (EPSG:32630) 
## source(s)   : memory
## varname     : raster1 
## names       : raster1_1, raster1_2, raster1_3 
## min values  :      8881,      7609,      6846 
## max values  :     65535,     65535,     65535
plotRGB(mosaico_promedio, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Mosaico creado con la función promedio")

Mosaico aplicando la función min() a la región superpuesta

Para aplicar esta función, se asigna el argumento fun = min, que toma el valor mínimo de los píxeles de las regiones superpuestas de las imágenes ráster que se van a unir y lo asigna a la imagen de salida.

mosaico_minimo <- mosaic(raster1, raster2, raster3, raster4, 
                      fun = min)
mosaico_minimo 
## class       : SpatRaster 
## dimensions  : 2595, 1994, 3  (nrow, ncol, nlyr)
## resolution  : 30, 30  (x, y)
## extent      : 414855, 474675, 4438845, 4516695  (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 30N (EPSG:32630) 
## source(s)   : memory
## varname     : raster1 
## names       : raster1_1, raster1_2, raster1_3 
## min values  :      8881,      7609,      6846 
## max values  :     65535,     65535,     65535
plotRGB(mosaico_minimo, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Mosaico creado con la función mínimo")

Mosaico aplicando la función max() a la región superpuesta

En este caso se asigana el valor máximo de píxel de las regiones superpuestas, pasando el argumento fun = max. Crea la mejor imagen sin pérdida de datos, y por eso es la elegida habitualmente.

mosaico_maximo <- mosaic(raster1, raster2, raster3, raster4, 
                      fun = max)
mosaico_maximo 
## class       : SpatRaster 
## dimensions  : 2595, 1994, 3  (nrow, ncol, nlyr)
## resolution  : 30, 30  (x, y)
## extent      : 414855, 474675, 4438845, 4516695  (xmin, xmax, ymin, ymax)
## coord. ref. : WGS 84 / UTM zone 30N (EPSG:32630) 
## source(s)   : memory
## varname     : raster1 
## names       : raster1_1, raster1_2, raster1_3 
## min values  :      8881,      7609,      6846 
## max values  :     65535,     65535,     65535
plotRGB(mosaico_maximo, 
        r=3, g=2, b=1, 
        stretch ='lin', 
        axes =TRUE, 
        main = "Mosaico creado con la función máximo")

GRABAR EN DISCO

writeRaster(mosaico_promedio, "D:/G174_2025/LABORATORIO_3_Preprocesamiento_imagenes_satelite/datos/mosaico_imagenes/originales/mosaico.tif")