检测栅格是否在 SpatialPolygons 对象之内、之外或相交

Question

我有许多栅格，我想检查它们是否完全包含在空间多边形内、完全没有空间多边形或与空间多边形相交（这可能意味着多边形完全在栅格内，或者多边形和栅格重叠）。我正在做这项检查，以便尽可能避免时间密集型掩蔽。

这是一个例子：

                                      # create 3 example rasters
r <- raster()
r[] <- rnorm(n = ncell(r))
e1 <- extent(c(45,55,45,50))
r1 <- crop(r,e1)
e2 <- extent(c(20,25,25,30))
r2 <- crop(r,e2)
e3 <- extent(c(38,55,57,65))
r3 <- crop(r,e3)


                                        #create SpatialPolygons

x <- c(40,60)
y <- c(40,60)
m <- expand.grid(x,y)
m <- m[c(1,2,4,3),]
p1 <- Polygon(m)
p1 <- Polygons(list(p1),1)

x <- c(10,15)
y <- c(10,15)
m <- expand.grid(x,y)
m <- m[c(1,2,4,3),]
p2 <- Polygon(m)
p2 <- Polygons(list(p2),2)

x <- c(30,45)
y <- c(70,80)
m <- expand.grid(x,y)
m <- m[c(1,2,4,3),]
p3 <- Polygon(m)
p3 <- Polygons(list(p3),3)


poly <- SpatialPolygons(list(p1,p2,p3))

绘制这些：

我将分别读取每个栅格并检查它是否在 SpatialPolygons 之内、之外或与 SpatialPolygons 相交。

您认为在 R 中最有效的方法是什么？我计划并行屏蔽数千个 4mb 栅格，并希望通过此检查加快处理速度。

注意，还有这个问题：https://gis.stackexchange.com/questions/34535/detect-whether-there-is-a-spatial-polygon-in-a-spatial-extent

但是，我认为它没有提供我正在寻找的细节。例如，所有栅格都在空间多边形的范围内，但并非所有栅格都在空间多边形内。

类似 rgeos 中的函数（gIntersects、gContains）可能会很方便。我不确定这些是否最有效，或者我应该如何将栅格（或其范围）转换为 sp 对象。

谢谢！

Answer 1

您也可以为此使用gRelate。它返回一个DE-9IM 代码，描述两个几何图形的内部、边界和外部组件之间的关系。

library(rgeos)
x <- sapply(rlist, function(x) 
  gsub('[^F]', 'T', gRelate(as(extent(x), 'SpatialPolygons'), poly)))

然后您可以将字符串与感兴趣的关系进行比较。例如，我们可以定义within、disjoint 和overlaps 如下（但请注意，对于给定的关系，其他一些交集是可选的 - “内部”是 defined by GEOS 为 T*F**F***，“不相交”为 @987654331 @，并且“重叠”为T*T***T**)：

pat <- c(TFFTFFTTT='within', FFTFFTTTT='disjoint', TTTTTTTTT='overlaps')
pat[x]

##  TFFTFFTTT  FFTFFTTTT  TTTTTTTTT 
##   "within" "disjoint"  "overlaps" 

它似乎比the gContainsProperly/gIntersects approach 稍快，但@Tedward 的post 更易于理解，并且更符合GEOS 定义（尽管可能需要创建特定关系定义的能力）。

---

DE-9IM 字符串的元素按顺序表示：

1. 几何体 A 的内部是否与几何体 B 的内部相交？
2. A 的边界是否与 B 的内部相交？
3. A 的外部是否与 B 的内部相交？
4. 几何体 A 的内部是否与几何体 B 的边界相交？
5. A 的边界是否与 B 的边界相交？
6. A 的外部是否与 B 的边界相交？
7. 几何体 A 的内部是否与几何体 B 的外部相交？
8. A 的边界是否与 B 的外部相交？
9. A 的外部是否与 B 的外部相交？

Answer 2

这是我为解决问题所做的：

library(rgeos)
rlist <- list(r1,r2,r3)

lapply(rlist, function(raster) {
  ei <- as(extent(raster), "SpatialPolygons")
  if (gContainsProperly(poly, ei)) {
    print ("fully within")
  } else if (gIntersects(poly, ei)) {
    print ("intersects")
  } else {
    print ("fully without")
  }
})

如果您知道更有效的解决方案，请告诉我。