我正在使用 Landsat 图像来计算净辐射。 在我转换的反射率栅格图层中有非常低的负值（例如：-0.000003）。我想确保我的反射率为 0 - 1，以减少我未来计算中的错误。
如何在 R 中将小于 0 的栅格值替换为“NA”。类似于栅格计算函数。我不确定如何提供一个例子，但我敢肯定你们中的一个可以帮助我，是吗？
这是我从 Bastiaanssen 等人得出的反射率方程。 (1998)
假设，pb1 表示 Landsat 波段 1 的反射率，pi = 3.14...，lb1 = 波段 1 的辐射率，ESUN = 波段 1 的外大气层值，dr = 一年中一天中地球太阳的相对距离.
#Calculate reflectivity per band. QC: Always 0 to 1
pb1 = (pi * lb1)/(ESUN1 * cos(solzen) * dr)
创建此栅格后，我要做的就是将小于 0 的 pb1 值设置为 NA。
帮助？
【问题讨论】：
使用raster 包，内存安全的方法是使用reclassify 或clamp
library(raster)
r <- raster(ncol=10, nrow=10)
values(r) <- rnorm(100)
x <- clamp(r, useValues=FALSE)
y <- reclassify(r, cbind(-Inf, 0, NA), right=FALSE)
请注意right=FALSE，不要将零单元格设置为NA。
这些方法是内存安全的，您还可以提供文件名参数，这样您之后就不需要调用writeRaster。
【讨论】：
right 不应该是 TRUE 以将 0 作为 NA 吗？ r
< 0，而不是<= 0
library(raster)
values(pb1)[values(pb1) < 0] = NA
或者，正如@jbaums 所建议的那样：
pb1[pb1 < 0] <- NA
如果要保留原始栅格对象，请记住在运行上述代码之前将原始栅格分配给新的对象名称。
【讨论】：
pb1[pb1 < 0] <- NA。 （还有@MaeAntoinette，这将修改原来的pb1。如果你也想保留原来的，那么你可以这样做，例如pb1c <- pb1; pb1c[pb1c < 0] <- NA）。
raster::clamp 是一种简单灵活的方法。可以将高于和/或低于阈值的所有内容设置为该阈值，或者通过设置useValues=FALSE 然后将高于/低于阈值的值设置为NA。例如，仅较低的值：
r <- raster(ncol=3, nrow=3)
values(r) <- seq_len(9)
x <- clamp(r, lower=3, useValues=FALSE)
values(x)
# [1] NA NA 3 4 5 6 7 8 9
仅限上限值：
x <- clamp(r, upper=6, useValues=FALSE)
values(x)
# [1] 1 2 3 4 5 6 NA NA NA
以及上下两个值：
x <- clamp(r, lower=3, upper=6, useValues=FALSE)
values(x)
# [1] NA NA 3 4 5 6 NA NA NA
注意useValues=TRUE（默认）的区别：
x <- clamp(r, lower=3, upper=6)
values(x)
# [1] 3 3 3 4 5 6 6 6 6
在本例中使用raster_2.8-19。
【讨论】：
使用terra 包，内存安全的方法是使用classify clamp 或ifel
library(terra)
r <- rast(ncol=10, nrow=10)
values(r) <- rnorm(100)
x <- clamp(r, 0, values=FALSE)
y <- classify(r, cbind(-Inf, 0, NA), right=FALSE)
z <- ifel(r < 0, NA, r)
这些方法是内存安全的，您还可以提供文件名参数，这样您之后就不需要调用writeRaster。
您可以使用app，但这有点笨拙且效率低得多
a <- app(r, function(x) { x[x<0] <- NA; x })
这样的方法非常适合使用小型数据集进行交互式探索，但它们通常不应用于“生产”代码。
r[r < 0] <- NA
【讨论】：
另一种选择是
pb1 <- raster::calc(pb1, function(x){x[x<0]<-NA; return(x)})
【讨论】：