Skip to content

Fix encoding of Long Primitive Point Index List in ANN #467

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We鈥檒l occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
fedorov wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Fix encoding of Long Primitive Point Index List in ANN #467

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
22 changes: 22 additions & 0 deletions docs/release_notes.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,28 @@ Brief release notes may be found on `on Github
<https://github.com/MGHComputationalPathology/highdicom/releases>`_. This page
contains migration notes for major breaking changes to the library's API.

.. _ann-point-index-list:

Encoding of Long Primitive Point Index List in Bulk Annotations
---------------------------------------------------------------

Up to version 0.28.1, highdicom encoded the Long Primitive Point Index List
(0066,0040) of POLYGON and POLYLINE annotation groups in
:class:`highdicom.ann.AnnotationGroup` as the (one-based) index of the first
*coordinate value* of each annotation within the Point Coordinates Data. This
is not consistent with the standard, which requires the index of the first
*point* (coordinate tuple). For example, for two 2D polygons where the first
has 4 points, the correct value is ``[1, 5]``, whereas earlier versions of
highdicom wrote ``[1, 9]``.

Newly created annotation groups now use the correct encoding. When reading,
:meth:`highdicom.ann.AnnotationGroup.get_graphic_data()` detects the legacy
encoding where the values cannot be valid point indices, decodes it
accordingly and issues a warning. Other applications that read files created
with earlier versions of highdicom (or that produce files read by highdicom)
may need to be updated accordingly. See `issue #466
<https://github.com/ImagingDataCommons/highdicom/issues/466>`_.

.. _add-segments-deprecation:

Deprecation of `add_segments` method
Expand Down
118 changes: 62 additions & 56 deletions src/highdicom/ann/content.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
"""Content that is specific to Annotation IODs."""
from copy import deepcopy
from typing import cast
import warnings
from collections.abc import Sequence
from typing_extensions import Self

Expand Down Expand Up @@ -391,13 +392,10 @@ def __init__(
if len(unique_z_values) == 1:
self.CommonZCoordinateValue = unique_z_values.item()
coordinates_data = coordinates[:, 0:2].flatten()
dimensionality = 2
else:
coordinates_data = coordinates.flatten()
dimensionality = 3
else:
coordinates_data = coordinates.flatten()
dimensionality = 2

if coordinates.dtype == np.double:
self.DoublePointCoordinatesData = coordinates_data.tobytes()
Expand All @@ -410,7 +408,9 @@ def __init__(
GraphicTypeValues.POLYGON,
GraphicTypeValues.POLYLINE,
):
spans = [item.shape[0] * dimensionality for item in graphic_data]
# One-based index of the first point (coordinate tuple, not
# individual coordinate value) of each annotation
spans = [item.shape[0] for item in graphic_data]
point_indices = np.cumsum(spans, dtype=np.int32) + 1
point_indices = np.concatenate([
np.array([1], dtype=np.int32),
Expand Down Expand Up @@ -609,12 +609,11 @@ def get_graphic_data(
GraphicTypeValues.POLYGON,
):
# Variable number of coordinates per point
point_indices = np.frombuffer(
self.LongPrimitivePointIndexList,
dtype=np.int32
) - 1
split_param = (
point_indices // stored_coordinate_dimensionality
split_param = self._get_point_offsets(
number_of_points=len(decoded_coordinates_data),
stored_coordinate_dimensionality=(
stored_coordinate_dimensionality
),
)[1:]
else:
raise ValueError(
Expand Down Expand Up @@ -718,68 +717,75 @@ def get_measurements(
units = []
return (names, value_array, units)

def _get_coordinate_index(
def _get_point_offsets(
self,
annotation_number: int,
coordinate_dimensionality: int,
number_of_coordinates: int
number_of_points: int,
stored_coordinate_dimensionality: int,
) -> np.ndarray:
"""Get coordinate index.
"""Get offsets of the first point of each annotation.

Values of Long Primitive Point Index List index points (coordinate
tuples), not individual coordinate values. However, highdicom
versions up to 0.28.1 (and possibly other implementations) encoded
the index of the first coordinate value instead. Such legacy
encodings are detected and decoded with a warning, where possible.
Legacy values that are also valid point indices cannot be
distinguished and are interpreted as point indices.

Parameters
----------
annotation_number: int
One-based identification number of the annotation
coordinate_dimensionality: int
Dimensionality of coordinate points
number_of_coordinates: int
Total number of coordinate points
number_of_points: int
Total number of points (coordinate tuples) in the group
stored_coordinate_dimensionality: int
Number of coordinate values stored per point

Returns
-------
numpy.ndarray
One-dimensional array of zero-based index values to obtain the
coordinate points for a given annotation
One-dimensional array of zero-based offsets of the first point of
each annotation

""" # noqa: E501
annotation_index = annotation_number - 1
graphic_type = self.graphic_type
if graphic_type in (
GraphicTypeValues.POLYGON,
GraphicTypeValues.POLYLINE,
):
point_indices = np.frombuffer(
self.LongPrimitivePointIndexList,
dtype=np.int32
"""
offsets = np.frombuffer(
self.LongPrimitivePointIndexList,
dtype=np.int32
).astype(np.int64) - 1

if len(offsets) == 0 or offsets[0] != 0:
raise ValueError(
'Invalid Long Primitive Point Index List: first value must '
'be 1.'
)
start = point_indices[annotation_index] - 1
try:
end = point_indices[annotation_index + 1] - 1
except IndexError:
end = number_of_coordinates
else:
if hasattr(self, 'CommonZCoordinateValue'):
stored_coordinate_dimensionality = 2
else:
stored_coordinate_dimensionality = coordinate_dimensionality
if graphic_type in (
GraphicTypeValues.ELLIPSE,
GraphicTypeValues.RECTANGLE,
if np.any(np.diff(offsets) <= 0):
raise ValueError(
'Invalid Long Primitive Point Index List: values must be '
'strictly increasing.'
)

if offsets[-1] >= number_of_points:
# Not a valid point index. Check whether the values index
# coordinate values instead (legacy encoding).
d = stored_coordinate_dimensionality
if (
np.all(offsets % d == 0) and
offsets[-1] < number_of_points * d
):
length = 4 * stored_coordinate_dimensionality
elif graphic_type == GraphicTypeValues.POINT:
length = stored_coordinate_dimensionality
warnings.warn(
'Values of Long Primitive Point Index List appear to '
'index individual coordinate values rather than points. '
'This non-standard encoding was used by highdicom '
'versions up to 0.28.1. The values will be '
'interpreted accordingly.',
UserWarning,
)
offsets = offsets // d
else:
raise ValueError(
'Encountered unexpected graphic type '
f'"{graphic_type.value}".'
'Invalid Long Primitive Point Index List: values exceed '
'the number of points.'
)
start = annotation_index * length
end = start + length

coordinate_index = np.arange(start, end)

return coordinate_index
return offsets

@classmethod
def from_dataset(
Expand Down
130 changes: 130 additions & 0 deletions tests/test_ann.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -486,6 +486,136 @@ def test_construction_with_wrong_graphic_data_rectangle(self):
)


class TestAnnotationGroupPointIndexList:

# First annotation has 4 points, second has 3 points
_graphic_data_2d = [
np.array([[0.0, 0.0], [10.0, 0.0], [10.0, 10.0], [0.0, 10.0]]),
np.array([[20.0, 20.0], [30.0, 20.0], [25.0, 30.0]]),
]
_graphic_data_3d_common_z = [
np.array([
[0.0, 0.0, 5.0],
[10.0, 0.0, 5.0],
[10.0, 10.0, 5.0],
[0.0, 10.0, 5.0],
]),
np.array([[20.0, 20.0, 5.0], [30.0, 20.0, 5.0], [25.0, 30.0, 5.0]]),
]
_graphic_data_3d = [
np.array([
[0.0, 0.0, 1.0],
[10.0, 0.0, 2.0],
[10.0, 10.0, 3.0],
[0.0, 10.0, 4.0],
]),
np.array([[20.0, 20.0, 1.0], [30.0, 20.0, 1.0], [25.0, 30.0, 1.0]]),
]

@staticmethod
def _create_group(graphic_data, graphic_type=GraphicTypeValues.POLYGON):
return AnnotationGroup(
number=1,
uid=UID(),
label='foo',
annotated_property_category=codes.SCT.MorphologicallyAbnormalStructure, # noqa: E501
annotated_property_type=codes.SCT.Neoplasm,
graphic_type=graphic_type,
graphic_data=graphic_data,
algorithm_type=AnnotationGroupGenerationTypeValues.MANUAL,
)

@staticmethod
def _reread(group):
return AnnotationGroup.from_dataset(Dataset(group))

@pytest.mark.parametrize(
'graphic_data,coordinate_type,stored_dimensionality',
[
(_graphic_data_2d, '2D', 2),
(_graphic_data_3d_common_z, '3D', 2),
(_graphic_data_3d, '3D', 3),
]
)
@pytest.mark.parametrize(
'graphic_type',
[GraphicTypeValues.POLYGON, GraphicTypeValues.POLYLINE]
)
def test_encoding(
self,
graphic_data,
coordinate_type,
stored_dimensionality,
graphic_type,
):
group = self._create_group(graphic_data, graphic_type)

# Values index points (coordinate tuples), not coordinate values,
# irrespective of the number of values stored per point
point_indices = np.frombuffer(
group.LongPrimitivePointIndexList,
dtype=np.int32
)
np.testing.assert_array_equal(point_indices, [1, 5])

coordinates = np.frombuffer(
group.DoublePointCoordinatesData,
dtype=np.float64
)
assert coordinates.size == 7 * stored_dimensionality

decoded = self._reread(group).get_graphic_data(coordinate_type)
assert len(decoded) == len(graphic_data)
for retrieved, expected in zip(decoded, graphic_data):
np.testing.assert_array_equal(retrieved, expected)

@pytest.mark.parametrize(
'graphic_data,coordinate_type,legacy_point_indices',
[
(_graphic_data_2d, '2D', [1, 9]),
(_graphic_data_3d_common_z, '3D', [1, 9]),
(_graphic_data_3d, '3D', [1, 13]),
]
)
def test_decoding_legacy_encoding(
self,
graphic_data,
coordinate_type,
legacy_point_indices,
):
group = self._create_group(graphic_data)
group.LongPrimitivePointIndexList = np.array(
legacy_point_indices,
dtype=np.int32
).tobytes()

with pytest.warns(UserWarning, match='Long Primitive Point Index'):
decoded = self._reread(group).get_graphic_data(coordinate_type)

assert len(decoded) == len(graphic_data)
for retrieved, expected in zip(decoded, graphic_data):
np.testing.assert_array_equal(retrieved, expected)

@pytest.mark.parametrize(
'point_indices',
[
[2, 5], # does not start at 1
[1, 1], # not strictly increasing
[1, 8], # neither a valid point index nor value index
[1, 15], # exceeds number of values
]
)
def test_decoding_invalid(self, point_indices):
group = self._create_group(self._graphic_data_2d)
group.LongPrimitivePointIndexList = np.array(
point_indices,
dtype=np.int32
).tobytes()

with pytest.raises(ValueError, match='Long Primitive Point Index'):
self._reread(group).get_graphic_data('2D')


class TestMicroscopyBulkSimpleAnnotations:

@pytest.fixture(autouse=True)
Expand Down
Loading