Table Of ContentWavelet Multi-scale Edge Detection
for Extraction of Geographic
Features to Improve Vector Map
Databases
Marlin L. Gendron
(US Naval Research Laboratory)
Dr. Juliette W. Ioup
(University of New Orleans)
Althoughnumerousatsmallergeographicscales,vectordatabasesoftendonotexistatthe
moredetailed,largerscales.Apossiblesolutionistheuseofimageprocessingtechniquesto
detect edges in high-resolution satellite imagery. Features such as roads and airports are
formedfromtheedgesandmatchedupwithsimilarfeaturesinexistinglow-resolutionvector
map databases. By replacing the old features with the new more accurate features, the
resolution of the existing map database is improved. To accomplish this, a robust edge
detectionalgorithmisneededthatwillperformwellinnoisyconditions.Thispaperstudies
andtestsonesuchmethod,theWaveletMulti-scaleEdgeDetector.Thewavelettransform
breaksdownasignalintofrequencybandsatdifferentlevels.Noisepresentatlowerscales
smoothes out at higher levels. It is demonstrated that this property can be used to detect
edgesinnoisysatelliteimagery.Onceedgesarelocated,anewmethodwillbeproposedfor
storingtheseedgesgeographicallyso thatfeaturescanbeformed and paired with existing
features inavectormap database.
KEY WORDS
1. Maps}Charts. 2. Database. 3. Wavelets.
1. INTRODUCTION. While ancient navigators struggled to locate features
crudelynotedontheirhand-drawnmaps,moderntravellersenjoytheluxuryofhigh-
speedcomputers,globalpositioningsystems(GPS),andaerospacetechnology.Early
adventurers relied on the wind to carry their ships and on their eyes to spot land.
Today,oureyesextendaroundtheglobeintheformofmulti-milliondollarsatellites,
and our current location is tracked by a symbol moving across a computer screen
(Lohrenz et al., 1997).
From tourists using moving map GPS systems in some rental cars to on-board
computers that guide smart weapons to their targets, the demand for accurate and
reliable maps is increasing. There are many digital map databases available to the
military, commercial sector, and the public. Off-the-shelf software offers worldwide
coverage,andmapsofalmosteverymajorcityareavailableviatheInternet.Without
a doubt, maps are important, but maintaining and updating these map databases is
costly and time consuming. In many cases, humans are still required to examine
79
Report Documentation Page Form Approved
OMB No. 0704-0188
Public reporting burden for the collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and
maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information,
including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 1204, Arlington
VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to a penalty for failing to comply with a collection of information if it
does not display a currently valid OMB control number.
1. REPORT DATE 3. DATES COVERED
JAN 2000 2. REPORT TYPE 00-00-2000 to 00-00-2000
4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER
Wavelet Multi-scale Edge Detection for Extraction of Geographic
5b. GRANT NUMBER
Features to Improve Vector Map Databases
5c. PROGRAM ELEMENT NUMBER
6. AUTHOR(S) 5d. PROJECT NUMBER
5e. TASK NUMBER
5f. WORK UNIT NUMBER
7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION
Naval Research Laboratory,Code 7440.1,Stennis Space Center,MS,39529 REPORT NUMBER
9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S)
11. SPONSOR/MONITOR’S REPORT
NUMBER(S)
12. DISTRIBUTION/AVAILABILITY STATEMENT
Approved for public release; distribution unlimited
13. SUPPLEMENTARY NOTES
Journal of Navigation 53:1, pp. 79-92. January
14. ABSTRACT
15. SUBJECT TERMS
16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 18. NUMBER 19a. NAME OF
ABSTRACT OF PAGES RESPONSIBLE PERSON
a. REPORT b. ABSTRACT c. THIS PAGE Same as 14
unclassified unclassified unclassified Report (SAR)
Standard Form 298 (Rev. 8-98)
Prescribed by ANSI Std Z39-18
80 M.L.GENDRON AND J.W.IOUP VOL.53
satelliteimageryanddigitisefeaturessuchasroadsandairportsintopointsthatform
vector representations. These points are stored in new map databases or used to
update existing ones (Edmonds, 1998).
This paper explores one possible method of improving the resolution of existing
vectormapdatabasesbyusingsemi-automatedimageprocessingtechniques.Features
from high-resolution satellite images are extracted via a Wavelet Multi-scale Edge
Detector(WMED)andplacedinageographicbitmap(Gendronetal.,1997).From
thebitmap,vectorsareformedandidenticalfeaturesarelocatedandupdatedinthe
map database.
2. TRANSFORMS. Manyproblemsinsciencecannotbesolveddirectly.Tools
are required to alter these problems in such a way that proven mathematical
techniques yield results. For example, certain operations performed on equations in
thetimedomainaredifficultandoftenimpossible.AFouriertransformation,which
transforms a function from the time domain to the frequency domain, makes many
oftheseoperationsthatareimpossibleinthetimedomainmanageablebyperforming
them in the frequency domain. After operations are completed in the frequency
domain, an inverse transform provides the means to return to the time domain
(Bracewell, 1986).
The Fourier transform operates on a signal in the time domain and generates
coefficients in the frequency domain. These coefficients describe the composition of
thesignalintermsofinfinite-lengthsineandcosinefunctions.Analternativetousing
suchlongwavesistouseshortwavescalledwavelets.Portionsoftheinputsignalthat
mostresemblethewaveletproducethehighestcoefficientsintheoutput(Strangand
Nguyen,1996).LiketheFouriertransform,thewavelettransformanditsinverseare
defined for both continuous and discrete signals. The discrete wavelet transform
(DWT)liesattheheartoftheedgedetectorusedonthesatelliteimagesstudiedinthis
paper.Theelementsintheinputseriesarepixelswhosemagnitudesarebetween0and
255.Thepositionoftheelementdoesnotdenotetime,butratherindicatesitslocation
in the image.
Thewavelettransformisexplainedmorefullylaterinthispaper,butthetime-scale
(spatial-scale) property that wavelets exhibit is worth mentioning now. By using a
pairofshortwaveletfiltersandapplyingthemtoaninputsignal,levelonecoefficients
are produced whose position in time relates directly back to the input. The wavelet
filters are then rescaled, shifted, and applied to the output to produce level two
coefficients. These coefficients are lower in resolution, but maintain their position
relative to the original time series (Misiti and Misiti, 1996).
The time-scale (spatial-scale) properties of the DWT can be used to extract true
featuresfromanoisyimage,becauserandomwhitenoisepresentathighresolutions
tends to smooth out at lower resolutions. Peaks detected in lower level coefficients,
whicharenotpresentinupperlevels,arediscardedasnoise.Peaksspatiallycloseat
all or most levels are retained as true edges. This is the basic concept behind the
WMED (Hajj et al., 1996a).
3. DIGITAL MAPS. Amapcanbestoredonacomputerindigitalformsimply
as a picture; for example, a digital raster map image file. Most map files originate
frompaperchartsscannedintothecomputerbyascanner.Whenpixelsintheimage
correspondtopointsonasurface,suchastheearth,theimageisgeo-rectifiedorgeo-
NO.1 WAVELET MULTI-SCALE EDGE DETECTION 81
registered.Severalmapfilesstoredtogetherinastructuredwayformamapdatabase.
Digital raster map databases exist at many different resolutions, or map scales, and
define how individual map files can be pieced together to form complete maps
(Lohrenz, 1991).
Frequently, when these maps are viewed, the current position of the cursor in
latitude and longitude is displayed on the screen. Cursor readout implies that the
computer knows the location of every pixel in the image. This is true and
accomplished by geo-rectifying the image. This geo-rectifying process occurs by
taking control points, at known latitude and longitude points on the earth, to
generatecoefficients.Thesecoefficientsareusedtocomputenewxandypositionsin
theimage.Thenumberofcontrolpointsandtheirprecisiondeterminestheaccuracy
of the raster map. The geographic extent of the image and its size in pixel space
determines the map scale of the image (Hager, et al., 1990; Landrum, 1989).
Raster maps make a lot of sense to the human brain. People are able to extract
features–such as roads and highways–intuitively and determine their location
relativetootherfeatures.Althoughthecomputer‘knows’thelocationofeachpixel
at a given resolution, the computer is unable to extract, identify, and determine the
location of a specific feature (such as a road) easily.
Analternativetostoringarasterpictureistostoreavectordatabase.Avectormap
database stores individual features, such as airports, roads and highways as points.
Aroadfeatureisstoredinthedatabaseasasetofpoints,includingabeginningand
endingpositionandintermediatepointsalongitspathwhereturnsoccur.Theroad’s
position, magnitude and direction are reconstructed from these points. A digital
vector map database contains many feature types and is considered accurate to a
given map scale (Edmonds, 1998). Today, most features in a vector database are
extractedfromrastermapsorpapercharts,whichisacostlyprocedure.Whenmore
accurateinformation(e.g.,higherresolutionsatelliteimagery)becomesavailable,the
vectoring process must start afresh to create a new vector database or correct an
existing one. A common method of correcting an existing database is to update
individual features or feature layers. This paper will explore a method of updating
existing features by using a WMED.
4. WAVELETS. Fourieranalysistransformsasignalbasedintimetoonebasedin
frequency. Although Fourier analysis is an extremely useful tool, it has one serious
disadvantage. Time information disappears when a signal is transformed to the
frequency domain. When analysing the frequency component, it is impossible to
know when an event occurred in time (Bracewell, 1986). The advent of the Short-
Time (or Finite) Fourier Transform helped to solve the problem by windowing the
signal and mapping it into a function oftime and frequency. Although astep in the
right direction, the chosen window size remains the same for all frequencies. This is
a serious drawback. A new windowing technique solved the problem by allowing
smaller windows for higher frequencies and larger windows for lower frequencies.
Thisnewtechnique,calledwaveletanalysis,usesthewavelettransformtobreakdown
signals into frequency bands at different scales or levels (Misiti and Misiti, 1996).
Awavelettransformcontainsafilterbankthatseparatesasignalintohighandlow
frequency bands by passing the signal through a highpass and lowpass filter. The
highpass filter is a moving difference filter capable of detecting sharp changes in the
data, and the lowpass filter is a moving average filter that smoothes the data. The
82 M.L.GENDRON AND J.W.IOUP VOL.53
w(t) w(2t-1)
w(2t)
1
1 1
t t
0 0
Haar wavelet Haar rescaled
and shifted
Figure1.HaarWaveletrescaledandshifted.
L 2 Approximation
S
H 2 Detail
Figure2.OneLevelDWT.
output of the lowpass and highpass filters are coefficients referred to as the
‘approximation’ and ‘detail’ coefficients, respectively (Strang and Nguyen, 1996).
Awaveletisashortwavewithanaveragevalueofzero.Byscalingthewaveletwith
a scale factor and shifting its position, approximation and detail coefficients can be
generatedatdifferentlevels.FortheDWT,thepositionandscalefactorarebasedon
powersoftwo.Thisiscalledadyadicapproach(StrangandNguyen,1996).Thescale
factor, derived from the scaling function, is associated with the lowpass filter. The
actual wavelet corresponds to the highpass filter. The formula that describes this
relationshipandproducesthescalingfunctionandthewaveletequationisknownas
the dilation equation. The dilation equation is a two-scale difference function that
governs the wavelet’s ability to scale and shift. For example, dilating and shifting a
simple boxcar function creates the Haar wavelet (Morgan, 1997). Both the Haar
wavelet and its rescaled version, shown in Figure 1, exhibit an extraordinary ability
to find edges in noisy signals (Hajj et al., 1996a).
Figure 2 shows a one-level DWT. The lowpass filter produces approximation
coefficientsbyconvolvingtheinputsignalwithlowpassfiltercoefficientsderivedfrom
the dilation equation. Reversing the lowpass filter and alternating the signs of the
filter coefficients produces the highpass filter. Convolving this filter with the input
signal results in the detail coefficients (Misiti and Misiti, 1996).
Thenumberofapproximationcoefficientsplusthenumberofdetailcoefficientsis
twice that of the input. To avoid doubling the number of terms after each level, the
coefficientsaredown-sampledbyafactoroftwoattheoutputofeachfiltertoremove
redundant information. Down-sampling by a factor of two means every other term
is removed. In any given filter output, the lost terms normally cannot be recovered,
NO.1 WAVELET MULTI-SCALE EDGE DETECTION 83
L 2 A
L 2
H 2 D2
S
H 2 D1
Figure3.TwoLevelDWTTree.
scale 1
DWT scale 2 Local
S filter scale 3 maxima Binary map
bank scale 4 detector
scale 5
Figure4.FiveLevel1-DWMED.
but remaining components from both filters can always be used to recover fully the
original signal (Strang and Nguyen, 1996).
Morelevelsareobtainedbyiteration.Successiveapproximationsaredecomposed
sothatonesignalisbrokendownintomanylowerresolutioncomponentsorlevels.
Figure 3 shows this wavelet decomposition tree.
By up-sampling the coefficients and reversing the filtering process, the signal is
reconstructed. Up-sampling is accomplished by adding zeros between every term,
thustakingtheplaceofthecoefficientremovedbydown-sampling.Althoughhalfof
thecoefficientsarereplacedbyzerosateachlevel,prefectreconstruction(PR)occurs
ifthecoefficientsarenotmanipulatedinanyotherway,andthesamewaveletisused.
PRmeanstheoriginaloutputsignalisreturnedunaltered(StrangandNguyen,1996).
Reconstruction, via the Inverse Discrete Wavelet Transform (IDWT), is not needed
fortheWMEDstudiedinthispaper.Onlythedetailcoefficientsthatcontainthehigh
frequency information are used to find the edges. The magnitude of the detail
coefficientsrelatesdirectlytohowcloselythatpartofthesignalresemblesthewavelet
atthegivenscale.Atanygivenlevel,thenumberoflargecoefficients(thosewithhigh
energy) is small (Strang and Nguyen, 1996).
Groupsoflargemagnitudedetailcoefficients,calledlocalmaxima,possiblydefine
edges (Hwang and Chang, 1996). Small coefficients are likely to represent noise and
should be removed using wavelet shrinkage or thresholding. There are two types of
thresholding: hard and soft. Hard thresholding is accomplished by choosing a
thresholdvalueandsettingallcoefficientsthatfallbelowthethresholdtozero(Laine
andZong,1996).Alternatively,insoftthresholding,coefficientsbelowthethreshold
are set to zero, but the larger coefficients are manipulated for good compression
performanceornoisereduction(LaineandZong,1996).Althoughsoftthresholding
84 M.L.GENDRON AND J.W.IOUP VOL.53
is preferable for most applications, hard thresholding is used for the WMED
implementedinthispaperbecausethefinaloutputofthedetectorisbinary(Gendron
et al., 1997).
Forthispaper,coefficientsbelowthethresholdaresetto0,andcoefficientsabove
arechangedtoavalueof1.Ifthenewoutputvalueis1,anedgeexists,otherwiseit
does not. The edge information encoded in this manner can easily be stored in a
bitmap. Later in this paper, a novel approach is introduced to equate an edge bit
stored in a bitmap to a geographic location on the earth (Gendron et al., 1997).
5. ONE-DIMENSIONAL SIGNALS. Thefirst step in designinga WMED
istodeterminethefilterlengthofthelevelonehighpassfilter(Figure4).Ifthefilteris
tooshortinlength,theresolutionwillbetoohigh,andsomeedgeswillbemissed.The
Haar wavelet, used in the DWT for this paper, has only two coefficients. By
lengthening the wavelet, or adding more coefficients, the filter overlap with the data
is increased, and the resolution reduced (Hajj et al., 1996a).
To determine the optimal length of the Haar highpass filter for this research, a
simpleone-dimensional(1-D)blocksignal,containing512samples,wasrunthrough
a one-dimensional DWT using a Haar wavelet of length 4. At level one, the high
frequency coefficients should reflect all the edges for the noiseless case via local
maxima. An examination of the output detail coefficients revealed that some edges
weremissing.Theresolutionwastoohigh.Thelengthofthefilterwasmodifiedusing
upsampling, which increases the filter length and doubles the frequency response.
(Strickland, 1997). Adding zeros between each filter coefficient increases the filter
length from two to four terms. The original signal was run through the DWT with
the new length. Now all edges were reflected in the level one output.
With the length problem addressed, the block signal was run through a five level
WMED. Figure 5 shows the block function and detail coefficients produced by the
Signal
d1
d2
d3
d4
d5
Figure5.FiveLevelDWTofNoiselessSignal.
NO.1 WAVELET MULTI-SCALE EDGE DETECTION 85
Signal
d1
d2
d3
d4
d5
Figure6.LocalMaximaofDetailCoefficients.
Signal
Edges of
signal
Figure7.Outputof1-DWMED.
DWTateachlevel.Atthehighestresolution,orlowestlevel,d1,thepeaksaresharp
and well defined, but become broader and wider at higher levels. However, even at
the lowest resolution, d5, the edges are easily detected by local maxima.
Next, local maxima are found at each level by hard thresholding to a value of 0
or1.Figure6showsaplotofthebinarysequenceofindicatorsateachlevel.Notice
most indicators match up with the step edges at all levels.
Therulefortracingedgesacrosslevelsistoconnectanindicatoratthelowestlevel
to the spatially closest indicator at the higher levels. If indicators are found at most
or all levels, it is considered an edge. Indicators at the lowest level that do not have
spatially close matches at upper levels are ignored (Hajj et al., 1996a). The smallest
stepintheinputsignal,theseventhedgefromtheleftmostpointattfl0,isdetected
atlevelone,butlostinsubsequentlevels.Itispossibletodetectthisedgebylowering
thethresholdvalue.Thisisacceptableinthiscase,buttherelaxationofthethreshold
can have profound implications in noisy conditions.
Figure7showsthefinaloutputoftheone-dimensionalWMEDafterindicatorsare
matched up. Notice that the seventh edge is not detected.
86 M.L.GENDRON AND J.W.IOUP VOL.53
Noisy signal
d1
d2
d3
d4
d5
Figure8.FiveLevelDWTofSignalwithNoise.
Noisy signal
d1
d2
d3
d4
d5
Figure9.WMEDIndicatorsforSignalwithNoise.
6. NOISE. To determine how well the WMED would perform in noisy
conditions, Gaussian noise is added to the same one-dimensional (1-D) block
function. Figure 8 shows the detail coefficients from the output of the DWT at five
levels. Notice the positions of the edges are far more difficult to discern.
NO.1 WAVELET MULTI-SCALE EDGE DETECTION 87
Noisy signal
Edges of
noisy signal
Figure10.Outputfrom1-DWMED.
1-D Wavelet
Col scanline multiscale
edge detector
Binary
Image feature
map
1-D Wavelet
Row scanline multiscale
edge detector
Figure11.2-DWMED.
Applying the DWT to the noisy signal, detail coefficients at level one (d1) reflect
edges and noise. It is not until higher levels, around d3, that the noise subsides and
trueedgesemerge.Figure9showsthelocalmaximaatfivelevels.Theindicatorsare
found by thresholding the maxima to values of 1 or 0.
Next, indicators, or possible edges, at d1 are matched up with spatially close
indicatorsathigherlevels(d2-d4).Therearetwoapproachestothematchingprocess,
andboth aread-hoc.Inonemethod,featuresare reportedif indicatorsare foundat
mostlevels.Anotherapproachdetectsafeaturewhenindicatorsarefoundatallthe
levels (Hajj et al., 1996a). Because a higher threshold value has been used to reduce
noise,onlythefirstfourlevelsareexamined.Figure9addscredibilitytothisdecision,
since the lowest resolution level (d5) shows some edge dropout. The resulting edges
are shown in Figure 10. The results are identical to the noiseless case, including the
absent seventh step edge.
7. TWO-DIMENSIONAL IMAGES. The WMED can be extended to two
dimensions (2-D) by scanning individual rows and columns of an image. Each
scanline passes through the one-dimensional WMED to form a binary map. Row
scanlines detect vertical edges, and column scanlines detect horizontal edges. By
combiningthetwosetsofedgesintoatwo-dimensionalbinarymap,allpossibleedges
emerge including diagonal ones. Figure 11 shows the 2-D WMED (Hajj et al.,
1996b).
Ifthetwo-dimensionalinputimagehappenstobeageo-rectifiedrastermapimage,
edges in the output equate to geographical locations in the image. By knowing the
