#r #spatial #r-raster
#r #пространственный #r-растр
Вопрос:
Я довольно новичок в работе с пространственными фреймами данных, и у меня есть то, что я считал относительно простой задачей: возьмите фрейм данных из 6 точек, со столбцами x и y, представляющими позиции широты / длины этих точек, и спроецируйте их так, чтобы их можно было использовать в пространственном фрейме данных, который у меня естьсделано.
Вот как я закодировал в 6 точках:
d1 <- structure(list(latitude = c(37.427733, 37.565759, 37.580956, 37.429285, 37.424270, 37.502496), longitude = c(-108.011061, -107.814039, -107.676662, -107.677166, -108.898826, -108.586042)))
d2 <- as.data.frame(d1)
d3 <- SpatialPointsDataFrame(c(d2[,c('longitude','latitude')]), data = d2)
И я попытался изменить / назначить проекцию для них (эти данные широты / длины были взяты из Google Maps), но, похоже, я не могу заставить это работать. Проекция для данных, на которые я хочу наложить эти точки, выглядит следующим образом:
proj=aea lat_0=23 lon_0=-96 lat_1=29.5 lat_2=45.5 x_0=0 y_0=0 datum=WGS84 units=m no_defs
Итак, в основном мой вопрос: как я могу преобразовать эти широты / длины в формат для x / y, который использует эта проекция? Вот объем набора данных, который я хочу наложить на него для справки, показывая, что он явно не в простой широте / длине:
class : Extent
xmin : -1145835
xmax : -1011345
ymin : 1613205
ymax : 1704855
Заранее большое вам всем спасибо!
Ответ №1:
Вам нужно (повторно) спроецировать свои пространственные точки в ту же проекцию, что и ваш другой источник данных. Я больше знаком с sf
пакетом для работы с пространственной векторной информацией в R, но похоже, что вы используете sp
packcage.
Первым шагом является присвоение вашему фрейму данных правильной проекции. Широта и долгота обычно указаны в WGS84 или epsg: 4326, поэтому:
library(sp)
d1 <- structure(list(latitude = c(37.427733, 37.565759, 37.580956, 37.429285, 37.424270, 37.502496), longitude = c(-108.011061, -107.814039, -107.676662, -107.677166, -108.898826, -108.586042)))
d2 <- as.data.frame(d1)
d3 <- SpatialPointsDataFrame(c(d2[,c('longitude','latitude')]), data = d2)
proj4string(d3) <- CRS(" init=epsg:4326")
sf::st_bbox(d3)
# xmin ymin xmax ymax
# -108.89883 37.42427 -107.67666 37.58096
Глядя на сводку или этот экстент, вы можете видеть, что широта и длина соответствуют оригиналу. Теперь мы можем перепроектировать, используя предоставленную вами proj4string
target_crs = CRS(" proj=aea lat_0=23 lon_0=-96 lat_1=29.5 lat_2=45.5 x_0=0 y_0=0 datum=WGS84 units=m no_defs") # This is your string assigned to an object
d4 = spTransform(d3, target_crs) ## object used to transform your data frame
sf::st_bbox(d4)
# xmin ymin xmax ymax
# -1127246 1661667 -1018854 1679942
Ваши координаты теперь представлены в целевом пространстве и должны быть нанесены поверх ваших фоновых данных.
Если вы хотите использовать альтернативный пакет, sf::st_transform()
это немного проще в использовании, и sf
пакет, как правило, более удобен для пользователя.
Комментарии:
1. Большое вам спасибо! Это отлично работает. Я думаю, вы правы, что формат sp несколько сбил меня с толку. Большое спасибо за добрые и полезные отзывы.