2

Using the LabelShapeStatisticFilter, I can extract oriented regions of interest from the original image correctly. I want to plot those oriented bounding boxes over the original image.

When I try to view the output of the GetOrientedBoundingBoxVertices() method, it's not clear to me in what coordinate system these vertices are defined. They do not appear to be in the original image coordinate system.

I am confident I am using the LabelShapeStatisticFilter class as intended (see below), following this excellent notebook: http://insightsoftwareconsortium.github.io/SimpleITK-Notebooks/Python_html/35_Segmentation_Shape_Analysis.html

bacteria_labels = shape_stats.GetLabels()
bacteria_volumes = [shape_stats.GetPhysicalSize(label) for label in bacteria_labels] 
num_images = 5 # number of bacteria images we want to display

bacteria_labels_volume_sorted = [label for _,label in sorted(zip(bacteria_volumes, bacteria_labels))]

resampler = sitk.ResampleImageFilter()
aligned_image_spacing = [10,10,10] #in nanometers

for label in bacteria_labels_volume_sorted[0:num_images]:
    aligned_image_size = [ int(ceil(shape_stats.GetOrientedBoundingBoxSize(label)[i]/aligned_image_spacing[i])) for i in range(3) ]
    direction_mat = shape_stats.GetOrientedBoundingBoxDirection(label)
    aligned_image_direction = [direction_mat[0], direction_mat[3], direction_mat[6], 
                               direction_mat[1], direction_mat[4], direction_mat[7],
                               direction_mat[2], direction_mat[5], direction_mat[8] ] 
    resampler.SetOutputDirection(aligned_image_direction)
    resampler.SetOutputOrigin(shape_stats.GetOrientedBoundingBoxOrigin(label))
    resampler.SetOutputSpacing(aligned_image_spacing)
    resampler.SetSize(aligned_image_size)

    obb_img = resampler.Execute(img)
    # Change the image axes order so that we have a nice display.
    obb_img = sitk.PermuteAxes(obb_img,[2,1,0])
    gui.MultiImageDisplay(image_list = [obb_img],                   
                          title_list = ["OBB_{0}".format(label)])

I expect to be able to draw these bounding boxes over the original image, but I'm not sure how.

UPDATE

Perhaps this can illustrate what I mean better. Resampled Oriented Bounding Box, output as expected:

resampled

However, after using original_label_image.TransformPhysicalPointToContinousIndex(), the oriented bounding box points in the original image space appear incorrect (shape_stats.OrientedBoundingBoxVertices() in original index space):

original label image

UPDATE 2

Using shape_stats.GetCentroid(), I can correctly get real-world coordinates of the centroids of each label and plot them correctly:

centroids

It also appears that output of shape_stats.GetOrientedBoundingBoxOrigin() is plausibly in real-world coordinates. One element of shape_stats.OrientedBoundingBoxVertices() corresponds to shape_stats.GetOrientedBoundingBoxOrigin().

enter image description here

Kevin
  • 31
  • 5

2 Answers2

1

The vertices are defined in physical spaces and not index space. You may need to use the Image class’s TransformPhysicslPointToIndex.

blowekamp
  • 1,401
  • 7
  • 7
  • From which image should I call TransformPhysicalPointToIndex? The vertices do not appear to be defined in the original image's real-world coordinates. – Kevin Jul 25 '19 at 17:45
  • SimpleITK uses i,j,k indexes but numpy uses k,j,i indexes. Was this conversion done? – blowekamp Jul 26 '19 at 06:32
  • Yes, the conversion was done. I can correctly plot other real-world points of interest. I will update with more images, but it appears that shape_stats.GetOrientedBoundingBoxOrigin(), which provides one vertex, is in its proper location. – Kevin Jul 26 '19 at 12:49
  • I can correctly plot other real-world points of interest generated from LabelShapeStatisticsFilter (see plot of the centroids in update 2). It appears that shape_stats.GetOrientedBoundingBoxOrigin(), which provides one vertex, are plausibly in world coordinates (see update 2). – Kevin Jul 26 '19 at 12:58
1

I believe I have figured it out: the oriented bounding box vertices are neither entirely in original image coordinates, or in the coordinates of the bounding box.

The origin of the oriented bounding box, returned by shape_stats.GetOrientedBoundingBoxOrigin(), is in original image world-coordinates. This origin also corresponds to one vertex of the oriented bounding box.

Each vertex of the oriented bounding box, returned by shape_stats.OrientedBoundingBoxVertices(), can be recovered in the real-world coordinates by a rotation about the origin using shape_stats.GetOrientedBoundingBoxDirection().

I don't know if this representation of vertices was intentional, but it was confusing to me at first (though I am a relative newcomer to sitk).

Kevin
  • 31
  • 5
  • This looks like a bug to me. Please report it to the ITK GitHub issues: https://github.com/InsightSoftwareConsortium/ITK/issues – blowekamp Jul 30 '19 at 14:40