将 lat/long 或 URM 坐标添加到 R 中的 shapefile

Question

我是 R 中的 GIS 新手，我正在尝试将 lat/long 或 UTM 坐标添加到 shapefile。我从这里下载了芝加哥市 (City_Boundaries.shp) 的边界：http://www.cityofchicago.org/city/en/depts/doit/supp_info/gis_data.html

我加载了 maptools 和 rgeos 库：

library(maptools)
library(rgeos)
library(rgdal)

我将数据导入 R 并尝试为 16T 区添加 UTM 代码：

city<- readShapeSpatial("C:/Users/Luke.Gerdes/Documents/GIS Files/Chicago/City_Boundary.shp")
proj4string(city) <- CRS("+proj=utm +north +zone=16T + datum=WGS84")

但是，生成的数据对我来说没有意义。当我查看 "city" 中的 "coords" 槽时，坐标值例如为 X=1092925 和 Y=1944820。我使用外部工具 (http://home.hiwaay.net/~taylorc/toolbox/geography/geoutm.html) 查找 GPS 坐标。结果是：lat=17.511，long=-81.42。这是牙买加和洪都拉斯海岸之间的某个地方。在我看来，shapefile 坐标存在于它们自己的宇宙中；正如名称 shapefile 所暗示的那样，坐标提供了城市形状的准确描述，但这些坐标不会自动映射到地球上。

我的最终目标是确定是否在芝加哥发生了一些带有地理标记的事件。我很乐意将 lat/long 列出的这些点转换为 UTM，如下所示：

  SP <- SpatialPoints(cbind(-87.63044, 41.79625), proj4string=CRS("+proj=longlat +datum=WGS84"))
  SP <- spTransform(SP, CRS("+proj=utm +north +zone=16T +datum=WGS84"))

如果这样更有意义，我也愿意尝试以原始格式处理事件数据。最终，我需要帮助将芝加哥城市边界变成一个多边形，以便我可以根据我的事件数据进行检查。像这样：

  gContains(city, SP)

如何将我的 shapefile 和事件数据放入一个通用（且准确）的参考系中，以便我可以检查城市中是否有点？如果您能提供任何帮助，我们将不胜感激。

Answer 1

一般建议使用 rgdal::readOGR 读取 shapefile，如

library(rgdal)
city = readOGR(".", "City_Boundary")

这不仅会给 city 适当的 CRS:

proj4string(city)
[1] "+proj=tmerc +lat_0=36.66666666666666 +lon_0=-88.33333333333333 +k=0.9999749999999999 +x_0=300000 +y_0=0 +datum=NAD83 +units=us-ft +no_defs +ellps=GRS80 +towgs84=0,0,0"

但如果您想在不重新投影的情况下更改它，也会警告您：

proj4string(city) <- CRS("+proj=utm +north +zone=16T + datum=WGS84")
Warning message:
In `proj4string<-`(`*tmp*`, value = <S4 object of class "CRS">) :
  A new CRS was assigned to an object with an existing CRS:
+proj=tmerc +lat_0=36.66666666666666 +lon_0=-88.33333333333333 +k=0.9999749999999999 +x_0=300000 +y_0=0 +datum=NAD83 +units=us-ft +no_defs +ellps=GRS80 +towgs84=0,0,0
without reprojecting.
For reprojection, use function spTransform in package rgdal

意思是这样做之后，city是错误的。您使用 spTransform 重新投影 city，如

city = spTransform(city, CRS("+proj=utm +north +zone=16T + datum=WGS84"))

您的 rgeos::gContains 调用可能会失败，因为 utm 的两个 CRS 字符串不同；删除前者 + datum=WGS84 中的 space 就可以了 (TRUE).

将 lat/long 或 URM 坐标添加到 R 中的 shapefile

Add lat/long or URM coordinates to shapefiles in R

r

gis

rgdal

maptools