`st_buffer`函数中默认的`dist`参数单位是什么?

11
我有以下的墨西哥地图。它展示了所有的市镇和大约400个气象站。

Mexico

我希望创建一个以每个站点为中心的10公里缓冲区,并最终将每个缓冲区内的市政当局与一个站点关联起来。
地图和站点分别存储在不同的sf对象中。我尝试了以下方法:
buffers <- st_buffer(stations, dist = 1)

我认为dist参数被设置为公里,所以我尝试dist = 10。不幸的是,这会为每个站点返回非常大的缓冲区。这就是我使用dist = 1的原因,但即使这些缓冲区也很大,就像一个州一样!这个问题建议我将我的站点转换为Irish Grid,但我无法复制已接受的答案。我现在想知道dist参数设置的单位是什么。

根据前面提到的问题,我假设它被设置为度数。如何在每个站周围设置10公里的缓冲区?

额外信息:

我的CRS在两个对象(墨西哥地图和车站)上都设置为4326。

这是我的stations数据:

> dput(head(stations))
structure(list(station_number = c(1004L, 1005L, 1008L, 1012L, 
1017L, 1018L), station_alt = c(1925, 1844, 2323, 1589, 2172, 
2053), month = c(9L, 9L, 9L, 9L, 9L, 9L), Mean_min = c(11.6, 
12.75, 12.25, 13.9666666666667, 12.9, 12.6833333333333), Mean_max = c(26.9333333333333, 
26.85, 24.0833333333333, 29.0333333333333, 24.8666666666667, 
26.1333333333333), months_observed = c(5L, 5L, 5L, 5L, 5L, 5L
), geometry = structure(list(structure(c(-102.199, 22.001), class = c("XY", 
"POINT", "sfg")), structure(c(-102.372, 21.781), class = c("XY", 
"POINT", "sfg")), structure(c(-102.135, 22.203), class = c("XY", 
"POINT", "sfg")), structure(c(-102.802, 21.794), class = c("XY", 
"POINT", "sfg")), structure(c(-102.444, 22.233), class = c("XY", 
"POINT", "sfg")), structure(c(-102.415, 22.141), class = c("XY", 
"POINT", "sfg"))), class = c("sfc_POINT", "sfc"), precision = 0, bbox = structure(c(xmin = -102.802, 
ymin = 21.781, xmax = -102.135, ymax = 22.233), class = "bbox"), crs = structure(list(
    epsg = NA_integer_, proj4string = NA_character_), class = "crs"), n_empty = 0L)), sf_column = "geometry", agr = structure(c(station_number = NA_integer_, 
station_alt = NA_integer_, month = NA_integer_, Mean_min = NA_integer_, 
Mean_max = NA_integer_, months_observed = NA_integer_), .Label = c("constant", 
"aggregate", "identity"), class = "factor"), row.names = c(NA, 
6L), class = c("sf", "data.frame"))

文档似乎说st_crs(stations)$units的值会被使用。这不是这种情况吗?如果您包含一个简单的可重现示例,其中包含可以用于测试和验证可能的解决方案的样本输入和期望输出,那么帮助您会更容易。 - MrFlick
@MrFlick 我更新了这篇帖子(添加了使用 dputstations)。至于你的问题,我不知道单位是否设置为与“站点”相同的单位。我该如何验证呢? - Arturo Sbr
1个回答

18

你的坐标是经纬度,所以距离将以度为单位。你应该先将其投影到以米为单位的空间参考,然后取10,000米。

st_buffer 的手册对 dist 参数说:

如果 dist 是一个单位对象,则它应该可转换为弧度,如果 x 具有地理坐标,则应可转换为 st_crs(x)$units。

如果您将坐标保留在 4326 中,您应该能够输入大约 0.1 这样的值,这应该是墨西哥 11 公里左右,但是您将看到一个警告信息:

在 st_buffer.sfc(st_geometry(x), dist, nQuadSegs, endCapStyle = endCapStyle, : st_buffer does not correctly buffer longitude/latitude data

因此,首先将其转换为另一种(以米为单位)投影并以米为单位输入距离。这应该有效,使用 EPSG 7801

library(sf)

pois <- st_as_sf(stations)
st_crs(pois) <- 4326
pois <- st_transform(pois, crs = 7801)
plot(st_geometry(pois))

buff <- st_buffer(pois, dist = 10000)
plot(st_geometry(buff), add = TRUE)

使用Leaflet和测量工具控制:

buff <- st_transform(buff, crs = 4326)

library(leaflet)

leaflet() %>% 
  addTiles() %>% 
  addMeasure(primaryLengthUnit = "meters") %>% 
  addMarkers(data = pois) %>% 
  addPolygons(data = buff)

谢谢!我按照您的方法做了,然后将“buff”更改回EPSG 4326,以便我可以像以前一样继续使用墨西哥地图、站点和缓冲区。您认为这样做有什么不方便吗? - Arturo Sbr
1
不,那完全正确。 :) 虽然我不确定 6829 是否正确。7801 对我来说更正确。 - SeGa
你是指 pois <- st_transform(pois, crs = 6829) 中的 6829 吗? - Arturo Sbr
1
是的,我已将其更改为7801。我会添加一些leaflet代码,以便您可以测量缓冲区。 - SeGa
你怎么知道7801或6829是适当的并且以米为单位? - Arthur Yip

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接