将shapefile添加到ggmap

Question

有点类似于here所描述的问题我遇到了对齐shapefile和ggmap对象的麻烦。

1. 我的shapefile由澳大利亚维多利亚州的局部区域边界组成，我试图将它们覆盖在州（维多利亚州）的谷歌地图上。 source shapefile具有以下PROJ4字符串（从prj文件中提取） [+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs] 这对应于EPSG:4283。 这是我的shapefile对象sp的摘要：

   > summary(sp)
   Object of class SpatialPolygonDataFrame
   Coordinates:
           min        max
   x  96.81677 159.109219
   y -43.74051  -9.142176
   Is projected: FALSE
   proj4string:
   [+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs]

2. 我转换shapefile以匹配谷歌的伪墨卡托投影（至少那是我认为我在做的） sp <- spTransform(sp, CRS("+proj=longlat +init=epsg:3857")) 并将sp变成强化数据框df.sp。
3. 然后我用 qazxsw poi 得到维多利亚的谷歌地形图，和（ggplot主题 qazxsw poi

map <- get_map("Victoria", zoom = 7, maptype = "terrain", source = "google")

从结果图中可以清楚地看出shapefile坐标和googlemap坐标不重叠。我是否在坐标转换方面做错了什么？如何正确匹配shapefile和googlemap坐标？我很感激对此事的任何帮助/见解。

Answer 1

这是我的第一篇文章，请原谅我，如果造型不是很棒。我发现无论出于何种原因，ggmap(map) + geom_polygon(data = df.sp, aes(x = long, y = lat, group = group)) + coord_equal() + theme_map() 包中的默认都没有在ESPG:3857函数中正确转换投影。

• 在GIS StackExchange中使用rgdal，我被引导到spTransform网站和gene's wonderful answer网站。我仔细检查了ABS的Prj2EPSG字符串（我假设这是数据的来源），就像你说的那样（这是Spatial Reference）。
• 在空间参考中输入EPSG：3857时，我得到了这个字符串（使用Proj4链接）proj4。
• 我用EPSG:4283取代+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs来投射经度和纬度坐标。
• 使用+proj=merc替换的字符串似乎比仅使用+proj=longlat给我更好的结果，这表明rgdal包处理它的方式有些不同。

我很想知道如何进一步改进这一点，但我希望也许这可以帮助你在此期间获得更好的结果。

---

编辑：我还发现，当使用merc绘图时添加+init=epsg:3857似乎可以在近距离时修复形状文件的对齐方式，但在进一步缩小时会让它们略微偏离。在缩放比例7及以上，我似乎得到完美的对齐，但6和更低，我没有。似乎无论出于何种原因，不同的缩放级别会略微改变投影。

所以试试吧

ggplot2