Я пытаюсь определить, к какому многоугольнику (ZCTA ... он же аналог почтового индекса) принадлежит данная точка, учитывая набор точек и шейп-файл. Хотя есть несколько вопросов такого типа, почти все, кажется, отсылают меня к QGIS. Пока я пойду и изучу другой инструмент, если понадобится, есть ли простой способ сделать это в R? Я имею опыт работы в среде R ... не так много в области ГИС.
Шейп-файл, который я использую, расположен здесь: //ftp.gisdata.mn.gov/pub/gdrs/data/pub/us_mn_state_mngeo/bdry_zip_code_tabulation_areas/shp_bdry_zip_code_tabulation_areas.zip
Моя первая попытка заключалась в том, чтобы загрузить шейп-файл как SpatialPolygonsDataFrame, точки как SpatialPointsDataFrame, а затем использовать «over ()», чтобы получить признаки совпадающих многоугольников:
library(maptools)
library(maps)
library(sp)
mn.zip.map <- readShapePoly("zip_code_tabulation_areas.shp")
# The shapefile is the one referenced in the link above
latlon <- data.frame(matrix(0,nrow=2,ncol=1))
latlon$lat <- c(44.730178, 44.784711)
latlon$lon <- c(-93.235381, -93.476415)
latlon[1] <- NULL
coordinates(latlon) = ~lon+lat
indices <- over(latlon, mn.zip.map)
С результатами:
> indices
ZCTA5CE10 GEOID10 CLASSFP10 MTFCC10 FUNCSTAT10 ALAND10 AWATER10 INTPTLAT10 INTPTLON10
1 <NA> <NA> <NA> <NA> <NA> NA NA <NA> <NA>
2 <NA> <NA> <NA> <NA> <NA> NA NA <NA> <NA>
Shape_Leng Shape_Area
1 NA NA
2 NA NA
Я надеялся, что первая строка будет выводить ZCTA5CE10 == 55124, а вторая строка - ZCTA5CE10 == 55379. Однако явно этого не происходит.
Кажется, что системы координат не выровнены ... но они обе должны быть Широта / Долгота, верно?
Что мне здесь не хватает? Заранее спасибо.