Преобразование широты и долготы в пространственные точки и проецирование их на другую карту в 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

Если вы хотите повторить это упражнение, вы можете скачатькарта здесь .

Это образец моего 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 - это огромные числа, а координаты станций выглядят как реальные геопространственные привязки.Я полагаю, что не правильно проецировал их на Мекса.Я ожидал, что 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-объект и наносить их на карту.

#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 ) )