R 中包含栅格计算 returns 错误的函数

Function containing raster calculation returns error in R

我正在尝试使用其中的栅格计算创建我自己的函数。但是当我调用函数时 returns 出现了一些错误。代码如下

BT <- function(B_1, B_2){
  tir_10 <- raster::raster(B_1)
  tir_11 <- raster::raster(B_2)
  
  l_lambda_10 <- 3.3420*10^-4*tir_10 + 0.1
  l_lambda_11 <- 3.3420*10^-4*tir_11 + 0.1
  
  BT_10 <- (1321.0789/(log(774.8853/l_lambda_10 + 1)))
  BT_11 <- (1201.1442/(log(480.8883/l_lambda_11 + 1)))
  
  return(list(BT_10 = BT_10, BT_11 = BT_11))
}

现在我尝试使用以下代码调用该函数

#Read the bands
library(raster)

a <- raster(ncol=100, nrow=100)
set.seed(2)     
values(a) = runif(10000, min=27791, max=30878)

b <- raster(ncol=100, nrow=100)
set.seed(2)     
values(b) = runif(10000, min=25686, max=28069)

BT(B_1 = a, B_2 = b)

但是它returns以下错误

Error in 3.342 * 10^-4 * tir_10 : hasValues(e2) is not TRUE

如何解决?

您必须 运行 runif 在函数内部的输入中使用来自 RasterLayers B_1 和 B_2 的最小值和最大值,这样:

BT <- function(B_1, B_2){
  tir_10 <- raster::raster(B_1)
  set.seed(2)     
  values(tir_10) = runif(10000, min=round(min(getValues(B_1)),0),max=round(max(getValues(B_1)),0))
  
  tir_11 <- raster::raster(B_2)
  set.seed(2)     
  values(tir_11) = runif(10000, min=round(min(getValues(B_2)),0), max=round(max(getValues(B_2)),0))
  
  l_lambda_10 <- 3.3420*10^-4*tir_10 + 0.1
  l_lambda_11 <- 3.3420*10^-4*tir_11 + 0.1
  
  BT_10 <- (1321.0789/(log(774.8853/l_lambda_10 + 1)))
  BT_11 <- (1201.1442/(log(480.8883/l_lambda_11 + 1)))
  
  return(list(BT_10 = BT_10, BT_11 = BT_11))
}

这两行是错误的,因为 tir_10tir_11 没有值。

tir_10 <- raster::raster(B_1)
tir_11 <- raster::raster(B_2)

您可以改为:

BT <- function(tir_10, tir_11){ 
  l_lambda_10 <- 3.3420*10^-4*tir_10 + 0.1
  l_lambda_11 <- 3.3420*10^-4*tir_11 + 0.1
  
  BT_10 <- (1321.0789/(log(774.8853/l_lambda_10 + 1)))
  BT_11 <- (1201.1442/(log(480.8883/l_lambda_11 + 1)))
  stack(list(BT_10=BT_10, BT_11=BT_11))
}

library(raster)
a <- raster(ncols=100, nrows=100)
set.seed(2)     
values(a) = runif(10000, min=27791, max=30878)
b <- setValues(a, runif(10000, min=25686, max=28069))

x <- BT(a, b)
x    
#class      : RasterStack 
#dimensions : 100, 100, 10000, 2  (nrow, ncol, ncell, nlayers)
#resolution : 3.6, 1.8  (x, y)
#extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#crs        : +proj=longlat +datum=WGS84 +no_defs 
#names      :    BT_10,    BT_11 
#min values : 298.5254, 297.9016 
#max values : 305.6353, 304.4019