空间操作:如何按面积而不是距离(从边缘)缓冲多边形

spatial operations: how to buffer a polygon by area, not distance (from edge)

我想缓冲区域 A 的空间多边形 P 以便缓冲的要素 P_buffered达到定义的区域A_buffered。函数 sf::st_buffer 按距离 d(从边缘)而不是按面积增长特征。

到目前为止,我尝试了:

  1. P(封闭圆的半径,bbox 的对角线)的一些径向测量近似 d,假设径向增加因子 x 将使该区域膨胀 x²(大致取决于 P 的旋转对称性)。
  2. 通过间隔减半
  3. 迭代地逼近d

(1) 的准确性随特征的形状而变化,而 (2) 太慢(至少我的实现是这样)

对于具有相应功能 (pseudo::buffer_area()) 或我似乎忽略的代码解决方案的某些包的提示,我将不胜感激。

回复。从 R 访问其他 GIS 可执行文件,请注意代码必须 运行 在只有 R 包的可用性被认为是理所当然的机器上。

示例数据dput 转储示例多边形 class sfc与 {sf} 一起使用:https://gist.github.com/1O/bc3798468b48f19ab2533f16c99c2268

我仍然建议使用 while 循环(您的第二种方法)...也许有人可以调整 growth/decline 算法...

library(sf)

points = matrix(c(0,0,20,0,10,10,0,10,0,0),ncol=2, byrow=TRUE)
pts = list(points)
pl1 = st_polygon(pts)

# x = polygon
# y = desired buffer area
# z = allowed delte from area 
# example y = 100, z = 0.01, the bufferarea aloowed is bewteen 99 and 101
buffer_distance <- function(x, y, z) {
  buff_dist <- 1
  buffer_area <- st_area(st_buffer(x, buff_dist)) - st_area(x)
  while (!data.table::between(buffer_area, 
                              y - z * y,
                              y + z * y)) {
    if (buffer_area > y) {
      buff_dist <- buff_dist * (1 - y / st_area(st_buffer(x, buff_dist) ))
    } else {
      buff_dist <- buff_dist * (1 + y / st_area(st_buffer(x, buff_dist) ))
    }
    buffer_area <- st_area(st_buffer(x, buff_dist)) - st_area(x)
  }
  return(buff_dist)
}

buffer_distance(pl1, 100, 0.01)
#[1] 1.688372

st_area(st_buffer(pl1, buffer_distance(pl1, 100, 0.01))) - st_area(pl1)
# [1] 100.3634

plot(st_buffer(pl1, buffer_distance(pl1, 100, 0.01)))
plot(pl1, add = TRUE)   

实际上有一个简单的解决方案:sf 可以通过将几何乘以标量 来线性缩放几何,例如:

## scale feature to triple linear extent (e.g. height):
zoom = 3
scaled_feature = zoom * feature

## scale feature to triple area:
scaled_feature = sqrt(zoom) * feature

必要的缓冲距离计算,例如,作为特征的内切圆半径与缩放特征之间的差异。

示例:


library(spData) ## contains sample data, e.g. London areas
library(sf)
library(dplyr)

## original feature (example: Kingston upon Thames in London,
## reprojected area-preserving Lambert (LCC Europe, EPSG 3034)
feature <- spData::lnd[1,] %>% st_transform(3034) %>% st_geometry 
zoom  <- sqrt(1.5) ## to scale area by 1.5
feature_zoomed <- feature * zoom
## get radial buffer distance as half the vertical difference
## between bounding boxes of inscribed circles):
buffer_distance <- 0.5 * (
    (feature_zoomed %>% st_inscribed_circle %>% st_bbox %>% .[c(4,2)] %>% dist) -
    (feature %>% st_inscribed_circle %>% st_bbox %>% .[c(4,2)] %>% dist)
)
feature_buffered <- st_buffer(feature, buffer_distance)

精度会随着顶点数和 st_buffer 设置(线段、连接样式、端盖样式)而变化,将多边形分解为多边形 (st_cast) 和 st_simplifying 可以提供帮助以及将初始 buffer_distance 送入近似循环,如@Wimpel 的回答所示。