The Scale Space Skeletonization
Luciano da Fontoura Costa and Leandro F. Estrozi
This page presents a scale space skeletonization technique introduced by Costa and Estrozi [13]. Exact dilations  which are related to the Eikonal Equation \({\partial \vec{C} \over \partial t} = F \cdot \hat{n} \) (i.e. a curve \(\vec{C}\) propagating with speed \(F\) in the normal direction) [2]  are used in order to calculate the exact Euclidean distance transform and to propagate labels associated to each of the pixels in a 2D contour. In this page you will also find the complete distance data structure (the SEDR) used in the exact dilation process, animations illustrating the described methodology, pseudocodes for multiresolution skeletonization, and an executable code (Windows) illustrating the involved concepts and techniques.
Exact Dilations, Exact Distance Transform And Label Propagation
The concept of Exact Dilations , introduced in [2], allows all the possible dilations of a binary 2D shape in the orthogonal lattice to be obtained. This technique, which is extremely simple, is based on a sorted list containing all the sorted distances that are found inside a ball of limited radius in the orthogonal lattice. This data structure, which is henceforth abbreviated as SEDR  standing for Sorted Exact Distances Representation [2], also includes relative vectors to all the possible points exhibiting the respective distance. Click here for an animated illustration of this data structure and to download a file containing this structure.
The exact dilations of a binary shape can be obtained by scanning the shape contour elements, in sequential fashion (requires a contour following algorithm), for each of the distinct distances in the SEDR structure. At each of such interaction, the relative vectors for the current distance are added to the coordinates of the contour element and, in case such positions are empty, they are marked as 1 (the original image is initialized with zeros). If the distance value, instead of the value 1, is assigned to the positions, the exact Euclidean distance transform is obtained. If the value of the label of the respective contour element is assigned, the process of label propagation is implemented. These processes are illustrated in the pseudo code to be found here. The main feature of such an label propagation implementation is that it provides the exact solution of the Eikonal equation for constant speed.
The figure above illustrates a original binary shape (left) and the respective label propagation (right) obtained by the above described label propagation algorithm.
Scale Space Shape Representation Without Border Shifting
Once the labels are propagated, a simple difference scheme involving the pixels in a 3x3 mask in the label image is used to produce a difference image (illustrated by the gray level traces in the figure below, left).
By thresholding this difference image at successive integer values (click here for a pseudocode describing this process), a family of scale space skeletons is obtained for the original contour that has remarkable properties, such as the fact of being connectedof8 and that the previous portions of the skeleton are not shifted along the image space as details are removed (by increasing the threshold). In addition, the obtained skeletons have unitary width and are completely connected.
The above animation (left) presents a scale space family of skeletons. The evolution of the blue lines corresponds to the skeletons obtained for subsequent threshold values.
These scale space skeletons can be straightforwardly used to define a scale space representation of the original shape presenting the important property that the frontiers of large scale portions of the object are not shifted when smaller scale details are removed. Such scale space reconstructions can be easily obtained by combining (by logical or) filled circles centered at each of the skeleton elements, having as respective radius the value of the exact distance transform of the original contour at those points. The above animation (right) presents the scale space reconstruction of the original shape.
The SEDR Data Structure
i 
j

A compressed ASCII file containing the SEDR structure can be found here.
Applications
The proposed multiresolution skeletonization method, which can be extended to 3D [2], allies the advantages of scale space analysis to the power of skeletons as shape representation. Therefore, the proposed technique has many interesting properties and presents great potential for applications, including neuronal shape analysis [3], particle system analysis (forthocoming paper), pattern recognition, among many others.
Algorithms (Scilab)
Distance Transform and Label Propagation  Scale Space Skeletonization 
// Calculates the distance transform and label propagation // according to exact Euclidean distances SEDR is a matrix // containing the subsequent Euclidean distances and relative // position vectors (see main page for explanation) // // Luciano da F Costa and Leandro F Estrozi // Cybernetic Vision Research Group, USP, Brazil // 25th Jan 2000 === Copyright CVRG // // Obs : list_x and list_y are the lists which store the // original shape contour coordinates obtained using // a stardard contour following algorithm. for x = 1 to N {Initializes img_dist and img_lbl} for y = 1 to M img_dist(x,y) = 1 img_lbl(x,y) = 1 end end for i = 1 to Nd {Scans the subsequent Euclidean dist.} d = SEDR(i,1); for j = 1 to SEDR(i,2) {Scans all respective pos. vectors} dx = SEDR(i,2j+1); dy = SEDR(i,2j+2); for k = 1 to L {Scans the contour elements} x = list_x(k)+dx; y = list_y(k)+dy; if img_dist(x,y) = 1 then img_dist(x,y) = d; img_lbl(x,y) = k; end end end end 
// Obtains scale space skeleton for the spatial scale // corresponding to a threshold T. The result is the // NxM image img_skl. // // Luciano da F. Costa and Leandro F. Estrozi // Cybernetic Vision Research Group, USP, Brazil // January 2002 === Copyright CVRG for x = 1 to N { Initializes img_skl} for y = 1 to M img_skl(x,y) = 0 end end K = number of contour elements for x=2 to N1 for y=2 to M1 max = 0 for every pixel (i,j) 4neighbour of (x,y) dif = img_lbl(i,j)  img_lbl(x,y) dif = MIN(dif,ndif) max = MAX(max, dif) end if max > T img_skl(x,y) = 1 end end end 
Demonstration Software
This demo program (for Windows) illustrates the potential of the multiresolution methodology. This program also includes high accuracy estimation of the curvature along the obtained multiresolution representation of the shapes, which uses the technique proposed in [4,5]. A fully operational version of this program can be obtained free of charge by contacting us.

Related Manuscripts:
[1]
Costa, L.da F.
and
Estrozi, L.F.
Multiresolution shape representation without border shifting.
Electronics Letters, Volume 35, Issue 21, Pages 18291830.
(Oct, 1999)
[2]
Costa, L.da F.
Multidimensional Scale Space Shape Analysis.
In Proc.International Workshop on SyntheticNatural Hybrid Coding and Three Dimensional Imaging,
Santorini, Greece, 1517, Pages 214217.
(Sep, 1999)
[3]
Costa, L.da F.
Campos, A.G.; Estrozi, L.F.; Rios, L.G. and BOSCO, A.
A BiologicallyMotivated Approach to Image Representation and Its Applications to Neuromorphology.
In IEEE International Workshop on Biologically Motivated Computer Vision,
Santorini, Seoul, Korea.
(2000)
[4]
Costa, L.da F.
Particle Systems Analysis by Using Skeletonization and Exact Dilations.
Part. Part. Syst. Charact.,
Volume 16, Issue 6, Pages 273277.
(1999)
[5]
CesarJr., R.M.
and
Costa, L.da F.
Piecewise linear segmentation of digital contours in \(O(n\log(n))\) through a technique based on effective curvature estimation.
Journal of RealTime Imaging,
Volume 1, Issue 6, Pages 409417.
(1995)
[6]
CesarJr., R.M.
and
Costa, L.da F.
Towards effective planar shape representation with multiscale digital curvature analysis based on signal processing techniques.
Pattern Recognition,
Volume 29, Issue 9, Pages 15591569.
(1996)
[7]
Costa, L.da F.
Robust Skeletonization through Exact Euclidean Distance Transform and its Application to Neuromorphometry.
RealTime Imaging,,
Volume 6, Pages 415431.
(2000)