0

我正在对一个非常大的栅格（437760000 个像元）进行一些后期处理，使用相同范围/crs 的其他栅格层进行约束。代码大部分都在工作，但我遇到了一个问题。

r1[r2== 6 & r3>= 40 & r3<= 60] <- sample(2:4, length(r1[r2== 6 & r3>= 40 & r3 <= 60]), replace = T)

其中 r1、r2 和 r3 是唯一的栅格图层。r1 正在根据约束进行更新，目的是改进地图。

此代码执行没有问题，但在完成时会引发以下警告：

Warning message:
In .local(x, i, j = j, ..., value) :
  the first replacement value is used for all cells

我想确保所有三个值都是随机选择的（最终我想使用 sample 中的 prob 参数来加权其中一个值）。我尝试了许多修复程序，它们都抛出了相同的警告消息，我认为这意味着只有三个值中的一个被应用于整个栅格。我正在为此工作。

有什么想法吗？谢谢！

4

2 回答 2

1

这是 Bastien 的两种替代（但相似）解决方案。他们避免使用values对于非常大的数据集可能会出现问题的方法。

library(terra)
set.seed(123)
r1 <- rast(matrix(round(runif(400, 0, 100)), 20, 20))
r2 <- rast(matrix(round(runif(400, 0, 10)), 20, 20))
r3 <- rast(matrix(round(runif(400, 30, 70)), 20, 20))

#1，使用lapp

x <- c(r1, r2, r3)
z <- lapp(x, function(r1, r2, r3) {
    i <- r2== 6 & r3>= 40 & r3<= 60
    r1[i] <- sample(2:4, sum(i), replace = T)
    r1
})

#2，使用global（而不是长度，取栅格的全局总和与TRUE/FALSE值）

i <- r2== 6 & r3>= 40 & r3<= 60
n <- unlist(global(i, "sum"))
r1[i] <- sample(2:4, n, replace = T)
于 2021-01-11T09:34:40.280 回答
1

这是您的问题的可重现示例：

library(terra)

set.seed(123)
r1 <- rast(matrix(round(runif(400, 0, 100)), 20, 20))
plot(r1)

r2 <- rast(matrix(round(runif(400, 0, 10)), 20, 20))
r3 <- rast(matrix(round(runif(400, 30, 70)), 20, 20))

即使我无法重现您的警告代码，我认为您的问题在于您对这个调用的解释：r2== 6 & r3>= 40 & r3 <= 60. 这条线产生一个栅格：

r2== 6 & r3>= 40 & r3 <= 60
class       : SpatRaster 
dimensions  : 20, 20, 1  (nrow, ncol, nlyr)
resolution  : 0.05, 0.05  (x, y)
extent      : 0, 1, 0, 1  (xmin, xmax, ymin, ymax)
coord. ref. :  
source      : memory 
name        : lyr.1 
min value   :     0 
max value   :     1

因此，进行此调用r1[r2== 6 & r3>= 40 & r3 <= 60]会产生一个数据框：

str(r1[r2== 6 & r3>= 40 & r3 <= 60])
'data.frame':   17 obs. of  1 variable:
 $ lyr.1: num  2 2 2 2 2 2 2 2 2 2 ...

您不希望这样，因为length1 列 data.frame = 1 并且因为您不能使用 data.frame 进行值替换。

试试这个：

pixel_to_change <- values(r2== 6 & r3>= 40 & r3<= 60) == 1
r1[pixel_to_change] <- sample(2:4, sum(pixel_to_change), replace = T)

这可能是您正在寻找的东西。

于 2021-01-07T13:45:01.240 回答