由于地理单位不同,我无法将shapefile添加到我的ggmap中 [英] I'm having trouble adding a shapefile to my ggmap due to differing geographic units

查看:192
本文介绍了由于地理单位不同,我无法将shapefile添加到我的ggmap中的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在尝试将ESRI shapefile(.shp)添加到北卡罗莱纳州的ggmap图中,该图具有以下代码:

x < - get_map(location =North Carolina,zoom = 6,maptype =terrain)



ggmap(x)+
geom_polygon(data = citylim83.df,aes(x = long,y = lat),fill =red,alpha = 0.2)



我已经加载并强化到 citylim83.df
这里是用于将shapefile加载到ggplot中的代码:

  epsgs<  -  make_EPSG()
citylim < - readOGR(dsn =。,layer =MunicipalBoundaries_polys)`

在完成EPSG搜索后,MunicipalBoundaries投影的单位是国家飞机系统的ft-US。即使这个.shp有一个NAD83的地理坐标系统,我也想把它投影到NAD83以摆脱状态平面系统(NAD83(UTM-17N)的EPSG代码是26917):

  citylim83 < -  spTransform(citylim,CRS(+ init = epsg:26917))
汇总(citylim83)

类的对象SpatialPolygonsDataFrame
预计:TRUE
[+ init = epsg:26917 + proj = utm + zone = 17 + ellps = GRS80 + datum = NAD83 + units = m

citylim83.df< - fortify(citylim83)

然后使用此数据框在上面显示的ggmap代码中。



即使它现在已经投影到NAD83中,它仍然不会显示在基本的ggmap上。什么是我导入的get_map对象的基础投影?有没有一个命令来找到这个,所以我可以匹配我的地图和我想要显示在其上的shapefile?我是否需要取消我的citylim对象?如果不明确,shapefile文件是北卡罗莱纳州每个城市的城市边界。任何帮助将非常感谢,因为我对ggplot2 / ggmap社区非常陌生。

解决方案

答案:

  library(ggmap)
library(maptools)
shapefile< - readShapeSpatial(' MunicipalBoundaries_polys.shp',proj4string = CRS(+ proj = lcc + lat_1 = 34.33333333333334 + lat_2 = 36.16666666666666 + lat_0 = 33.75 + lon_0 = -79 + x_0 = 609601.2199999997 + y_0 = 0 + ellps = GRS80 + towgs84 = 0,0, 0,0,0,0,0 + units = us-ft + no_defs))
shp< - spTransform(shapefile,CRS(+ proj = longlat + datum = WGS84))
数据< - fortify(shp)

nc< - get_map(North Carolina,zoom = 6,maptype ='terrain')
ncmap< - ggmap(nc,extent =device)
ncmap +
geom_polygon(aes(x = long,y = lat,group = group),data = data,
color ='gray',fill ='black ',alpha = .4,size = .1)

坐标需要转换为长条投影一个它正在工作。


I am trying to add an ESRI shapefile (.shp) to a ggmap plot of the state of North Carolina, which has the following code:

x <- get_map(location="North Carolina", zoom=6, maptype="terrain")

ggmap(x) + geom_polygon(data=citylim83.df, aes(x=long, y=lat), fill="red", alpha=0.2)

The shapefile I've loaded and fortified into citylim83.df . Here is the code used to load the shapefile into ggplot:

    epsgs <- make_EPSG()
    citylim <- readOGR(dsn=".", layer="MunicipalBoundaries_polys")`

The units of the projection of MunicipalBoundaries, after doing a search in EPSG, are ft-US for the state-plane system. Even though this .shp has a geographic coordinate system of NAD83, I want to also project it to NAD83 to get rid of the state plane system (the EPSG code for NAD83 (UTM-17N) is 26917):

    citylim83 <- spTransform(citylim, CRS("+init=epsg:26917"))
    summary(citylim83)

    Object of class SpatialPolygonsDataFrame
    Is projected: TRUE 
    [+init=epsg:26917 +proj=utm +zone=17 +ellps=GRS80 +datum=NAD83 +units=m

    citylim83.df <- fortify(citylim83)

This data frame was then used in the ggmap code shown above.

Even though it is now definitely projected into NAD83, it still will not show up on the base ggmap. What is the base projection of the get_map object I've imported? Is there a command to find this out so I can match my map with the shapefile I want displayed on top of it? Do I have to "unproject" my citylim object? FYI the shapefile is the city limit boundary of each municipality in North Carolina, if it wasn't clear. Any help would be much appreciated on this, as I'm very new to the ggplot2/ggmap community.

解决方案

Here You have an answer:

library(ggmap)
library(maptools)
shapefile <- readShapeSpatial('MunicipalBoundaries_polys.shp', proj4string = CRS("+proj=lcc +lat_1=34.33333333333334 +lat_2=36.16666666666666 +lat_0=33.75 +lon_0=-79 +x_0=609601.2199999997 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"))
shp <- spTransform(shapefile, CRS("+proj=longlat +datum=WGS84"))
data <- fortify(shp)

nc <- get_map("North Carolina", zoom = 6, maptype = 'terrain')
ncmap <- ggmap(nc,  extent = "device")
ncmap +
  geom_polygon(aes(x = long, y = lat, group = group), data = data,
               colour = 'grey', fill = 'black', alpha = .4, size = .1)

The coordinates needs to be transform to longlat projection and it's working.

这篇关于由于地理单位不同,我无法将shapefile添加到我的ggmap中的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆