R-如何在一个点周围绘制半径并使用该结果来过滤其他点

Question

我希望在lat long点上绘制半径，然后使用该缓冲区过滤其他适合其中的点。例如：

#stores datasets
stores = data.frame(store_id = 1:3,
                    lat = c("40.7505","40.7502","40.6045"),
                    long = c("-73.8456","-73.8453","-73.8012")
                    )
#my location
me  = data.frame(lat = "40.7504", long = "-73.8456")
#draw a 100 meter radius around me 

#use the above result to check which points in dataset stores are within that buffer

不确定如何处理。我之前曾与over合作以相交点和多边形，但不确定如何在孤独点上运行类似的情况。

Answer 1

您可以假设尝试在球体或椭球表面上进行几何计算，但通常在执行几何映射操作的情况进入平坦的表面。

这是如何使用 sf 软件包完成此操作。首先，在LON-LAT坐标中创建要点：

library(sf)
lat <- c(40.7505, 40.7502, 40.6045)
lon <- c(-73.8456, -73.8453, -73.8012)
stores <- st_sfc(st_multipoint(cbind(lon, lat)), crs = 4326)
me <- st_sfc(st_point(c(-73.8456, 40.7504)), crs = 4326)

crs = 4326参数指定LON-LAT坐标系的EPSG代码。接下来，我们需要选择一个地图Projecton。对于此示例，我将使用UTM区域18，其中包含以上几点：

stores_utm <- st_transform(stores, "+proj=utm +zone=18")
me_utm     <- st_transform(me, "+proj=utm +zone=18")

现在，我们可以缓冲代表自己的点100米以产生半径为100米的圆：

circle <- st_buffer(me_utm, 100)

现在，我们几乎可以使用几何谓词来测试哪些点在圆圈中。但是，stores_utm当前是MULTIPOINT，因此几何谓词将其视为一个几何实体。我们可以通过将stores_utm铸造为POINT来解决此问题，这将为我们提供三个单独的观点的集合：

stores_utm_column <- st_cast(stores_utm, "POINT")
stores_utm_column
# Geometry set for 3 features 
# geometry type:  POINT
# dimension:      XY
# bbox:           xmin: 597453 ymin: 4495545 xmax: 601422.3 ymax: 4511702
# epsg (SRID):    32618
# proj4string:    +proj=utm +zone=18 +ellps=WGS84 +units=m +no_defs
# POINT (597453 4511702)
# POINT (597478.7 4511669)
# POINT (601422.3 4495545)

现在我们可以测试哪些点在圆圈中：

> st_contains(circle, stores_utm_column, sparse = FALSE)
#      [,1] [,2]  [,3]
# [1,] TRUE TRUE FALSE

表明前两个点在圆圈中，第三点不是。

当然，每个地图投影都会引入一些扭曲。您选择投影将取决于问题的性质。

Answer 2

farking points_in_circle((来自空道级包装处理。

例如，使用您的数据：

library(spatialrisk)
# Stores 
stores <- data.frame(store_id = 1:3,
                     lat = c(40.7505, 40.7502, 40.6045),
                     long = c(-73.8456, -73.8453, -73.8012))
# My location
me <- data.frame(lat = 40.7504, long = -73.8456)
> spatialrisk::points_in_circle(stores, me$long[1], me$lat[1], radius = 100, lon = long)
# store_id     lat     long distance_m
#        1 40.7505 -73.8456   11.13195
#        2 40.7502 -73.8453   33.70076