Растворить перекрывающиеся полигоны, используя разность и объединение в R

Иметь файл формы, который имеет несколько полигонов с логическим разделением зон и графиков. Участки перекрывают зоны. Задача состоит в том, чтобы растворить / объединить участки с зонами без наложения.

Вот spplot файла формы. Здесь участки находятся над полевыми зонами. Также вот шейп-файл с перекрывающимися полигонами (Зоны и участки): Шейп-файл

В QGIS то же самое было достигнуто с помощью извлечения зон и участков, нахождения разницы и последующего растворения с помощью Union.Now необходимо запрограммировать то же самое в 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')
Источник

2 ответа

Решение

Прочитать шейп-файл

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
Источник

Чтобы решение работало на всех полях, к вышеуказанному решению добавлена ​​дополнительная строка кода для добавления геометрии буфера.

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")
Источник
Другие вопросы по тегам