`st_buffer` 中的 `dist` 参数默认设置为啥单位？

Posted 2023-02-16

【中文标题】`st_buffer` 中的 `dist` 参数默认设置为啥单位？

【英文标题】：What unit is the `dist` argument in `st_buffer` set to by default?`st_buffer` 中的 `dist` 参数默认设置为什么单位？

【发布时间】：2019-07-12 06:00:21

【问题描述】：

我有以下墨西哥地图。它显示了其所有城市和大约 400 个气象站。

我想在每个车站周围创建一个 10 公里的缓冲区，并最终将每个市政当局与位于每个半径内的车站相关联。

地图和车站存储在单独的 sf 对象中。我厌倦了以下内容：

buffers <- st_buffer(stations, dist = 1)

我认为dist 参数设置为公里，所以我尝试了dist = 10。不幸的是，这为每个站点返回了 HUGE 缓冲区。这就是我使用dist = 1 的原因，但即使这些缓冲区也与状态一样大！ This question，建议我将我的站点转换为 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 的值将被使用。不是这样吗？如果您包含一个简单的reproducible example，其中包含可用于测试和验证可能解决方案的示例输入和所需输出，则更容易为您提供帮助。

@MrFlick 我更新了帖子（使用dput 添加了stations）。至于你的问题，我不知道单位是否设置为与stations 相同的单位。我该如何验证？

【参考方案1】：

您的坐标是长/纬度，因此距离将以度为单位。您应该首先投影到以米为单位的空间参考，然后取 10 000 米。

st_buffer 的手册是这样描述 dist 参数的：

如果 dist 是一个单位对象，它应该可以转换为 arc_degree 如果 x 有地理坐标，并且到 st_crs(x)$units 否则

如果您将坐标保留在 4326 中，您应该能够采取类似 0.1 的坐标，对于墨西哥来说应该是大约 11 公里，但您会看到一条警告消息：

在 st_buffer.sfc(st_geometry(x), dist, nQuadSegs, endCapStyle = endCapStyle, : st_buffer 没有正确缓冲 经纬度数据

所以首先转换为另一个投影（以米为单位）并输入以米为单位的距离。这应该可以工作，它使用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)

使用传单和测量工具进行控制：

buff <- st_transform(buff, crs = 4326)

library(leaflet)

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

【讨论】：

谢谢！我做了你所做的，然后将buff改回EPSG 4326，这样我就可以像以前一样继续使用墨西哥地图、车站和缓冲区。你觉得这样做有什么不便吗？

不，那是完全正确的。 :) 虽然我不确定 6829 是否正确。 7801 对我来说看起来更正确。

你是指pois <- st_transform(pois, crs = 6829)中的6829吗？

是的，我将其更改为 7801。我将添加一些传单代码，以便您可以测量缓冲区。

您怎么知道 7801 或 6829 是否合适以及以米为单位？