轉換rasterstack到csv中的R

Question

並行處理I具有大的柵格（> 20 GB），並希望每個轉換成CSV文件的特殊格式如下：

unique_key_column 
x_coordinate 
y_coordinate 
layer1_values 
layer2_values 

等

library(raster) 
r <- raster(nrows=10,ncols=10) 
r[] <- rnorm(10) 
stack <- stack(r,r,r,r,r) 

    #create function to convert coordinate to special format 
    # -34.9 will be 1034900000 
    # sxxxdddddd, where s= sign (-=1, +=2), x=degrees (34=034), 
    # and d = decimal (.9=900000) 

formatCoordinate <- function(x){ 
    first_part <- ifelse(x < 0 , "1","2") 
    second_part <- abs(as.integer(x)) 
    #make sure 3 part has 6 decimal places, then convert it to string 
    third_part <- substr(gsub(".+\\.","",as.character(format(round(x, 2), 
          nsmall = 6))),1,6) 
    result <- sprintf("%s%03d%s",first_part,second_part,third_part) 
    result 
} 

    #the actual processing 

stack =readAll(stack) 
names(stack) <-c("l1", "l2", "l3", "l4", "l5") 
#convert rasterStack to dataframe 
stackPoints <- as.data.frame(rasterToPoints(stack)) 
#format x and x coordinates 
colX <- formatCoordinate(stackPoints$x) 
colY <- formatCoordinate(stackPoints$y) 
#combine formatted x and y coordinates to compose a unique key 
pK <- paste0(colX, colY) 
stackPoints["key"] <- pK 
col_idx <- grep("key", names(stackPoints)) 
stackPoints <- stackPoints[, c(col_idx, (1:ncol(stackPoints))[-col_idx])] 
#write results to a csv file 
write.table(stackPoints, "r.csv", row.names=F, sep=";", dec=",", append=F) 

上面的代碼適用於小型光柵，但對於大型光盤，我無法將堆疊加載到RAM中。 有沒有辦法將我的代碼轉換爲使用並行處理？即讀取光柵並使用多核寫入csv，而不需要將光柵加載到RAM（Mac OSX 10.11和Ubuntu 14.04，每個8核）。 最好的，

Answer 1

首先你想弄清楚如何在單線程上編寫一個循環，因爲從for()到foreach()將是非常簡單的。我對RasterStack對象不熟悉，但看起來它具有可與nlayers(x)一起計數的圖層，可用x[[i]]提取。

所以首先我將編寫和調試是這樣的：

for(i in 1:nlayers(stack)){ 
    #convert layer of rasterStack to dataframe 
    layer_pts <- as.data.frame(rasterToPoints(stack[[i]])) 

    #write layer_pts to a csv file 
} 

然後foreach()容易。請記住，您需要使用raster軟件包啓動每個線程。爲了更快地合併，我推薦data.table。

library(foreach) 
library(doMC) 
library(data.table) 
registerDoMC(detectCores() - 2) # for me this is 40 - 2 = 38 
layer_list <- 
    foreach(i = 1:nlayers(stack), .packages = c('raster', 'data.table')) %dopar% { 
    #convert layer of rasterStack to data.table 
    layer_pts <- as.data.table(rasterToPoints(stack[[i]])) 
    setkey(layer_pts, x, y) # data.table can key on x and y, no synthetic key needed 
    layer_pts 
    } 

tbl_out <- Reduce(merge, layer_list) # uses keys from setkey 

# if you wanted the "key" column (but not essential) 
tbl_out[, key:= paste0(formatCoordinate(x), formatCoordinate(y)) ] 

write.csv(tbl_out, 'r.csv') 

請注意，如果內存不足，則可能必須減少使用的內核數量。例如，基於試驗和錯誤的registerDoMC(4)。