小编lbu*_*ett的帖子

从包含相当大数量(约20000)可能部分重叠的多边形的shapefile开始,我需要提取所有通过交叉它们的不同"边界"而产生的子多边形.

在实践中,从一些模拟数据开始:

library(tibble)
library(dplyr)
library(sf)

ncircles <- 9
rmax     <- 120
x_limits <- c(-70,70)
y_limits <- c(-30,30)
set.seed(100) 
xy <- data.frame(
  id = paste0("id_", 1:ncircles), 
  x = runif(ncircles, min(x_limits), max(x_limits)),
  y = runif(ncircles, min(y_limits), max(y_limits))) %>% 
  as_tibble()

polys <- st_as_sf(xy, coords = c(2,3)) %>% 
  st_buffer(runif(ncircles, min = 1, max = 20)) 

plot(polys[1])  

我需要派生一个包含ALL的多边形sf或sp多边形,并且只需要交叉点生成的多边形,例如:

(请注意,颜色仅用于举例说明预期结果,其中每个"不同颜色"区域是一个单独的多边形,不覆盖任何其他多边形)

我知道我可以通过一次分析一个多边形,识别并保存所有交叉点然后"擦除"那些区域形成完整的多边形并继续循环来解决问题,但这很慢.

我觉得应该有一个更有效的解决方案,但我无法弄清楚,所以任何帮助将不胜感激!(这两个sf和sp基础的解决方案,欢迎)

更新:

最后,我发现即使是"一次一个多边形",任务也远非简单!我真的在努力解决这个显而易见的"简单"问题!任何提示？即使是一个缓慢的解决方案或提示在正确的道路上开始将不胜感激!

更新2:

也许这会澄清一些事情:所需的功能类似于这里描述的功能:

https://it.mathworks.com/matlabcentral/fileexchange/18173-polygon-intersection?requestedDomain=www.mathworks.com

更新3:

我将赏金奖给了@ shuiping-chen(谢谢!),他的回答正确地解决了所提供的示例数据集上的问题.然而,"方法"一般被推广到"四倍"或"n-uple"交叉点可能的情况.如果我管理的话,我会在未来几天尝试解决这个问题并发布一个更通用的解决方案!

大家好,

我正在努力解决这个问题,希望有人能提出一个简单的解决方案.

我的目标是在多边形的范围内创建一个规则的多边形网格,但是按用户定义的角度旋转.

我知道我可以在sf使用中轻松创建一个北/南多边形网格,例如:

library(sf)
#> Linking to GEOS 3.6.2, GDAL 2.2.3, proj.4 4.9.3
inpoly <- st_read(system.file("shape/nc.shp", package="sf"))[1,] %>% 
  sf::st_transform(3857) %>% 
  sf::st_geometry()
grd <- sf::st_make_grid(inpoly, cellsize = 3000)
plot(inpoly, col = "blue")
plot(grd, add = TRUE)

我也知道我可以使用以下方法轻松旋转给定的角度:

rotang = 20
rot = function(a) matrix(c(cos(a), sin(a), -sin(a), cos(a)), 2, 2)
grd_rot <- (grd - st_centroid(st_union(grd))) * rot(rotang * pi / 180) +
  st_centroid(st_union(grd))
plot(inpoly, col = "blue")
plot(grd_rot, add = TRUE)

我的问题是,根据旋转角度,输入多边形的一般"方向"和单元格大小,旋转的网格可能不再覆盖多边形的整个范围,如下所示:

rotang = 45 …

图片比文字好,所以请看看

我拥有的是什么

• 一个RasterLayer对象(此处仅填充随机值仅用于说明目的,实际值无关紧要)
• 一个SpatialPolygons对象,里面有很多很多的多边形

您可以使用以下代码重新创建我用于图像的示例数据:

library(sp)
library(raster)
library(rgeos)

# create example raster
r <- raster(nrows=10, ncol=15, xmn=0, ymn=0)
values(r) <- sample(x=1:1000, size=150)

# create example (Spatial) Polygons
p1 <- Polygon(coords=matrix(c(50, 100, 100, 50, 50, 15, 15, 35, 35, 15), nrow=5, ncol=2), hole=FALSE)
p2 <- Polygon(coords=matrix(c(77, 123, 111, 77, 43, 57, 66, 43), nrow=4, ncol=2), hole=FALSE)
p3 <- Polygon(coords=matrix(c(110, 125, 125, 110, 67, 75, 80, 67), nrow=4, ncol=2), hole=FALSE)

lots.of.polygons <- SpatialPolygons(list(Polygons(list(p1, p2, p3), 1)))
crs(lots.of.polygons) <- crs(r) # copy crs …

我有一个从互联网上传的形状.根据计算出的质心,我想从它绘制一条50度线,并找到与轮廓相交的坐标.不知道怎么办呢？

谢谢.

脚本:

library(ggplot2)

df = read.table("E:/cloud1.txt")  #stored at https://ufile.io/zq679
colnames(df)<- c("x","y")

meanX <- mean(df$x)
meanY <- mean(df$y)

ggplot(df, aes(x, y))+ geom_point()+ geom_point(aes(meanX, meanY),colour="green",size=2)