Elsevier

Pattern Recognition

Volume 36, Issue 7, July 2003, Pages 1529-1544
Pattern Recognition

A topology-preserving parallel 3D thinning algorithm for extracting the curve skeleton

https://doi.org/10.1016/S0031-3203(02)00348-5Get rights and content

Abstract

We introduce a new topology-preserving 3D thinning procedure for deriving the curve voxel skeleton from 3D binary digital images. Based on a rigorously defined classification procedure, the algorithm consists of sequential thinning iterations each characterized by six parallel directional sub-iterations followed by a set of sequential sub-iterations. The algorithm is shown to produce concise and geometrically accurate 3D curve skeletons. The thinning algorithm is also insensitive to object rotation and only moderately sensitive to noise. Although this thinning procedure is valid for curve skeleton extraction of general elongated objects, in this paper, we specifically discuss its application to the orientation modeling of trabecular biological tissues.

Introduction

In this paper, we introduce a new 3D topology-preserving thinning algorithm specifically designed for the reduction of complex 3D digital images of elongated objects, such as trabeculated biological tissues, to 1D voxel representations (curve skeletons). The algorithm is based on a rigorously defined classification procedure and follows a directional sub-iteration approach. Central to the algorithm is a set of simple albeit rigorous and efficient conditions for identifying ambiguous simple sets to be operated upon in a parallel thinning procedure. The algorithm is shown to yield skeletons only moderately sensitive to noise and to the object rotation.

Thinning, defined as an iterative reduction process that transforms a discrete structure to its lower dimension, is an important operation in both 2D and 3D digital image processing. Earlier work on thinning has concentrated mostly on 2D images, see for example the extensive literature survey in Ref. [1]. However, with 3D images becoming more readily available, research has increasingly focused on the problem of 3D thinning, especially in biomedical image reduction [2], [3], [4], [5]. Lobregt et al. [6] first presented a 3D thinning algorithm based on the preservation of Euler characteristics; Ma and Sonka [4] introduced a curve thinning algorithm using deleting templates; Saha et al. [7] developed a parallel algorithm that preserves the numbers of object components, cavities, and tunnels; recently, Borgefors et al. [8] introduced a partially parallel algorithm based on the D6 distance labeling. However, most of the existing algorithms are difficult to implement due to either complex checking conditions or lengthy matching templates involved.

A number of useful applications of 3D thinning, such as 3D image condensation and recovery [8], topological segmentation [9], and computer-aided manufacturing [10], have been discussed in the literature. In this paper, we use thinning to determine the local structural orientation of trabeculated biological tissues, such as found in the embryonic heart [11] and in trabecular bone [12], [13]. According to experimental observations, the mechanical properties of trabeculated tissues depend greatly on the structural orientation. In particular, the accurate measurement of structural orientation is critical for the finite element modeling of the trabeculated myocardial tissue for biomechanical studies of early development [14]. Denslow et al. [15] developed an eigenvalue decomposition method for numerical characterization of the directionality of trabeculated tissues from 2D images. Their approach, however, does not solve the problem of computing the heterogeneous structural orientation of trabeculated tissues in 3D. In fact, structural orientation can be computed only on the 1D representation derived from the original 3D image. To solve this problem, we develop a topology-preserving 3D thinning algorithm for computing directly the curve skeleton of a trabeculated tissue from its 3D image, based on the definition of 3D simple (removable) voxels [7], [16].

In the following section, we first review the basic notions and definitions of a 3D digital image; then we introduce the notion of simple voxel and establish a directional approach for solving the problem of ambiguous simple sets in a parallel thinning cycle; finally, we present a simple partially parallel thinning algorithm and discuss several test results.

Section snippets

Basic definitions

A 3D digital image is a finite array of solid cubic elements, or voxels. In the binary 3D image considered in this work, a voxel has a value of either 1 or 0, called black or object voxel, and white or background voxel, respectively. We denote this image by the triple (P,B,W), where P is a set consisting of all voxels within the digital image, B (BP) is the set of all black voxels, and W(W=P−B) denotes the set of all white voxels. In turn, a voxel p is uniquely defined by an integer triple (x1,

Theoretical considerations

The basic notion of thinning the 3D binary image (P,B,W) is to progressively change black voxels to white voxels (i.e., to progressively delete black voxels) until the minimal set of black voxels, or the skeleton set, is reached. A thinning process is said to be topology-preserving if it does not change the numbers of components, cavities, and tunnels in B. As shown in Ref. [6], the checking of topological changes following the deletion of a black voxel p can be restricted to N26(p). The

Preliminary

The basic operation in the examination and removal of a black voxel p is the α-component counting in various subsets of Nα(p). We adopt the following simple recursive scanning procedure (α-scan) for counting α-components of the voxel set S.

  • 1.

    Set the initial scanning front S′:S′⊂S.

  • 2.

    Set initial scanned set A:A=S′.

  • 3.

    Find the α-envelope of S′ within S:E=Nαe(S′)∩S.

  • 4.

    Update first the scanning front S′ and then the scanned set A:S′=E−A, A=AE.

  • 5.

    Repeat from step 3 to step 4 until S′=∅.

The final scanned set A is

Pre-processing

The proposed thinning procedure yields a curve skeleton which is relatively insensitive to the presence of noise (Fig. 8). However, an excessive amount of noise may cause the creation of spurious noisy branches and significant local skeleton distortions (Fig. 9a).

In order to contain noise effects, we preprocess the original image with a noise-reduction procedure. The maximum and minimum filters are commonly adopted to reduce noise in digital images [25]. However, application of a high-order

Examples

We first apply the proposed thinning procedure to 2D and 3D binary images extracted from a stack of grayscale confocal images of the trabeculated heart of a HH stage 21 chick embryo. The original grayscale image resolution is 768×512, at 3.3μm per pixel and the distance between consecutive section is 18μm. The typical trabecular structure of the myocardial wall is clearly visible in the 3D image of the stack, viewed frontally at the apex of the heart—Fig. 10.

From the first image of the stack,

Discussion

We generally follow Saha et al.'s [7] classification of simple voxels. However, we achieve the parallel deletion of simple voxels with a directional sub-iteration approach different to Saha et al.'s sub-field approach. Our parallel thinning algorithm, which is based on a set of additional checking condition, allows for the derivation of curve skeleton from the 3D image with a single thinning pass.

Recognizing that our checking conditions for parallel thinning are overly restrictive, and thus are

Conclusion

In this paper, we present a topology-preserving thinning algorithm for reducing a 3D binary image to its 1D (curve) representation. The algorithm consists of sequential thinning iterations each characterized by six parallel directional sub-iterations followed by a set of sequential sub-iterations. Examples indicate that the algorithm yields concise and geometrically sound curve skeletons for 3D images. Unlike other multi-pass thinning algorithms [8], [20], our curve thinning algorithm produces

Acknowledgements

We wish to thank Larry Taber for continuous advice and support throughout this research. We also thank Scott Hollister for useful discussions and for providing the images of the trabecular bone. Finally, we gratefully acknowledge the support of NIH through grants R01 46367 and R01 HL64347-02.

About the AuthorWENJIE XIE received his M.Sc. degree (Mechanics, 1995) from Beijing University and Ph.D. degree (Biomedical Engineering, 2002) from University of Rochester. He is currently a senior development engineer at ANSYS, Inc. (Canonsburg, PA). His research interests include computational biomechanics, volume image based geometrical modeling, and motion tracking.

References (28)

  • P.K Saha et al.

    3d digital topology under binary transformation with applications

    Comput. Vision Image Understand.

    (1996)
  • Y.F Tsao et al.

    A parallel thinning algorithm for 3-d structure

    Comput. Graph. Image Process.

    (1981)
  • K Palagyi et al.

    A 3d 6-subiteration thinning algorithm for extracting medial lines

    Pattern Recognition Lett.

    (1998)
  • L Latecki et al.

    Well-composed sets

    Comput. Vision Image Understand.

    (1995)
  • Cited by (84)

    • Improved skeleton extraction method considering surface feature of natural micro fractures in unconventional shale/tight reservoirs

      2018, Journal of Petroleum Science and Engineering
      Citation Excerpt :

      By looping from the edges of the pore space to the interior, this method can continuously identify the simple points in the pore space and delete them so that the medial axis skeleton model can be extracted. In reality, this algorithm is realized by mapping the weight topology parameters of the direct space neighboring domains (Han et al., 2003; Dolean et al., 2008; Xie et al., 2003; Lohou and Bertrand, 2002; Zhou et al., 2007). Traditional medial surface extraction methods include the maximum ball algorithm (Arand and Hesser, 2017; Silin and Patzek, 2006), simulated annealing algorithm (Zhao et al., 2007), and reduction method (Toriwaki and Yoshida, 2009).

    • 2D skeleton extraction based on heat equation

      2018, Computers and Graphics (Pergamon)
      Citation Excerpt :

      Thinning/peeling based methods. Most of the algorithms of skeleton generation are based on the thinning technique or the wavefront/grassfire transform [26,29,30]. These algorithms shrink the boundary iteratively to determine an object’s skeleton.

    View all citing articles on Scopus

    About the AuthorWENJIE XIE received his M.Sc. degree (Mechanics, 1995) from Beijing University and Ph.D. degree (Biomedical Engineering, 2002) from University of Rochester. He is currently a senior development engineer at ANSYS, Inc. (Canonsburg, PA). His research interests include computational biomechanics, volume image based geometrical modeling, and motion tracking.

    About the AuthorROBERT P. THOMPSON schooled at Rice University (B.A. 1968, Biology/Chemistry) and Colorado State University (Ph.D. 1974, Biophysics) and is now Professor of Cell Biology and Anatomy at the Medical University of South Carolina. His research focuses upon supra-cellular organization and morphogenesis of the embryonic heart.

    About the AuthorRENATO PERUCCHIO holds a Doctorate in Aeronautical Engineering (University of Pisa, Italy, 1977) and a Ph.D. in Civil Engineering: Structural (Cornell University, 1983), and is Professor of Mechanical Engineering and Biomedical Engineering at the University of Rochester. His research interests include computational biomechanics, solid and geometric modeling, and structural mechanics. He is currently working on an interdisciplinary investigation of growth and morphogenesis of the embryonic heart.

    View full text