R语言sfnetworks包, 使用空间谓词查询边
# Sat Jul 31 00:42:27 2021 edit
# 字符编码:UTF-8
# R 版本:R 4.1 x64 for window 11
# cgh163email@163.com
# 个人笔记不负责任,拎了个梨🍐🍈
#.rs.restartR()
require()
rm(list = ls());gc()
# 使用空间谓词查询边
# edge_intersects() edge_is_disjoint() edge_touches() edge_crosses() edge_is_within() edge_contains() edge_contains_properly() edge_overlaps() edge_equals() edge_covers() edge_is_covered_by() edge_is_within_distance()
library(sf, quietly = TRUE)
library(tidygraph, quietly = TRUE)
# Create a network.
net = as_sfnetwork(roxel) %>%
st_transform(3035)
# Create a geometry to test against.
p1 = st_point(c(4151358, 3208045))
p2 = st_point(c(4151340, 3207520))
p3 = st_point(c(4151756, 3207506))
p4 = st_point(c(4151774, 3208031))
poly = st_multipoint(c(p1, p2, p3, p4)) %>%
st_cast('POLYGON') %>%
st_sfc(crs = 3035)
# Use predicate query function in a filter call.
intersects = net %>%
activate(edges) %>%
filter(edge_intersects(poly))
oldpar = par(no.readonly = TRUE)
par(mar = c(1,1,1,1))
plot(st_geometry(net, "edges"))
plot(st_geometry(intersects, "edges"), col = "red", lwd = 2, add = TRUE)
par(oldpar)
dev.copy(png, "1.png");dev.off()
# Use predicate query function in a mutate call.
net %>%
activate(edges) %>%
mutate(disjoint = edge_is_disjoint(poly)) %>%
select(disjoint)
#> # A sfnetwork with 701 nodes and 851 edges
#> #
#> # CRS: EPSG:3035
#> #
#> # A directed multigraph with 14 components with spatially explicit edges
#> #
#> # Edge Data: 851 x 4 (active)
#> # Geometry type: LINESTRING
#> # Dimension: XY
#> # Bounding box: xmin: 4150707 ymin: 3206375 xmax: 4152367 ymax: 3208565
#> from to disjoint geometry
#> <int> <int> <lgl> <LINESTRING [m]>
#> 1 1 2 FALSE (4151491 3207923, 4151474 3207946)
#> 2 3 4 FALSE (4151398 3207777, 4151390 3207727, 4151370 3207673)
#> 3 5 6 FALSE (4151408 3207539, 4151417 3207573, 4151421 3207592)
#> 4 7 8 TRUE (4151885 3206698, 4151861 3206711, 4151845 3206725, 4151…
#> 5 9 10 TRUE (4151732 3207017, 4151721 3206809)
#> 6 11 12 TRUE (4152152 3206984, 4152143 3206932, 4152147 3206923)
#> # … with 845 more rows
#> #
#> # Node Data: 701 x 1
#> # Geometry type: POINT
#> # Dimension: XY
#> # Bounding box: xmin: 4150707 ymin: 3206375 xmax: 4152367 ymax: 3208565
#> geometry
#> <POINT [m]>
#> 1 (4151491 3207923)
#> 2 (4151474 3207946)
#> 3 (4151398 3207777)
#> # … with 698 more rows
image.png
网友评论