使用自定义条件从栅格砖/堆栈对每一层进行栅格重新分类

问题描述 投票:0回答:0

我有一个名为“wnt”的 rasterBrick/stack,它有 8 层。我想根据相应层的均值 (mn) 和标准差 (sd) 对每一层进行重新分类。分类方案就像(如果图层的单元格值 = x)那么

x>mn + sd = 1, mn +0.5sd < x < mn+ sd = 2, mn+ 0.5sd < x < mn + 0.5sd = 3, mn+ sd < x < mn+ 0.5sd = 4, x < mn + sd = 5.

可以找到 rasterBrick/stack 'wnt' here

我期待有 8 个重新分类层的最终栅格堆栈。

以下我试过的代码。但没有得到预期的结果。

# create an empty raster brick to store the reclassified data
wnt_reclass <- brick(nrows=nlayers(wnt), ncols=ncell(wnt), 
                     xmn=extent(wnt)[1], xmx=extent(wnt)[2], 
                     ymn=extent(wnt)[3], ymx=extent(wnt)[4], 
                     crs=projection(wnt))

# loop through each layer of the raster brick, calculate the mean and standard deviation 
# for the current layer, and reclassify the data based on the current layer's mean and sd
for (i in 1:nlayers(wnt)) {
  # extract the current layer
  layer <- wnt[[i]]
  
  # calculate the mean and sd for the current layer
  mn <- mean(layer[], na.rm = TRUE)
  sd <- sd(layer[], na.rm = TRUE)
  
  # reclassify the data based on the current layer's mean and sd
  x <- layer[]
  x[x > (mn + sd)] <- 1
  x[(mn + 0.5 * sd) < x & x < (mn + sd)] <- 2
  x[(mn + 0.5 * sd) == x] <- 3
  x[(mn + sd) < x & x < (mn + 0.5 * sd)] <- 4
  x[x < (mn + sd)] <- 5
  
  # set the reclassified data for the current layer in the output raster brick
  wnt_reclass[[i]] <- setValues(layer, x)
}

# Set the names of the reclassified layers
names(wnt_reclass) <- paste0("wnt_", 1:nlayers(wnt_reclass))

# Set the extent and projection of the reclassified rasterBrick to match the original 'wnt' rasterBrick
extent(wnt_reclass) <- extent(wnt)
projection(wnt_reclass) <- projection(wnt)
r r-raster
© www.soinside.com 2019 - 2024. All rights reserved.