情報地質 第26巻 第2号 060-061頁 2015年
Geoinformatics, vol.26, no.2, pp.060-061, 2015

Estimating the Optimal Parameter Values for BS-Horizon DEM
Generation Algorithm in Flat Lowland Area
Thi An TRAN *, Shinji MASUMOTO**, Venkatesh RAGHAVAN *, Susumu NONOGAKI ***Tatsuya NEMOTO** and Go YONEZAWA** Graduate School for Creative Cities, Osaka City University, 3-3-138 Sugimoto, Sumiyoshi-ku,Osaka 558-8585, Japan. E-mail: tranthian.gis@gmail.com.
** Graduate School of Science, Osaka City University, 3-3-138 Sugimoto, Sumiyoshi-ku, Osaka558-8585, Japan.
*** Geological Survey of Japan, National Institute of Advanced Industrial Science and Technology,
Central 7, 1-1-1 Higashi, Tsukuba, Ibaraki 305-8567, Japan.

Key words : BS-Horizon, DEM interpolation, Bi-cubic spline, equality-inequality
constraint

1. Introduction
This study aims to generate high resolution DEM offlat lowland area in Danang city, Vietnam for thepurpose of flood inundation and flood hazard mapping.
Topography plays an important role in flood modelling,hence it is necessary to create reliable DEM. An
interpolation method based bi-cubic spline algorithm
was used for generate 5m resolution DEM. The
generated DEM was used as input data for further
extraction of geomorphology and flood risk zoning of thisstudy area.

2. Methodology
2.1 Study area and data set
Study area comprises of 98km2 of lowland area
located in the south of Danang City with elevation
ranging from 0m to 10m. The data used for DEM
generation are in-situ elevation points surveyed byDepartment of Natural Resource and Environment,
Danang city in 2009, include 79,600 elevation points.
This alluvial lowland area prone to experienced several
serious flood events in past.
2.2. DEM generation process
In order to generate high resolution DEM of thestudyarea,
anequality-inequality
constrainedinterpolation method implemented as BS-Horizon
program (Nonogaki et al., 2012) was applied. The in-situ
spot height data that includes 79,600 elevation points
was taken into account as equality constraint tointerpolate topographic surface. Further, the inequality
constraint was characterized by the elevation range from0m to 10m of the study area. In general, the possible
constraint of an elevation zi that is obtained at a point(xi, yi ) can be expressed as follow:
f (xi, yi)  zi = 0(1)(2a)f (xi, yi)  zi < 0f (xi, yi)  zi > 0(2b)Equality constraint (1) is used in case that the surface
passes through the known elevation point. Inequality
constraint (2a) is applied in cases that the surface passes
under the point, and the inequality constraint (2b) isused in cases that the surface passes above the point(Nonogaki et al., 2012).
This equality-inequality constrained data was used
as input data for BS-Horizon interpolation algorithm.
BS-Horizon is a bi-cubic spline approximation that relies
on an exterior penalty function to determine the optimalsurface. Exterior penalty function controls the balance
between the smoothness and the violation of theconstraint using a penalty parameter (alpha). The
function can be shown in the equation 3:
Q(x, α) = J(x) + α R(x)(3)

Where Q(x) is called an exterior penalty function, J(x)evaluates the smoothness of the surface, R(x) evaluates
the degree of violation of constraints and α control aweight balance between J(x) and R(x).
In order to run BS-Horizon the equality-inequality
constrained data needs to be prepared. Firstly, wecolored the study area in a certain color ramp. A color
table was prepared includes corresponding color ramp
and the elevation constraint (0m to 10m). This colored
map was saved as bitmap 24 bit image and imported toan application namely BMP2DAT. In this application thecalculation region, the output cell size and the referencecolor table need to be specified. The output data is a text
file of data points with long/lat location and inequality
elevation constraints. Next, the spot height data wasconverted to equality constraint in the same format as
the inequality data. Finally, both of these files weremerged into an equality-inequality constraint data (incsv format) and input to BS-Horizon program. The
process of BS-Horizon DEM generation can be describedin Figure 1.
2.3 Estimating optimal parameters for BS-Horizon
According to BS-Horizon algorithm implemented byNonogaki et al. (2012), it is observed that two main
parameters effecting on the output DEM surface are M(Mx or My, the number of sub-regions that constitute thesurface domain, in which bi-cubic spline function works)
and penalty parameter (alpha). In order to obtain thebest parameter values for BS-Horizon, batch processing
was carried out by setting J(x) and R(x) values according
to M and α and examine the results are shown in a graph(Figure 2). Three cases of Mx and My were used in this61

Table 1. Comparison of DEM with different parameter settings.
Alpha
2.07E+06
6.95E+06
2.34E+07
2.07E+06
6.95E+06
2.34E+07
2.07E+06
6.95E+06
2.34E+07

Mx189
189189643
643
643
1072
1072
1072

Min
1.9
3.3
15.4
0.7
2.0
2.5
0.3
0.4
0.3

My123
123
123
438
438
438
730
730
730

Max
13.5
16.1
31.5
10.6
10.9
10.6
10.2
10.3
10.1

Mean
3.05
3.06
3.07
3.03
3.03
3.03
3.03
3.03
3.03

STD
1.98
2.0
2.1
1.99
2.0
2.0
1.98
1.99
1.99

RMSE
0.40

0.16

0.11

Table 2. BS-Horizon DEM compared with other methods.
Figure 1. DEM generation workflow.

DEM
BS-Horizon DEM
IDW DEM
RST DEMMin
0.32.471099.0Max
10.1
10.6
10.3

RMSE
0.11
0.29
0.21(a) Mx = 189, My = 123

R(x) – case1
J(x) – case1

R(x) – case2
J(x) – case2

R(x) – case3
J(x) – case3

Figure 2. Representation of J(x) and R(x) according todifferent M and alpha.

experiment, namely: case 1: Mx = 189, My = 123; case 2:
Mx = 643, My = 438 and case 3: Mx = 1072, My = 730. Alphawas used as values between 1 to 1010.
Figure 2 shows the transition tendency of J(x) andR(x) with the progress of iterative computation (Noumi,2003). It can be observed that there are three differentzones. Zone A shows the steady increase of J and steadydecrease of R, while in zone B, J increases continuously
and R goes to minimum values or the generated surface
shows smaller errors. In zone C, J continuously increases
and R shows spikes. It is suggested that B zone is thesuitable for generating optimal surface. Several possible
DEMs in this B zone were selected to calculate somestatistical parameters pertaining to error estimation.
Based on the given constraint with elevation range from0 to 10m and the root mean square error (RMSE) ofDEMs, the optimal surface was selected (Table 1). Also,
it is understandable to see from Figure 3 that DEM
surface in case 3 shows more details of topography thanthe one in case 1.

3. Result
Considering the given inequality constraint, it isclear from Table 1 that the best parameters for
BS-Horizon is α = 2.34107, Mx = 1072 and My = 730. As aresult, the output DEM has elevation from 0.3-10.06m,
the RMSE is only 0.11m. Compared to otherinterpolation methods, such as Inverse distanceweighting (IDW) or Regularized spline with tension(RST), Horizon-DEM provides a better surface
estimation and the vertical errors are also minimized asshown in (Table 2).

(b) Mx = 1072, My = 730

Figure 3. DEM surface two case of Mx, My.

4. Conclusion
Field survey elevation point data combined withinequality constraint used in BS-Horizon interpolation
method can provide reliable DEM. Experiment withdifferent conditions Mx, My and alpha parameter inBS-Horizon reveals that the most suitable for DEM
generation for this flat lowland area are Mx = 1072, My=730 and alpha = 2.34107. Horizon DEM has better
RMSE compare to other DEM interpolation methods
such as IDW or RST.

References
Nonogaki, S., Masumoto, S. and Shiono, K. (2012) Griddingof geological surfaces based on equality-inequality
constraints from elevation data and trend data.
International Journal of Geoinformatics, vol.8, no.4,
pp.49-60.
Noumi, Y. (2003) Generation of DEM using inter-contour
height information on topographic map. Journal ofGeosciences, Osaka City University, vol.46, pp.217-230.