Download SRTM for an entire country

Austria SRTM

Happy new year 2017!

It’s January 1st and I made the new years resolution to post more this year. 🙂 So here is a short tutorial on how to download SRTM (Shuttle Radar Topography Mission) 90m resolution data (3 arc seconds) for an entire country using the {raster} and {rgeos} R packages.

Basics (Download single SRTM tile)

As described in a previous post it’s possible to download the 90m SRTM using the getData() function as follows:

srtm <- getData('SRTM', lon=16, lat=48)
  1. Select Dataset: The first argument specifies the dataset. ‘SRTM’ returns the SRTM 90 elevation data.
  2. Specify Lon: The second argument specifies the lon of the SRTM tile.
  3. Specify Lat:  The second argument specifies the lat of the SRTM tile.

The code above will return a sginle SRTM Tile somewhere around Vienna:SRTM Austria

Advanced Download (entire country)

Specifying the longitude and latitude for the data download is often not very handy. Especially when we want to download SRTM data for an entire country. Usually a country consists of more than one tile and looking up the coordinates of each tile is a pain.

Here is a short script where you only have to specify the ISO country code and the tiles are automatically being mosaiced and plotted. You can find the code and all the data needed on github.

Preparation: Before we can start, download the shapefile which contains the the SRTM tile grid. It can be found here or on github. Without this file, the code below won’t work.

Once you have downloaded the shapefile, we can start the script. First load the necessary libraries which are:

library(raster)
library(rgeos)
library(rasterVis)

Then specify the country for your SRTM data download:

#Specify target ISO country code and path to downloaded shapefile
country_name <- "AUT"                        #Austria
shp          <- shapefile("srtm/tiles.shp")  #Path to downloaded SRTM Tiles

After these initial steps, you can run the code from below and lean back. The script automatically determines all the tiles that are needed for the download, mosaics and crops them into a single file:

#Get country geometry first
country <- getData("GADM", 
                   country = country_name, 
                   level=0)

#Intersect country geometry and tile grid
intersects <- gIntersects(country, shp, byid=T)
tiles      <- shp[intersects[,1],]

#Download tiles
srtm_list  <- list()
for(i in 1:length(tiles)) {
  lon <- extent(tiles[i,])[1]  + (extent(tiles[i,])[2] - extent(tiles[i,])[1]) / 2
  lat <- extent(tiles[i,])[3]  + (extent(tiles[i,])[4] - extent(tiles[i,])[3]) / 2
  
  tile <- getData('SRTM', 
                  lon=lon, 
                  lat=lat)
  
  srtm_list[[i]] <- tile
}

#Mosaic tiles
srtm_list$fun <- mean 
srtm_mosaic   <- do.call(mosaic, srtm_list)

#Crop tiles to country borders
srtm_crop     <- mask(srtm_mosaic, country)

#Plot
p <- levelplot(srtm_mosaic)
p + layer(sp.lines(country, 
                   lwd=0.8, 
                   col='darkgray'))

Here the result for Austria plotted with {rasterVis}:

Austria SRTM

Austria SRTM

If you have any questions, feel free to comment below.

Cheers!

Martin

About This Author

Martin was born in Czech Republic and studied at the University of Natural Resources and Life Sciences, Vienna. He is currently working at GeoVille - an Earth Observation Company based in Austria, specialised in Land Monitoring. His main interests are: Open-source applications like R, (geospatial) statistics and data-management, web-mapping and visualization. He loves travelling, geocaching, photography and sports.

10 Comments

You can post comments in this post.


  • Hi Martin,

    great blog.

    I have a troubles getting the country-wide SRTM download, mosaicking and cropping to work. The download stops after the first itereation because the GeoTiff file can’t be found. I am working on Ubuntu 16.04 and I guess the problem lies in the fact that the getData() function searches for pattern = “.TIF” and my unzipped images end on lower case pattern = “.tif”. I tried to source the getData() function but unfortunately it didn’t work. Do you have any work around suggestions?

    Many thanks,
    Valentin

    Valentin Louis 6 years ago Reply


    • Hey Valentin! This is interesting,… I never tried the function on Linux.
      Try setting your working directory manually before the download. I.e. `setwd(“your/custom/dir/for/srtm/”)`.
      Then wrap the `getData()` call inside a try-clause: `try(getData(“…”))`. This should make sure you go through all the iterations.
      After everything is downloaded, call `rasters <- list.files(full.names=T, pattern=".tif$")` and stack them using `stack <- stack(rasters)`. Good luck!

      Martin 6 years ago Reply


      • Thanks, Martin.

        The try-clause worked, but the stack didn’t because of the different extent. lapply() on the rasters did the job.

        Best wishes,
        Valentin

        Valentin Louis 6 years ago Reply


        • Happy to hear that you made it work! Cheers Martin

          Martin 6 years ago Reply


  • Awesome tutorial! I’ve been wrestling with various topographic datasets recently, and this method seems super handy. Any tips for how I might incorporate multiple countries? Or better yet, a particular region of the world (e.g. Central America and the Caribbean)? I’m trying to make a topographic map with multiple countries that’s still cropped to the coastline.

    -Hollis

    Hollis Dahn 6 years ago Reply


  • I keep getting this error. Please guide me:

    country_name  shp <- shapefile("srtm/tiles.shp")  #Path to downloaded SRTM Tiles
    Error: file.exists(extension(x, ".shp")) is not TRUE

    Muhammad Waqas 5 years ago Reply


    • Hey Muhammad,
      it seems that the path to the shapefile is wrong. Follow these steps:
      1. Make sure you downloaded and unzipped the shapefile: http://www.gis-blog.com/wp-content/uploads/2017/01/srtm.zip
      2. Find the location of your shapefile on your computer. For example: “C:/Downloads/srtm/tiles.shp” and copy it into your script instead of the path provided here: country_name shp <- shapefile("srtm/tiles.shp") Hope this works! C

      Martin 5 years ago Reply


  • I have been attempting to carry out this workflow using sf objects, I am struggling to understand how to extract the lon/lat from sf objects to work within the function.

    Would you advise the operation of this in sf objects over sp?

    My goal is to download SRTM DEM for 121 countries coastlines. I would also like to have this in 1 arc second but trying this at the moment.

    Amy 3 years ago Reply


  • Nice post, I have a question: how to merge 2 countries when the second is within the first bigger country? if so how

    Abou 3 years ago Reply


  • Hey Martin,
    thank you for your code and the all the handy information.
    I tried your code for Norway and everything seemed to work, but the plottet map doesn’t seem to exceed 60°N so I just got the most southern tip of Norway. Is there any way to get Data above 60°N using the same method in R?

    Thank you very much!
    Josephine

    Josephine 3 years ago Reply


Post A Reply

*