在 Google 地图图块上绘制 shapefile [关闭]答案

【问题标题】：Plotting shapefiles on top of Google map tiles [closed]在 Google 地图图块上绘制 shapefile [关闭]
【发布时间】：2010-01-16 18:32:07
【问题描述】：

我有一些要在 Google 地图图块上绘制的 shapefile。最有效的方法是什么？一种方法可能是使用 pkg RgoogleMaps，但是，我仍然不清楚如何执行此操作。我假设使用 PlotonStaticMap 与重新格式化 shapefile 数据的某种组合

【问题讨论】：

标签： google-maps r mapping shapefile

【解决方案1】：

从开发人员那里，它似乎工作得很好。

  shpFile <- system.file("shapes/sids.shp", package="maptools");  
  shp<-importShapefile(shpFile,projection="LL");  
  bb <- qbbox(lat = shp[,"Y"], lon = shp[,"X"]);  
  MyMap <- GetMap.bbox(bb$lonR, bb$latR, destfile = "SIDS.jpg");  
  #compute regularized SID rate  
  sid <- 100*attr(shp, "PolyData")$SID74/(attr(shp, "PolyData")$BIR74+500)  
  b <- as.integer(cut(sid, quantile(sid, seq(0,1,length=8)) ));  
  b[is.na(b)] <- 1;  
  opal <- col2rgb(grey.colors(7), alpha=TRUE)/255;  opal["alpha",] <- 0.2;  
  shp[,"col"] <- rgb(0.1,0.1,0.1,0.2);  
  for (i in 1:length(b)) shp[shp[,"PID"] == i,"col"] <- rgb(opal[1,b[i]],opal[2,b[i]],opal[3,b[i]],opal[4,b[i]]);  
  PlotPolysOnStaticMap(MyMap, shp, lwd=.5, col = shp[,"col"], add = F);

【讨论】：

很好的解决方案，但你应该添加一些上下文： install.packages('PBSmapping', 'RgoogleMaps', 'maptools', 'rgdal', 'ReadImages');图书馆（PBS映射）；图书馆（RgoogleMaps）；图书馆（地图工具）；另请注意，如果 SIDS.jpg 不起作用，您可以尝试 SIDS.png

【解决方案2】：

数据来自data.gov.sg。加载绘制地图所需的包

library(rgdal)  #for shapefiles
library(ggmap)  #plotting maps
library(ggplot2) #general plots

使用 readOGR 读取 shapefile

 shpfile <- readOGR(dsn = 'data/street-and-places/','StreetsandPlaces',verbose = FALSE)
head(coordinates(shpfile),5)

     coords.x1 coords.x2
[1,]  28640.40  29320.77
[2,]  29429.83  28548.69
[3,]  29224.01  28360.38
[4,]  29827.20  29664.26
[5,]  28451.00  29451.78

我们观察到坐标位于不同的投影系统中。所以我们必须把它转换成网络墨卡托投影系统。我们使用 spTransform 转换 shapefile。

crsobj <- CRS("+proj=longlat +datum=WGS84")   #Web Mercator projection system
shpfile_t <- spTransform(shpfile,crsobj)      #Applying projection transformation
df <- as.data.frame(coordinates(shpfile_t))   #Converting to a data frame
head(df,5)

coords.x1<dbl>coords.x2<dbl>
1   103.8391    1.281441        
2   103.8462    1.274459        
3   103.8443    1.272756        
4   103.8497    1.284547        
5   103.8374    1.282626

我们注意到转换后的投影系统。让我们将它绘制在一张底图上。

sgmap <- get_map(location="Singapore", zoom=11,   
             maptype="roadmap", source="google") #Using Osm base map of Singapore
p <- ggmap(sgmap) +
    geom_point(data = df,
    aes(x = coords.x1, y = coords.x2),
    color = 'orange',size = 1)
p

【讨论】：