Skip to content

Commit

Permalink
Fix plane ordering bug; split spatial functions
Browse files Browse the repository at this point in the history
  • Loading branch information
@CPBridge
CPBridge committed Feb 2, 2024
1 parent ec26339 commit f598555
Show file tree
Hide file tree
Showing 3 changed files with 217 additions and 83 deletions.
58 changes: 44 additions & 14 deletions src/highdicom/seg/content.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,11 @@
)
from highdicom.enum import CoordinateSystemNames
from highdicom.seg.enum import SegmentAlgorithmTypeValues
from highdicom.spatial import map_pixel_into_coordinate_system
from highdicom.spatial import (
_get_slice_distances,
get_normal_vector,
map_pixel_into_coordinate_system,
)
from highdicom.sr.coding import CodedConcept
from highdicom.uid import UID
from highdicom.utils import compute_plane_position_slide_per_frame
Expand Down Expand Up @@ -605,7 +609,8 @@ def get_index_position(self, pointer: str) -> int:

def get_index_values(
self,
plane_positions: Sequence[PlanePositionSequence]
plane_positions: Sequence[PlanePositionSequence],
image_orientation: Optional[Sequence[float]] = None,
) -> Tuple[np.ndarray, np.ndarray]:
"""Get values of indexed attributes that specify position of planes.
Expand All @@ -626,6 +631,15 @@ def get_index_values(
plane_indices: numpy.ndarray
1D array of planes indices for sorting frames according to their
spatial position specified by the dimension index
image_orientation: Union[Sequence[float], None], optional
An image orientation to use to order frames within a 3D coordinate
system. By default (if ``image_orientation`` is ``None``), the
plane positions are ordered using their raw numerical values and
not along any particular spatial vector. If ``image_orientation``
is provided, planes are ordered along the positive direction of the
vector normal to the specified. Should be a sequence of 6 floats.
This is only valid when plane position inputs contain only the
ImagePositionPatient.
Note
----
Expand Down Expand Up @@ -659,21 +673,37 @@ def get_index_values(
for p in plane_positions
])

# Build an array that can be used to sort planes according to the
# Dimension Index Value based on the order of the items in the
# Dimension Index Sequence.
_, plane_sort_indices = np.unique(
plane_position_values,
axis=0,
return_index=True
)
if image_orientation is not None:
if not hasattr(plane_positions[0][0], 'ImagePositionPatient'):
raise ValueError(
'Provided "image_orientation" is only valid when '
'plane_positions contain the ImagePositionPatient.'
)
normal_vector = get_normal_vector(image_orientation)
origin_distances = _get_slice_distances(
plane_position_values[:, 0, :],
normal_vector,
)
_, plane_sort_indices = np.unique(
origin_distances,
return_index=True,
)
else:
# Build an array that can be used to sort planes according to the
# Dimension Index Value based on the order of the items in the
# Dimension Index Sequence.
_, plane_sort_indices = np.unique(
plane_position_values,
axis=0,
return_index=True
)

if len(plane_sort_indices) != len(plane_positions):
raise ValueError(
"Input image/frame positions are not unique according to the "
"Dimension Index Pointers. The generated segmentation would be "
"ambiguous. Ensure that source images/frames have distinct "
"locations."
'Input image/frame positions are not unique according to the '
'Dimension Index Pointers. The generated segmentation would be '
'ambiguous. Ensure that source images/frames have distinct '
'locations.'
)

return (plane_position_values, plane_sort_indices)
Expand Down
96 changes: 52 additions & 44 deletions src/highdicom/seg/sop.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1460,9 +1460,15 @@ def __init__(
# number) plane_sort_index is a list of indices into the input
# planes giving the order in which they should be arranged to
# correctly sort them for inclusion into the segmentation
sort_orientation = (
plane_orientation[0].ImageOrientationPatient
if self._coordinate_system == CoordinateSystemNames.PATIENT
else None
)
plane_position_values, plane_sort_index = \
self.DimensionIndexSequence.get_index_values(
plane_positions
plane_positions,
image_orientation=sort_orientation,
)

are_spatial_locations_preserved = (
Expand Down Expand Up @@ -1501,49 +1507,6 @@ def __init__(
"the source image."
)

# Dimension Organization Type
dimension_organization_type = self._check_tiled_dimension_organization_type(
dimension_organization_type=dimension_organization_type,
is_tiled=is_tiled,
omit_empty_frames=omit_empty_frames,
plane_positions=plane_positions,
rows=self.Rows,
columns=self.Columns,
)
if self._coordinate_system == CoordinateSystemNames.PATIENT:
spacing = get_regular_slice_spacing(
image_positions=np.array(plane_position_values[:, 0, :]),
image_orientation=np.array(
plane_orientation[0].ImageOrientationPatient
),
sort=False,
enforce_postive=True,
)

if spacing is not None and spacing > 1.0:
# The image is a regular volume, so we should record this
dimension_organization_type = (
DimensionOrganizationTypeValues.THREE_DIMENSIONAL
)
# Also add the slice spacing to the pixel measures
(
self.SharedFunctionalGroupsSequence[0]
.PixelMeasuresSequence[0]
.SpacingBetweenSlices
) = spacing
else:
if (
dimension_organization_type ==
DimensionOrganizationTypeValues.THREE_DIMENSIONAL
):
raise ValueError(
'Dimension organization "3D" has been specified, '
'but the source image is not a regularly-spaced 3D '
'volume.'
)
if dimension_organization_type is not None:
self.DimensionOrganizationType = dimension_organization_type.value

# Find indices such that empty planes are removed
if omit_empty_frames:
if tile_pixel_array:
Expand Down Expand Up @@ -1589,6 +1552,51 @@ def __init__(
else:
unique_dimension_values = [None]

# Dimension Organization Type
dimension_organization_type = self._check_tiled_dimension_organization_type(
dimension_organization_type=dimension_organization_type,
is_tiled=is_tiled,
omit_empty_frames=omit_empty_frames,
plane_positions=plane_positions,
rows=self.Rows,
columns=self.Columns,
)
if self._coordinate_system == CoordinateSystemNames.PATIENT:
spacing = get_regular_slice_spacing(
image_positions=np.array(
plane_position_values[plane_sort_index, 0, :]
),
image_orientation=np.array(
plane_orientation[0].ImageOrientationPatient
),
sort=False,
enforce_positive=True,
)

if spacing is not None and spacing > 0.0:
# The image is a regular volume, so we should record this
dimension_organization_type = (
DimensionOrganizationTypeValues.THREE_DIMENSIONAL
)
# Also add the slice spacing to the pixel measures
(
self.SharedFunctionalGroupsSequence[0]
.PixelMeasuresSequence[0]
.SpacingBetweenSlices
) = spacing
else:
if (
dimension_organization_type ==
DimensionOrganizationTypeValues.THREE_DIMENSIONAL
):
raise ValueError(
'Dimension organization "3D" has been specified, '
'but the source image is not a regularly-spaced 3D '
'volume.'
)
if dimension_organization_type is not None:
self.DimensionOrganizationType = dimension_organization_type.value

if (
has_ref_frame_uid and
self._coordinate_system == CoordinateSystemNames.SLIDE
Expand Down

0 comments on commit f598555

Please sign in to comment.