3d skeleton

[thanasi]

Hi,

 I am trying to skeletonize a 3d dataset, and run into errors when using the mahotas skeletonize with a 3d structure element. I saw there was an old thread about this that seemed to have died. If you don't see this as a feature that you plan to implement soon, then can you point me to a paper or other resource that would guide me on how to use mahotas/ndimage functionality to build a 3d-skeleton creator? I am concerned mainly about efficiency and while I have some experience with c/c++, I have no experience wrapping it into python.

Thanks!

Thanasi

[Luis Pedro Coelho]

Yes, the skeleton functionality in mahotas is 2D-only.

I think you'd have to do it from scratch. Googling brings up too many
alternatives to list here.

*

For implementation, you have 2 options

1) Cython, like:

while thin_more:
for i in xrange(im.shape[0]):
for j in xrange(im.shape[1]):
for k in xrange(im.shape[2]):
thin_here(im, i, j, k)

This would be the easiest to implement. Look at skimage code for examples.

2) C/C++. mahotas has a very efficient method to the thinning in 2D. It peeks at
the ndarray structure to get arrays of pointer offsets and then does multiple
passes as fast at possible.

This is harder to get right (at least it was for me—a whole day just
transforming this function into this form), but significantly faster (I
obsessively measured it while coding it).

So, it depends on your computer time/programmer time tradeoffs.

HTH
Luis

--
Luis Pedro Coelho | Institute for Molecular Medicine | http://luispedro.org

[Sprindzuk Matvey]

Does anybody know where to find a detailed description of 13 vectors
for Haralick features computation?

Thank you

2012/3/27 Luis Pedro Coelho <lu...@luispedro.org>:

[Thouis (Ray) Jones]

There's an implementation here, with pointers at the top of the file
for the original reference. It's not the most straightforward to read
code, but it's not that hard to understand:

https://github.com/thouis/CellProfiler/blob/master/cellprofiler/cpmath/haralick.py#L110

[thanasi]

I am still a little unfamiliar with the details of thinning, and would like your opinion on this idea.

I've seen that thinning can be achieved by subtracting a HitOrMiss transform from the original binary image (see  http://homepages.inf.ed.ac.uk/rbf/HIPR2/thin.htm). This is also how thinning is achieved in 1-3D in the IDL library.  Why is this approach not taken/what are the benefits of other thinning methods over this?  Thanks!

[Luis Pedro Coelho]

On Thursday, March 29, 2012 11:59:13 AM thanasi wrote:
> I am still a little unfamiliar with the details of thinning, and would like
> your opinion on this idea.
> I've seen that thinning can be achieved by subtracting a HitOrMiss
> transform from the original binary image (see
> http://homepages.inf.ed.ac.uk/rbf/HIPR2/thin.htm). This is also how
> thinning is achieved in 1-3D in the IDL library. Why is this approach not
> taken/what are the benefits of other thinning methods over this?

This is what mahotas does for 2D.

If you have the right templates for 3D thinning, you can just have a wrapper
around hitmiss() and call it. This works well and my first implementation of
thinning for 2D worked like that.

It's just not especially fast (you have to loop in Python for the multiple
hitmiss() calls and each of those computes filter offsets).

HTH,

--
Luis Pedro Coelho | Institute for Molecular Medicine | http://luispedro.org

LxMLS 2012: Lisbon Machine Learning School
http://lxmls.it.pt

[thanasi]

Thanks - I didn't completely understand the first implementation then.  I'll try this - will I gain noticeable speed by typecasting variables using Cython?

[Luis Pedro Coelho]

On Thursday, March 29, 2012 12:18:44 PM thanasi wrote:
> Thanks - I didn't completely understand the first implementation then.
> I'll try this - will I gain noticeable speed by typecasting variables
> using Cython?

You need to measure it.

hitmiss() [at least from mahotas] is C++, so it's fast. What is slow is that
you won't be calling it optimally and doing some extra computation in every
loop.

Still, you might not need the fastest possible code, if you have less than
hundreds of thousands of images.