【问题标题】:Rasters not written to list when parallel procesing并行处理时未写入列表的栅格
【发布时间】:2021-09-30 08:34:52
【问题描述】:

我可能做错了,但这就是我想要做的...... 我有两个需要纠正的大栅格。由于 R 只在一个核心上工作,我想使用并行处理。因此,我将两个栅格拆分为较小的栅格,我想在这些栅格上执行我的任务。之后我会将较小的栅格合并为一个。

我的问题是,在foreach 循环中,应该创建新对象 - 我将它们保存到列表中。但是列表仍然是空的。

library(raster)
library(rgdal)
library(SpaDES)
library(doParallel)


dmp <- raster("exampe1.tif"))
luc <- raster("example2.sdat"))

UseCores <- parallel::detectCores()-1  #Define how many cores you want to use
cl       <- makeCluster(UseCores)      #Register CoreCluster
registerDoParallel(cl)

dmp_s <- splitRaster(dmp, nx= 10, ny = 10)   #Split raster to smaller rasters
luc_s <- splitRaster(luc, nx= 10, ny = 10)   #Split raster to smaller rasters

dmp_s_copy <- dmp_s
luc_s_copy <- luc_s


dmp1 <- vector(mode = 'list', length = length(dmp_s))
m <- vector(mode = 'list', length = length(dmp_s))
dmp2 <- vector(mode = 'list', length = length(dmp_s))
dmp_list <- vector(mode = 'list', length = length(dmp_s))
dmp_stack <- vector(mode = 'list', length = length(dmp_s))
dmp_final1 <- vector(mode = 'list', length = length(dmp_s))

foreach(i=1:length(dmp_s)) %dopar% {
  dmp1[i] <- dmp_s[[i]] * luc_s[[i]]
  
  m[[i]] <- matrix(1, ncol=5, nrow=5)
  dmp2[[i]] <- raster::focal(dmp1[[i]], w=m[[i]], fun=max, na.rm=TRUE, pad=TRUE, NAonly=TRUE)
  
  dmp_list[[i]] <- list(dmp_s[[i]], dmp2[[i]])
  dmp_stack[[i]] <- stack(dmp_list[[i]])
  
  dmp_final1[[i]] <- min(dmp_stack[[i]])
  
}

dmp_final <- mergeRaster(dmp_final1)

如果我不使用并行处理,代码可以工作,但我需要它。如何让它发挥作用?

【问题讨论】:

  • 你需要研究一下 foreach 包的小插曲。 foreach 循环与for 循环根本不同,它更类似于lapply 循环。在并行处理期间,工作人员无法访问主 R 进程的全局环境。您需要重写循环体,使其具有可以传递给主 R 进程并由其组合的返回值。
  • 你也可以试试terra::focal。这可能会加快速度(而且您不必担心边缘效应)。

标签: r parallel-processing raster r-raster doparallel


【解决方案1】:

我认为这可能是因为并行集群中启动了多个列表位置。根据您的操作系统,这可能会导致问题。在这种情况下,我使用“fork”类型(仅在 linux 中)来节省一些内存,但如果您使用 Windows,只需将“FORK”替换为“PSOCK”。无论如何,我实现了两个thigs:

  1. 更快的栅格分割功能。
  2. 一个整洁的并行进程(避免可能导致内存问题的中间对象)。

代码如下:

library(raster)
library(rgdal)
library(SpaDES)
library(doParallel)
library(mapview)

# read data
dmp <- raster("./data/example1.tif")
luc <- raster("./data/example2.tif")

# split rasters (this is slow)
# dmp_s <- splitRaster(dmp, nx= 100, ny = 100)   #Split raster to smaller rasters
# luc_s <- splitRaster(luc, nx= 100, ny = 100)   #Split raster to smaller rasters

# Function to split rasters (adapted from https://stackoverflow.com/questions/29784829/r-raster-package-split-image-into-multiples)
SplitRas <- function(raster, split){
        h <- ceiling(ncol(raster)/split)
        v <- ceiling(nrow(raster)/split)
        agg <- aggregate(raster,fact=c(h,v))
        agg[] <- 1:ncell(agg)
        agg_poly <- rasterToPolygons(agg)
        names(agg_poly) <- "polis"
        r_list <- list()
        for(i in 1:ncell(agg)){
                e1 <- extent(agg_poly[agg_poly$polis==i,])
                r_list[[i]] <- crop(raster,e1)
        }
        return(r_list)
}

# Apply function to split rasters
dmp_s <- SplitRas(dmp, 10)
luc_s <- SplitRas(luc, 10)

# start cluster & process
UseCores <- parallel::detectCores()-1  #Define how many cores you want to use
cl       <- makeCluster(UseCores, type="FORK")      #Register CoreCluster
doParallel::registerDoParallel(cl)

dmp_final1 <- foreach(i=1:length(dmp_s)) %dopar% {
        x1 <- dmp_s[[i]] * luc_s[[i]]
        m <- matrix(1, ncol=5, nrow=5)
        x2 <- raster::focal(x1, w=m, fun=max, na.rm=TRUE, pad=TRUE, NAonly=TRUE)
        min(stack(list(x1, x2)))
}

# stop cluster
stopCluster(cl)

# merge products
dmp_final <- mergeRaster(dmp_final1)

# showmap
mapview(dmp_final)

【讨论】:

  • 谢谢。我在你的代码中有一个错误。我不得不用x1 &lt;- dmp_s[[i]] x1@data@values &lt;- dmp_s[[i]]@data@values * luc_s[[i]]@data@values 纠正forach 循环
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2021-03-31
  • 1970-01-01
  • 2015-06-15
  • 1970-01-01
  • 2015-10-10
相关资源
最近更新 更多