Geometric Medial Axis Representation from Binary Thinning


(Andrew Wilson) #1

How can I turn the medial axis produced by ITKs binary thinning into a mesh? Isosurface extraction methods would not be able to generate such a non manifold mesh. I can sort of imagine a way to do it but was hoping there was something in ITK or VTK to do it.

Also how does binary thinning compare to other methods like powercrust and q-mat (quadratic error minimization not quadtree).


(Dženan Zukić) #2

You could use @phcerdan’s Spatial Graph Extractor.

If you want to do it yourself, one way would be to iterate through the thinned image and for each pixel:

  1. add non-zero pixels to a set of nodes
  2. add “previous” neighbors from the neighborhood to the undirected list of neighbors

Finally turn this into a line mesh data structure of your favorite visualization library.


(grothausmann.roman) #3

I once created a simple CLI that basically does what @dzenanz suggested (https://github.com/romangrothausmann/ITK-CLIs/blob/d029ee0838fe4e4865c1b8bfb57ac01dd87f8d50/vo2ve.cxx)
I used it for converting the output of itkBinaryThinningImageFilter3D into a VTP file to load e.g. in ParaView and to apply boost graph algorithms on it.
Note, that it has a limitation
https://github.com/romangrothausmann/ITK-CLIs/blob/d029ee0838fe4e4865c1b8bfb57ac01dd87f8d50/vo2ve.cxx#L6-L17
that I was not able to resolve.
Please let me know if you know a way to remove that limitation.
@phcerdan do you have any idea how to solve that?


(Pablo Hernandez-Cerdan) #4

Spatial Graph Extractor (SGext) uses a graph from the thin image directly (using the boost graph library via the DGtal library). I refer to it as the raw graph, where all the non-zero voxels are nodes, and each one is connected to other voxels in the neighborhood with a 26 connectivity.

Then, it uses a depth-first search (DFS) visitor to merge/reduce this raw graph into a more useful spatial graph, where nodes with degree 2 are stored in edge points of a new edge connecting end-points (degree 1) or junctions (degree >2).

I used a different approach, so not sure from that comment in https://github.com/romangrothausmann/ITK-CLIs/blob/d029ee0838fe4e4865c1b8bfb57ac01dd87f8d50/vo2ve.cxx#L6-L17, what is your problem or how to solve it.

This was motivated because I implemented a 3D thinning algorithm in DGtal based on Couprie and Bertrand: Asymmetric parallel 3D thinning scheme and algorithms based on isthmuses. See the Voxel Complex docs for an intro.

A great addition to ITK would be to create a bridge/wrap with the boost graph library for binary images: See https://dgtal.org/doc/stable/moduleBoostGraphWrapping.html for the DGtal approach.


(grothausmann.roman) #5

Many thanks for your reply.
Instead of DFS I tried to solve the multiple connections of nodes with degree > 2 by labelling these voxel cluster and taking their center as node position to avoid any directional preference depending on a root node (which I guess DFS needs?).

Yes, I want to try your 3D thinning algorithm in DGtal in relation to this problem but did not find the time yet, sorry for that. Is your algorithm bound to DGtal because of the missing ITK bridge to the boost graph library?

I would also very much welcome an ITK bridge/wrap to the boost graph library. One use case would be to reduce a skeleton with BoostBiconnectedComponents, i.e. to remove all “dangeling” side-branches of centerline-loops (or the other way round). I’ve once create an ITK CLI for that as well (prune_ends). I can’t remember the full logic any more, but the above comment in question might be a stray comment from that CLI, which contains it as well.


(Pablo Hernandez-Cerdan) #6

The 3D thinning was developed and it is bound to DGtal, but SGExt is not. SGExt would just need the aforementioned bridge with boost to work without DGtal. However DGtal is pretty lightweight and CMakefied, if you just want to go ahead and give it a try. Its license is LGPL. Don’t hesitate to ask, or open an issue in SGext, it is in a beta state, but I have used it quite a bit in research projects, so hopefully you won’t find any outstanding bug. However, tutorials are missing. Any feedback is appreciated! :slight_smile:

If you want to give SGext a try, build DGtal with their ITK bridge. And try the thin script.

The nice thing about the Couprie and Betrand algorithm is that you can remove dangling branches automatically with the persistence parameter. I developed it because I was not very satisfied with the 3Dthinning in ITK for my use case. If you can give it a try and compare it with your cases, it would be great to know. It is great to know that you have also dealt with thinning and graphs @grothausmann.roman!