Search code examples
rmergerastermaskterra

Simplify R scripts to load multiple SpatVectors, SpatRaster and then mask the spatial rasters by using SpatVectors


I have an R script for masking several rasters by using corresponding spatial vectors, as the rasters have a 5 km buffer area need to be cropped before I can merge all of them together. I load the data manually one by one from folders, which look a lot of repetiton.

I would like to know if it is possible to simplify my scripts to avoid loading date one by one and do each process one by one?

I have a lot of raster files (more than 100 .tif files) need to be masked by corresponding spatial vectors and merged to a final large map.

library(terra)

setwd("I:/France_grid_1e5")

##Load data
A10 <- vect("./fr_grid_A10.shp")
A100 <- vect("./fr_grid_A100.shp")
A101 <- vect("./fr_grid_A101.shp")
##more than 80 .shp files need to be loaded here and follow the same procedures below.####

country_shp <- vect("I:/VoCC/Data/France_shapefiles/France.shp")

#A10
raster10 <- rast("I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/10FVoCC_fr_Euclidean_25m_ssp370_1e5grid_5kmBuffer.tif")
plot(raster10)
raster10
plot(A10, add = TRUE)
raster10 <- mask(raster10, A10)
writeRaster(raster10, filename = "I:/VoCC/Data/Output/RemoveBuffer_fr/10FVoCC_fr_Euclidean_25m_ssp370_1e5grid.tif", overwrite=TRUE)

#A100
raster100 <- rast("I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/100FVoCC_fr_Euclidean_25m_ssp370_1e5grid_5kmBuffer.tif")
plot(raster100)
raster100
plot(A100, add = TRUE)
raster100 <- mask(raster100, A100)
writeRaster(raster100, filename = "I:/VoCC/Data/Output/RemoveBuffer_fr/100FVoCC_fr_Euclidean_25m_ssp370_1e5grid.tif", overwrite=TRUE)

#A101
raster101 <- rast("I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/101FVoCC_fr_Euclidean_25m_ssp370_1e5grid_5kmBuffer.tif")
plot(raster101)
raster101
plot(A101, add = TRUE)
raster101 <- mask(raster101, A101)
writeRaster(raster101, filename = "I:/VoCC/Data/Output/RemoveBuffer_fr/101FVoCC_fr_Euclidean_25m_ssp370_1e5grid.tif", overwrite=TRUE)


####Merge the raster tiles####

raster_files <- list.files("I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/", pattern = ".tif$", full.names = TRUE)

rasters <- lapply(raster_files, terra::rast)

raster_sprc <- terra::sprc(rasters)

#raster_mosaic<- terra::mosaic(raster_sprc) 

raster_merge <- terra::merge(raster_sprc)
##In areas where the SpatRasters overlap, the values of the SpatRaster that is 
##first in the sequence of arguments (or in the SpatRasterCollection) will be retained 
plot(raster_merge)
raster_merge
final_raster <- mask(raster_merge, country_shp)
plot(final_raster)

writeRaster(final_raster, filename = "I:/VoCC/Data/Output/RemoveBuffer_BE/Final_FVoCC_BE_Euclidean_25m_ssp370_1e5grid.tif", overwrite=TRUE)

Does anyone have some ideas how I could simplify these R scripts?

I would like to get some ideas about how I can avoid loading my data line by line manually. The number assigned to my raster tiles is not in a sequence, you can see the number for the raster tile can be 'A10' and then skip to 'A100', some numbers are being removed because of the raster tiles has no data values.

Based on Till's answer, I modified my codes but met some problems when I read in the data. In my output, the 'vect_list' and 'rast_list' are both empty and shown as 'List of 0'. Anyone has a clue why this happened?

   library(purrr)
library(terra)

# read data ---------------------------------------------------------------

# these patterns will be used to select files in your directory you are
# looking to use. the patterns might not be specific enough to
# avoid overselecting in your folders

vect_files_to_keep <- c("A10", "A100", "A101")
rast_files_to_keep <- c("/10", "/100", "/101")

vect_list <-
  dir("I:\\France_grid_1e5", pattern = "\\.shp") |>
  keep(\(x) x %in% vect_files_to_keep) |>
  map(vect)

rast_list <-
  dir("I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/", pattern = "\\.tif") |>
  keep(\(x) x %in% rast_files_to_keep) |>
  map(rast)  

I modified part of the scripts again and now it can work for me. Here is a small part of my codes below:

library(terra)

setwd("I:/France_grid_1e5/")
# Define the filenames to keep
vect_files_to_keep <- c("A10.shp", "A100.shp", "A101.shp")

# Get the list of shapefiles
shapefile_list <- dir("I:/France_grid_1e5/", pattern = "\\.shp")

# Filter the list of shapefiles based on vect_files_to_keep
filtered_shapefile_list <- shapefile_list[grepl(paste(vect_files_to_keep, collapse = "|"), shapefile_list)]

# Read the shapefiles using the vect function from terra
vect_list <- lapply(filtered_shapefile_list, vect)

Solution

  • You can serialize the repeated tasks your script executes. The purrr package is great for that.

    Here is some code, that should give you a starting point. I don't know your exact folder structure, so you might have to do some adjustments.

    library(purrr)
    library(terra)
    
    # read data ---------------------------------------------------------------
    
    # these patterns will be used to select files in your directory you are
    # looking to use. the patterns might not be specific enough to
    # avoid overselecting in your folders
    
    vect_files_to_keep <- c("A10", "A100", "A101")
    rast_files_to_keep <- c("/10", "/100", "/101")
    
    vect_list <-
      dir("path/to/data", pattern = "\\.shp") |>
      keep(\(x) x %in% files_to_keep) |>
      map(vect)
    
    rast_list <-
      dir("path/to/data", pattern = "\\.tif") |>
      keep(\(x) x %in% rast_files_to_keep) |>
      map(rast)
    
    # plot --------------------------------------------------------------------
    
    map(vect_list, plot)
    map(rast_list, plot)
    
    # mask --------------------------------------------------------------------
    
    masked_list <-
      map2(vect_list,
           rast_list,
           mask)
    
    # merge -------------------------------------------------------------------
    
    merged <-
      masked_list |>
      raster_sprce() |>
      terra::merge()
    
    plot(merged)
    
    # final -------------------------------------------------------------------
    
    country_shp <- vect("path/to/data/France.shp")
    final_raster <- mask(merged, country_shp)
    plot(final_raster)