R:点周围的方形缓冲区

sqyvllje  于 2023-02-10  发布在  其他
关注(0)|答案(3)|浏览(140)

我一直试图弄清楚如何在点周围创建方形缓冲区,但我最接近的是使用terra::buffer和quadsegs = 1生成一个菱形缓冲区。下面是可复制的代码。任何建议都非常感谢!
附言:上传剧情的时候出了点问题,不过我相信是堆栈溢出的问题

library(terra)
library(geosphere)
创建数据
lon <- seq(from = 10, by = 3/3600, length.out = 4)
lat <- rep(0, 4)
lon.lat <- cbind(lon, lat)
crs.lon.lat <- "epsg:4326"
grid <- terra::vect(lon.lat, crs = crs.lon.lat)
grid$id <- 1:length(grid)
以米为单位设置缓冲区大小并创建缓冲区
res.7as <- geosphere::distGeo(c(0, 0), c(1, 0))*7/3600
grid.buf <- terra::buffer(grid,
                          width = res.7as,
                          quadsegs = 1)
plot(grid.buf)
plot(grid, add = T)
wkyowqbh

wkyowqbh1#

可以使用terra::spin旋转矢量几何。
示例数据

library(terra)
lon <- seq(from = 10, by = 3/3600, length.out = 4)
lon.lat <- cbind(lon, 0)
pts <- terra::vect(lon.lat, crs = "epsg:4326")
buf <- buffer(pts, width=50, quadsegs = 1)

旋转中心是矢量化的,所以你可以

xy <- crds(pts)
s <- spin(buf, 45, xy[,1], xy[,2])
plot(s); points(pts)
hpcdzsge

hpcdzsge2#

例如,像这样,从您的数据开始:

library(sf)
library(dplyr) ## for convenient dataframe manipulation

angle = pi/4 ## desired angle of rotation

## function to rotate a feature (of class sfc from package sf)
## around its centroid 
rot = function(feature, a = 0){
    cnt  <- st_centroid(feature)
    (feature - cnt) *
        matrix(c(cos(a), sin(a), -sin(a), cos(a)), 2, 2) +
        cnt
}

SpatVect类转换为sfc类以应用{sf}方法:

grid.buf.sf <- st_as_sf(grid.buf)

根据需要旋转:

grid.buf.rotated <- 
    grid.buf.sf %>%
    rowwise %>%
    mutate(geometry = rot(geometry, angle))

如果需要,重新转换为SpatVect以与{terra}一起使用:

grid.buf.rotated <- vect(grid.buf.rotated)

旋转函数取自https://r-spatial.github.io/sf/articles/sf3.html#affine-transformations

vyu0f0g1

vyu0f0g13#

下面是一个完全基于sf的替代解决方案-st_buffer(..., endCapStyle = "SQUARE")为我们完成了这项工作。

library("sf")

example_point_1 = sf::st_point(c(0, 10)) 
example_point_2 = sf::st_point(c(110, 40))

example_points_sf <- sf::st_sfc(example_point_1, example_point_2) %>% 
                                sf::st_sf()

sf::st_buffer(example_points_sf,
              dist = 50,
              endCapStyle = "SQUARE") %>% plot() 
plot(example_points_sf, add = T)

相关问题