当点位于国家海岸限制（在 R 中）时，坐标中的国家/省名称问题[重复]

Question

我有一组大部分位于伊比利亚半岛内或附近的坐标。

我有兴趣找出大量坐标的国家名称和省份名称。但是，在国家或省内未检测到的国家边界坐标（沿海地区附近的坐标）存在一些问题。首先，我试图找出省份名称和国家名称之后，我发现这两种不同方法存在这个问题。

我知道有些坐标可能是拼写错误，但我已经手动检查了一些坐标是否在国家/地区范围内。

关于如何克服这个问题的任何建议？

在这里创建一个最小的可重现示例：

首先，我用有问题的坐标创建了一个数据框。最后一个坐标 (36.76353, -4.425162) 是“正确的”，以表明函数有效。

lat <- c(36.81973, 43.69739, 36.51103, 38.50000, 42.25167, 42.25177, 39.31316, 36.76353)
long <- c(-2.411557, -5.919138, -4.635675, -0.100000, -8.804174, -8.790300, 2.995276, -4.425162)
coords <-data.frame(long, lat)

然后通过添加 library(mapSpain) 重新调整同事的答案 Latitude Longitude Coordinates to State Code in R我创建了这个函数，应该告诉省份不同点在哪里。

#I create a function that should tell the provinces of Spain from coordinates
library(sf)
library(spData)
library(mapSpain)

## pointsDF: A data.frame whose first column contains longitudes and
##           whose second column contains latitudes.
##
## states:   An sf MULTIPOLYGON object with 50 states plus DC.
##
## name_col: Name of a column in `states` that supplies the states'
##           names.
lonlat_to_state <- function(pointsDF,
                            states = mapSpain::esp_get_prov(),
                            name_col = "ine.prov.name") {
  ## Convert points data.frame to an sf POINTS object
  pts <- st_as_sf(pointsDF, coords = 1:2, crs = 4326)
  
  ## Transform spatial data to some planar coordinate system
  ## (e.g. Web Mercator) as required for geometric operations
  states <- st_transform(states, crs = 3857)
  pts <- st_transform(pts, crs = 3857)
  
  ## Find names of state (if any) intersected by each point
  state_names <- states[[name_col]]
  ii <- as.integer(st_intersects(pts, states))
  state_names[ii]
}

lonlat_to_state(coords)

[1] NA       NA       NA       NA       NA       NA       NA       "Málaga"

如果我使用，我也会遇到同样的问题

library(maps)
map.where(x = coords$long, y = coords$lat)
[1] NA      NA      NA      NA      NA      NA      NA      "Spain"

欢迎任何建议，感谢您的宝贵时间！

Answer 1

st_intersects 按预期工作，因为点（除了一个）位于多边形之外。我们可以通过绘制它们来看到这一点（放大每个）

library(ggplot2)
for (i in seq_len(nrow(pts))) {
  xlim = unlist(st_geometry(pts[i,]))[1] + c(-5000, 5000)
  ylim = unlist(st_geometry(pts[i,]))[2] + c(-5000, 5000)
  g = ggplot() +
    geom_sf(data=states) +
    geom_sf(data=pts) +
    coord_sf(xlim = xlim, ylim = ylim)
  print(g)
}

地块的几个例子：