File Exchange

image thumbnail

Skeleton3D

version 1.12.0.1 (16.6 KB) by Philip Kollmannsberger
Calculates the 3D skeleton of an arbitrary binary volume using parallel medial axis thinning.

44 Downloads

Updated 05 Dec 2018

View License on GitHub

Editor's Note: This file was selected as MATLAB Central Pick of the Week

Skeleton3D: Parallel medial axis thinning of a 3D binary volume
This code calculates the 3D medial axis skeleton of an arbitrary 3d binary volume. It is an optimized MATLAB implementation of the homotopic thinning algorithm described in [1]. We developed it to quantify the network of cell processes in bone [2], but it should work on images of any tubular or filamentous structures. An example volume (testvol.mat) is included, along with an example script (Test_Skeleton3D.m). Any comments, corrections or suggestions are highly welcome.
Usage:
skel = Skeleton3D(bin)
where "bin" is a 3D binary image, and "skel" the resulting image containing only the skeleton voxels, or
skel = Skeleton3D(bin,mask)

to mask all foreground voxels in 'mask' from skeletonization, e.g. to preserve certain structures in a image volume.

For additional cleanup, e.g. pruning of short branches, please use my Skel2Graph3D package on MATLAB File Exchange.

This code is inspired by the ITK implementation by Hanno Homann [3] and the Fiji/ImageJ plugin by Ignacio Arganda-Carreras [4]. If you include this in your own work, please cite our original publicaton [2].

Philip Kollmannsberger 09/2013
philipk@gmx.net

References:

[1] Ta-Chih Lee, Rangasami L. Kashyap and Chong-Nam Chu
"Building skeleton models via 3-D medial surface/axis thinning algorithms."
Computer Vision, Graphics, and Image Processing, 56(6):462–478, 1994.

[2] Kollmannsberger, Kerschnitzki et al.,
"The small world of osteocytes: connectomics of the lacuno-canalicular network in bone."
New Journal of Physics 19:073019, 2017.

[3] http://hdl.handle.net/1926/1292

[4] http://fiji.sc/wiki/index.php/Skeletonize3D

Cite As

Kollmannsberger, Kerschnitzki et al., "The small world of osteocytes: connectomics of the lacuno-canalicular network in bone." New Journal of Physics 19:073019, 2017.

Comments and Ratings (23)

Cody Crosby

alphin j

Hi Philip,
How can I use it for getting curves from a point cloud(set of dense curves)?

Martin Denk

Hello Philip,

in the file p_oct_label.m
you checked the adjacency of octant 8 according to 17 twice.
""""
idx = cube(:,17) == 1;
if any(idx)
cube(idx,17) = label(idx);
cube(idx,:) = p_oct_label(4,label(idx),cube(idx,:));
end
"""
I think you can delete one of them.

How much effort will it cost to generate surface skeletons like in [1]?

Thanks for answering me Philip, I got it! Now, I wanna do something more. I have another file with the diameter size value of the vessel in it (in bmp format) and which is of the same matrix size as the 3D-Skeleton that I obtain. My idea is to color the points of the skeleton with a sort of colormap where hot values means big vessel size and blue small vessel size. Any idea to help me in this challenge would be greatfull for me. Thanks!

Alex Cocco

Thank you, Phillip!

@Alex: the easiest way is to compute the distance transform on the original volume using bwdist() and to multiply the result with the final skeleton. This should give you the skeleton with a local thickness estimate stored as pixel intensity.

@Julien: this sounds like "skel" is not of type logical - try "skel = logical(skel)", or replace line 42 by "cands=zeros(width,height,depth);".

I have the following error appearing using a 3D matrix as input and full of 0 and 1... Do you have any idea on how I can fix that? I really not see why this message appear since my matrix is already supposed to be logical and comes from slices of a CT-scan...
"Error using false
Input following 'like' is not a logical array.

Error in Skeleton3D (line 42)
cands=false(width,height,depth,'like',skel);"

Alex Cocco

Is there a way to have the code output an estimate of local thickness at each point along the skeleton?

Xiaoyu Fu

great function, thanks!

I updated the code by merging some bug fixes and performance improvements made by David Dreher some time ago.

@Eric Chadwick: can you send me a small test dataset or code how you generate the data?

Eric Chadwick

Works well, except that I have some strange artifacts when I run this on a volume with intersecting vessels. I am working on processing a much larger vasculature scan, but have created 3 intersecting cylinders to practice on. At the intersection, the function preserves the outline of the cylinders as well as the centreline. Do you know why this might be?

Eric Chadwick

Sarun

Sarun (view profile)

Lucy Robinson

This seems like a very useful function! Can it deal with shapes like a torus (doughnut), that have holes?

Thanks

Lucy

To remove short branches, please try my Skel2Graph3D script also here on FEX.

A slight directional bias can result from the fact that the medial axis is not uniquely defined in a voxel structure, but the error should never be more than +/- 1 voxel.

In the last update, a H1 line has been added and the code was cleaned up a bit.

Philip

William

Be useful if minimum branch length (noise removal) option was included

Alec Jacobson

Seems somehow biased to the left direction. That is, Skeleton3D(BW) and fliplr(Skeleton3D(fliplr(BW))) do not produce the same result

William

I have an MRA vessel that I'm trying to find the centre-line for, running it through the Skeleton3D code doesn't seem to isolate anything from the matrix. Any thoughts on whats going wrong?

maryam

maryam (view profile)

Sven

Sven (view profile)

Almost 5 stars. Wonderfully useful tool, but the help is a little sparse (it would be useful to have a H1 line so we know what the parameters and outputs are), and it looks like the function once had global variables which are now (thankfully) replaced but their global designation remains.

Since this is a port from other code I bet it can be vectorised in places for speed but I've found it useful so far - thanks.

zq

zq (view profile)

Updates

1.12.0.1

updated reference

1.12.0.0

Fixed title

1.11.0.0

Linked to github

1.1.0.0

Removed designation of global variables, and added H1 line for help. For cleaning noise such as short artifacts, please use my Skel2Graph3D package.

MATLAB Release Compatibility
Created with R2010b
Compatible with any release
Platform Compatibility
Windows macOS Linux
Acknowledgements

Inspired: Skel2Graph 3D

Discover Live Editor

Create scripts with code, output, and formatted text in a single executable document.


Learn About Live Editor