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

huangapple go评论74阅读模式
英文:

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 repetition.

I would like to know if it is possible to simplify my scripts to avoid loading data 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)
英文:

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(&quot;I:/France_grid_1e5&quot;)

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

country_shp &lt;- vect(&quot;I:/VoCC/Data/France_shapefiles/France.shp&quot;)

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

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

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


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

raster_files &lt;- list.files(&quot;I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/&quot;, pattern = &quot;.tif$&quot;, full.names = TRUE)

rasters &lt;- lapply(raster_files, terra::rast)

raster_sprc &lt;- terra::sprc(rasters)

#raster_mosaic&lt;- terra::mosaic(raster_sprc) 

raster_merge &lt;- 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 &lt;- mask(raster_merge, country_shp)
plot(final_raster)

writeRaster(final_raster, filename = &quot;I:/VoCC/Data/Output/RemoveBuffer_BE/Final_FVoCC_BE_Euclidean_25m_ssp370_1e5grid.tif&quot;, 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 &lt;- c(&quot;A10&quot;, &quot;A100&quot;, &quot;A101&quot;)
rast_files_to_keep &lt;- c(&quot;/10&quot;, &quot;/100&quot;, &quot;/101&quot;)

vect_list &lt;-
  dir(&quot;I:\\France_grid_1e5&quot;, pattern = &quot;\\.shp&quot;) |&gt;
  keep(\(x) x %in% vect_files_to_keep) |&gt;
  map(vect)

rast_list &lt;-
  dir(&quot;I:/VoCC/Data/Output/France_1e5grid_5kmbuffer/&quot;, pattern = &quot;\\.tif&quot;) |&gt;
  keep(\(x) x %in% rast_files_to_keep) |&gt;
  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(&quot;I:/France_grid_1e5/&quot;)
# Define the filenames to keep
vect_files_to_keep &lt;- c(&quot;A10.shp&quot;, &quot;A100.shp&quot;, &quot;A101.shp&quot;)

# Get the list of shapefiles
shapefile_list &lt;- dir(&quot;I:/France_grid_1e5/&quot;, pattern = &quot;\\.shp&quot;)

# Filter the list of shapefiles based on vect_files_to_keep
filtered_shapefile_list &lt;- shapefile_list[grepl(paste(vect_files_to_keep, collapse = &quot;|&quot;), shapefile_list)]

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

答案1

得分: 0

你可以序列化你的脚本执行的重复任务。purrr 包非常适合这个任务。

以下是一些代码,应该可以为你提供一个起点。我不知道你的确切文件夹结构,所以你可能需要做一些调整。

library(purrr)
library(terra)

# 读取数据 ---------------------------------------------------------------

# 这些模式将用于选择你的目录中要使用的文件。这些模式可能不够具体,可能会在你的文件夹中选择过多的文件。

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)

# 绘图 --------------------------------------------------------------------

map(vect_list, plot)
map(rast_list, plot)

# 掩膜 --------------------------------------------------------------------

masked_list <-
  map2(vect_list,
       rast_list,
       mask)

# 合并 -------------------------------------------------------------------

merged <-
  masked_list |>
  raster_sprce() |>
  terra::merge()

plot(merged)

# 最终 -------------------------------------------------------------------

country_shp <- vect("path/to/data/France.shp")
final_raster <- mask(merged, country_shp)
plot(final_raster)

希望这可以帮助你开始处理你的任务。

英文:

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 &lt;- c(&quot;A10&quot;, &quot;A100&quot;, &quot;A101&quot;)
rast_files_to_keep &lt;- c(&quot;/10&quot;, &quot;/100&quot;, &quot;/101&quot;)

vect_list &lt;-
  dir(&quot;path/to/data&quot;, pattern = &quot;\\.shp&quot;) |&gt;
  keep(\(x) x %in% files_to_keep) |&gt;
  map(vect)

rast_list &lt;-
  dir(&quot;path/to/data&quot;, pattern = &quot;\\.tif&quot;) |&gt;
  keep(\(x) x %in% rast_files_to_keep) |&gt;
  map(rast)

# plot --------------------------------------------------------------------

map(vect_list, plot)
map(rast_list, plot)

# mask --------------------------------------------------------------------

masked_list &lt;-
  map2(vect_list,
       rast_list,
       mask)

# merge -------------------------------------------------------------------

merged &lt;-
  masked_list |&gt;
  raster_sprce() |&gt;
  terra::merge()

plot(merged)

# final -------------------------------------------------------------------

country_shp &lt;- vect(&quot;path/to/data/France.shp&quot;)
final_raster &lt;- mask(merged, country_shp)
plot(final_raster)

huangapple
  • 本文由 发表于 2023年7月17日 23:32:28
  • 转载请务必保留本文链接:https://go.coder-hub.com/76706029.html
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定