r - 在 R 或 Grass GIS 中填充栅格孔

Question

示例数据

x <- raster(x=matrix(data=1:36, nrow=6), xmn=-1000, xmx=1000, ymn=-100, ymx=900)
x[c(8, 15, 16, 17, 22, 25, 26, 30, 31)] <- NA
plot(x)

问题

我如何(从算法上)区分栅格中的孔洞，即单元格 c(15:17, 22) 与其他非孔洞间隙(即其余空单元格)界定的区域？

这样就可以只对光栅的孔/岛区域进行操作，用自定义值填充孔等。

实际光栅有大约 30000 个孔，因此速度很重要。我对 R 和 Grass GIS 解决方案都感兴趣。非常感谢您的帮助，非常感谢!

Answer 1

这是一个没有多边形化的解决方案:(它不优雅，但它有效)。 但是您必须将洞/岛重新分类为值(即 999)，并将所有其他非岛重新分类为 NA。像这样:

x <- raster(x=matrix(rep(NA,36), nrow=6), xmn=-1000, xmx=1000, ymn=-100, ymx=900)
x[c(8, 15, 16, 17, 22, 25, 26, 30, 31)] <- 999

plot(x)

现在我使用 clump() 函数来检查是否有岛屿，该函数最酷的地方在于，它还返回这些岛屿的 ID:

#Get Islands with IDs
cl <- clump(x,directions=8)
plot(cl)

然后我根据岛屿的频率创建一个数据框(这只是为了获取每个岛屿的 ID)

freqCl <- as.data.frame(freq(cl))

#remove the (row) which corresponds to the NA values (this is important for the last step)
freqCl <- freqCl[-which(is.na(freqCl$value)),]

检查是否有岛屿接触边界:

#Check if the island touches any border and therefore isn't a "real island" (first and last column or row)

noIslandID <- c()
#First row
if(any(rownames(freqCl) %in% cl[1,])){
  eliminate   <- rownames(freqCl)[rownames(freqCl) %in% cl[1,]]
  noIslandID  <- append(noIslandID, eliminate)
}
#Last row
if(any(rownames(freqCl) %in% cl[nrow(cl),])){
  eliminate  <- rownames(freqCl)[rownames(freqCl) %in% cl[nrow(cl),]]
  noIslandID <- append(noIslandID, eliminate)
}
#First col
if(any(rownames(freqCl) %in% cl[,1])){
  eliminate   <- rownames(freqCl)[rownames(freqCl) %in% cl[,1]]
  noIslandID  <- append(noIslandID, eliminate)
}
#Last col
if(any(rownames(freqCl) %in% cl[,ncol(cl)])){
  eliminate   <- rownames(freqCl)[rownames(freqCl) %in% cl[,ncol(cl)]]
  noIslandID  <- append(noIslandID, eliminate)
}

消除那些接触边界的岛屿:

noIslandID <- unique(noIslandID)
IslandID   <- setdiff(rownames(freqCl), noIslandID)

在最后一步中，为初始栅格中的每个“真实岛屿”分配一个 1:

for(i in 1:length(IslandID)) {
  x[cl[]==as.numeric(IslandID[i])] <- 1
}

plot(x)