判断一系列点和一系列多边形是否相交的函数

Question

我有一个点列表和一个多边形列表，我想找到每个点所在的多边形。虽然使用 sf 包中的命令 st_intersects 很简单，我发现这里不太好，因为点和多边形都太多了，结果矩阵太大，我的电脑无法工作。

我尝试通过制作一个函数来解决这个问题，该函数可以查看所有多边形以找到每个点的正确多边形。我做了一个，但效果不佳，我不知道它有什么问题。所以我需要一些帮助！谢谢！

以下是我创建点列表的方法

pts = st_sfc(st_point(c(.5,.5)), st_point(c(1.5, 1.5)), st_point(c(2.5, 
2.5))) %>% st_sf() %>% mutate(namepoint=c(1,2,3))

Simple feature collection with 3 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.5 ymin: 0.5 xmax: 2.5 ymax: 2.5
epsg (SRID):    NA
proj4string:    NA
   namepoint        geometry
1         1 POINT (0.5 0.5)
2         2 POINT (1.5 1.5)
3         3 POINT (2.5 2.5)

以下是我如何创建多边形列表

pol_down = st_polygon(list(rbind(c(0,0), c(2,0), c(2,2), c(0,2), c(0,0))))
pol_up=pol_down+c(2,1) 
pol_add <- list(pol_down,pol_up) %>% st_sfc() %>% st_sf %>% 
mutate(namepgn=c('A','B'))

Simple feature collection with 2 features and 1 field
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 4 ymax: 3
epsg (SRID):    NA
proj4string:    NA
   namepgn                       geometry
1       A POLYGON ((0 0, 2 0, 2 2, 0 ...
2       B POLYGON ((2 1, 4 1, 4 3, 2 ...

我做的函数是

 fun_pgn <- function(x,y){
                location_m<- st_intersects(x,
                                           y,
                                           sparse = F)
                name <- pull(y,1) %>% .[location_m]
                return(name)
             }

使用命令st_intersects获取他们的关系

st_intersects(pts,pol,sparse = F)

      [,1]  [,2]
[1,]  TRUE FALSE
[2,]  TRUE FALSE
[3,] FALSE  TRUE

然后我想要的是在数据pts中添加一个变量，比如pgnname来记录每个点的多边形，

        namepoint    pgnname    geometry
1               1       "A"     POINT (0.5 0.5)
2               2       "A"     POINT (1.5 1.5)
3               3       "B"     POINT (2.5 2.5)

但是我的功能不起作用。在运行宁以下代码

pts %>% mutate(pgnname=fun_pgn(geometry,pol))

接着是

Simple feature collection with 3 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.5 ymin: 0.5 xmax: 2.5 ymax: 2.5
epsg (SRID):    NA
proj4string:    NA
  namepoint pgnname        geometry
1         1       A POINT (0.5 0.5)
2         2       B POINT (1.5 1.5)
3         3    <NA> POINT (2.5 2.5)

多边形的名称不正确！但是我的功能的关键部分运行良好。例如

> pull(pol,1) %>% .[st_intersects(pts[3,],pol,sparse = F)]
[1] "B"

所以我很困惑。当我运行我的函数时会发生什么？

Answer 1

你说你不能用st_intersects因为"the resulting matrix is too big for my computer"。最简单的解决方案是不生成密集矩阵，您可以使用 sparse = F 强制生成。没有这个参数的 st_intersects 的默认行为是生成一个稀疏列表，这就是你试图用你的函数做的。

x = st_intersects(pts, pol)
pts %>% mutate(pgnname = pol$namepgn[unlist(x)])
# Simple feature collection with 3 features and 2 fields
# geometry type:  POINT
# dimension:      XY
# bbox:           xmin: 0.5 ymin: 0.5 xmax: 2.5 ymax: 2.5
# epsg (SRID):    NA
# proj4string:    NA
#   namepoint pgnname        geometry
# 1         1       A POINT (0.5 0.5)
# 2         2       A POINT (1.5 1.5)
# 3         3       B POINT (2.5 2.5)

请注意，如果一个点可能落在多个多边形中，而您只想为该点分配一个多边形，则可以改用以下方法：

pgn = sapply(x, `[`, 1)
pts %>% mutate(pgnname = pol$namepgn[pgn])

判断一系列点和一系列多边形是否相交的函数

A function for determining whether a series of points and a series polygons intersect

r

function

sf