使用等高线在多边形层下方切割多边形

问题描述

| 我想根据高度将多边形层切成两部分（上部和下部）。多边形可能会凸出或凹入，并且切割位置可能会彼此不同。等高线的间隔为5m，这意味着我可能需要生成一个轮廓线非常密集的等高线，例如1m的间隔。关于如何做到这一点的任何想法，是在ArcGIS中还是在R中更好？ 以下是Q的运行示例：

library(sp)
library(raster)
r<-raster(ncol=100,nrow=100)
values(r)<-rep(1:100,100)
plot(r)   ### I have no idea why half of the value is negative...
p1<-cbind(c(-100,-90,-50,-100),c(60,70,30,60))
p2<-cbind(c(0,50,100,0),c(0,-25,10,0))
p1p<-polygons(list(polygon(p1,hole=T)),\"p1\")
p2p<-polygons(list(polygon(p2,\"p2\")
p<-Spatialpolygons(list(p1p,p2p),1:2)
plot(p,add=T)
segments(-90,80,20)  ##where the polygon Could be devided
segments(50,20,-30)  ##

提前谢谢〜 马可

解决方法

如果我理解正确，则可以在R中使用ѭ1包和相关的Spatial工具。 我用技巧来缓冲一条相交的线，然后从该站点生成\“ difference \”多边形： http://www.chopshopgeo.com/blog/?p=89 生成示例栅格和一个上覆的多边形。

vdata <- list(x = 1:nrow(volcano),y = 1:ncol(volcano),z = volcano)

## raw polygon data created using image(vdata); xy <- locator()

xy <- structure(list(x = c(43.4965355534823,41.7658494766076,36.2591210501883,25.560334393145,13.7602020508178,18.7949251835441,29.179041644792,40.6645037913237,44.2832110429707,47.272577903027,47.5872480988224
),y = c(30.0641086410103,34.1278207016757,37.6989616034726,40.900674136118,32.7732500147872,27.4781100569505,22.5523984682652,22.7986840476995,24.5226831037393,29.3252519027075,33.8815351222414
)),.Names = c(\"x\",\"y\"))

## close the polygon
coords <- cbind(xy$x,xy$y)
coords <- rbind(coords,coords[1,])

library(sp)

## create a Spatial polygons object
poly <- SpatialPolygons(list(Polygons(list(Polygon(coords,hole = FALSE)),\"1\")))


## create a contour line that cuts the polygon at height 171
cl <- contourLines(vdata,levels = 171)

## for ContourLines2SLDF
library(maptools)

clines <- ContourLines2SLDF(cl)

现在，将多边形与线相交，然后稍微缓冲线，然后再次与多边形相区别以给出多部分多边形。

library(rgeos)
lpi <- gIntersection(poly,clines)

blpi <- gBuffer(lpi,width = 0.000001)

dpi <- gDifference(poly,blpi)

绘制原始数据，然后从“空间”对象中手动提取半个多边形。

par(mfrow = c(2,1))

image(vdata)
plot(poly,add = TRUE)

plot(SpatialPolygons(list(Polygons(list(dpi@polygons[[1]]@Polygons[[1]]),\"1\"))),add = TRUE,col = \"lightblue\")

image(vdata)
plot(poly,add = TRUE)
cl <- contourLines(vdata,levels = 171)

plot(SpatialPolygons(list(Polygons(list(dpi@polygons[[1]]@Polygons[[2]]),\"2\"))),col = \"lightgreen\")

这适用于这种相当简单的情况，可能对您的方案有用。