ANN: Long Primitive Point Index List (0066,0040) encoded/decoded as value index instead of point index #466

Note on how this report was produced: The analysis and the text of this issue were generated by Claude Code (Anthropic's Claude Opus 5.5 model), working in a local clone of this repository at the request of @fedorov, who prompted it with the WG-26 discussion summarized below. Claude did all of the following:

• inspected the highdicom source and git history
• wrote and ran the reproduction script against master
• queried the IDC index
• downloaded and decoded the example IDC file

The outputs below are copied verbatim from those runs. The proposed fix is a suggestion and has not been implemented or tested.

Summary

For POLYGON and POLYLINE annotation groups in Microscopy Bulk Simple Annotations, highdicom encodes and decodes Long Primitive Point Index List (0066,0040) as an index into the flat array of coordinate values. The standard intends an index into the list of points (coordinate tuples).

Take two 2D polygons, the first with 4 points (8 values) and the second with 3 points:

Writing and reading make the same mistake, so highdicom round-trips its own files without error. Files that follow the standard are read incorrectly and silently. Files highdicom writes are not standard-conformant. This has been the behavior since the original ANN implementation (#114, b5bcb62).

Background: the standard's interpretation

The question came up on the DICOM WG-26 mailing list (thread "MBSA - How to interpret Long Primitive Point Index List (0066,0040)", 2 Oct 2026). David Clunie's answer (paraphrased):

• The list references points, not individual values. For multiple 2D polygons where the first has 4 points (8 values), the list is [1, 5], not [1, 9].
• "A list of point indices" in the attribute description means exactly that. The "index used" text only says that indexing starts at 1 and that the first entry refers to the first point of the first annotation.
• The wording is ambiguous enough that a CP to clarify the text and improve the Part 17 examples is likely. Any clarification has to stay consistent with the mechanisms this one descends from, namely C.27.4 Surface Mesh Primitives Macro and Track Point Index List (0066,0129) in C.8.33.2 Tractography Results Module. Both index points (tuples), not values.

Where in the code

Encoding, src/highdicom/ann/content.py L409-L419:

spans = [item.shape[0] * dimensionality for item in graphic_data]  # number of *values* per annotation
point_indices = np.cumsum(spans, dtype=np.int32) + 1

dimensionality is the number of values stored per point: 2 for 2D or for 3D with Common Z Coordinate Value, and 3 for explicit 3D.

Decoding, AnnotationGroup.get_graphic_data(), L607-L618:

point_indices = np.frombuffer(self.LongPrimitivePointIndexList, dtype=np.int32) - 1
split_param = (point_indices // stored_coordinate_dimensionality)[1:]

This divides by the stored dimensionality, which assumes value-based indices. It does not validate the indices, so a standard-conformant file is split in the wrong places without any error or warning.

Unused helper, AnnotationGroup._get_coordinate_index(), L721-L782: it is never called anywhere in the package, and it also treats the list entries as value offsets.

Tests: no test in tests/test_ann.py checks the raw encoded value of LongPrimitivePointIndexList. Every test is an encode/decode round trip, which passes because both directions share the same interpretation.

How the analysis was done

1. Code inspection

grep -rn "LongPrimitivePointIndexList" src tests finds the three places above and no tests. git log -S LongPrimitivePointIndexList shows the encoding has not changed since b5bcb62.

2. Reproduction against current master (57ce37d)

The script below was run with PYTHONPATH=src uv run --no-project --with pydicom --with numpy --with pillow --with pillow-jpls --with typing-extensions python repro.py. For each case it does two things:

• encodes two polygons (4 points + 3 points) and prints the resulting index list
• replaces the list with the standard-conformant point-based list [1, 1 + n_points_first] and decodes it with get_graphic_data()

import numpy as np
from pydicom.sr.codedict import codes
import highdicom as hd

cases = [
    # 2D
    [np.array([[0, 0], [10, 0], [10, 10], [0, 10]], float),
     np.array([[20, 20], [30, 20], [25, 30]], float)],
    # 3D, single z plane -> CommonZCoordinateValue, 2 stored values per point
    [np.array([[0, 0, 5], [10, 0, 5], [10, 10, 5], [0, 10, 5]], float),
     np.array([[20, 20, 5], [30, 20, 5], [25, 30, 5]], float)],
    # 3D, explicit z -> 3 stored values per point
    [np.array([[0, 0, 1], [10, 0, 2], [10, 10, 3], [0, 10, 4]], float),
     np.array([[20, 20, 1], [30, 20, 1], [25, 30, 1]], float)],
]

for gd in cases:
    g = hd.ann.AnnotationGroup(
        number=1, uid=hd.UID(), label='x',
        annotated_property_category=codes.SCT.MorphologicallyAbnormalStructure,
        annotated_property_type=codes.SCT.Neoplasm,
        graphic_type='POLYGON', graphic_data=gd, algorithm_type='MANUAL',
    )
    print(gd[0].shape[1], np.frombuffer(g.LongPrimitivePointIndexList, np.int32))

    # Simulate a standard-conformant file with point-based indices
    n = [len(x) for x in gd]
    g2 = hd.ann.AnnotationGroup.from_dataset(g)
    g2.LongPrimitivePointIndexList = np.array([1, 1 + n[0]], np.int32).tobytes()
    g2._graphic_data = {}
    ct = '2D' if gd[0].shape[1] == 2 else '3D'
    print('  point-based read ->', [a.shape for a in g2.get_graphic_data(ct)])

Output:

2 [1 9]
  point-based read -> [(2, 2), (5, 2)]
3 [1 9]
  point-based read -> [(2, 3), (5, 3)]
3 [ 1 13]
  point-based read -> [(1, 3), (6, 3)]

The expected point-based list is [1, 5] in all three cases, and the expected decoded shapes are [(4, d), (3, d)]. Notice that highdicom's encoded value even changes between "3D with common Z" and "3D explicit" for the same geometry. A point index would not change.

3. Effect on existing data in NCI Imaging Data Commons (IDC)

The IDC index was queried for POLYGON/POLYLINE annotation groups:

SELECT i.collection_id, i.analysis_result_id, g.GraphicType,
       COUNT(DISTINCT g.SeriesInstanceUID) AS n_series,
       SUM(g.NumberOfAnnotations) AS n_annotations
FROM ann_group_index g JOIN index i USING (SeriesInstanceUID)
WHERE g.GraphicType IN ('POLYGON', 'POLYLINE')
GROUP BY ALL ORDER BY n_series DESC

Result: 6,075 series containing about 6.17 billion POLYGON annotations, all in the pan_cancer_nuclei_seg_dicom analysis result, across 14 TCGA collections (BRCA, GBM, LUAD, UCEC, COAD, SKCM, LUSC, PRAD, BLCA, STAD, CESC, PAAD, READ, UVM). There are no POLYLINE groups.

To confirm how these files are actually encoded, the smallest one was downloaded and its raw attributes decoded with pydicom:

• SeriesInstanceUID 1.2.826.0.1.3680043.10.511.3.61063851916899676364195955814616690, 11 KB, s3://idc-open-data/64dbbcfb-0934-4f06-a733-ce325b462d56/
• Software Versions: github.com/ImagingDataCommons/idc-sm-annotations-conversion, which uses highdicom

g = ds.AnnotationGroupSequence[0]
coords = np.frombuffer(g.PointCoordinatesData, np.float32)
idx = np.frombuffer(g.LongPrimitivePointIndexList, np.int32)

POLYGON, 8 annotations, 2D, 276 values = 138 points
LongPrimitivePointIndexList = [1 31 57 97 117 155 199 243]

Every entry is odd, which fits offsets into (x, y) value pairs. The last entry, 243, is larger than the number of points (138), so the list cannot hold point indices. These files use the value-based encoding.

Not verified: how dicom-microscopy-viewer / Slim interpret this attribute. Because they display these IDC files correctly, they very likely use the value-based interpretation as well, and would need to change at the same time.

Proposed fix

1. Encoder: use point counts, spans = [item.shape[0] for item in graphic_data].
2. Decoder: use the indices directly as point offsets, split_param = (point_indices - 1)[1:]. To avoid breaking the large body of existing legacy files, detect value-based encoding and warn. With N = total number of points and d = stored values per point:
   • any (idx - 1) % d != 0 → point-based, which is the only possibility
   • max(idx) - 1 >= N → value-based legacy file (warn and divide by d)
   • otherwise → ambiguous; follow the standard (point-based). For a legacy value-based file this only happens when the last annotation holds more than (1 - 1/d) of all points, i.e. more than half for 2D and more than two-thirds for explicit 3D.
   • Groups with a single annotation give [1] under both interpretations and are unaffected.
3. Tests: assert the raw LongPrimitivePointIndexList against hand-built, standard-conformant expected values for 2D, 3D with common Z, and explicit 3D, plus decoding of a legacy value-based list.
4. Cleanup: remove _get_coordinate_index() or fix it.
5. Coordination: release together with dicom-microscopy-viewer, and decide with IDC how to handle the 6,075 existing series (re-conversion vs. relying on tolerant readers). The examples from the expected WG-26 CP would be good test fixtures.