Adding elevation point data and lines to topography plot using ggplot2

Question

I have pieced together some code to make a topography plot using elevation data and ggplot2.

However, now I would like to overlay onto this plot the country outlines, along with points to show the location / elevation of points of interest e.g. cities and also other polygons and lines to show for example routes.

How do I do this?

My code is as follows:

library(dplyr)
library(ggplot2)
library(raster)
library(scales)
library(rgeos)

#get elevation data
dem1<- getData("SRTM",lat=34.5,lon=33)
dem2<- getData("SRTM",lat=34.5,lon=36)
dem <- merge(dem1,dem2)

#create polygon to crop the elevation data file   
library(sp)
coords = matrix(c(34, 31,
                  34, 32.5,
                  36, 32.5,
                  36, 31.0,
                  34, 31), 
                ncol = 2, byrow = TRUE)
P1 = Polygon(coords)
Ps1 = SpatialPolygons(list(Polygons(list(P1), ID = "a")), proj4string=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))

#crop the elevation data using the polygon
dem <- crop(dem,Ps1)

#lower the reolution to enable faster plotting
dem_lower_res <- aggregate(dem, fact=10)
dem.p  <-  rasterToPoints(dem_lower_res )
df <-  data.frame(dem.p)
colnames(df) = c("lon", "lat", "alt")

#plot the data
ggplot(df, aes(lon,lat)) +
  geom_raster(aes(fill = alt))+
  scale_fill_gradientn(colours = terrain.colors(10))

Answer 1

Just as an example, I used GADM country outline and two city locations (notional).

• GADM dataset was downloaded as sp object, and plotted by geom_polygon()
• A dataframe of two points (my_points) was converted to sf object, and plotted by geom_sf. Also geom_sf_label to print the label (city name) may be of interest for you.

It highly depends what kind of point/polygon/line you have... hope this helps.

---

library(raster)
#get elevation data
dem1<- getData("SRTM",lat=34.5,lon=33)
dem2<- getData("SRTM",lat=34.5,lon=36)
dem <- merge(dem1,dem2)

# get country border ################################################
gadm <- getData('GADM', country = "ISR", level = 0)
# use st_read() if you have a shapefile

# Add points  #######################################################
my_points = data.frame(
  city = c('Jerusalem', 'Tel Aviv-Yafo'),
  lon =  c(35.1, 34.9),
  lat =  c(31.8, 32.1)
)
library(sf)
my_points = st_as_sf(my_points, coords = c("lon", "lat"), crs = 4326)

# simplified your code a little, using extent() function ############
#create polygon to crop the elevation data file   
Ps1 = as(extent(34, 36, 31, 32.5), 'SpatialPolygons')
crs(Ps1) = "+proj=longlat +datum=WGS84 +no_defs"

#crop the elevation data using the polygon
dem = crop(dem, Ps1, snap= 'out')

#lower the reolution to enable faster plotting
dem_lower_res <- aggregate(dem, fact=10)
dem.p  <-  rasterToPoints(dem_lower_res )
df <-  data.frame(dem.p)
colnames(df) = c("lon", "lat", "alt")

#plot the data
library(ggplot2)
ggplot(df) +
  geom_raster(aes(lon, lat, fill = alt))+
  scale_fill_gradientn(colours = terrain.colors(10))+
  geom_sf(data= my_points, color= 'red')+                                         # point
  geom_sf_label(data= my_points, aes(label= city), nudge_y = 0.1,  color= 'red')+ # label
  geom_polygon(data= gadm, aes(long, lat), color= 'blue', fill= NA) +             # polygon
  coord_sf(xlim = c(34, 36), ylim = c(31, 33))+
  theme_bw()