% preamble
\documentclass[12pt]{article}
\topmargin -0.5 in
\oddsidemargin 0.25 in
\evensidemargin 0.25 in
\textwidth 6.0 in
\textheight 9.0 in

\author{}
\date{May 1, 2026}

\usepackage{graphicx}

\begin{document}

\bibliographystyle{acm} 

\title{Augmenting Grid-Based Contours to Improve Thin Plate DEM Generation}

\maketitle
put your name here!!
Department of Mathematics \& Computer Science
Wheaton College
Norton, MA \  02766
Wm. Randolph Franklin
Electrical, Computer, and Systems Engineering Department
Rensselaer Polytechnic Institute
Troy, NY \  12180

\thispagestyle{empty}
\newpage

We present two new pre-processing techniques that improve thin plate
Digital Elevation Model (DEM) approximations from grid-based contour data.
One method computes gradients from an initial interpolated or approximated
surface.  The aspects are used to create gradient paths that are
interpolated using Catmull-Rom splines.   The computed elevations are
added to the initial contour data set.  Thin plate methods are applied to
all of the data.
The splines allow information
to flow across contours, improving the final surface.
The second method successively 
computes new, intermediate contours in between existing
isolines, which provide additional data for subsequent thin plate processing.
Both methods alleviate artifacts visible
in previous thin plate methods.
The surfaces are tested with published methods to show qualitative
and quantitative
improvements over previous methods. 


Introduction

Digital Elevation Models (DEM) are often used to store three-dimensional
elevation data via a regular grid.  Because DEMs are not available for
many areas and/or because they are storage intensive, 
they are often interpolated
or approximated from sparse data.   We have chosen isoline data from which 
to compute DEMs because contour maps are readily available for many geographic
locations in the
form of Digital Line Graphs (DLG), a standard of the United States Geological
Survey (USGS).  We use a grid-based approach because such methods often produce
DEMs that preserve terrain morphology better than other methods, such as
those using a Triangulated Irregular Network (TIN).
Examples of systems that generate DEMs from contours are TOPOGRID,
available
in ArcInfo, TAPES-C, and TOPOG.

Thin Plate Splines

The notion of minimizing a thin plate to interpolate or approximate a
surface is an old (i.e., ) and trusted technique.
Given $N$ data points, where i 1..N, the differential equation
that models a thin plate is given by:

equation here

Although the thin plate equation has been used in the surface 
reconstruction problem, one problem manifests itself more often 
when using contour line
data as opposed to scattered data.  In simplest terms, a common solution to
the thin plate equation at a particular node
can be stated as the weighted average of the node's neighbors.
Consider contour data depicting hilly or mountainous terrain.
Furthermore, consider a contour line $A$ with a certain elevation, and
a second contour line $B$ which is at the next higher elevation (see top
of Figure~

Improving Thin Plate Methods

Our approach to improving any of the thin plate methods is to add more,
accurate elevation points into the initial contour data set,
without the need for operator intervention (such as adding additional
peak elevation data points).  This results
in a two-stage method where the contours are processed first, creating
a richer data set.   We create additional data through the use of
gradient paths and intermediate contours.
The second stage applies any of the thin plate
methods to yield the final surface.

The algorithm for finding all of the gradient paths in a grid of contours is as
follows:
Compute initial thin plate surface				\\
Compute gradient at each grid point 		\\
For each point $P_{i,j}$ on the grid not visited                        \\
For \\
   create empty path                                                 \\
   repeat                                                            \\
   For \\
   mark $P_{i,j}$ as visited                                       \\
   add $P_{i,j}$ to path                                           \\
   if $P_{i,j}$ contains a valid gradient direction                \\
   For \\
   move to neighboring point $P_{k, l}$ following gradient direction        \\
   until there is no valid neighbor                                  \\
   apply Catmull-Rom spline to path, using contour elevations as knots \\
   copy new computed elevations from path back to grid               \\

Results

Evaluation criteria

The criteria used to assess the quality of a computed DEM are as follows:

The total squared curvature must be as low as possible.  Although natural
surfaces exhibit some curvature, artifacts such as the aforementioned
Gibbs' phenomena contribute greatly to the total curvature.
For N = n2 total points, 
this is found by comparing each computed elevation value to its four
neighbors:

equation here

where each $u$ represents the elevation at the grid location indexed by $i$
and $j$.
The lower the squared curvature, the smoother the surface. This is
useful for direct comparisons of results from different algorithms working on
the same data. 

DEM elevations falling on the original contour lines must have values
equal to (interpolation) or almost equal to (approximation) the contour
labels.
This is measured by the root mean square error ($RMSE$)
of the surface: 

equation here

where something...

Tests

The test case is shown in Figure ,
taken from a USGS DLG of
Crater Lake, Oregon.  The contours were rasterized into a $900 \times 900$ grid.
Elevations are given in feet and the grid spacing is 
in meters; the contour interval is 40 feet.  
As is evident in the contours, the Crater Lake
data has both steep
sections (rising from the lake in the lower left) and flatter sections,
yielding a good test for reconstruction techniques.  

The results of the quantitative
tests two through four are shown in Table~1.
The thin plate with springs 
approximation yields a total squared curvature of 72,678,
which indicates a globally smooth surface in relation to the other surfaces.  
The average curvature is among the highest, however, suggesting that there
must be large areas of high curvature which may be evidence of Gibbs'
phenomena.  
The $RMSE$ of 1.29 is 3.2\% of the contour interval, which 
falls within the standard of five per cent.
Adding tension to the thin plate creates a true
interpolation ($RMSE$ = 0), but at the cost of overall curvature.  
Visually, the gradient paths method produced a less-terraced surface, but
with a few artifacts.  

table here

Finally, to evaluate criteria six, plots were made of the relative heights
of the DEMs of Crater Lake produced by each of the procedures.  The
plot can be seen in Figure .
The frequency of the first height class for thin plate approximation is actually
235479; similarly, the frequency of the first height class for the thin plate
under tension is 320682.  Both of these indicate that the surfaces change
rapidly right at the contour lines.  The overall pattern of the graphs shows
the terracing effect.  
The TOPOGRID procedure shows a very regular pattern.

Conclusions

The thin plate interpolation, approximation, or the addition of
tension may compute smooth DEMs from contour input,
but often Gibbs' phenomena and especially terracing effects are visible.
The problem worsens as contour spacing increases, and is readily apparent
in shaded relief maps.
Automatically adding additional data 
in a pre-processing step through either gradient paths
or intermediate contours visually improve the 
surface created by subsequent thin plate processing, compared to
thin plate methods alone.  
In all cases, the new methods produce smooth surfaces as
shown by the total squared curvature and average curvature, while still
being faithful to the original contour data, as measured by the $RMSE$.  
The profiles and height class plots
show that the new methods create better surfaces in between
contours than previous thin plate procedures alone.
The surfaces compare favorably to ArcInfo's TOPOGRID procedure,
the latter of which exhibited a much higher $RMSE$.


\end{document}