Morphological snakes

[grahamross123]

Morphological snakes function based on scikit image morphological_chan_vese:

To do:

• Fix smoothing to be the same as in skimage implementation
• Insert hardcoded input image and reference into test
• Clean up level_set_segmentation notebook

[StRigaud]

Digging into the smoothing operation done in scikit-image, the smoothing is actually not doing an opening and closing.

They define 2 operations inferior_superior and superior_inferior, which respectively correspond to the minimum value of a dilation operation along all the direction and the maximum value of an erosion along all the direction.

In a more codish description:

the superior_inferior

erode the image respectively with all the elements

1 0 0  0 0 1  0 1 0  0 0 0
0 1 0  0 1 0  0 1 0  1 1 1
0 0 1  1 0 0  0 1 0  0 0 0

and take the maximum pixels across the 4 different results dimensions

the inferior_superior

dilate the image respectively with all the elements

1 0 0  0 0 1  0 1 0  0 0 0
0 1 0  0 1 0  0 1 0  1 1 1
0 0 1  1 0 0  0 1 0  0 0 0

and take the minimum pixels across the 4 different results dimensions

Code is here

I am not sure how doable this is with the current morphological. operator available in clesperanto nor if it is GPU compatible.

We may need an alternative operation ...

[StRigaud]

Digging into the smoothing operation done in scikit-image, the smoothing is actually not doing an opening and closing.

They define 2 operations inferior_superior and superior_inferior, which respectively correspond to the minimum value of a dilation operation along all the direction and the maximum value of an erosion along all the direction.

In a more codish description:

the superior_inferior

erode the image respectively with all the elements

1 0 0  0 0 1  0 1 0  0 0 0
0 1 0  0 1 0  0 1 0  1 1 1
0 0 1  1 0 0  0 1 0  0 0 0

and take the maximum pixels across the 4 different results dimensions

the inferior_superior

dilate the image respectively with all the elements

1 0 0  0 0 1  0 1 0  0 0 0
0 1 0  0 1 0  0 1 0  1 1 1
0 0 1  1 0 0  0 1 0  0 0 0

and take the minimum pixels across the 4 different results dimensions

Code is here

I am not sure how doable this is with the current morphological. operator available in clesperanto nor if it is GPU compatible.

We may need an alternative operation ...

[haesleinhuepf]

We may need an alternative operation ...

Agreed. If possible, I would pack this functionality in one or two new opencl-kernels.

[grahamross123]

It looks possible to do this on the GPU but we would need a few more opencl kernels. If we use the same structure as the other pyclesperanto erode / dilate functions, we'd need one kernel for each 2D matrix (8 kernels total) unless we can make a general kernel which can take a 2D matrix as input.

So you could loop through erosions / dilations with each 2D matrix and then take the max / min of the resulting array using maximum_z_projection.

As for 3D, we'd need much more kernels for each 3D matrix (18 total, unless using a more general implementation) and we might need a new function to obtain the max / min along an axis of the resulting 4D image after running each erosion / dilation.

If that makes sense and I'm not missing anything, I'm happy to work on this.

[StRigaud]

we'd need one kernel for each 2D matrix (8 kernels total) unless we can make a general kernel which can take a 2D matrix as input.

That's one way, though it would make too many 'kernel.ocl'. A general kernel which compute the correct kernel based on direction and dimensions (2d or 3d) on the fly is tricky imho. I would push for a hard coded version of each directional kernel (2d and 3d), its repetitive and ugly but works well while remaining efficient.

we might need a new function to obtain the max / min along an axis of the resulting 4D image

This, I'm not sure its possible... thus we may consider the computation of the min/max value directly during kernel execution and not rely on the maximum_z_projection.

That would mean 2 kernels,: supinf and infsup. Each of them having 8 + 18 hard coded directional 3x3(x3) kernels and would, for each pixel, apply an erosion (or dilation) followed by a maximum (minimum).

[grahamross123]

That would mean 2 kernels,: supinf and infsup. Each of them having 8 + 18 hard coded directional 3x3(x3) kernels and would, for each pixel, apply an erosion (or dilation) followed by a maximum (minimum).

I agree that sounds cleaner overall and (probably) easier to implement.