OSM

下載OSM助教課投影片

library(GISTools);library(rgdal);library(OpenStreetMap)
library(ggplot2);library(gridExtra);library(ggspatial)
setwd("D:/1082SA/Data")
TW=readOGR("Popn_TWN2.shp")
TP=TW[TW$COUNTY=="臺北市",]

par(mar=c(0,0,0,0))
gplot=function(epsg) ggplot()+geom_polygon(data=fortify(spTransform(TP,CRS(paste0("+init=epsg:",epsg))),region="TOWN"),aes(x=long,y=lat,group=group),color='white')+coord_fixed()+labs(title=paste0("EPSG:",epsg),x=NULL,y=NULL)
grid.arrange(gplot(3826),gplot(4326),gplot(3857),nrow=1)

map=openmap(c(21.8,119),c(25.5,123),zoom = 7,type="osm")
map84=openproj(map,CRS("+init=epsg:4326"))
grid.arrange(autoplot(map),autoplot(map84),nrow=1)

par(mfrow=c(1,3))
TP = spTransform(TP,CRS("+init=epsg:4326"))
#
ul = c(TP@bbox[2,2],TP@bbox[1,1])
lr = c(TP@bbox[2,1],TP@bbox[1,2])
map = openmap(ul,lr,12, "osm")
plot(map)
plot(spTransform(TP,osm()), add=T)
#
map=openmap(rev(TP@bbox[,1]),rev(TP@bbox[,2]),11,"osm")
plot(map)
plot(spTransform(TP,osm()), add=T)
#
map=openproj(map,CRS("+init=epsg:4326"))
plot(map)
plot(TP, add=T)

#ggplot畫osm方法
map=openmap(c(lat= 24.95, lon= 121.4),c(lat= 25.25, lon= 121.8),zoom = 11,type="osm")
autoplot(openproj(map))+layer_spatial(TP,fill='transparent',col='red')

## Coordinate system already present. Adding new coordinate system, which will replace the existing one.