R 中带孔的 SpatialPolygons 栅格化
Rasterize SpatialPolygons with holes in R
我正在尝试以某种方式栅格化 SpatialPolygons 对象,以获得与栅格的每个像元重叠的多边形数量。我的多边形有洞。我遇到的问题是,如果两个孔重叠,它会在计数中包含此重叠。
我创建了一个玩具示例:
library(sp)
library(raster)
p1 <- rbind(c(-180,-20), c(-140,55), c(10, 0), c(-140,-60), c(-180,-20))
hole <- rbind(c(-150,-20), c(-100,-10), c(-110,20), c(-150,-20))
p1 <- list(p1, hole)
p2 <- rbind(c(-10,0), c(140,60), c(160,0), c(140,-55), c(-10,0))
p3 <- rbind(c(-180,-20), c(-130,55), c(0,60), c(40,5), c(15,-45), c(-180,-20))
p3 <- list(p3, hole)
pols <- spPolygons(p1, p2, p3)
pols@polygons[[1]]@Polygons[[2]]@hole
pols@polygons[[3]]@Polygons[[2]]@hole
r1 <- raster(ncol=180, nrow=180)
r <- rasterize(pols, r1, fun="count")
plot(pols, col=rgb(1,0,0,0.5))
plot(r)
我认为这可能是导致与问题相关的问题的原因 use R to rasterize shapefile with holes? , but there is no final solutions. In question rasterize ESRI shapefile with holes but FALSE hole slots 问题是孔的格式不正确,但我的情况不是这样。如果您查看孔,它们都是 TRUE,请参阅:
pols@polygons[[1]]@Polygons[[2]]@hole
pols@polygons[[3]]@Polygons[[2]]@hole
如有任何帮助,我们将不胜感激!
这是一个解决方法。首先栅格化每个多边形,然后将它们相加。通过这样做,孔将为零。
library(sp)
library(raster)
# Create a list of Spatial Polygons
pols_list <- lapply(list(p1, p2, p3), spPolygons)
# Create a blank raster
r1 <- raster(ncol = 180, nrow = 180)
# Rasterize each polygon
r_list <- lapply(pols_list, rasterize, y = r1, fun = "count")
# Create a raster stack
s <- stack(r_list)
# Calculate the sum
r2 <- sum(s, na.rm = TRUE)
# Plot r2
plot(r2)
使用包 fasterize
进行光栅化(除了速度更快)似乎可以避免这个问题:
library(fasterize)
library(sf)
r1 <- raster(ncol=180, nrow=180)
r <- fasterize::fasterize(sf::st_as_sf(pols),
r1, fun = "count")
plot(r, col = rainbow(3))
HTH!
我正在尝试以某种方式栅格化 SpatialPolygons 对象,以获得与栅格的每个像元重叠的多边形数量。我的多边形有洞。我遇到的问题是,如果两个孔重叠,它会在计数中包含此重叠。
我创建了一个玩具示例:
library(sp)
library(raster)
p1 <- rbind(c(-180,-20), c(-140,55), c(10, 0), c(-140,-60), c(-180,-20))
hole <- rbind(c(-150,-20), c(-100,-10), c(-110,20), c(-150,-20))
p1 <- list(p1, hole)
p2 <- rbind(c(-10,0), c(140,60), c(160,0), c(140,-55), c(-10,0))
p3 <- rbind(c(-180,-20), c(-130,55), c(0,60), c(40,5), c(15,-45), c(-180,-20))
p3 <- list(p3, hole)
pols <- spPolygons(p1, p2, p3)
pols@polygons[[1]]@Polygons[[2]]@hole
pols@polygons[[3]]@Polygons[[2]]@hole
r1 <- raster(ncol=180, nrow=180)
r <- rasterize(pols, r1, fun="count")
plot(pols, col=rgb(1,0,0,0.5))
plot(r)
我认为这可能是导致与问题相关的问题的原因 use R to rasterize shapefile with holes? , but there is no final solutions. In question rasterize ESRI shapefile with holes but FALSE hole slots 问题是孔的格式不正确,但我的情况不是这样。如果您查看孔,它们都是 TRUE,请参阅:
pols@polygons[[1]]@Polygons[[2]]@hole
pols@polygons[[3]]@Polygons[[2]]@hole
如有任何帮助,我们将不胜感激!
这是一个解决方法。首先栅格化每个多边形,然后将它们相加。通过这样做,孔将为零。
library(sp)
library(raster)
# Create a list of Spatial Polygons
pols_list <- lapply(list(p1, p2, p3), spPolygons)
# Create a blank raster
r1 <- raster(ncol = 180, nrow = 180)
# Rasterize each polygon
r_list <- lapply(pols_list, rasterize, y = r1, fun = "count")
# Create a raster stack
s <- stack(r_list)
# Calculate the sum
r2 <- sum(s, na.rm = TRUE)
# Plot r2
plot(r2)
使用包 fasterize
进行光栅化(除了速度更快)似乎可以避免这个问题:
library(fasterize)
library(sf)
r1 <- raster(ncol=180, nrow=180)
r <- fasterize::fasterize(sf::st_as_sf(pols),
r1, fun = "count")
plot(r, col = rainbow(3))
HTH!