删除重叠的多边形,点太多
Removing overlapping polygons, too many points
我有很多数据点。其实分太多了。 None 个点重叠,但有些点彼此非常接近。我想少一点点,但不移动任何位置。
我最终会得到尽可能多的点,但只有与任何其他点相距至少 ~5.7 公里的点。 (如果有一点点重叠也没关系--0.5km的误差是可以接受的)
我尝试 write an algorithm in R 完成此操作,但出现了很多意想不到的结果。我有一些数据覆盖地球大约 300,000 个点。我还有一些其他数据是几百万。当我执行算法时,我可以按国家/地区对数据进行细分,这可以将这些数字减少到 20,000 到 100,000 的范围内。如果点的位置无关紧要,那么我可能只会制作一个插值栅格并将其称为好但对于这个问题,我需要保持特定位置不变。
我尝试的另一件事是制作一个 0.028 度的规则网格和 运行 NNJoin 以找到最近的数据点。这比我的 R 代码好一点,但结果有点可笑,正如您可能想象的那样。
我的另一个想法是缓冲点,计算有多少点与缓冲层相交。我还在研究这个
是否有一个已经确定的方法来得出这个结果?如果有可以执行此操作的包或库,我很乐意使用 PostGIS、QGIS、Python、R。
tl;dr 如何减少密集的点簇,同时用减少的点集保持覆盖?
这是一个方法。
示例数据
x <- runif(10000, -180, 180)
y <- runif(10000, -90, 90)
pts <- cbind(x, y)
解决方案
library(raster)
# you will want a lower resolution than this
r <- raster(nrow=18, ncol=36, vals=1)
# get cell numbers
cells <- cellFromXY(r, pts)
# pick one point per cell
sel <- aggregate(pts, list(cells), function(i)i[1]) # or sample
让我们看看
plot(r)
points(pts, cex=.1)
points(sel[,2:3], pch=20, col="red")
请注意,这里使用了 lon/lat,因此跨纬度的距离并不相同。不确定这是否重要;但如果是这样,你可以改变。
之后:
有多种方法可以通过更改范围或在创建 RasterLayer 时创建偏移变化。有关更多信息,请参见 ?raster 和 ?extent。您也可以使用 shift
#add a row and a column
r1 <- raster(nrow=19, ncol=37, xmx=190, ymn=-100)
r2 <- shift(r1, -.5*xres(r1), -.5*yres(r1))
plot(as(r1, "SpatialPolygons"))
lines(as(r2, "SpatialPolygons"), col="red")
我有很多数据点。其实分太多了。 None 个点重叠,但有些点彼此非常接近。我想少一点点,但不移动任何位置。
我最终会得到尽可能多的点,但只有与任何其他点相距至少 ~5.7 公里的点。 (如果有一点点重叠也没关系--0.5km的误差是可以接受的)
我尝试 write an algorithm in R 完成此操作,但出现了很多意想不到的结果。我有一些数据覆盖地球大约 300,000 个点。我还有一些其他数据是几百万。当我执行算法时,我可以按国家/地区对数据进行细分,这可以将这些数字减少到 20,000 到 100,000 的范围内。如果点的位置无关紧要,那么我可能只会制作一个插值栅格并将其称为好但对于这个问题,我需要保持特定位置不变。
我尝试的另一件事是制作一个 0.028 度的规则网格和 运行 NNJoin 以找到最近的数据点。这比我的 R 代码好一点,但结果有点可笑,正如您可能想象的那样。
我的另一个想法是缓冲点,计算有多少点与缓冲层相交。我还在研究这个
是否有一个已经确定的方法来得出这个结果?如果有可以执行此操作的包或库,我很乐意使用 PostGIS、QGIS、Python、R。
tl;dr 如何减少密集的点簇,同时用减少的点集保持覆盖?
这是一个方法。
示例数据
x <- runif(10000, -180, 180)
y <- runif(10000, -90, 90)
pts <- cbind(x, y)
解决方案
library(raster)
# you will want a lower resolution than this
r <- raster(nrow=18, ncol=36, vals=1)
# get cell numbers
cells <- cellFromXY(r, pts)
# pick one point per cell
sel <- aggregate(pts, list(cells), function(i)i[1]) # or sample
让我们看看
plot(r)
points(pts, cex=.1)
points(sel[,2:3], pch=20, col="red")
请注意,这里使用了 lon/lat,因此跨纬度的距离并不相同。不确定这是否重要;但如果是这样,你可以改变。
之后:
有多种方法可以通过更改范围或在创建 RasterLayer 时创建偏移变化。有关更多信息,请参见 ?raster 和 ?extent。您也可以使用 shift
#add a row and a column
r1 <- raster(nrow=19, ncol=37, xmx=190, ymn=-100)
r2 <- shift(r1, -.5*xres(r1), -.5*yres(r1))
plot(as(r1, "SpatialPolygons"))
lines(as(r2, "SpatialPolygons"), col="red")