RStudio中对rasterstack执行遇0重置的逐像素条件求和问题

错误原因 #

• 你调用overlay时传入的Fun_Duration不符合函数要求：raster包的overlay/calc等栅格计算函数的fun参数要求传入接收向量输入、返回等长向量的向量化函数，你写的函数是直接操作RasterStack对象的，完全不匹配参数要求，所以触发not vectorized报错。
• 循环内部的overlay逻辑错误：overlay传入匿名函数的参数是对应位置的像素值向量，不存在层的索引结构，你在匿名函数里用g[[i-1]]、ff[[i]]对栅格层做索引的逻辑完全无效。
• 最后一行冗余的b <- overlay(jj, fun = Fun_Duration)是直接触发报错的代码，你已经提前调用了Fun_Duration(jj)，不需要再套一层overlay。

正确实现方案 #

你的需求本质是逐像素计算连续1的累计长度，遇到0就重置，直接用calc函数即可实现，calc会自动遍历每个像素的时间序列向量，传入你定义的处理函数，完整可运行代码如下：

library(raster)
# 测试数据生成逻辑不变
r1 <- raster(ncol=100,nrow=100)
r2 <- raster(ncol=100,nrow=100)
r3 <- raster(ncol=100,nrow=100)
r4 <- raster(ncol=100,nrow=100)
r5 <- raster(ncol=100,nrow=100)

r1[] <- round(runif(n=ncell(r1), min=0, max = 1))
r2[] <- round(runif(n=ncell(r1), min=0, max = 1))
r3[] <- round(runif(n=ncell(r1), min=0, max = 1))
r4[] <- round(runif(n=ncell(r1), min=0, max = 1))
r5[] <- round(runif(n=ncell(r1), min=0, max = 1))

jj <- stack(r1, r2, r3, r4, r5)

# 定义单个像素时间序列的处理函数
consec_1_sum <- function(x) {
  current_count <- 0
  out <- vector(length = length(x))
  for (i in seq_along(x)) {
    if (x[i] == 0) {
      current_count <- 0
    } else {
      current_count <- current_count + 1
    }
    out[i] <- current_count
  }
  return(out)
}

# 直接用calc计算，返回结果是和输入层数一致的RasterBrick
b <- calc(jj, fun = consec_1_sum)

结果验证 #

可以随机抽取像素验证计算逻辑是否符合要求：

# 抽取第10个像素的原始序列
raw_vals <- as.vector(jj[10])
# 抽取第10个像素的计算后序列
calc_vals <- as.vector(b[10])
print(raw_vals)
print(calc_vals)

内容的提问来源于stack exchange，提问作者Johnny_Basar