使用 R 中的 terra 包在多边形中点

89_*_*ple 3 r terra

我想使用 terra 包进行多边形分析中的一个点。我有一组点,我想提取它们落在哪些多边形中。下面的示例数据显示了单个多边形

library(terra)
crdref <- "+proj=longlat +datum=WGS84"
  
longitude <- c(-116.7, -120.4, -116.7, -113.5, -115.5, -120.8, -119.5, -113.7, -113.7, -110.7)
latitude <- c(45.3, 42.6, 38.9, 42.1, 35.7, 38.9, 36.2, 39, 41.6, 36.9)
lonlat <- cbind(longitude, latitude)
pts <- vect(lonlat, crs=crdref)
  
lon <- c(-116.8, -114.2, -112.9, -111.9, -114.2, -115.4, -117.7)
lat <- c(41.3, 42.9, 42.4, 39.8, 37.6, 38.3, 37.6)
lonlat <- cbind(id=1, part=1, lon, lat)
pols <- vect(lonlat, type="polygons", crs=crdref)
  
plot(pols, border='blue', col='yellow', lwd=3)
points(pts, col='red', pch=20, cex=3)
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述

terra::extract(pts, pols)
       id.y id.x
[1,]    1   NA
Run Code Online (Sandbox Code Playgroud)

但输出对我来说并不清楚。我只需要每个点属于哪个多边形。

Rob*_*ans 7

示例数据

library(terra)
longitude <- c(-116.7, -120.4, -116.7, -113.5, -115.5, -120.8, -119.5, -113.7, -113.7, -110.7)
latitude <- c(45.3, 42.6, 38.9, 42.1, 35.7, 38.9, 36.2, 39, 41.6, 36.9)
pts <- vect(cbind(longitude, latitude), crs="+proj=longlat")
  
lon <- c(-116.8, -114.2, -112.9, -111.9, -114.2, -115.4, -117.7)
lat <- c(41.3, 42.9, 42.4, 39.8, 37.6, 38.3, 37.6)
lonlat <- cbind(id=1, part=1, lon, lat)
pols <- vect(lonlat, type="polygons", crs="+proj=longlat")
Run Code Online (Sandbox Code Playgroud)

这里有四种方法

## 1
e <- extract(pols, pts) 
e[!is.na(e[,2]), 1]
#[1] 3 4 8 9
 
## 2
relate(pols, pts, "contains") |> which()
#[1] 3 4 8 9

## 3 
is.related(pts, pols, "intersects") |> which()
#[1] 3 4 8 9

## 4
pts$id <- 1:nrow(pts)
intersect(pts, pols)$id 
#[1] 3 4 8 9
Run Code Online (Sandbox Code Playgroud)

所以你的思路是正确的,extract但你的论点顺序是错误的。extract如果您有多个多边形,这可能是最简单的。

更多例子在这里