Photo-realistic shaded relief using Blender
The past days, whilst procrastinating from my thesis work, I was working on a map for a presentation I’ll be giving at the Dutch Ornithological Union congress in a little more than 2 weeks. I had in mind to use Daniel Huffman’s Blender tutorial to make a map of the Caucasus region with sort-of a photo-realistic relief, to emphasise the barrier the Greater and Lesser Caucasus form, without it looking unrealistic.
Daniel’s tutorial is an excellent, easy to follow example on how Blender can be incorporated in the map making GIS workflow. Although I rarely make maps that show topography, it is a technique I surely want to use more often. As I made a few modifications to the process, e.g. by adding cloudless satellite imagery, the rest of the post is mostly a reference for in the future of elements not described in the tutorial.
Why use Blender instead of the default Hillshade generators in GIS software?
Let’s take an extreme example. The image below shows the Caucasus region as rendered by a hillshade algorithm. It is a very harsh representation of topography. In some cases, this can of course be very useful as it exaggerates topographical features, but it is not necessarily very realistic.
Hillshade rendering of the Caucasus region. Very harsh and overwhelming.
The primary factor used to ‘light up’ elements through a hillshade is the direction towards a virtual sun. In reality, light scatters and bounces back and forth between topographical elements, but that is not incorporated at all in the hillshade algorithm. Blender, on the other hand, is perfectly capable of realistically scattering and bouncing light around in a scene and the resulting map (see below) is much more realistic and pleasing to the eye. I’d argue it also is easier to interpret to the untrained eye. For example: the Lesser Caucasus appears quite similar to the Greater Caucasus in the hillshade rendering, but looks substantially lower in Blender’s rendering, as it probably should.
Blender’s more realistic rendering of Caucasus topography shows a much clearer contrast between the higher Greater Caucasus and the lower Lesser Caucasus.
In the end the purpose of the map determines whether an exaggerated hillshade rendering is a better choice than Blender’s rendering, but if the purpose is photo-realism, the latter surely is a better choice.
Making a shaded relief of the Caucasus region
I can add nothing to what’s already described in the tutorial. I used the SRTM DEM, which I downloaded using this tool by Derek Watkins. At the end of the tutorial you will end up with a rendering that looks similar to Blender’s rendering above.
Adding satellite imagery to the map
Although the greyscale image is pleasing to look at already, I wanted to add satellite imagery from the region to add a little bit of context. Furthermore, large parts of the map are very flat, making it hard to understand what you’re looking at right away: are these plains wedged in between mountains or are these lakes/seas? Where in the world is this actually? Surely satellite imagery can help.
Cloudless Sentinel Imagery
I wasn’t particularly enthused by the idea of having to pick cloudless satellite photos from one of the myriad of providers, but luckily there is a nice cloudless Sentinel 2 version (by EOX IT Services GmbH) that’s suitable.
What should be easy turned out to be quite complicated. Although it may just be because of my lack of experience, exporting a high-resolution version of a WMS map in QGIS turned out to be problematic. For some reason I kept getting empty maps or no exports at all, so I decided to try offline downloaders for WMS maps. After trying many different bits of code and pieces of software, I ended up settling with this Python wms-downloader, which worked excellently after it was configured properly.
The following config.yml did the trick:
service: version: 1.0.0 url: https://tiles.maps.eox.at/wms?service=wms&request=getcapabilities srs: EPSG:4326 format: jpeg layer: s2cloudless-2018 bbox: west: 35.999860000 south: 38.999860000 east: 51.000140000 north: 45.000140000 size: 0.70312500000000000000 resolution: 1024 timeout: 300 projection: EPSG:4326 directory: cloudless-images vrtfile: tiles.vrt tmpfile: /tmp/wms.xml
I am not 100% sure how the size config was supposed to work, but I ended up taking the value from the resolution of the s2cloudless-2018 layer in the WMS XML file (see url in .yml bit above and the code snippet below). These values determine the zoom level, with smaller resolution values resulting in more zoomed in tiles and consequently a higher number of tiles needed to cover your whole area.
<TileSet> <SRS>EPSG:4326</SRS> <BoundingBox SRS="EPSG:4326" minx="-180.000000" miny="-90.000000" maxx="180.000000" maxy="90.000000"/> <Resolutions> 0.70312500000000000000 0.35156250000000000000 0.17578125000000000000 0.08789062500000000000 0.04394531250000000000 0.02197265625000000000 0.01098632812500000000 0.00549316406250000000 0.00274658203125000000 0.00137329101562500000 0.00068664550781250000 0.00034332275390625000 0.00017166137695312500 0.00008583068847656250 0.00004291534423828120 0.00002145767211914060 0.00001072883605957030 0.00000536441802978516 </Resolutions> <Width>256</Width> <Height>256</Height> <Format>image/jpeg</Format> <Layers>s2cloudless-2018</Layers> <Styles/> </TileSet>
The resolution of images is the width and height of each of the tiles you download from the WMS server. bbox contains coordinates for the bounding box from which to download the tiles, for which I recommend to take the extent of the DEM.
Once all the cloudless Sentinel 2 tiles are downloaded, which can take quite a while, you can load the in QGIS, combine them to a virtual raster and export them as a GeoTIFF. This is also the moment to change projection of both the DEM and cloudless imagery if necessary (I used EPSG:900913). Make sure both files (DEM and cloudless imagery) are of the same resolution.
In Blender’s Node Editor view, which has been covered extensively in Daniel’s tutorial, it’s time to add an image texture. You can do this through Add > Texture > Image Texture.
Add an image texture node in the Node Editor view.
Subsequently you can configure the image texture node with your generated Cloudless GeoTIFF, and use the same settings you used before for adding an image texture with the DEM. Instead, rather than connecting it to the Displacement option, we will feed it into the Color setting of the Diffuse BSDF. You should end up with something that looks like the image below.
Overview of the used nodes and connections. This should roughly be your final ‘setup’ before rendering starts.
Before you start to render, make sure you’ve been through the Repeating the Process section here.
And voilà, after a very long render, this is the final outcome.
The Caucasus region. A smaller version of a very large 14,000+ pixels wide original rendering, which you can find below.
There are a few artifacts in the image resulting from the processing done by the creators of the cloudless Sentinel 2 map (see the red spots and purple bars west of the Greater Caucasus), but other than that it looks pretty nice. There’s also a full-resolution, 14,000+ pixels wide (18MB), version available for those of you interested. Feel free to use both versions, but make sure to credit both the maker and the data sources (see bottom of this post).
Moving the rendered image back into GIS software
Since your original files contain georeferencing data, georeferencing the rendered file from Blender is as simple as copying the georeferencing information from the original file to the render. You can do that using some utilities included with a QGIS installation.
listgeo original.tif > original.gtf geotifcp -g original.gtf retouched.tif output.tif
The first line of the above snippet stores georeferencing metadata in a separate .gtf file, the second line copies that metadata to output.tif.
Data used
Sentinel-2 cloudless 2018 by EOX IT Services GmbH is licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License. (Contains modified Copernicus Sentinel data 2017 & 2018)"
NASA JPL (2013). NASA Shuttle Radar Topography Mission Global 1 arc second [Data set]. NASA EOSDIS Land Processes DAAC. doi: 10.5067/MEaSUREs/SRTM/SRTMGL1.003. Downloaded using 30-Meter SRTM Tile Downloader by Derek Watkins.