将纬度和经度转换为空间点并将它们投影到 R 中的另一张地图上

Question

我从政府官方网站下载了墨西哥所有州的 SHP 地图。我还有一个包含 300 个气象站及其坐标的数据集。 我想将这 300 个站点添加到 SHP 对象中。 但是，我的坐标似乎不在同一个地理空间区域上。

我阅读并绘制了墨西哥地图。然后，我把所有300个站点的经纬度都变成了空间点。不幸的是，当我将它们添加到地图时，我无法找到它们。

这是我的代码：

# Read official SHP map
Mex <- readOGR(dsn = "/Desktop",
               layer = "areas_geoestadisticas_estatales")
# Turn `stations` into Spatial Points and project them onto `Mex`
Geo_stations <- SpatialPoints(stations[, c("stations_lat", "stations_long")],
                              proj4string = CRS(proj4string(mex)))
# Plot `Mex`, then add `Geo_stations`
plot(Mex) 
plot(Geo_stations, col = "red", add = TRUE) # Nowhere to be found

如果您想重复此练习，you can download the map here。

这是我的 stations 数据集的样本：

    station_number station_lat station_long
 1:          10003      25.100     -106.567
 2:          10018      24.944     -106.259
 3:          10031      24.523     -105.952
 4:          10038      23.554     -105.411
 5:          10042      24.174     -105.967
 6:           1004      22.001     -102.199
 7:          10050      25.063     -106.531
 8:           1005      21.781     -102.372
 9:          10064      24.148     -105.953
10:          10087      25.129     -106.363

比较了coordinates(Mex)和coordinates(Geo_stations)的结果，我发现Mex的"coordinates"是一个巨大的数字，而站里的那些看起来像是实际的地理空间参考。我想我没有将它们正确投射到 Mex 上。本来希望plot(Geo_stations, add = TRUE)在全国加一层站

Answer 1

library(sf)
library(mapview)

#load the shapefile
sh = st_read( "./702825217341_s/conjunto_de_datos/areas_geoestadisticas_estatales.shp",
              stringsAsFactors = FALSE )

#set encoding of second column, because of the special characters used
Encoding(sh$NOM_ENT) <- "UTF-16"
#set crs of the current map
st_set_crs( sh, 5332 )
#transform to WGS84-coordinates
st_transform( sh, 4326 )

#let's see what we have so far
mapview::mapview(sh)

看起来不错。

现在您可以将站点添加到 sf-object，然后绘制到地图中。

#sample stations
stations <- data.frame( station_number = c(10003),
                        station_alt = c(1754),
                        station_long = c(-106.567),
                        station_lat = c(26)
)

#make a sf-object out of them
stations <- st_as_sf( stations, coords = c( "station_long", "station_lat" ) )

#and plot
mapview::mapview( list( sh, stations ) )

将纬度和经度转换为空间点并将它们投影到 R 中的另一张地图上

Convert latitude and longitude to Spatial Points and project them over another map in R

r

rgdal