From 04dedb0c06428bfd88fd89d1b9b7055d5a734f6e Mon Sep 17 00:00:00 2001 From: Andrey Fedorov Date: Tue, 6 Oct 2026 13:29:39 -0400 Subject: [PATCH] Fix encoding of Long Primitive Point Index List in ANN Long Primitive Point Index List (0066,0040) of POLYGON and POLYLINE annotation groups was encoded and decoded as the one-based index of the first coordinate value of each annotation in Point Coordinates Data, rather than the index of the first point (coordinate tuple), as intended by the standard (clarified on DICOM WG-26). For two 2D polygons where the first has 4 points, highdicom wrote [1, 9] instead of [1, 5]. - Encode point indices in AnnotationGroup - Decode point indices in get_graphic_data(), validating the values and detecting the legacy value-based encoding (with a warning) where the values cannot be valid point indices - Replace unused _get_coordinate_index() helper - Add tests asserting the raw encoded values, legacy decoding and rejection of invalid values - Add migration note to release notes Addresses #466 Co-Authored-By: Claude Opus 5.5 --- docs/release_notes.rst | 22 ++++++ src/highdicom/ann/content.py | 118 ++++++++++++++++--------------- tests/test_ann.py | 130 +++++++++++++++++++++++++++++++++++ 3 files changed, 214 insertions(+), 56 deletions(-) diff --git a/docs/release_notes.rst b/docs/release_notes.rst index bc06a2d1b..151a873f1 100644 --- a/docs/release_notes.rst +++ b/docs/release_notes.rst @@ -7,6 +7,28 @@ Brief release notes may be found on `on Github `_. 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 +`_. + .. _add-segments-deprecation: Deprecation of `add_segments` method diff --git a/src/highdicom/ann/content.py b/src/highdicom/ann/content.py index 3eb50d1b8..c1860d2d9 100644 --- a/src/highdicom/ann/content.py +++ b/src/highdicom/ann/content.py @@ -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 @@ -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() @@ -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), @@ -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( @@ -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( diff --git a/tests/test_ann.py b/tests/test_ann.py index 432d17f21..eccf5fc03 100644 --- a/tests/test_ann.py +++ b/tests/test_ann.py @@ -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)