查找多线串和多边形Sf r之间的交点 [英] find intersection between multilinestring and polygons sf r
本文介绍了查找多线串和多边形Sf r之间的交点的处理方法，对大家解决问题具有一定的参考价值，需要的朋友们下面随着小编来一起学习吧！
问题描述
我在地图上有一组空间坐标，并有一个穿过它们的多行字符串。我需要计算出这条线在每种颜色的多边形中花费的长度。
注意事项：
- 我的最小表示法不是很好，我的实际数据是一个数据。Frame有一个颜色列(在真实数据中是组)，一些标识列，以及一列多边形几何图形。我真的不知道如何制作，所以我的Reprex只有3个单独的多边形。
- 有时两种颜色的多边形会重叠，在这种情况下，我希望将这条线与重叠形状相交的长度与该线只与一个形状相交时的长度区分开来。
最小Reprex：
poly1 <-
# create list of matrices and the first point same as last point
list(
matrix(
c(0, 0,
4, 0,
5, 1,
4, 2,
3, 2,
1, 1,
0, 0),
ncol=2, byrow=T
)
)
poly1 <- sf::st_polygon(poly1)
poly2 <-
# create list of matrices and the first point same as last point
list(
matrix(
c(4, 1,
7, 0,
7, 1,
6, 3,
3, 2,
2, 1,
4, 1),
ncol=2, byrow=T
)
)
poly2 <- sf::st_polygon(poly2)
poly3 <-
# create list of matrices and the first point same as last point
list(
matrix(
c(7, 1,
10, 1,
12, 2,
11, 4,
8, 3,
7, 2,
7, 1),
ncol=2, byrow=T
)
)
poly3 <- sf::st_polygon(poly3)
line <-
# create list of matrices and the first point same as last point
list(
matrix(
c(0, 1,
2, 0,
4, 1,
6, 3,
8, 2,
10, 1,
12, 1),
ncol=2, byrow=T
)
)
line <-
sf::st_multilinestring(line)
ggplot() +
geom_sf(data = poly1, fill = "green",alpha=.5) +
geom_sf(data = poly2, fill = "blue",alpha=.5) +
geom_sf(data = poly3, fill = "green", alpha=.5)+
geom_sf(data = line, color="black", size=2) +
ggthemes::theme_map()
所需输出：
- 只裁剪到直线的st多边形的可视化表示(我使用的是st_cross())，如下所示：
poly1_cropped <- st_intersection( line, poly1)
poly2_cropped <- st_intersection( line, poly2)
poly3_cropped <- st_intersection( line, poly3)
ggplot() +
geom_sf(data = poly1, fill = "green",alpha=.5) +
geom_sf(data = poly2, fill = "blue",alpha=.5) +
geom_sf(data = poly3, fill = "green", alpha=.5)+
geom_sf(data = line, color="black", size=2) +
ggthemes::theme_map() +
geom_sf(data = poly1_cropped, color="green", size=2) +
geom_sf(data = poly2_cropped, color="blue", size=2) +
geom_sf(data = poly3_cropped, color="green", size=2)
然后是量化线条穿过每个形状的数据框，类似于：
shape | color | shapes_overlapping | length
poly1 green 0 3
poly1 green 1 .5
poly2 blue 1 .5
poly2 blue 0 2
poly3 green 0 2.5
推荐答案
我有一个适用于您的示例的解决方案，您可以根据实际数据对其进行修改。
我做的第一件事是将所有3个POLYGON对象合并为一个Simple feature collection，最初使用一个MULTIPOLYGON，但然后使用st_cast()拆分成每个单独的POLYGON。
# combine polygons into a single simple feature collection
c(poly1, poly2, poly3) %>%
# make into simple feature
st_sfc %>%
st_sf %>%
# split into individual polygons
st_cast('POLYGON') %>%
{. ->> int1}
int1
# Simple feature collection with 3 features and 0 fields
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 0 ymin: 0 xmax: 12 ymax: 4
# CRS: NA
# geometry
# 1 POLYGON ((0 0, 4 0, 5 1, 4 ...
# 2 POLYGON ((4 1, 7 0, 7 1, 6 ...
# 3 POLYGON ((7 1, 10 1, 12 2, ...
st_intersection()，但与您的示例中找到两个几何图形之间的交集不同，我们将其应用于单个sf对象，该对象返回其自交，包括n.overlaps(层数)和origins(重叠区域中的原始多边形)。有关更多详细信息，请参阅sf页here。我们还创建了一个新的id列，该列对于每个多边形和层数都是唯一的。
# now, we need to find if any polygons overlap, and how many layers there are
int1 %>%
st_intersection %>%
filter(
st_geometry_type(.) == 'POLYGON'
) %>%
mutate(
id = row_number()
) %>%
{. ->> int2}
int2
# Simple feature collection with 4 features and 3 fields
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 0 ymin: 0 xmax: 12 ymax: 4
# CRS: NA
# n.overlaps origins geometry id
# 1 1 1 POLYGON ((4 0, 0 0, 1 1, 3 ... 1
# 1.1 2 1, 2 POLYGON ((4 2, 5 1, 4.75 0.... 2
# 2 1 2 POLYGON ((7 1, 7 0, 4.75 0.... 3
# 3 1 3 POLYGON ((7 2, 8 3, 11 4, 1... 4
为了演示这一点，我们可以根据id列绘制int2和颜色(fill)每个多边形。
int2_plot <- int2 %>%
ggplot()+
geom_sf(aes(fill = factor(id)))+
geom_sf(data = line)
int2_plot
ggplot2::ggplot_build()提取绘图的组成部分-特别是对应于每个多边形的颜色id。我们使用plotrix::color.id()将十六进制颜色代码转换为颜色名称。这可能是必需的，也可能不是必需的，具体取决于您如何解释颜色。
# now, get colours of each polygon from ggplot
int2_plot %>%
ggplot_build %>%
.$data %>%
.[[1]] %>%
tibble %>%
select(fill, group, id = group) %>%
rowwise %>%
mutate(
colour = plotrix::color.id(fill),
) %>%
select(id, colour) %>%
{. ->> plot_colours}
plot_colours
# # A tibble: 4 x 2
# # Rowwise:
# id colour
# <int> <chr>
# 1 1 salmon
# 2 2 chartreuse3
# 3 3 turquoise3
# 4 4 mediumpurple1
现在我们求出line与每个多边形相交的长度。我们使用st_intersection()找出与多边形重叠的LINESTRING对象，然后st_length()计算其长度。因为这个示例数据没有CRS，所以这个长度是无单位的，并作为numeric向量返回。如果您的数据包含CRS，st_length()将返回一个units向量，其中包含数字和距离单位，如3.726 [m]；可以使用as.numeric()检索数字。
# now, what's the length of line intersecting each polygon (colour)?
int2 %>%
mutate(
line_overlap = st_intersection(int2, line) %>% st_length
) %>%
{. ->> int3}
int3
# Simple feature collection with 4 features and 4 fields
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 0 ymin: 0 xmax: 12 ymax: 4
# CRS: NA
# n.overlaps origins geometry id line_overlap
# 1 1 1 POLYGON ((4 0, 0 0, 1 1, 3 ... 1 3.7267800
# 1.1 2 1, 2 POLYGON ((4 2, 5 1, 4.75 0.... 2 0.7071068
# 2 1 2 POLYGON ((7 1, 7 0, 4.75 0.... 3 2.1213203
# 3 1 3 POLYGON ((7 2, 8 3, 11 4, 1... 4 2.9814240
origins)。要为这些面中的每个面创建单独的行，我们将origins转换为list，然后可以使用unnest()将其拆分为每个列表元素(原始面)的单独行。
# now, duplicate rows where n.overlaps > 1 - so there's a row for each polygon in the overlap
int3 %>%
rowwise %>%
mutate(
origins = list(origins)
) %>%
tidyr::unnest(cols = c('origins')) %>%
{. ->> int4}
int4
# Simple feature collection with 5 features and 4 fields
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 0 ymin: 0 xmax: 12 ymax: 4
# CRS: NA
# # A tibble: 5 x 5
# n.overlaps origins geometry id line_overlap
# <int> <int> <POLYGON> <int> <dbl>
# 1 1 1 ((4 0, 0 0, 1 1, 3 2, 2 1, 4 1, 4.75 0.75, 4 0)) 1 3.73
# 2 2 1 ((4 2, 5 1, 4.75 0.75, 4 1, 2 1, 3 2, 4 2)) 2 0.707
# 3 2 2 ((4 2, 5 1, 4.75 0.75, 4 1, 2 1, 3 2, 4 2)) 2 0.707
# 4 1 2 ((7 1, 7 0, 4.75 0.75, 5 1, 4 2, 3 2, 6 3, 7 1)) 3 2.12
# 5 1 3 ((7 2, 8 3, 11 4, 12 2, 10 1, 7 1, 7 2)) 4 2.98
最后，我们加入前面拉出的每个多边形的颜色。
# now, join in the colours of each polygon
int4 %>%
left_join(plot_colours) %>%
{. ->> int5}
int5
# Simple feature collection with 5 features and 5 fields
# Geometry type: POLYGON
# Dimension: XY
# Bounding box: xmin: 0 ymin: 0 xmax: 12 ymax: 4
# CRS: NA
# # A tibble: 5 x 6
# n.overlaps origins geometry id line_overlap colour
# <int> <int> <POLYGON> <int> <dbl> <chr>
# 1 1 1 ((4 0, 0 0, 1 1, 3 2, 2 1, 4 1, 4.75 0.75, 4 0)) 1 3.73 salmon
# 2 2 1 ((4 2, 5 1, 4.75 0.75, 4 1, 2 1, 3 2, 4 2)) 2 0.707 chartreuse3
# 3 2 2 ((4 2, 5 1, 4.75 0.75, 4 1, 2 1, 3 2, 4 2)) 2 0.707 chartreuse3
# 4 1 2 ((7 1, 7 0, 4.75 0.75, 5 1, 4 2, 3 2, 6 3, 7 1)) 3 2.12 turquoise3
# 5 1 3 ((7 2, 8 3, 11 4, 12 2, 10 1, 7 1, 7 2)) 4 2.98 mediumpurple1
如果您不想将结果作为sf对象，我们可以使用st_drop_geometry()仅返回属性组件。
# keep the attribute data only
int5 %>%
st_drop_geometry
# # A tibble: 5 x 5
# n.overlaps origins id line_overlap colour
# * <int> <int> <int> <dbl> <chr>
# 1 1 1 1 3.73 salmon
# 2 2 1 2 0.707 chartreuse3
# 3 2 2 2 0.707 chartreuse3
# 4 1 2 3 2.12 turquoise3
# 5 1 3 4 2.98 mediumpurple1
这篇关于查找多线串和多边形Sf r之间的交点的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
查看全文
