TWELVE DIFFERENT INTERPOLATION METHODS: A CASE STUDY

OF SURFER 8.0
Chin-Shung Yang* Szu-Pyng Kao* Fen-Bin Lee** Pen-Shan Hung**
*National Chung Hsing University, Taichung, Taiwan, ROC yangchin@ms1.emax.net.tw

**Feng Chia University, Taichung, Taiwan, ROC

majorlfp@mail.apol.com.tw

KEY WORDS: DTM (digital terrain model), Interpolation

ABSTRACT
SURFER is a contouring and 3D surface mapping program, which quickly and easily transforms random surveying

data, using interpolation, into continuous curved face contours. In particular, the new version, SURFER 8.0,
provides over twelve interpolation methods, each having specific functions and related parameters. In this study, the
5 meter DTM was used as test data to compare the various interpolation results; the accuracy of these results was
then discussed and evaluated.

1. INTRODUCTION following paper, we will introduce every

How to adequately use exist numerous interpolation method in SURFER.
wide-distributed height points has been an important 2. SURFER INTERPOLATION METHODS
topic in the field of spatial information. Normally, 2.1 The Inverse Distance to a Power method
contouring is the way to accurately describe the The Inverse Distance to a Power method is a
terrain relief by means of Scenography, Shading, weighted average interpolator, which can be either
Hachure and Layer Tinting in a way which is best fit exact or smoothing. With Inverse Distance to a
to the habit of human vision. Power, data are weighted during interpolation, so
Presently, discretely collected height points have to that the influence of one point, relative to another,
be interpolated to form curved faces, the selection of declines with distance from the grid node. Weighting
spatial interpolation methods decide the quality, is assigned to data through the use of a weighting
accuracy and follow-up analysis applications. power, which controls how the weighting factors
Interpolation methods are used here to calculated the drop off as distance from the grid node increases.
unknown heights of interested points by referring to The greater the weighting power, the less effect the
the elevation information of neighboring points. points, far removed from the grid node, have during
There are a great many commercial interpolation interpolation. As the power increases, the grid node
software, however, most of them are tiny and value approaches the value of the nearest point. For a

designed to solve specific problems with limited smaller power, the weights are more evenly
versatility. The SURFER is a software developed by distributed among the neighboring data points.

US GOLDEN company, and the newest version 8.0 Normally, Inverse Distance to a Power behaves as an
contains up to 12 interpolation methods to been free exact interpolator. When calculating a grid node, the
chosen for various needs. Users are suggested to first weights assigned to the data points are fractions, the
have the basic understanding of every interpolation sum of all the weights being equal to 1.0. When a
methods before he or she can effectively select particular observation is coincident with a grid node,

parameters in every interpolation methods. In the the distance between that observation and the grid
node is 0.0, that observation is given a weight of 1.0; while attempting to honor your data as closely as

all other observations are given weights of 0.0. Thus, possible. Minimum Curvature is not an exact

the grid node is assigned the value of the coincident interpolator, however. This means that your data are

observation. The smoothing parameter is a not always honored exactly.

mechanism for buffering this behavior. When you

2.4 The Modified Shepard's Method
assign a non-zero smoothing parameter, no point is
The Modified Shepard's Method uses an inverse
given an overwhelming weight, meaning that no
distance weighted least squares method. As such,
point is given a weighting factor equal to 1.0. One of
Modified Shepard's Method is similar to the Inverse
the characteristics of Inverse Distance to a Power is
Distance to a Power interpolator, but the use of local
the generation of "bull's-eyes" surrounding the
least squares eliminates or reduces the "bull's-eye"
observation position within the grid area. A
appearance of the generated contours. Modified
smoothing parameter can be assigned during Inverse
Shepard's Method can be either an exact or a
Distance to a Power to reduce the "bull's-eye" effect
smoothing interpolator. The Surfer algorithm
by smoothing the interpolated grid.
implements Franke and Nielson's (1980) Modified

2.2 The Kriging Method Quadratic Shepard's Method with a full sector search

Kriging is a geostatistical gridding method that has as described in Renka (1988).

proven useful and popular in many fields. This

method produces visually appealing maps from 2.5 The Natural Neighbor Method
irregularly spaced data. Kriging attempts to express The Natural Neighbor method is quite popular in

trends suggested in your data, so that, for example, some fields. What is the Natural Neighbor

high points might be connected along a ridge rather interpolation? Consider a set of Thiessen polygons

than isolated by bull's-eye type contours. Kriging is a (the dual of a Delaunay triangulation). If a new point

very flexible gridding method. The Kriging defaults (target) were added to the data set, these Thiessen

can be accepted to produce an accurate grid of your polygons would be modified. In fact, some of the

data, or Kriging can be custom-fit to a data set, by polygons would shrink in size, while none would

specifying the appropriate variogram model. Within increase in size. The area associated with the target's

SURFER, Kriging can be either an exact or a Thiessen polygon that was taken from an existing

smoothing interpolator, depending on the polygon is called the "borrowed area." The Natural

user-specified parameters. It incorporates anisotropy Neighbor interpolation algorithm uses a weighted

and underlying trends in an efficient and natural average of the neighboring observations, where the

manner. weights are proportional to the "borrowed area". The

Natural Neighbor method does not extrapolate
2.3 The Minimum Curvature Method
Minimum Curvature is widely used in the earth contours beyond the convex hull of the data

sciences. The interpolated surface generated by locations (i.e. the outline of the Thiessen polygons).

Minimum Curvature is analogous to a thin, linearly

elastic plate passing through each of the data values, 2.6 The Nearest Neighbor Method

with a minimum amount of bending. Minimum The Nearest Neighbor method assigns the value of

Curvature generates the smoothest possible surface the nearest point to each grid node. This method is
useful when data are already evenly spaced, but need between data points. The original points are

to be converted to a SURFER grid file. Alternatively, connected in such a way that no triangle edges are

in cases where the data are close to being on a grid, intersected by other triangles. The result is a

with only a few missing values, this method is patchwork of triangular faces over the extent of the

effective for filling in the holes in the data. grid. This method is an exact interpolator. Each

Sometimes with nearly complete grids of data, there triangle defines a plane over the grid nodes lying

are areas of missing data that you want to exclude within the triangle, with the tilt and elevation of the

from the grid file. In this case, you can set the Search triangle determined by the three original data points

Ellipse to a certain value, so the areas of no data are defining the triangle. All grid nodes within a given

assigned the blanking value in the grid file. By triangle are defined by the triangular surface.

setting the search ellipse radii to values less than the Because the original data are used to define the

distance between data values in your file, the triangles, the data are honored very closely.

blanking value is assigned at all grid nodes where Triangulation with Linear Interpolation works best

data values do not exist. when your data are evenly distributed over the grid

2.7 The Polynomial Regression Method area. Data sets containing sparse areas result in

Polynomial Regression is used to define large-scale distinct triangular facets on the map.

trends and patterns in your data. Polynomial

Regression is not really an interpolator because it 2.10 The Moving Average Method
does not attempt to predict unknown Z values. There The Moving Average method assigns values to grid

are several options you can use to define the type of nodes by averaging the data within the grid node's

trend surface. search ellipse. To use Moving Average, a search

ellipse must be defined and the minimum number of
2.8 The Radial Basis Function Interpolation
data to use, specified. For each grid node, the
Method
neighboring data are identified by centering the
Radial Basis Function interpolation is a diverse
search ellipse on the node. The output grid node
group of data interpolation methods. In terms of the
value is set equal to the arithmetic average of the
ability to fit your data and produce a smooth surface,
identified neighboring data. If there are fewer, than
the Multiquadric method is considered by many to
the specified minimum number of data within the
be the best. All of the Radial Basis Function methods
neighborhood, the grid node is blanked.
are exact interpolators, so they attempt to honor your
data. You can introduce a smoothing factor to all the
2.11 The Data Metrics Methods
methods in an attempt to produce a smoother
The collection of data metrics methods creates grids
surface.
of information about the data on a node-by-node

2.9 The Triangulation with Linear Interpolation basis. The data metrics methods are not, in general,

Method weighted average interpolators of the Z-values. For

The Triangulation with Linear Interpolation method example, you can obtain information such as:

in SURFER uses the optimal Delaunay triangulation. a) The number of data points used to interpolate each

This algorithm creates triangles by drawing lines grid node.

If the number of data points used are fairly equal at
each grid node, then the quality of the grid at each

grid node can be interpreted. 3. EXPERIMENT DESIGN

b) The standard deviation, variance, coefficient of The purpose of this experiment was primarily to take
variation, and median absolute deviation of the data a spatial interpolation method, with the assistance of
at each grid node. SURFER software, to move the spatial interpolation

These are measures of the variability in space of the of DTM from 40 m to 5 m, and to report
grid, which is important information for statistical comparisons and results. The actual operation made
analysis. use of the Chen-Yu-Lan river region the ground 5 m

c) The distance to the nearest data point. DTM of the results, and not pass by the mathematics
For example, if the XY values of a data set are to calculate but the direct to take the 40 m DTM ,
sampling locations, the Distance to the nearest data again then this 40 meter manuscript input the
metric can be used to determine new sampling interpolation of SURFER software put to compute 5

locations. A contour map of the distance to the m DTM , and 5 meter manuscript ratio then right
nearest data point can quantify where higher acquire its interpolation to put the analysis of error
sampling density may be desired. margin of the calculation. The flow charts and sketch

maps of the experimental district were designed as

2.12 The Local Polynomial Method follows:

The Local Polynomial method assigns values to grid

nodes by using a weighted least squares fit, with data

within the grid node's search ellipse.

 Figure 1:Flowchart of experimental design

Figure 2:Sketch map of experiment district

4. COMPARISONS OF ACCURACY AND experiment the scope of district is big, and not easily

ANALYSIS present the difference of the result of calculation on

This experiment tries to make the district to on the the screen of computer, so only pick part of districts

spot measure 5 m DTM to manuscripts with the and establish the result of Layer Tinting and shadow
region of Chen-Yu-Lan river, and the total area is to display as follows :

roughly 215 square kilometer, because of the

Source Inverse Distance to a Power Kriging

Minimum Curvature Modified Shepard's Method Natural Neighbor

 Nearest Neighbor Polynomial Regression Radial Basis Function

Triangulation with Linear Moving Average Data Metrics

Interpolation

Local Polynomial

The personal computer used in this experiment was actions, the reference of the performance. Below

a Pentium 4,2GHz memory 768MB. record every then for since 40 m DTM put to 5 m, and reduce the

kind of time that interpolation method use when acquisition the covariance of result that get with 5 m
putting calculation, with the accuracy of conduct and manuscript.

Interpolation method Use time Min. Max. Mean STD. Dev.

Inverse Distance to a Power 00:23:06 -110.34 126.19 -0.023 8.227

Kriging 01:55:29 -104.93 128.63 -0.017 3.602

Minimum Curvature 00:02:55 -472.88 510.79 -0.023 10.641

Modified Shepard's Method 00:01:08 -119.6 144.34 -0.016 3.586

Natural Neighbor 00:11:57 -106.65 126.61 -0.019 3.880

Nearest Neighbor 00:01:21 -143.5 175.39 -0.068 8.495

Polynomial Regression 00:00:02 -697.88 956.95 0.304 304.533

Radial Basis Function 02:22:57 -114.57 132.4 -0.016 3.472

Triangulation with Linear 00:00:06 -106.31 126.76 -0.018 4.061

Interpolation

Moving Average 00:00:10 -83.86 113.16 0.015 13.612

Data Metrics 00:18:08 -138.04 154.02 -1.069 20.638

Local Polynomial 00:25:24 -161.73 149.95 -0.113 20.133

From statistics that above the form can get, all the aspect of mean error, are impracticable margin
methods use the time most longer is Radial Basis for error of Polynomial Regression biger outside, the
Function, and the shortest Polynomial Regression rest of worth and all very small.
discrepancy very big, especially at the standard error With the scope of error margin worth for horizontal

aspect, every kind of method have the obvious sit the mark, and the quantity is for vertical sit to
margin as well, its inside than is interesting of is, at mark as follows of histogram:

Local Polynomial Inverse Distance to a Power Kriging

Minimum Curvature Modified Shepard's Method Natural Neighbor

Nearest Neighbor Polynomial Regression Radial Basis Function

Triangulation with Linear Moving Average Data Metrics

Interpolation

Top each statistical chart of form in, put the Regression, can then obviously find its distributing

method due to each interpolation horizontal sit to the appearance not good, the rest sketch all make the
mark the scope all different, therefore can't directly rule symmetry to distribute.
judge mutually the density of its distribute the 5. CONCLUSION
difference, but in the statistical chart of Polynomial In this paper we compare 12 different interpolation
methods. For each method, we analyze its 2. Barnett V. Interpreting multivariate data [M].
applicability, algorithm, efficiency and advantage. NewYork, 1981:21.

There is no absolutely best method but only the 3. Franke R. Scattered Data Interpolation: Test of
optimal choice under certain circumstances. One Some Methods [J]. Mathematics of Computations,
should first review the characteristic and theorem of 1982,33(157):181.
each method as well as the property and spatial 4. Franke R, Nielson G. Smooth Interpolation of

analysis of data before he or she can successfully Large Sets of Scattered Data [J]. International
select a spatial interpolation method which is Journal for Numerical Methods in Engineering,
relatively best in certain situation. However, the 1980,15(2):1691.

outcome should be evaluated by conscientious 5. Lee DT, Schachter BJ. Two Algorithms for
experiences. Constructing a Delaunay Triangulation,

References: International Journal of Computer and Information

1. Briggs IC. Machine Contouring Using Minimum Sciences [J].1980,9(3):219.

Curvature [J], Geophysics, 1974,39(1):39. 6. SURFER on-line manual