diff --git a/docs/release_notes.rst b/docs/release_notes.rst index bc06a2d1..151a873f 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 3eb50d1b..c1860d2d 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 432d17f2..eccf5fc0 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)