使用 R 中的差异和联合溶解重叠多边形答案

【问题标题】：Dissolve Overlapping Polygons using difference and union in R使用 R 中的差异和联合溶解重叠多边形
【发布时间】：2019-04-03 06:49:00
【问题描述】：

拥有包含多个多边形的形状文件，这些多边形具有区域和地块的逻辑划分。图在区域上重叠。任务是解散/合并带有无重叠区域的地块。

这是形状文件的 spplot。这里的地块位于 Field Zones 的顶部。这里还有重叠多边形（区域和绘图）的 shapefile：Shapefile

在 QGIS 中，使用提取区域和绘图，找到差异然后使用 Union 溶解来实现相同的目的。现在需要在 R 中进行相同的编程。

在 R 中尝试了以下步骤，但无法获得正确的结果，正在寻找如何在 R 中溶解这种类型的重叠多边形：

library(sp);
library(raster);
library(rgeos)

#Importing the shape files

field_boundary_fp <- "Database/gadenstedt2_outer_field 3 -26_0_zoned-plotm.shp"
poly_map_proj_str <- "+proj=longlat +datum=WGS84 +no_defs";
utm_target_proj   <- "+init=epsg:32632";

field_boundary_sdf <- maptools::readShapePoly(fn = field_boundary_fp,
                                          proj4string =  CRS(poly_map_proj_str),
                                          repair = T,
                                          delete_null_obj = T,
                                          verbose = T);
spplot(
field_boundary_sdf,"Rx"
)

# Extracting the Zones and Plots#

Zone_sdf <- field_boundary_sdf[field_boundary_sdf@data$Type == "Zone", ]
Plot_sdf <- field_boundary_sdf[field_boundary_sdf@data$Type == "Plot", ]
plot(Plot_sdf)
plot(Zone_sdf)

#Finding the Intersection Part between the both
test <- gIntersection(Zone_sdf, Plot_sdf,id="ZoneIdx")
plot(test)
plot(test, add = T, col = 'blue')

# Finding the difference

test2 <- gDifference(Zone_sdf,Plot_sdf,id="ZoneIdx")
plot(test2)
plot(test2, add = T, col = 'red')

#Trying for Union then
polygon3 <- gUnion(test2, Plot_sdf,id="ZoneIdx")
plot(polygon3)
plot(polygon3, add = T, col = 'yellow')

【问题讨论】：

我们没有你的数据，所以我们不能运行你的代码，所以我们看不出哪里出了问题。请更新您的问题。
感谢@Spacedman 的回复。已上传包含多边形和相关数据的 shapefile。链接：filedropper.com/shapefile

标签： r r-raster sp rgeo-shapefile

【解决方案1】：

读取 shapefile

library(raster)
fields <- shapefile("gadenstedt2_outer_field 3 -26_0_zoned-plotm.shp")

首先将区域和字段分开。

zone <- fields[fields$Type == "Zone", ]
plot <- fields[fields$Type == "Plot", ]

从区域中删除绘图

d <- erase(zone, plot)

然后将plot 附加到d

r <- bind(plot, d)

现在聚合

rd <- aggregate(r, "Rx")
spplot(rd, "Rx")

---- 现在用一个可重现的例子，让其他人也可以受益；你不应该问依赖于需要下载的文件的问题----

示例数据（两个 SpatialPolygonDataFrame 对象）

library(raster)
p <- shapefile(system.file("external/lux.shp", package="raster"))
p <- aggregate(p, "NAME_1")
p$zone <- 10 + (1:length(p))
r <- raster(ncol=2, nrow=2, vals=1:4, ext=extent(6, 6.4, 49.75, 50), crs=crs(p))
names(r) <- "zone"
b <- as(r, 'SpatialPolygonsDataFrame')

擦除和追加

e <- erase(p, b)
pb <- bind(e, b)

data.frame(pb)
#        NAME_1 zone
#1     Diekirch   11
#2 Grevenmacher   12
#3   Luxembourg   13
#4         <NA>    1
#5         <NA>    2
#6         <NA>    3
#7         <NA>    4

【讨论】：

谢谢罗伯特的回答。 Intersect 仅给出地块和区域相交地块。还需要剩余的区域部分，以便解散地块和区域
再次感谢罗伯特，但添加的区域没有像绘图那样的数据值 - 现场生物质区域识别需要 Rx 值
再次感谢罗伯特的更新。在将答案设为正确之前，我继续分析数据，发现绘图的“Rx”值现在与区域相同。在解散地块和区域后，地块“Rx”值应该保持不变，就像我们现在对区域一样。你能帮忙吗？
我简化了，现在它做到了。很难准确理解您的问题（您对“溶解”一词的使用不正确）。
哦，所以你确实想要一个dissolve。目前的解决方案如何？

【解决方案2】：

为确保解决方案适用于所有字段，在上述解决方案中添加了额外的代码行以添加缓冲区几何。

fields <- gBuffer(fields, byid=TRUE, width=0) # Expands the given geometry to include 
the area within the specified width 

zone <- fields[fields$Type == "Zone", ]
plot <- fields[fields$Type == "Plot", ]

d <- erase(zone, plot)
spplot(d, "Rx")

r <- bind(plot, d)

rd <- aggregate(r, "Rx")

spplot(rd, "Rx")

【讨论】：