Я хотел бы сделать настоящий ГИС-клип в R из полигонов почвы, используя серию одиночных граничных полигонов, но я не могу найти функцию R, чтобы правильно это сделать. Он должен работать так же, как clip
функция в ArcMap ESRI. Я попробовал over
метод в sp
пакете, но он не работает для polys поверх polys.
Одно из предложений заключалось в том, чтобы использовать пакет gIntersection
in rgeos
в качестве клипа, используя следующий код:
#------------------------------------
library(rgeos)
library(maptools)
#Read layers as SpatialPolygonsDataFrame (both the same Albers projection)
Soils_poly = readShapePoly("Soils_polygons") #Note - Has 400 polygons
clipper_poly = readShapePoly("clipper_polygon") #Note - Has 1 polygon
#Try gintersection as clip
Clipped_polys = gIntersection(Clipper_Tile_poly, Soils_poly)
#-----------------------------------
Это займет 5 минут (слишком медленно) и ошибки с этим:
Ошибка в RGEOSBinTopoFunc (spgeom1, spgeom2, byid, id, drop_not_poly, "rgeos_intersection"): исключение TopologyException: не найден исходящий dirEdge в -721459.77681285271 2009506.5980877089
Я также попробовал этот код для проверки на совпадение:
gIntersects(Clipper_Tile_poly, Soils_poly)
и результат был ИСТИННЫМ. clip
Функция в ESRI ArcMap отлично работает для этих данных.
Кто-нибудь знает функцию R, чтобы правильно сделать клип на пространственных многоугольниках, используя пространственные многоугольники?
Ответы:
Подсказка @mdsummer об использовании
byid=TRUE
работает точно.Смотрите воспроизводимый пример ниже:
источник
Вы также можете использовать растровый пакет
raster::intersect(spdf1, spdf2)
. Преимущество заключается в поддержании атрибутов в случае, если у вас есть SpatialPolygonsDataFrame.Спасибо за этот вопрос, чтобы указать на это и для примера кода.
источник