Calculate volume of largest and minimum ellipse around a segmentation

ral_san · January 31, 2020, 3:39pm

Is there any function in SimpleITK similar to the OpenCV function FitEllipse --> https://docs.opencv.org/3.4.9/de/d62/tutorial_bounding_rotated_ellipses.html ?

What I am trying to do is to calculate the volume of the largest and minimal ellipsoid volume related to a segmented tumor.
Something like this figure, where 1 is my tumor segmentation and 2 and 3 would be the ellipsoid volumes that I need to calculate.

I am not looking for code but just for some direction.

Thank you for your help.

dzenanz · January 31, 2020, 3:47pm

I don’t think there is an ellipse fitting class. What you could look into is LabelShapeStatisticsImageFilter and LabelGeometryImageFilter, and possibly their source code.

blowekamp · January 31, 2020, 4:21pm

I am not sure how well defined or unique the “minimal” ellipsoid would be.

In the LabelShapeStatisticsImageFilter, the principle axis are computed. There are several derived parameters from that you could look into. The original Insight Journal article may be the best reference.

Beyond those existing attributes you would need to create your own. A method needs to be added to LabelShapStatisticsImageFilter to get a list of indices or points for a given label. In the mean time something like np.where(sitk.GetArrayViewFromImage(lableImage)) could provide similar functionality.

blowekamp · January 31, 2020, 9:40pm

Here is the PR adding the new functionality to SimpleITK:

ral_san · February 14, 2020, 10:21am

I think I managed to solve this issue.

gist.github.com

https://gist.github.com/raluca-san/5d8bd4dcb3278b8e2f6dfa959614f9c8

MinVolEllipse.py

# -*- coding: utf-8 -*-
"""
@author: Raluca Sandu
"""

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import numpy as np
import numpy.linalg as la

This file has been truncated. show original

What is missing is extracting the surface points from a segmentation mask. That could be extracted using the SimpleITK.LabelContour(img, fullyConnected=False).
Also the spacing information is missing, so the ellipsoid should be adjusted by placing it on the correct grid.
This can be achieved scikit-image.draw.ellipsoid https://scikit-image.org/docs/dev/api/skimage.draw.html#skimage.draw.ellipsoid.

Is it possible to propose adding this functionality to the official SimpleITK library?

dzenanz · February 14, 2020, 1:35pm

Thanks for sharing your solution. But the problem is somewhat specialized. Moreover, ITK and SimpleITK are written in C++ so you would have to translate your solution into C++ for incorporation.

blowekamp · February 14, 2020, 2:56pm

Thanks for sharing your implementation of computing the optimal ellipsoid.

I believe you are asking to get the points from a contour in SimpleITK.

The above linked PR was merged and you can download the latest binaries to use it.

Here is sample usage:


In [4]: labels = sitk.ReadImage("/scratch/blowekamp/SimpleITK/SimpleITK-build/ExternalData/Testing/Data/Input/2th_cthead1.png")                                                                                                                             

In [5]: sitk.LabelContour(labels, fullyConnected=False)                                                                                                                                                                                                     
Out[5]: <SimpleITK.SimpleITK.Image; proxy of <Swig Object of type 'std::vector< itk::simple::Image >::value_type *' at 0x7f8052bdbc00> >

In [6]: contours = sitk.LabelContour(labels, fullyConnected=False)                                                                                                                                                                                          

In [7]: shape = sitk.LabelShapeStatisticsImageFilter()                                                                                                                                                                                                      

In [8]: shape.Execute(contours)                                                                                                                                                                                                                             

In [9]: shape.GetLabels()                                                                                                                                                                                                                                   
Out[9]: (100, 200)

In [10]: idx = shape.GetIndexes(100)