尽可能快地对具有超过100万条目的两个简单特征{sf}进行空间连接。

Question

尽可能快地对具有超过100万条目的两个简单特征{sf}进行空间连接。

6

我希望这不是太琐碎了，但我真的找不到答案，而且我对这个主题还太陌生，无法自己提出替代方案。所以问题在这里：

我有两个shapefile x和y，它们表示Sentinel2卫星图像的不同处理级别。 x包含约1,300,000个多边形/片段，完全覆盖图像范围，没有任何其他重要信息。 y包含约500个多边形，表示图像中无云区域（也覆盖大部分图像，除了一些“云洞”），以及有关使用的图像的4个列（传感器、时间...）的信息。

我正在尝试将图像信息添加到x中，以便在y覆盖x的地方。很简单吧？但我找不到不花费数天的方法。

我将x作为简单要素{sf}读入，因为使用shapefile / readOGR读取需要很长时间。我尝试了不同的y操作

当我尝试merge(x,y)时，我只能采用一个sf，因为merge不支持两个sf。将x（作为sf）和y（作为shp）合并会给我带来错误“无法分配13.0 Gb大小的向量”

所以我尝试了sf::st_join(x,y)，它支持变量都是sf，但现在已经28个小时了仍未完成。 sf::st_intersect(x,y)对于10,000个片段的子集大约需要9分钟，因此对于整个片段而言可能并不会更快。

将x分成几个小块是否可以解决整个问题，或者是否有另一种简单的解决方案？我是否可以通过我的工作空间使合并工作，或者是否没有连接那么多多边形的捷径？

非常感谢您的帮助，希望我的描述不太模糊！

我的小工作站： win 7 64位 8 GB RAM Intel i7-4790 @ 3.6 GHz

- Matthias_Stack

您可能想通过子集更新shapefile。对存在y的x进行子集操作，然后将所需信息保存到x中即可。但是，如果您展示样本数据和期望输出，那会更容易些。 - manotheshark

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- gonzalez.ivan90 · Accepted Answer

我经常遇到这种问题，正如@manotheshark2所说，我更喜欢在循环中对我的向量图层进行子集化处理。以下是我的建议：

加载您的数据

library(raster)
library(rgdal)
x <- readOGR('C:/', 'sentinelCovers')
y <- readOGR('C:/', 'cloudHoles')

为了确定哪些x多边形与y多边形相交，并在x表中创建该列，需要为其分配一个唯一的ID。

x$xyID <- NA # Answer col
y$yID <- 1:nrow(y@data) # ID col

运行一个循环，对 x 进行子集操作。

for (posX in 1:nrow(x@data)){
  pol.x <- x[posX, ]
  intX <- raster::intersect(pol.x, y)
  # x$xyID[posX] <- intX@data$yID ## Run this if there's unique y polygons
  # x$xyID[posX] <- paste0(intX@data$yID, collapse = ',') ## Run this if there's multiple y polygons
}

您可以检查是在x层还是y层运行循环更好。

x$xyID <- NA # Answer col
x$xID <- 1:nrow(x@data) # ID Col

for (posY in 1:nrow(y@data)){
  pol.y <- y[posY, ]
  intY <- tryCatch(raster::intersect(pol.y, x), finally = 'NULL')
  if (is.null(intY)) next
  x$xyID[x@data$xID %in% intY@data$xID] <- pol.y$yID
}