Skip to content

Commit 5d96ae9

Browse files
committed
fix(dicom): keep distinct series with same orientation separate
SeparateOnImageOrientation kept a single cosinesToID lookup across all input volumes. When two series shared an ImageOrientationPatient (e.g. T1 MP-RAGE and a sagittal T2 CUBE acquired in the same session) the second volume's files were reassigned to the first volume's bucket, collapsing them into one scan. Scope cosinesToID per input volume so the orientation-based split only happens within a series. Also fixes the Y-cosine row indexing in areCosinesAlmostEqual: the second loop iteration walked indices 1,2,3 — a sliding window straddling the X and Y rows — instead of 3,4,5, the actual Y cosine row. For axial IOP [1,0,0,0,1,0] this happened to dot to 0 and made same-orientation volumes look mismatched, which masked the cosinesToID leak for axial acquisitions. With the leak fixed the windowed comparison no longer affects outputs, but it's still wrong. While there, honor the epsilon argument instead of always using the EPSILON constant. Adds an e2e regression test that synthesizes two minimal DICOM series with distinct SeriesInstanceUIDs and identical oblique IOP, loads them through two manifest URLs, and asserts two volume cards. Closes #861.
1 parent 9a84cde commit 5d96ae9

5 files changed

Lines changed: 280 additions & 23 deletions

File tree

src/io/itk-dicom/dicom.cpp

Lines changed: 15 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -75,18 +75,6 @@ void replaceChars(std::string &str, char search, char replaceChar) {
7575
}
7676
}
7777

78-
// doesn't actually do any length checks, or overflow checks, or anything
79-
// really.
80-
template <int N>
81-
double dotProduct(const std::vector<double> &vec1,
82-
const std::vector<double> &vec2) {
83-
double result = 0;
84-
for (int i = 0; i < N; i++) {
85-
result += vec1.at(i) * vec2.at(i);
86-
}
87-
return result;
88-
}
89-
9078
std::vector<double> ReadImageOrientationValue(const std::string &filename) {
9179
gdcm::Reader reader;
9280
reader.SetFileName(filename.c_str());
@@ -98,16 +86,18 @@ std::vector<double> ReadImageOrientationValue(const std::string &filename) {
9886
return gdcm::ImageHelper::GetDirectionCosinesValue(file);
9987
}
10088

101-
bool areCosinesAlmostEqual(std::vector<double> cosines1,
102-
std::vector<double> cosines2,
89+
bool areCosinesAlmostEqual(const std::vector<double> &cosines1,
90+
const std::vector<double> &cosines2,
10391
double epsilon = EPSILON) {
104-
for (int i = 0; i <= 1; i++) {
105-
std::vector<double> vec1{cosines1.at(i), cosines1.at(i + 1),
106-
cosines1.at(i + 2)};
107-
std::vector<double> vec2{cosines2.at(i), cosines2.at(i + 1),
108-
cosines2.at(i + 2)};
109-
double dot = dotProduct<3>(vec1, vec2);
110-
if (dot < (1 - EPSILON)) {
92+
// ImageOrientationPatient is two row vectors: X cosines at [0..2] and
93+
// Y cosines at [3..5]. Compare each row.
94+
for (int row = 0; row < 2; row++) {
95+
const int offset = row * 3;
96+
double dot = 0;
97+
for (int i = 0; i < 3; i++) {
98+
dot += cosines1.at(offset + i) * cosines2.at(offset + i);
99+
}
100+
if (dot < (1 - epsilon)) {
111101
return false;
112102
}
113103
}
@@ -116,8 +106,6 @@ bool areCosinesAlmostEqual(std::vector<double> cosines1,
116106

117107
VolumeMapType SeparateOnImageOrientation(const VolumeMapType &volumeMap) {
118108
VolumeMapType newVolumeMap;
119-
// Vector< Pair< cosines, volumeID >>
120-
std::vector<std::pair<std::vector<double>, std::string>> cosinesToID;
121109

122110
// append unique ID part to the volume ID, based on cosines
123111
// The format replaces non-alphanumeric chars to be semi-consistent with DICOM
@@ -141,6 +129,10 @@ VolumeMapType SeparateOnImageOrientation(const VolumeMapType &volumeMap) {
141129
};
142130

143131
for (const auto &[volumeID, names] : volumeMap) {
132+
// Scope the orientation lookup to a single input volume so distinct
133+
// series with identical orientations stay separate.
134+
std::vector<std::pair<std::vector<double>, std::string>> cosinesToID;
135+
144136
for (const auto &filename : names) {
145137
std::vector<double> curCosines = ReadImageOrientationValue(filename);
146138

902 Bytes
Binary file not shown.
811 Bytes
Binary file not shown.
Lines changed: 69 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,69 @@
1+
// Regression for https://github.com/Kitware/VolView/issues/861
2+
//
3+
// Two DICOM series with distinct SeriesInstanceUIDs but identical
4+
// ImageOrientationPatient must remain two volume cards, not collapse
5+
// into one merged scan.
6+
//
7+
// Synthetic DICOMs are generated on the fly so the test has no external
8+
// dependencies.
9+
import * as path from 'path';
10+
import * as fs from 'fs';
11+
import { volViewPage } from '../pageobjects/volview.page';
12+
import { TEMP_DIR } from '../../wdio.shared.conf';
13+
import { buildSyntheticDicom, newUid } from './syntheticDicom';
14+
15+
// Oblique ImageOrientationPatient — same value for both series. With
16+
// this orientation the pre-fix areCosinesAlmostEqual sees a non-zero
17+
// second window and the cross-volume cosinesToID leak collapses both
18+
// series into one.
19+
const SHARED_IMAGE_ORIENTATION_PATIENT = [
20+
-0.00964, 0.99248, 0.12202, 0.06932, 0.12239, -0.99006,
21+
] as const;
22+
23+
const SLICE_COUNT = 5;
24+
const STUDY_UID = newUid();
25+
26+
function writeSeries(label: 'A' | 'B', dirName: string, manifestName: string) {
27+
const seriesUid = newUid();
28+
const dir = path.join(TEMP_DIR, dirName);
29+
fs.mkdirSync(dir, { recursive: true });
30+
31+
const resources: { url: string; name: string }[] = [];
32+
for (let i = 0; i < SLICE_COUNT; i++) {
33+
const filename = `slice-${i}.dcm`;
34+
const bytes = buildSyntheticDicom({
35+
studyUid: STUDY_UID,
36+
seriesUid,
37+
sopUid: newUid(),
38+
instanceNumber: i + 1,
39+
imageOrientationPatient: SHARED_IMAGE_ORIENTATION_PATIENT,
40+
// Step along the slice-normal so GDCM can sort within the series.
41+
imagePositionPatient: [0, 0, i],
42+
});
43+
fs.writeFileSync(path.join(dir, filename), bytes);
44+
resources.push({ url: `tmp/${dirName}/${filename}`, name: filename });
45+
}
46+
47+
fs.writeFileSync(
48+
path.join(TEMP_DIR, manifestName),
49+
JSON.stringify({ resources })
50+
);
51+
}
52+
53+
describe('Multi-series load: two series with identical ImageOrientationPatient', () => {
54+
before(() => {
55+
writeSeries('A', 'multi-series-A', 'multi-series-A.json');
56+
writeSeries('B', 'multi-series-B', 'multi-series-B.json');
57+
});
58+
59+
it('keeps the two series separate (two volume cards)', async () => {
60+
await volViewPage.open(
61+
'?urls=[tmp/multi-series-A.json,tmp/multi-series-B.json]'
62+
);
63+
await volViewPage.waitForViews();
64+
await browser.waitUntil(
65+
async () => (await $$('.volume-card').length) === 2,
66+
{ timeout: 30000, timeoutMsg: 'expected exactly two volume cards' }
67+
);
68+
});
69+
});

tests/specs/syntheticDicom.ts

Lines changed: 196 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,196 @@
1+
// Minimal synthetic DICOM (Explicit VR Little Endian) for tests.
2+
// Emits just enough tags for ITK/GDCM to categorize and load a series:
3+
// SOP Class/Instance UIDs, Study/SeriesInstanceUID, SeriesNumber, Modality,
4+
// Patient identifiers, ImageOrientationPatient, ImagePositionPatient,
5+
// PixelSpacing, SliceThickness, image geometry, and zeroed PixelData.
6+
7+
const SOP_CLASS_MR = '1.2.840.10008.5.1.4.1.1.4';
8+
const TS_EXPLICIT_VR_LE = '1.2.840.10008.1.2.1';
9+
10+
const enc = new TextEncoder();
11+
12+
const padUi = (s: string) => (s.length % 2 === 0 ? s : `${s}\0`);
13+
const padText = (s: string) => (s.length % 2 === 0 ? s : `${s} `);
14+
15+
const writeShort = (v: number) => {
16+
const b = new Uint8Array(2);
17+
new DataView(b.buffer).setUint16(0, v, true);
18+
return b;
19+
};
20+
21+
const writeLong = (v: number) => {
22+
const b = new Uint8Array(4);
23+
new DataView(b.buffer).setUint32(0, v, true);
24+
return b;
25+
};
26+
27+
const tagBytes = (group: number, element: number) => {
28+
const b = new Uint8Array(4);
29+
const dv = new DataView(b.buffer);
30+
dv.setUint16(0, group, true);
31+
dv.setUint16(2, element, true);
32+
return b;
33+
};
34+
35+
const combine = (...arr: Uint8Array[]) => {
36+
const total = arr.reduce((acc, a) => acc + a.length, 0);
37+
const out = new Uint8Array(total);
38+
let off = 0;
39+
for (const a of arr) {
40+
out.set(a, off);
41+
off += a.length;
42+
}
43+
return out;
44+
};
45+
46+
// Short-form explicit VR (2-byte length): UI, CS, DA, DS, IS, LO, PN, SH, US, UL, ...
47+
const elemShort = (
48+
group: number,
49+
element: number,
50+
vr: string,
51+
value: Uint8Array
52+
) =>
53+
combine(
54+
tagBytes(group, element),
55+
enc.encode(vr),
56+
writeShort(value.length),
57+
value
58+
);
59+
60+
// Long-form explicit VR (2-byte reserved + 4-byte length): OB, OW, UN, SQ, ...
61+
const elemLong = (
62+
group: number,
63+
element: number,
64+
vr: string,
65+
value: Uint8Array
66+
) =>
67+
combine(
68+
tagBytes(group, element),
69+
enc.encode(vr),
70+
new Uint8Array(2),
71+
writeLong(value.length),
72+
value
73+
);
74+
75+
const ui = (g: number, e: number, v: string) =>
76+
elemShort(g, e, 'UI', enc.encode(padUi(v)));
77+
const cs = (g: number, e: number, v: string) =>
78+
elemShort(g, e, 'CS', enc.encode(padText(v)));
79+
const pn = (g: number, e: number, v: string) =>
80+
elemShort(g, e, 'PN', enc.encode(padText(v)));
81+
const lo = (g: number, e: number, v: string) =>
82+
elemShort(g, e, 'LO', enc.encode(padText(v)));
83+
const sh = (g: number, e: number, v: string) =>
84+
elemShort(g, e, 'SH', enc.encode(padText(v)));
85+
const da = (g: number, e: number, v: string) =>
86+
elemShort(g, e, 'DA', enc.encode(padText(v)));
87+
const is = (g: number, e: number, v: string) =>
88+
elemShort(g, e, 'IS', enc.encode(padText(v)));
89+
const ds = (g: number, e: number, v: string) =>
90+
elemShort(g, e, 'DS', enc.encode(padText(v)));
91+
const us = (g: number, e: number, v: number) =>
92+
elemShort(g, e, 'US', writeShort(v));
93+
94+
export type SyntheticSliceOptions = {
95+
studyUid: string;
96+
seriesUid: string;
97+
sopUid: string;
98+
instanceNumber: number;
99+
imageOrientationPatient: readonly [
100+
number,
101+
number,
102+
number,
103+
number,
104+
number,
105+
number,
106+
];
107+
imagePositionPatient: readonly [number, number, number];
108+
rows?: number;
109+
cols?: number;
110+
pixelSpacing?: readonly [number, number];
111+
sliceThickness?: number;
112+
modality?: string;
113+
patientName?: string;
114+
patientId?: string;
115+
seriesNumber?: number;
116+
studyDate?: string;
117+
};
118+
119+
export function buildSyntheticDicom(opts: SyntheticSliceOptions): Uint8Array {
120+
const {
121+
studyUid,
122+
seriesUid,
123+
sopUid,
124+
instanceNumber,
125+
imageOrientationPatient,
126+
imagePositionPatient,
127+
rows = 4,
128+
cols = 4,
129+
pixelSpacing = [1, 1] as const,
130+
sliceThickness = 1,
131+
modality = 'MR',
132+
patientName = 'TEST',
133+
patientId = 'TEST001',
134+
seriesNumber = 1,
135+
studyDate = '20260101',
136+
} = opts;
137+
138+
const dataset = combine(
139+
ui(0x0008, 0x0016, SOP_CLASS_MR),
140+
ui(0x0008, 0x0018, sopUid),
141+
da(0x0008, 0x0020, studyDate),
142+
da(0x0008, 0x0021, studyDate),
143+
cs(0x0008, 0x0060, modality),
144+
pn(0x0010, 0x0010, patientName),
145+
lo(0x0010, 0x0020, patientId),
146+
da(0x0010, 0x0030, '19700101'),
147+
cs(0x0010, 0x0040, 'O'),
148+
ds(0x0018, 0x0050, String(sliceThickness)),
149+
ui(0x0020, 0x000d, studyUid),
150+
ui(0x0020, 0x000e, seriesUid),
151+
sh(0x0020, 0x0010, '1'),
152+
is(0x0020, 0x0011, String(seriesNumber)),
153+
is(0x0020, 0x0013, String(instanceNumber)),
154+
ds(
155+
0x0020,
156+
0x0032,
157+
imagePositionPatient.map((n) => n.toString()).join('\\')
158+
),
159+
ds(
160+
0x0020,
161+
0x0037,
162+
imageOrientationPatient.map((n) => n.toString()).join('\\')
163+
),
164+
us(0x0028, 0x0002, 1),
165+
cs(0x0028, 0x0004, 'MONOCHROME2'),
166+
us(0x0028, 0x0010, rows),
167+
us(0x0028, 0x0011, cols),
168+
ds(0x0028, 0x0030, pixelSpacing.map((n) => n.toString()).join('\\')),
169+
us(0x0028, 0x0100, 16),
170+
us(0x0028, 0x0101, 16),
171+
us(0x0028, 0x0102, 15),
172+
us(0x0028, 0x0103, 0),
173+
elemLong(0x7fe0, 0x0010, 'OW', new Uint8Array(rows * cols * 2))
174+
);
175+
176+
const fileMetaBody = combine(
177+
elemLong(0x0002, 0x0001, 'OB', new Uint8Array([0x00, 0x01])),
178+
ui(0x0002, 0x0002, SOP_CLASS_MR),
179+
ui(0x0002, 0x0003, sopUid),
180+
ui(0x0002, 0x0010, TS_EXPLICIT_VR_LE)
181+
);
182+
const fileMeta = combine(
183+
elemShort(0x0002, 0x0000, 'UL', writeLong(fileMetaBody.length)),
184+
fileMetaBody
185+
);
186+
187+
return combine(new Uint8Array(128), enc.encode('DICM'), fileMeta, dataset);
188+
}
189+
190+
// Pseudo-UID unique per call within a test process. All-numeric so it
191+
// stays within DICOM UID character constraints.
192+
let counter = 0;
193+
export const newUid = () => {
194+
counter += 1;
195+
return `1.2.826.0.1.3680043.10.999.${process.pid}.${Date.now()}.${counter}`;
196+
};

0 commit comments

Comments
 (0)