如何在图形上用直线分割多边形？

Question

我已经为此苦苦挣扎了一段时间。我在这里发现了一些类似的问题，但对于这个特定问题来说，所有这些问题似乎都是多余的。我有两个部分重叠的多边形，我想彼此相减（从base_pol 中减去sel_pol）。之后，生成的多边形 (diff) 将沿垂直直线 x=20 进一步分割。这就是我卡住的地方。

首先，我使用了sf 包。这有助于我将两个多边形相减并检索得到的 diff 多边形的顶点。但是，我在sf 包中找不到任何功能来分割x=20 处生成的多边形，所以我转向另一个包lwgeom。

不幸的是，我在这里有两个主要问题。首先，st_split 函数以一种奇怪的格式返回顶点：GEOMETRYCOLLECTION (POLYGON ((5 0.55, 5 0.6, 19 0.6, 23 0.5, 20.5 0.5, 17 0.55, 5 0.55)))，我不确定如何使用它以与使用 sf 相同的直观方式沿我的代码进一步传递顶点，并且，其次，得到的顶点是错误的。它没有切割diff 多边形。

很明显，我在这里做错了什么。我也不喜欢这样一个基本操作必须使用两个独立的库的事实。如果有人能给我一些关于如何使这个过程更整洁的提示，我将非常感激。

我的可重现示例如下：

library(ggplot2)
library(sf)
library(lwgeom)

rm(list = ls())

group_prc <- data.frame(x = c(5, 23, 19, 5), y = c(0.5, 0.5, 0.6, 0.6), group = "prc")
group_exp <- data.frame(x = c(5, 24, 17, 5), y = c(0.45, 0.45, 0.55, 0.55), group = "exp")
base_pol<- Polygons(list(Polygon(group_prc[, c("x", "y")])), "base")
sel_pol <- Polygons(list(Polygon(group_exp[, c("x", "y")])), "sel")
shape <- SpatialPolygons(list(base_pol, sel_pol))

diff <- shape["base"]-shape["sel"]
diffVert <- diff@polygons[[1]]@Polygons[[1]]@coords
diffVert

dat <- rbind(group_prc, group_exp, data.frame(diffVert, group = "diff"))
dat
ggplot(dat, aes(x = x, y = y, fill = group)) + geom_polygon(alpha = 0.5) + geom_vline(xintercept = 20)

splitRes <- st_split(st_polygon(list(diff@polygons[[1]]@Polygons[[1]]@coords)), 
                     st_linestring(rbind(c(20,20), c(0, 1))))
splitRes

这是一个视觉表示。紫色区域是我在代码中称为diff 的区域，我想从diff 多边形的x=20 中获取位于左侧的所有顶点（以及最终区域）。

Answer 1

我设法通过一个技巧解决了它。不再需要lwgeom 包。基本上，我将一个足够长的矩形定义为data.frame 并将其转换为sf 兼容的SpatialPolygon 对象。

然后，有了三个SpatialPolygon 对象，我现在可以将它们彼此相减，就像我之前获得diff 时所做的那样。我从下面的代码中删除了diff，因为它只是一个临时解决方案。

library(sf)
library(rgeos)

group_prc <- data.frame(x = c(5, 23, 19, 5), y = c(0.5, 0.5, 0.6, 0.6), group = "prc")
group_exp <- data.frame(x = c(5, 24, 17, 5), y = c(0.45, 0.45, 0.55, 0.55), group = "exp")
group_lim <- data.frame(x = c(20, 20, 500, 500), y = c(0, 1, 1, 0), group = "lim")
prc_pol <- Polygons(list(Polygon(group_prc[, c("x", "y")])), "prc")
exp_pol <- Polygons(list(Polygon(group_exp[, c("x", "y")])), "exp")
lim_pol <- Polygons(list(Polygon(group_prcmax[, c("x", "y")])), "lim")
shape <- SpatialPolygons(list(prc_pol, exp_pol, lim_pol))

out_pol <-  shape["prc"]-shape["exp"]-shape["lim"]

out_pol@polygons[[1]]@Polygons[[1]]@coords

输出是这个数据框：

      x         y
[1,]  5 0.5500000
[2,]  5 0.6000000
[3,] 19 0.6000000
[4,] 20 0.5750000
[5,] 20 0.5071429
[6,] 17 0.5500000
[7,]  5 0.5500000