In this post, let us explore the R package rgdal for map shape plotting. We shall attempt to plot the map of Singapore and display major road networks in Singapore.
Pre-requisites
To get the data you need, you can go to a GIS provider. In this post, we shall be using diva-gis.
Steps:
1. Go to: http://www.diva-gis.org/gdata
2. Select ‘Spatial Data Download’
3. Select Country = ‘Singapore’
4. Select Subject = ‘Administrative areas (GADM)’
You will see a download link, where you can obtain the direct download url that we will be using in this post.
library(httr)
setwd("./map")
url <- "http://biogeo.ucdavis.edu/data/diva/adm/SGP_adm.zip"
GET(url, write_disk("SGP_adm.zip", overwrite=TRUE))
## Response [http://data.biogeo.ucdavis.edu/data/diva/adm/SGP_adm.zip]
## Date: 2015-05-04 20:57
## Status: 200
## Content-Type: application/zip
## Size: 64.2 kB
## <ON DISK> SGP_adm.zip
unzip("SGP_adm.zip",exdir="./SGP_adm")
dir("./map/SGP_adm")
## character(0)
After downloading and extracting the map files, let us determine which ‘layer’ of the map to draw and proceed to read and plot it.
library(rgdal)
## Loading required package: sp
## rgdal: version: 0.9-2, (SVN revision 526)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 1.11.1, released 2014/09/24
## Path to GDAL shared files: C:/R/win-library/3.1/rgdal/gdal
## GDAL does not use iconv for recoding strings.
## Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]
## Path to PROJ.4 shared files: C:/R/win-library/3.1/rgdal/proj
ogrListLayers(dsn="./map/SGP_adm/SGP_adm0.shp")
## [1] "SGP_adm0"
## attr(,"driver")
## [1] "ESRI Shapefile"
## attr(,"nlayers")
## [1] 1
sgp <- readOGR(dsn="./map/SGP_adm/SGP_adm0.shp", layer="SGP_adm0")
## OGR data source with driver: ESRI Shapefile
## Source: "./map/SGP_adm/SGP_adm0.shp", layer: "SGP_adm0"
## with 1 features
## It has 70 fields
plot(sgp)
Alright, now we have the base map of Singapore. It looks a little plain now so let’s overlay it with the major roads in Singapore. To get the shape files for roads, refer back to Step 4 under pre-requisites and change the Subject to “Roads” to get the direct link.
url <- "http://biogeo.ucdavis.edu/data/diva/rds/SGP_rds.zip"
GET(url, write_disk("SGP_rds.zip", overwrite=TRUE))
## Response [http://data.biogeo.ucdavis.edu/data/diva/rds/SGP_rds.zip]
## Date: 2015-05-04 20:57
## Status: 200
## Content-Type: application/zip
## Size: 3.14 kB
## <ON DISK> SGP_rds.zip
unzip("SGP_rds.zip", exdir="./SGP_rds")
dir("./map/SGP_rds")
## [1] "SGP_roads.dbf" "SGP_roads.prj" "SGP_roads.shp" "SGP_roads.shx"
ogrListLayers(dsn="./map/SGP_rds/SGP_roads.shp")
## [1] "SGP_roads"
## attr(,"driver")
## [1] "ESRI Shapefile"
## attr(,"nlayers")
## [1] 1
roads <- readOGR(dsn="./map/SGP_rds/SGP_roads.shp", layer="SGP_roads")
## OGR data source with driver: ESRI Shapefile
## Source: "./map/SGP_rds/SGP_roads.shp", layer: "SGP_roads"
## with 22 features
## It has 5 fields
plot(sgp)
lines(roads,col="red")
