Improvements to DEM Merging with r.mblend

ıs Moreira de Sousa

and Jo

ao Paulo Leit

ISRIC - World Soil Information, Droevendaalsesteeg 3, Building 101, 6708 PB Wageningen, The Netherlands

Swiss Federal Institute of Aquatic Science and Technology (EAWAG),Urban Water Management Department (SWW),

Uberlandstrasse 133, CH-8600, Dubendorf, Switzerland

Keywords:

Digital Elevation Model (DEM), Terrain Analysis, Raster.

Abstract:

r.mblend is an implementation of the MBlend method for merging Digital Elevation Models (DEMs). This

method produces smooth transitions between contiguous DEMs of different spatial resolution, for instance,

when acquired by different sensors. r.mblend is coded on the Python API provided by the Geographic

Resources Analysis Support System (GRASS), being fully integrated in that GIS software. It introduces

improvements to the original method and provides the user with various parameters to ﬁne tune the merging

procedure. This article showcases the main differences between r.mblend and two conventional DEM merge

methods: Cover and Average.

1 INTRODUCTION

In the Geographic Information Systems (GIS) dom-

ain, the representation of terrain elevation has been

predominantly performed using the raster data format,

in what are called Digital Elevation Models (DEMs).

The discretisation of elevation by a regular grid is rat-

her useful in software development, with a direct cor-

respondence to a two dimensional array. This ease

of development has fostered the creation of numerous

spatial analysis methods (de Smith et al., 2015), ma-

king DEMs ever more convenient.

DEMs have traditionally been acquired by stereo-

scopic sensors on board of air-borne or space-borne

vehicles. For decades DEMs remained an expen-

sive and inaccessible type of data. The emergence of

technologies like Light Detection and Ranging (Li-

DAR) sensors, and small and easy to operate Unman-

ned Aerial Vehicles (UAVs) have made the acquisition

of high resolution DEMs considerably simpler and in-

expensive (K

ung et al., 2011).

With multiple DEMs obtained by different met-

hods and at different spatial resolutions available, spa-

tial analysts often face today the need to combine or

merge various of these data sets. However, there is no

obvious method for doing so; a direct merge of over-

lapping DEMs produces artefacts along borders, lea-

ding to inconsistent terrain aspects and slopes (Katzil

and Doytsher, 2005; Luedeling et al., 2007). Spatial

analysis conducted on such merged DEMs inevitably

results in ﬁckle results, be it in view-shed computa-

tion, overland water ﬂow, least cost path, etc.

The MBlend method (Leit

ao et al., 2016) pro-

poses to merge two overlapping DEMs by retaining

the highest spatial resolution DEM and introducing

a smooth transition into the lower resolution DEM.

Modiﬁcations are applied only to the lower resolu-

tion DEM, producing a single DEM that covers the

entire study area with the highest possible accuracy,

while also ensuring smooth transitions between the

original DEMs. r.mblend is an implementation of

the MBlend method, coded in the Python language as

an add-on to the Geographic Resources Analysis Sup-

port System (GRASS). It introduces an advanced and

ﬂexible computation of the transition between DEMs,

that the user may tune through various parameters.

This article compares results obtained using

r.mblend with those of two conventional DEM mer-

ging methods on two different test cases. Section 2

describes the methods used, Section 4 presents the

test cases used for comparison and Section 6 rounds

up the results.

2 RASTER MERGING

2.1 Conventional Methods

Common GIS programmes provide simple functi-

ons to merge raster data sets. They usually re-

quire the inputs to have the same cell size and be in

de Sousa, L. and Leitão, J.

Improvements to DEM Merging with r.mblend.

DOI: 10.5220/0006672500420049

In Proceedings of the 4th International Conference on Geographical Information Systems Theory, Applications and Management (GISTAM 2018), pages 42-49

ISBN: 978-989-758-294-3

the same coordinate system. These simple methods

can be classiﬁed in two different types: Cover and

Average (Eastman, 2012).

Cover type methods do not operate any adjust-

ments to the input DEMs, they are simply superim-

posed. The DEM resulting from this method assumes

cell values of the ﬁrst input across its entire extent and

values from the second input in areas not covered by

the ﬁrst. The resulting DEM can yield signiﬁcant ele-

vation discontinuities along the boundary between the

input DEMs, resulting in erroneous slope and aspect

values (Hickey, 2000).

Average methods assign to the merged DEM the

average elevation within areas where the input DEMs

overlap. Outside the overlapping area the resulting

DEM assumes the value of the existing input; only va-

lues within the overlapping area are changed. Some

of these methods try to tackle the discontinuities is-

sue using a weighted average, as is the case with

IDRISI (Eastman, 2012). A subset of these, usually

referred as Blend methods, go further, using a avera-

ging weighting function that may be linear, smoothed

(e.g. bicubic), or discontinuous; this way more weight

can be given to one of the inputs in certain areas, e.g.

closer to borders. It must be noted though, that these

averaging methods act as low pass ﬁlters, therefore

reducing the accuracy of the resulting DEM.

2.2 MBlend

The MBlend method differs from conventional met-

hods in two essential aspects: it is aware of the dif-

ferent spatial resolution of its inputs and modiﬁes

areas for which only low resolution data are availa-

ble (Leit

ao et al., 2016). This method works by identi-

fying two edges: the border between the low and high

resolution inputs (near edge) and the border around

the area of the low resolution input not overlapping

with the high resolution input (far edge). Points are

set along each of these two edges, those on the near

edge take the difference between the two inputs at the

location; those on the far edge take the value zero (see

Figure 1). From these points is computed a transition

surface, spatially restricted to the area of the low re-

solution input, not overlapping with the high resolu-

tion input. Finally, the transition surface is subtracted

from the low resolution input; the resulting DEM as-

sumes the values of the high resolution input within

its extent and outside the values of the low resolution

input minus the transition surface.

The MBlend method consist in seven essential

steps:

1. Obtain the low resolution extent - this is the ex-

tent of the study area that is only covered by the

low resolution DEM. It can be obtained by vecto-

rising the extent of each DEM and then applying

an intersection.

2. Compute differences - obtained by subtracting the

low resolution from the high resolution DEM.

3. Obtain the near edge - the differences map is vec-

torised into points. A buffer around the low re-

solution extent is then used to select from these

difference points those that lay along the border

between the two rasters (see Figure 1 (a)).

4. Obtain the far edge - the low resolution DEM is

vectorised to points and those along the border are

selected using an internal buffer to the low resolu-

tion extent.

5. Build interpolation points set - the value zero is

assigned to the points in the far edge; it is then

merged with the near edge into a single data set.

6. Interpolate smoothing surface - a new raster sur-

face is created by interpolation using the edges

points data set. The resulting surface smoothly

transitions from the full difference between the

two input DEMs along the near edge towards zero

along the far edge (see Figure 1 (b)).

7. Apply smoothing - the smoothing surface is added

to the low resolution raster. The result is then pa-

tched with the high resolution raster to obtain a

single data set covering the entire study area.

3 Implementation

3.1 The GRASS Add-on Development

Environment

GRASS is a Geographical Information System (GIS)

originally developed by the US Army Corps of En-

gineers with a focus on spatial data management and

analysis (Neteler et al., 2012). It is characterised by a

deep dataset management and archiving structure and

a vast roll of analysis operations, also known as mo-

dules. GRASS manages multiple data-set types: ras-

ter, vector, imagery and voxel (3D).

GRASS was originally written in C, with its mo-

dern structure now also coded in C++. Since 2012, in

the wake of version 6.4.2

, an API to the GRASS C

library was made avaialble for the Python program-

ming language (Sanner et al., 1999). This API greatly

simpliﬁed the development of new GRASS modules,

https://grass.osgeo.org/announces/announce grass642

.html

Improvements to DEM Merging with r.mblend

also facilitating the integration of popular Python li-

braries such as NumPy

or Pandas

A system to host new modules – called “add-

ons” – was also created, whereby third party develo-

pers commit their code to the GRASS repository, thus

making their module(s) automatically available to all

GRASS users. These “add-ons” can be added to every

GRASS installation with the module g.extension.

This module connects automatically to the GRASS

repository, downloads and installs the required bina-

ries or code.

Developed within this environment, r.mblend is

versioned and managed at GitHub

, and released un-

der the European Union Public Licence

3.2 Improvements Over the Original

MBlend Method

The main difference from the r.mblend implementa-

tion to the original method concerns the computation

of the far edge. r.mblend uses by default only those

points in the far edge that are farther away from the

near edge, with the aim of obtaining a geometrically

even transition in the smoothing surface. In detail,

this computation is performed by r.mblend as fol-

lows:

1. Compute a distance map within the low resolution

area relative to the high resolution raster.

2. Vectorise the distance raster into a points data set.

3. Normalise the distance values and select those

above a certain threshold (by default 95% of the

maximum distance).

4. Use an inner buffer to the interpolation area to se-

lect further those points only along the low reso-

lution raster border.

Figure 1 presents these differences to the original

proposal with a simple case. The user is able to adjust

the distance cut off thus dosing the weight of the far

edge on the smoothing surface interpolation.

The r.mblend implementation also provides the

user with the option to use the average difference be-

tween the two input DEMs as the value assigned to

the far edge interpolation points (instead only zero,

as in the original proposal). In this mode the resulting

DEM remains closer to the high resolution input. This

may be useful when the differences between the two

DEMs are spatially uncorrelated.

http://www.numpy.org/

http://pandas.pydata.org/

https://github.com/ldesousa/r.mblend

http://ec.europa.eu/idabc/eupl.html

(a) Edges in original method (b) Smooth surface

Figure 1: Interpolation points edges with original method

(a) and a 95% cut off to the maximum distance (c) and

the respective smoothing surfaces (b and d). Near edge in

green, far edge in yellow.

3.3 Model Parameters

The r.mblend module takes the following arguments:

• high - name of the high resolution DEM;

• low - name of the low resolution DEM (overlap-

ped by the high resolution input);

• output - name of the resulting blended DEM;

• far edge - percentage of the maximum distance

to the high resolution DEM used to determine the

far edge;

• inter points - number of points (from both ed-

ges) to use in interpolation;

• -a - optional ﬂag that indicates to assign the

average difference between the two input rasters

to the far edge (instead of zero).

GISTAM 2018 - 4th International Conference on Geographical Information Systems Theory, Applications and Management

The far edge argument is bounded between 0

and 100; by default is used a value of 95. Values clo-

ser to zero translate into a higher number of points in

the far edge, impacting the shape of the differences

raster.

The Inverse Distance Weighting (IDW) method is

used to interpolate the smoothing surface. This met-

hod takes as a parameter the number of points (from

the two edges) used to interpolate each new cell va-

lue. By default 50 points are used; inter points

provides the user a mean to tweak this value. The

higher the number of points used the smoother is the

resulting smoothing surface; however, it also means a

lengthier computation time.

4 TEST CASES

4.1 A - Lucern

For the ﬁrst test were employed two DEMs represen-

ting an urban catchment in the city of Lucerne in Swit-

zerland. This is a relatively smooth surface, but in-

cluding a number of detailed man-made features. The

lower resolution DEM was obtained with an air-borne

LiDAR sensor and provided to this study by the ofﬁ-

cial cadastral service of the Canton of Lucerne. It has

a cell side of 0.5 metres and a vertical accuracy of ap-

proximately 0.5 metres (Figure 2a). This dataset was

last updated in July of 2012 (Doe, 2014).

The high resolution DEM was obtained with a

conventional camera mounted aboard an electricity

powered, ﬁxed-wing UAV. This UAV made several

ﬂights at an altitude of 114 metres over the study area

in March of 2014. Overlapping images where acqui-

red from different angles allowing for stereoscopic

depth rendition. The resulting DEM has a spatial re-

solution of 0.5 metres and a vertical accuracy of 0.2

metres (Figure 2b).

4.2 B - North Carolina

The second test case was derived from the open spa-

tial data set from North Carolina distributed with

GRASS as sample data

. This data set includes a

10 metres cell side DEM representing relatively rug-

ged terrain with carved valleys and sparse man-made

features. A section of this DEM was cropped to be

used as high resolution input (Figure 3b). The origi-

nal DEM was then converted to a lower spatial resolu-

tion with 60 metres side cells, to which non spatially

correlated noise was added (Figure 3a).

https://grass.osgeo.org/download/sample-data/

(a) LiDAR

(b) UAV stereoscopic

Figure 2: Overlapping DEMs of different resolutions built

from the North Carolina data set.

5 RESULTS

To compare r.mblend with the Cover and Average

methods the Mosaic tool provided with the ArcGIS

software was used. This is tool is able to merge DEMs

using both types of conventional methods. All results

where then assessed with GRASS.

A high pass ﬁlter was used for a ﬁrst assessment

of the merged DEMs produced by each of the three

methods. A 5-by-5 cell ﬁlter was used in order to

highlight zones of transition, e.g. sharp edges, walls

and so forth. Figure 4 presents a detail of these re-

sults for the Lucerne test case. Immediately standing

out is the artiﬁcial step introduced by the Cover and

Average methods along the border between the two

input DEMs. The step is not so marked with Average,

but still present; contrariwise, at a closer inspection a

loss of detail is visible with this method, with various

Improvements to DEM Merging with r.mblend

(a) Low resolution

(b) High resolution

Figure 3: Overlapping DEMs of different resolutions built

from the North Carolina data set.

small transitions in the UAV DEM loosing magnitude.

As for r.mblend it shows no border at all, while pre-

serving the ﬁne detail in the high resolution DEM.

In the North Carolina test case the artiﬁcial step

introduced by the Cover and Average methods is also

present, even if less marked (Figure 5). In this case

the transitions within the larger cell areas stand out

considerably more. Since this is rugged terrain, the 60

metres cells introduce relevant ridges and cliffs. It is

also interesting to observe the effects of the Average

method on the high resolution area, introducing the

artiﬁcial ridges from the 60 meters DEM.

Taking the differences from the original to the out-

put DEMs provides another point of assessment. Fi-

gure 6 shows together the differences from the blen-

ded result to the high resolution DEM and the dif-

(a) Cover method

(b) Average method

Figure 4: Results of high pass ﬁlter applied on merged

DEMs in test case A.

GISTAM 2018 - 4th International Conference on Geographical Information Systems Theory, Applications and Management

(a) Cover method

(b) Average method

Figure 5: Results of a high pass ﬁlter applied on merged

rasters in test case B.

ferences from the result to the low resolution DEM

in the areas not covered by the high resolution input.

This analysis is not presented for the Cover method

since it does not change the inputs. The differen-

ces between Average and r.mblend are striking, ap-

pearing in opposite areas. r.mblend applies changes

only to the low resolution DEM, with a smooth tran-

sition surface; Average leaves the low resolution data

untouched, while applying irregular and many times

severe changes to the high resolution data.

A similar pattern in differences is patent in the

North Carolina test (Figure 7). r.mblend yields again

the smooth transition surface, applying changes solely

to the low resolution input. As before, the Average

method introduces broad changes to the area where

both inputs overlap, in this case coinciding with the

full extent of the high resolution DEM.

(a) Average method

(b) r.mblend

Figure 6: Differences from resulting blended DEM to inputs

in test case A.

Improvements to DEM Merging with r.mblend

(a) Average method

(b) r.mblend

Figure 7: Differences from resulting blended DEM to inputs

in test case B.

6 SUMMARY AND FUTURE

WORK

This article compared the r.mblend GRASS add-

on with conventional methods to merge overlapping

DEMs of different spatial resolution. Using two dif-

ferent test cases, it was possible to assess its advan-

tages on smooth and rugged terrain. r.mblend eli-

minates the steps introduced along the border of the

areas where the merging inputs overlap; these steps

are not so marked in rugged but still present. This

smooth transition is not achieved at the expense of

loss of detail, as the high resolution DEM is left un-

touched. This contrasts particularly with the Average

method, that visibly derides the information from the

high resolution input. r.mblend presents itself as cle-

arly superior alternative to the conventional methods

assessed.

Presently, r.mblend operates on a single execu-

tion thread. All operations conducted are relatively

straightforward, except for the interpolation of the

smoothing surface. For the Lucern case study pre-

sented above, this operation may take in the order

of dozens of minutes. However, it is possible to pa-

rallelise this operation, since there is no dependence

between cells of the resulting surface. The GRASS

Python API provides elementary tools for parallelisa-

tion, spawning GRASS commands as sub-processes.

Therefore, an obvious evolution to r.mblend is to

slice the interpolation area and run the interpolation

independently on each slice.

Other ways of improvement also concern the

smoothing surface interpolation. Alternative methods

beyond IDW can be made available to the user, as so

their respective parameters. This would provide the

user with further degrees of freedom to tune the mo-

dule output.

Finally, r.mblend can also be extended to auto-

matically apply an high pass ﬁlter on the resulting

DEM, providing it as a secondary output. This is a

useful asset to assess to quality of the resulting DEM,

either visually or in more elaborate analysis.

REFERENCES

de Smith, M. J., Goodchild, M. F., and Longley, P. A.

(2015). Geospatial Analysis: A Comprehensive Guide

to Principles, Techniques and Software Tools - Fifth

Edition, chapter Geocomputational methods and mo-

deling, pages 625–672. Winchelsea Press.

Doe, R. (2014). GIS Kanton Luzern. https://rawi.lu.

ch/themen/gis_kanton_luzern. Accessed: 30-10-

2014.

Eastman, J. (2012). IDRISI Selva. Clark University, MA,

USA.

Hickey, R. (2000). Slope angle and slope length solutions

for GIS. Cartography, 29(1):1–8.

Katzil, Y. and Doytsher, Y. (2005). Spatial rubber sheeting

of dtms. In Proceedings of the 6th Geomatic Week

Conference, Barcelona, Spain, volume 811.

ung, O., Strecha, C., Beyeler, A., Zufferey, J.-C., Flore-

ano, D., Fua, P., and Gervaix, F. (2011). The accuracy

of automatic photogrammetric techniques on ultra-

light uav imagery. In UAV-g 2011-Unmanned Aerial

Vehicle in Geomatics, number EPFL-CONF-168806.

Leit

ao, J., Prodanovi

c, D., and Maksimovi

C. (2016). Im-

proving merge methods for grid-based digital eleva-

tion models. Computers & Geosciences, 88:115 –

131.

Luedeling, E., Siebert, S., and Buerkert, A. (2007). Filling

the voids in the srtm elevation modela tin-based delta

GISTAM 2018 - 4th International Conference on Geographical Information Systems Theory, Applications and Management

surface approach. ISPRS Journal of Photogrammetry

and Remote Sensing, 62(4):283–294.

Neteler, M., Bowman, M. H., Landa, M., and Metz, M.

(2012). Grass gis: A multi-purpose open source gis.

Environmental Modelling & Software, 31:124–130.

Sanner, M. F. et al. (1999). Python: a programming lan-

guage for software integration and development. J

Mol Graph Model, 17(1):57–61.

Improvements to DEM Merging with r.mblend