2004 
  
data 
e the 
goes 
d M ; 
90 k 
ring 
?rate 
line 
4 THEINTERPOLATION PROCESS 
In order to interpolate the elevation at a point, the algorithm lo- 
cates an edge of the triangle of the dual of the Voronoi diagram 
(i.e. Delaunay triangulation) for a set of points and oriented line 
segments in which the given point lies. Then, it determines if the 
given point is a vertex or it lies on the line segment supporting a 
vertex (oriented half line segment) of the triangle in which it lies. 
If this is the case, it means that the given point is one of the data 
points, or it lies on a data line segment. In the first case, the ele- 
vation at that point is the same as the elevation at the data point. 
In the second case, we assume that the elevation at the point can 
take both of the elevations linearly interpolated from the eleva- 
tions at the extremities of the two oriented line segments. If it 
is not the case, the algorithm computes the list of objects (cor- 
responding to Delaunay triangulation vertices) from which the 
point would steal some area if it was inserted in the Delaunay 
triangulation. This is done without inserting the point in the De- 
launay triangulation. We are using the Quad-edge data structure 
(Guibas and Stolfi, 1985) for storing both the Delaunay triangu- 
lation and the Voronoi diagram. Starting from the located edge, 
and visiting the three edges of the enclosing triangle, the algo- 
rithm tests whether the given edge is safe with respect to the in- 
terpolation point (i.e. the edge would remain after the addition of 
the interpolation point in the Delaunay triangulation (Guibas and 
Stolfi, 1985)). If an edge is safe, then it is added to the circular list 
of safe edges enclosing the interpolated point. If an edge is not 
safe, the edge having the same origin and immediately after it in 
the clockwise orientation (the edge pointed by its Oprev operator 
(Guibas and Stolfi, 1985)), and the edge having the same desti- 
nation and immediately after it in the counterclockwise orienta- 
tion (the edge pointed by its Dnext operator (Guibas and Stolfi, 
1985)) are successively checked. The safe edges detected by this 
algorithm are inserted in a circular list in the counterclockwise 
order, and the “previous” pointer points to the previous edge in 
the counterclockwise orientation. The origin of the edge pointed 
by the Rot operator (Guibas and Stolfi, 1985) of each edge of this 
circular list is a natural neighbour of the point being interpolated. 
Once the list of enclosing safe edges is computed, the area stolen 
by the interpolated point to each associated neighbour and the in- 
terpolated elevation are computed and the total area and the sum 
of interpolated elevations are maintained. Once all the neigh- 
bours have been visited, the sum of the interpolated elevations is 
divided by the total area in order to get the interpolated elevation 
of the point being interpolated. The construction of the Quad- 
Edge data structure requires O (n log n) worst case time, where 
n is the number of sampled points. Each interpolation requires 
O (log n) amortized worst case time. 
5 EXPERIMENTAL RESULTS 
In Figures 5,6,7, and 8, we can see the results of the use of our 
natural neighbour interpolation based on the Voronoi diagram for 
a set of points and oriented line segments for simulating linear 
discontinuities like faults, bridges, brakes of slope, dams and 
faults in topographic modelling. Figures 5 and 6 show the mod- 
elling of a vertical fault as the result of the interpolation based 
on the Voronoi diagram for a set of points and one pair of ori- 
ented line segment. Figures 7 and 8 show the use of a smoothing 
function on top of the interpolation in order to get a smoother 
surface. The smoothing function we used is an Hermitian inter- 
polation function (Davis, 1975). It applies the following change 
on the local coordinates producing the new local coordinates: 
us (M) = Su, (M )! — 2u (M )* . Figure 7 has a different 
altimetric scale than Figures 5 and 8 in order to show the effect 
of the smoothing function on the top of the interpolated surface. 
1117 
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Vol XXXV, Part B3. Istanbul 2004 
  
Figure 5: An example of line Voronoi diagram based natural 
neighbour interpolation, the data objects are marked by plain 
disks 
  
Figure 6: A view of the same surface with illumination (the 
source is at 30 degrees above in the north east direction) 
6 DISCUSSION 
We have shown in this paper an extension of the natural neigh- 
bour interpolation from the ordinary Voronoi diagram to the Vo- 
ronoi diagram for a set of points and oriented line segments. We 
have presented an example of use of this interpolation technique 
for digital terrain modelling. 
7 ACKNOWLEDGMENTS 
This research work has received the financial support of NSERC 
Discovery Grant and University of Calgary Starter Grant to the 
first author and an Alberta Ingenuity Fund Fellowship to the sec- 
ond author. 
REFERENCES 
Anton, F., Gold, C. and Mioc, D., 1998. Local coordinates and 
interpolation in a Voronoi diagram for a set of points and line seg- 
ments. In: The Voronoi Conference on Analytic Number Theory 
and Space Tillings, pp. 9-12. 
Berger, M., 1977a. Géométrie. Vol. 1. CEDIC, Paris. Actions de 
groupes, espaces affines et projectifs. [Actions of groups, affine 
and projective spaces].