diff --git a/package-lock.json b/package-lock.json index c50b7d12a..8c1c0594f 100644 --- a/package-lock.json +++ b/package-lock.json @@ -23,6 +23,7 @@ "@sentry/vue": "^10.27.0", "@thi.ng/api": "^8.12.9", "@thi.ng/rasterize": "^1.0.171", + "@types/color-name": "^2.0.0", "@types/cors": "^2.8.19", "@types/deep-equal": "^1.0.4", "@types/express": "^5.0.5", @@ -39,6 +40,7 @@ "@wdio/spec-reporter": "^9.20.0", "@wdio/static-server-service": "^9.20.0", "@wdio/visual-service": "^10.1.0", + "color-name": "^1.1.4", "comlink": "^4.4.2", "concurrently": "^10.0.4", "core-js": "3.47.0", @@ -5142,6 +5144,13 @@ "assertion-error": "^2.0.1" } }, + "node_modules/@types/color-name": { + "version": "2.0.0", + "resolved": "https://registry.npmjs.org/@types/color-name/-/color-name-2.0.0.tgz", + "integrity": "sha512-63mTjolMJv75upGaUbT6J3lRDWl6pETPQsaWni9w3dMArhNBpgtHkX8ISb9zLV3YYLPA/SMk8ZGALa3k9WY/aQ==", + "dev": true, + "license": "MIT" + }, "node_modules/@types/connect": { "version": "3.4.38", "resolved": "https://registry.npmjs.org/@types/connect/-/connect-3.4.38.tgz", diff --git a/package.json b/package.json index 71cfe5912..93fe9da79 100644 --- a/package.json +++ b/package.json @@ -49,6 +49,7 @@ "@sentry/vue": "^10.27.0", "@thi.ng/api": "^8.12.9", "@thi.ng/rasterize": "^1.0.171", + "@types/color-name": "^2.0.0", "@types/cors": "^2.8.19", "@types/deep-equal": "^1.0.4", "@types/express": "^5.0.5", @@ -65,6 +66,7 @@ "@wdio/spec-reporter": "^9.20.0", "@wdio/static-server-service": "^9.20.0", "@wdio/visual-service": "^10.1.0", + "color-name": "^1.1.4", "comlink": "^4.4.2", "concurrently": "^10.0.4", "core-js": "3.47.0", diff --git a/src/store/__tests__/segmentImportPath.spec.ts b/src/io/resample/__tests__/ensureSameSpace.spec.ts similarity index 75% rename from src/store/__tests__/segmentImportPath.spec.ts rename to src/io/resample/__tests__/ensureSameSpace.spec.ts index 043e2e18e..910b1d77d 100644 --- a/src/store/__tests__/segmentImportPath.spec.ts +++ b/src/io/resample/__tests__/ensureSameSpace.spec.ts @@ -1,11 +1,8 @@ import { afterEach, describe, expect, it, vi } from 'vitest'; -import { createPinia, setActivePinia } from 'pinia'; import { resolve } from 'node:path'; import { InterfaceTypes, runPipelineNode, type Image } from 'itk-wasm'; import vtkDataArray from '@kitware/vtk.js/Common/Core/DataArray'; import vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; -import { useSegmentGroupStore } from '@/src/store/segmentGroups'; -import { useImageCacheStore } from '@/src/store/image-cache'; import { ensureSameSpace } from '@/src/io/resample/resample'; import * as wasm from '@/src/io/resample/itkWasmUtils'; @@ -53,26 +50,16 @@ describe('labelmap import index alignment', () => { expect(runWasm).not.toHaveBeenCalled(); }); - it('resamples a reflected child onto the parent grid with label interpolation', async () => { - setActivePinia(createPinia()); + it('reorients a reflected label child onto the parent grid without resampling', async () => { const runWasm = useNodeResampling(); const parent = makeImage(); const child = makeImage(); child.setOrigin([0, 2, 0]); child.setDirection(1, 0, 0, 0, -1, 0, 0, 0, 1); child.getPointData().getScalars().setComponent(0, 0, 7); - const cache = useImageCacheStore(); - cache.addVTKImageData(parent, 'CTA', { id: 'parent' }); - cache.addVTKImageData(child, 'cta-head-neck-total.seg.nii.gz', { - id: 'child', - }); - - const store = useSegmentGroupStore(); - const [id] = await store.convertImageToLabelmap('child', 'parent'); - const imported = store.dataIndex[id]; + const imported = await ensureSameSpace(parent, child, true); - expect(runWasm).toHaveBeenCalledOnce(); - expect(runWasm.mock.calls[0][1]).toContain('--label'); + expect(runWasm).not.toHaveBeenCalled(); expect(imported.getDirection()).toEqual(parent.getDirection()); expect(imported.getOrigin()).toEqual(parent.getOrigin()); expect(Array.from(imported.getPointData().getScalars().getData())).toEqual([ diff --git a/src/io/resample/__tests__/reorientLabelImage.spec.ts b/src/io/resample/__tests__/reorientLabelImage.spec.ts new file mode 100644 index 000000000..df79b6a03 --- /dev/null +++ b/src/io/resample/__tests__/reorientLabelImage.spec.ts @@ -0,0 +1,129 @@ +import { describe, expect, it } from 'vitest'; +import vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; +import vtkDataArray from '@kitware/vtk.js/Common/Core/DataArray'; +import { reorientLabelImage } from '../reorientLabelImage'; + +const image = (dimensions: [number, number, number]) => { + const result = vtkImageData.newInstance(); + result.setDimensions(dimensions); + result.getPointData().setScalars( + vtkDataArray.newInstance({ + values: Uint16Array.from( + { length: dimensions.reduce((a, b) => a * b, 1) }, + (_, i) => i + 1 + ), + }) + ); + return result; +}; + +describe('label-grid reorientation', () => { + it.each([0, 1, 2, 3, 4, 5, 6, 7])( + 'preserves every label under axis flips %i', + (flips) => { + const source = image([3, 4, 5]); + const target = image([3, 4, 5]); + const direction: [ + number, + number, + number, + number, + number, + number, + number, + number, + number, + ] = [1, 0, 0, 0, 1, 0, 0, 0, 1]; + const origin: [number, number, number] = [0, 0, 0]; + [3, 4, 5].forEach((size, axis) => { + if (flips & (1 << axis)) { + direction[axis * 4] = -1; + origin[axis] = size - 1; + } + }); + source.setDirection(direction); + source.setOrigin(origin); + const output = reorientLabelImage(target, source)!; + const values = output.getPointData().getScalars().getData(); + for (let k = 0; k < 5; k++) + for (let j = 0; j < 4; j++) + for (let i = 0; i < 3; i++) { + const x = flips & 1 ? 2 - i : i, + y = flips & 2 ? 3 - j : j, + z = flips & 4 ? 4 - k : k; + expect(values[i + 3 * (j + 4 * k)]).toBe(1 + x + 3 * (y + 4 * z)); + } + expect(output.getDirection()).toEqual(target.getDirection()); + expect(source.getPointData().getScalars().getData()[0]).toBe(1); + } + ); + + it('permutes unequal axes', () => { + const source = image([3, 4, 5]); + const target = image([4, 3, 5]); + target.setDirection([0, 1, 0, 1, 0, 0, 0, 0, 1]); + const output = reorientLabelImage(target, source)!; + expect(output.getExtent()).toEqual(target.getExtent()); + const values = output.getPointData().getScalars().getData(); + for (let k = 0; k < 5; k++) + for (let j = 0; j < 3; j++) + for (let i = 0; i < 4; i++) + expect(values[i + 4 * (j + 3 * k)]).toBe(1 + j + 3 * (i + 4 * k)); + }); + + it('defers nonzero extents to the general path', () => { + const source = image([3, 4, 5]); + source.setExtent(1, 3, 2, 5, 3, 7); + expect(reorientLabelImage(source, source)).toBeNull(); + }); + + it('accepts DICOM precision differences without changing labels', () => { + const source = image([3, 4, 5]); + const target = image([3, 4, 5]); + source.setOrigin([0.000004, 0.000004, 0.000012]); + expect([ + ...reorientLabelImage(target, source)! + .getPointData() + .getScalars() + .getData(), + ]).toEqual([...source.getPointData().getScalars().getData()]); + }); + + it('returns the source itself when it already sits on the target grid', () => { + const spacing: [number, number, number] = [0.7, 0.7, 3]; + const origin: [number, number, number] = [-120.1, -98.4, 33.7]; + const source = image([3, 4, 5]); + const target = image([3, 4, 5]); + [source, target].forEach((im) => { + im.setSpacing(spacing); + im.setOrigin(origin); + }); + expect(reorientLabelImage(target, source)).toBe(source); + + // The same geometry laid out along a flipped axis is a different grid and + // still has to go through the reslice. + const flipped = image([3, 4, 5]); + flipped.setSpacing(spacing); + flipped.setOrigin([origin[0] + spacing[0] * 2, origin[1], origin[2]]); + flipped.setDirection([-1, 0, 0, 0, 1, 0, 0, 0, 1]); + const output = reorientLabelImage(target, flipped)!; + expect(output).not.toBe(flipped); + const values = output.getPointData().getScalars().getData(); + for (let k = 0; k < 5; k++) + for (let j = 0; j < 4; j++) + for (let i = 0; i < 3; i++) + expect(values[i + 3 * (j + 4 * k)]).toBe(1 + (2 - i) + 3 * (j + 4 * k)); + }); + + it('defers fractional shifts, different sampling, and cropping to interpolation', () => { + const source = image([3, 4, 5]); + const target = image([3, 4, 5]); + source.setOrigin([0.25, 0, 0]); + expect(reorientLabelImage(target, source)).toBeNull(); + source.setOrigin([0, 0, 0]); + source.setSpacing([0.5, 1, 1]); + expect(reorientLabelImage(target, source)).toBeNull(); + source.setSpacing([1, 1, 1]); + expect(reorientLabelImage(image([1, 4, 5]), source)).toBeNull(); + }); +}); diff --git a/src/io/resample/reorientLabelImage.ts b/src/io/resample/reorientLabelImage.ts new file mode 100644 index 000000000..c7ef31375 --- /dev/null +++ b/src/io/resample/reorientLabelImage.ts @@ -0,0 +1,86 @@ +import vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; +import vtkImageReslice from '@kitware/vtk.js/Imaging/Core/ImageReslice'; +import { InterpolationMode } from '@kitware/vtk.js/Imaging/Core/AbstractImageInterpolator/Constants'; +import { mat4, vec3 } from 'gl-matrix'; + +/** Reorder an equivalent voxel grid without interpolating label boundaries. */ +export function reorientLabelImage(target: vtkImageData, source: vtkImageData) { + // ImageReslice produces zero-based output extents. + if ( + [target, source].some((image) => + [0, 2, 4].some((axis) => image.getExtent()[axis] !== 0) + ) + ) + return null; + const matrix = mat4.multiply( + mat4.create(), + source.getWorldToIndex(), + target.getIndexToWorld() + ); + const tolerance = 1e-3; + const axes = [0, 1, 2].map((column) => { + const values = [0, 1, 2].map((row) => matrix[4 * column + row]); + const axis = values.findIndex((value) => Math.abs(value) > 0.5); + if ( + axis < 0 || + values.some((value, row) => + row === axis + ? Math.abs(Math.abs(value) - 1) > tolerance + : Math.abs(value) > tolerance + ) + ) + return -1; + return axis; + }); + if (axes.includes(-1) || new Set(axes).size !== 3) return null; + const sourceSize = source.getDimensions(); + if ( + target.getDimensions().some((size, axis) => size !== sourceSize[axes[axis]]) + ) + return null; + + // Check the entire extent so rounding error cannot accumulate into a shift + // at the far edge. Matching physical bounds alone does not imply equal grids. + const from = target.getExtent(); + const to = source.getExtent(); + for (let corner = 0; corner < 8; corner++) { + const point = vec3.fromValues( + from[corner & 1], + from[2 + ((corner >> 1) & 1)], + from[4 + ((corner >> 2) & 1)] + ); + vec3.transformMat4(point, point, matrix); + if ( + [0, 1, 2].some( + (axis) => + Math.min( + Math.abs(point[axis] - to[axis * 2]), + Math.abs(point[axis] - to[axis * 2 + 1]) + ) > tolerance + ) + ) + return null; + } + + // The corner check just proved both grids coincide within `tolerance`, so an + // unpermuted, unflipped mapping means the source already sits on the target + // grid. Demanding a bit-exact identity here instead would reslice every real + // image, whose transforms never multiply back to exactly one. + if (axes.every((axis, column) => axis === column && matrix[5 * column] > 0)) + return source; + + const filter = vtkImageReslice.newInstance(); + filter.setOutputOrigin(target.getOrigin()); + filter.setOutputSpacing(target.getSpacing()); + filter.setOutputDirection(target.getDirection()); + filter.setOutputExtent(target.getExtent()); + filter.setOutputDimensionality(3); + filter.setTransformInputSampling(false); + filter.setInterpolationMode(InterpolationMode.NEAREST); + try { + filter.setInputData(source); + return filter.getOutputData() as vtkImageData; + } finally { + filter.delete(); + } +} diff --git a/src/io/resample/resample.ts b/src/io/resample/resample.ts index 41b302ec2..2d96fddaf 100644 --- a/src/io/resample/resample.ts +++ b/src/io/resample/resample.ts @@ -4,6 +4,7 @@ import vtkITKHelper from '@kitware/vtk.js/Common/DataModel/ITKHelper'; import { compareImageIndexGrids } from '@/src/utils/imageSpace'; import { shallowCopyImageData } from '@/src/utils/vtk-helpers'; import { runWasm } from './itkWasmUtils'; +import { reorientLabelImage } from './reorientLabelImage'; export async function resample(fixed: Image, moving: Image, label = false) { const labelFlag = label ? ['--label'] : []; @@ -30,7 +31,12 @@ export async function ensureSameSpace( ) { // Callers own what they get back and may hand it to something that disposes // it, so never return the candidate itself. - if (compareImageIndexGrids(target, resampleCandidate)) { + if (label) { + const reoriented = reorientLabelImage(target, resampleCandidate); + if (reoriented === resampleCandidate) + return shallowCopyImageData(resampleCandidate); + if (reoriented) return reoriented; + } else if (compareImageIndexGrids(target, resampleCandidate)) { return shallowCopyImageData(resampleCandidate); } const itkImage = await resample( diff --git a/src/segmentation/__tests__/model.spec.ts b/src/segmentation/__tests__/model.spec.ts new file mode 100644 index 000000000..f2fe69850 --- /dev/null +++ b/src/segmentation/__tests__/model.spec.ts @@ -0,0 +1,297 @@ +import { describe, expect, it } from 'vitest'; +import type { RGBAColor } from '@kitware/vtk.js/types'; + +import { TOOL_COLORS } from '@/src/config'; +import { + emptyExtent, + isEmptyExtent, + markedExtent, +} from '@/src/segmentation/geometry'; +import { + cssColorToRGBA, + tryCssColorToRGBA, + rgbaToCssColor, +} from '@/src/segmentation/color'; +import type { + LabelmapBinding, + SegmentMask, + Segmentation, +} from '@/src/segmentation/model'; +import type { Extent3D } from '@/src/segmentation/geometry'; +import { resolveSegmentAppearance } from '@/src/segmentation/segment'; +import vtkLabelMap from '@/src/vtk/LabelMap'; + +describe('emptyExtent', () => { + it('is the pinned empty sentinel', () => { + expect(emptyExtent()).toEqual([0, -1, 0, -1, 0, -1]); + }); + + it('returns a fresh extent per call', () => { + const first = emptyExtent(); + first[1] = 10; + + expect(emptyExtent()).toEqual([0, -1, 0, -1, 0, -1]); + }); +}); + +describe('isEmptyExtent', () => { + it('accepts the empty sentinel', () => { + expect(isEmptyExtent(emptyExtent())).toBe(true); + }); + + it('rejects an extent covering a whole image', () => { + expect(isEmptyExtent([0, 9, 0, 19, 0, 29])).toBe(false); + }); + + it('rejects a single voxel extent', () => { + // vtk.js extents are inclusive, so min === max is one voxel, not empty. + expect(isEmptyExtent([4, 4, 5, 5, 6, 6])).toBe(false); + }); + + it.each([ + ['i', [5, 4, 0, 9, 0, 9] as Extent3D], + ['j', [0, 9, 5, 4, 0, 9] as Extent3D], + ['k', [0, 9, 0, 9, 5, 4] as Extent3D], + ])('is empty when the %s axis is inverted', (_axis, extent) => { + expect(isEmptyExtent(extent)).toBe(true); + }); +}); + +describe('cssColorToRGBA', () => { + it('parses a tool color hex string', () => { + expect(cssColorToRGBA('#58f24c')).toEqual([88, 242, 76, 255]); + }); + + it('parses an alpha channel when present', () => { + expect(cssColorToRGBA('#58f24c80')).toEqual([88, 242, 76, 128]); + }); + + it('parses the named color used by the vector tool label defaults', () => { + expect(cssColorToRGBA('red')).toEqual([255, 0, 0, 255]); + }); + + it.each([ + ['orange', [255, 165, 0, 255]], + ['rebeccapurple', [102, 51, 153, 255]], + ])('parses the CSS color keyword %s', (name, expected) => { + expect(cssColorToRGBA(name)).toEqual(expected); + }); +}); + +describe('rgbaToCssColor', () => { + it.each([ + [[88, 242, 76, 255] as RGBAColor], + [[0, 0, 0, 255] as RGBAColor], + [[214, 0, 0, 128] as RGBAColor], + ])('round trips %j', (rgba) => { + expect(cssColorToRGBA(rgbaToCssColor(rgba))).toEqual(rgba); + }); + + it('emits a css color literal', () => { + expect(rgbaToCssColor([88, 242, 76, 255])).toMatch( + /^(#[0-9a-fA-F]{6,8}|rgba?\(.+\))$/ + ); + }); +}); + +describe('color conversion of existing label colors', () => { + it.each(TOOL_COLORS)('round trips %s', (css) => { + const rgba = cssColorToRGBA(css); + + expect(cssColorToRGBA(rgbaToCssColor(rgba))).toEqual(rgba); + }); +}); + +describe('mask model', () => { + it('holds no labelmap binding until voxels are allocated', () => { + const segment: SegmentMask = { + id: 'segment-1', + segmentId: 'segment-1', + representations: {}, + }; + + expect(segment.representations.labelmap).toBeUndefined(); + }); + + it('binds a segment to its own voxels and the box they cover', () => { + const image = vtkLabelMap.newInstance(); + const binding: LabelmapBinding = { + image, + extent: [0, 9, 0, 19, 0, 29], + name: 'Tumor', + }; + const segment: SegmentMask = { + id: 'segment-1', + segmentId: 'segment-1', + representations: { labelmap: binding }, + }; + + expect(segment.representations.labelmap?.image).toBe(image); + expect(isEmptyExtent(binding.extent)).toBe(false); + }); + + it('keeps record order separate from the records themselves', () => { + const makeMask = (id: string, segmentId: string): SegmentMask => ({ + id, + segmentId, + representations: {}, + }); + const segmentation: Segmentation = { + id: 'segmentation-1', + name: 'Segmentation', + parentImageId: 'image-1', + masks: { + 'segment-1': makeMask('segment-1', 'segment-1'), + 'segment-2': makeMask('segment-2', 'segment-2'), + }, + order: ['segment-2', 'segment-1'], + fillOpacity: 1, + outlineOpacity: 1, + outlineThickness: 2, + }; + + expect(segmentation.order).toEqual(['segment-2', 'segment-1']); + expect(Object.keys(segmentation.masks).sort()).toEqual([ + 'segment-1', + 'segment-2', + ]); + }); +}); + +describe('the appearance resolver', () => { + it('fills the app defaults for what a type leaves unset', () => { + const resolved = resolveSegmentAppearance({ + id: 'segment-1', + name: 'Tumor', + color: [255, 0, 0, 255], + visible: true, + locked: false, + }); + + expect(resolved).toMatchObject({ + name: 'Tumor', + cssColor: '#ff0000', + fillOpacity: 1, + outlineOpacity: 1, + }); + }); + + it('keeps what a type does state', () => { + const resolved = resolveSegmentAppearance({ + id: 'segment-1', + name: 'Tumor', + color: [255, 0, 0, 255], + visible: false, + locked: true, + fillOpacity: 0.25, + strokeWidth: 3, + }); + + expect(resolved.fillOpacity).toBe(0.25); + expect(resolved.strokeWidth).toBe(3); + expect(resolved.visible).toBe(false); + expect(resolved.locked).toBe(true); + }); + + it('answers for a type that is gone', () => { + const resolved = resolveSegmentAppearance(undefined); + + expect(resolved.name).toBe(''); + expect(resolved.fillOpacity).toBe(1); + expect(resolved.visible).toBe(true); + expect(resolved.locked).toBe(false); + }); + + describe('CSS colors', () => { + // Functional notation is not supported. A config using it is told so at the + // boundary that reads the file, rather than resolving to a plausible black. + it.each([ + 'rgb(0, 255, 0)', + 'rgb(0 255 0)', + 'rgba(255, 0, 0, 0.5)', + 'hsl(120, 100%, 50%)', + ])('does not parse %s', (css) => { + expect(tryCssColorToRGBA(css)).toBeUndefined(); + }); + + it('treats transparent as fully transparent, not black', () => { + expect(cssColorToRGBA('transparent')).toEqual([0, 0, 0, 0]); + }); + + it('reports unparseable input rather than silently blackening it', () => { + expect(tryCssColorToRGBA('not-a-color')).toBeUndefined(); + expect(cssColorToRGBA('not-a-color')).toEqual([0, 0, 0, 255]); + }); + + // A config.json label color and a 6.x state file both reach here unvalidated, + // and the migration that calls this is not wrapped in a try. + it.each(['constructor', '__proto__', 'toString', 'valueOf'])( + 'treats the inherited property name %s as unparseable', + (name) => { + expect(tryCssColorToRGBA(name)).toBeUndefined(); + expect(cssColorToRGBA(name)).toEqual([0, 0, 0, 255]); + } + ); + }); +}); + +// A binding's extent is the allocation: paint pads it on growth and an erase +// never shrinks it, so the box a segment occupies has to be read off the voxels. +describe('markedExtent', () => { + const VALUE = 3; + + // A mask bounded to `extent`, holding VALUE at each parent index in `marks`. + const scalarsOf = ( + extent: Extent3D, + marks: Array<[number, number, number]> + ) => { + const si = extent[1] - extent[0] + 1; + const sj = extent[3] - extent[2] + 1; + const sk = extent[5] - extent[4] + 1; + const scalars = new Uint8Array(si * sj * sk); + marks.forEach(([i, j, k]) => { + scalars[ + i - extent[0] + (j - extent[2]) * si + (k - extent[4]) * si * sj + ] = VALUE; + }); + return scalars; + }; + + it('bounds the marked voxels, not the allocation they sit in', () => { + const extent: Extent3D = [0, 5, 0, 5, 0, 5]; + const scalars = scalarsOf(extent, [ + [1, 2, 3], + [4, 2, 3], + [2, 5, 1], + ]); + + expect(markedExtent(scalars, extent, VALUE)).toEqual([1, 4, 2, 5, 1, 3]); + }); + + it('reads a single voxel as its own box', () => { + const extent: Extent3D = [2, 4, 2, 4, 2, 4]; + const scalars = scalarsOf(extent, [[3, 3, 3]]); + + expect(markedExtent(scalars, extent, VALUE)).toEqual([3, 3, 3, 3, 3, 3]); + }); + + it('reports an empty box for a mask holding nothing of that value', () => { + const extent: Extent3D = [0, 3, 0, 3, 0, 3]; + const scalars = scalarsOf(extent, [[1, 1, 1]]); + + expect(isEmptyExtent(markedExtent(scalars, extent, VALUE + 1))).toBe(true); + expect(isEmptyExtent(markedExtent(new Uint8Array(64), extent, VALUE))).toBe( + true + ); + }); + + it('reads a mask whose own origin is not the image origin', () => { + const extent: Extent3D = [4, 6, 7, 9, 1, 2]; + const scalars = scalarsOf(extent, [ + [5, 8, 1], + [6, 9, 2], + ]); + + expect(markedExtent(scalars, extent, VALUE)).toEqual([5, 6, 8, 9, 1, 2]); + }); +}); diff --git a/src/segmentation/__tests__/reframeMaskScalars.spec.ts b/src/segmentation/__tests__/reframeMaskScalars.spec.ts new file mode 100644 index 000000000..ba4554f60 --- /dev/null +++ b/src/segmentation/__tests__/reframeMaskScalars.spec.ts @@ -0,0 +1,44 @@ +import { describe, expect, it } from 'vitest'; +import { reframeMaskScalars } from '@/src/segmentation/masks/storage'; +import type { Extent3D } from '@/src/segmentation/geometry'; + +describe('reframing into a reusable mask buffer', () => { + it('overwrites old labels and padding when the source moves or is erased', () => { + const to: Extent3D = [4, 6, 5, 7, 6, 8]; + const output = new Uint8Array(27).fill(9); + const source = new Uint8Array([1]); + expect(reframeMaskScalars(source, [5, 5, 6, 6, 7, 7], to, output)).toBe( + output + ); + expect([...output]).toEqual( + Array.from({ length: 27 }, (_, i) => (i === 13 ? 1 : 0)) + ); + reframeMaskScalars(source, [6, 6, 7, 7, 8, 8], to, output); + expect([...output]).toEqual( + Array.from({ length: 27 }, (_, i) => (i === 26 ? 1 : 0)) + ); + source[0] = 0; + reframeMaskScalars(source, [6, 6, 7, 7, 8, 8], to, output); + expect([...output]).toEqual(Array(27).fill(0)); + }); + + it('clips source rows at the destination and clears disjoint copies', () => { + const from: Extent3D = [2, 5, 3, 4, 4, 4]; + const to: Extent3D = [3, 4, 4, 5, 4, 4]; + const source = new Uint8Array([1, 2, 3, 4, 5, 6, 7, 8]); + const output = new Uint8Array(4).fill(9); + reframeMaskScalars(source, from, to, output); + expect([...output]).toEqual([6, 7, 0, 0]); + expect([...source]).toEqual([1, 2, 3, 4, 5, 6, 7, 8]); + reframeMaskScalars(source, from, [3, 4, 9, 10, 4, 4], output); + expect([...output]).toEqual([0, 0, 0, 0]); + }); + + it('rejects a destination buffer of the wrong size before changing it', () => { + const output = new Uint8Array([9]); + expect(() => + reframeMaskScalars([1], [0, 0, 0, 0, 0, 0], [0, 1, 0, 0, 0, 0], output) + ).toThrow('Mask output size mismatch'); + expect([...output]).toEqual([9]); + }); +}); diff --git a/src/segmentation/color.ts b/src/segmentation/color.ts new file mode 100644 index 000000000..832a59fb5 --- /dev/null +++ b/src/segmentation/color.ts @@ -0,0 +1,41 @@ +import type { RGBAColor } from '@kitware/vtk.js/types'; +import colorNames from 'color-name'; +import { hexaToRGBA, rgbaToHexa } from '@/src/utils/color'; + +// A Map, not the package's plain object: prototype keys such as 'constructor' +// must not be mistaken for colors. +const NAMED_COLORS = new Map(Object.entries(colorNames)); + +const HEX_COLOR = /^#?([0-9a-f]{3,4}|[0-9a-f]{6}|[0-9a-f]{8})$/; + +const expandShorthandHex = (hex: string) => + hex.length <= 4 + ? hex + .split('') + .map((digit) => digit.repeat(2)) + .join('') + : hex; + +/** Parses `transparent`, the CSS colour keywords, and hex. */ +export function tryCssColorToRGBA(css: string): RGBAColor | undefined { + const value = css.trim().toLowerCase(); + if (value === 'transparent') return [0, 0, 0, 0]; + + const named = NAMED_COLORS.get(value); + if (named) return [...named, 255] as RGBAColor; + + const hex = HEX_COLOR.exec(value)?.[1]; + return hex ? hexaToRGBA(expandShorthandHex(hex)) : undefined; +} + +/** Falls back to opaque black for unparseable label colors. */ +export function cssColorToRGBA(css: string): RGBAColor { + return tryCssColorToRGBA(css) ?? [0, 0, 0, 255]; +} + +// Opaque colors keep the 6-digit form label colors are written in, so a color +// that round trips through a segment comes back byte-identical. +export function rgbaToCssColor(rgba: RGBAColor) { + const hexa = rgbaToHexa(rgba); + return rgba[3] === 255 ? hexa.slice(0, 7) : hexa; +} diff --git a/src/segmentation/editing/algorithms/__tests__/gaussianSmooth.spec.ts b/src/segmentation/editing/algorithms/__tests__/gaussianSmooth.spec.ts new file mode 100644 index 000000000..4a809e805 --- /dev/null +++ b/src/segmentation/editing/algorithms/__tests__/gaussianSmooth.spec.ts @@ -0,0 +1,122 @@ +import { describe, it, expect } from 'vitest'; +import { fullExtent, type Extent3D } from '@/src/segmentation/geometry'; +import { gaussianSmoothLabelMapWorker } from '@/src/segmentation/editing/algorithms/gaussianSmooth.worker'; + +const LABEL = 3; +type Dims = [number, number, number]; + +/** Reconstructs the whole parent, including growth outside the input mask. */ +function smooth( + data: Uint8Array, + dimensions: Dims, + maskExtent: Extent3D = fullExtent(dimensions), + parentDimensions: Dims = dimensions +) { + const { scalars, extent } = gaussianSmoothLabelMapWorker({ + data, + dimensions, + spacing: [1, 1, 1], + maskExtent, + parentDimensions, + params: { sigma: 1, label: LABEL }, + }); + const [pi, pj, pk] = parentDimensions; + const parent = new Uint8Array(pi * pj * pk); + let offset = 0; + for (let k = extent[4]; k <= extent[5]; k += 1) { + for (let j = extent[2]; j <= extent[3]; j += 1) { + for (let i = extent[0]; i <= extent[1]; i += 1) { + parent[i + pi * (j + pj * k)] = scalars[offset++]; + } + } + } + return parent; +} + +function boxInParent(extent: Extent3D, dimensions: Dims) { + const [pi, pj, pk] = dimensions; + const parent = new Uint8Array(pi * pj * pk); + for (let k = extent[4]; k <= extent[5]; k += 1) { + for (let j = extent[2]; j <= extent[3]; j += 1) { + for (let i = extent[0]; i <= extent[1]; i += 1) { + parent[i + pi * (j + pj * k)] = LABEL; + } + } + } + return parent; +} + +describe('gaussianSmoothLabelMapWorker', () => { + it('smooths an isolated voxel away even when the mask is that one voxel', () => { + const smoothed = smooth( + new Uint8Array([LABEL]), + [1, 1, 1], + [1, 1, 1, 1, 1, 1], + [3, 3, 3] + ); + expect(smoothed.every((value) => value === 0)).toBe(true); + }); + + it.each([ + [2, 4, 2, 4, 2, 4], + [1, 3, 1, 3, 1, 3], + [0, 2, 2, 4, 2, 4], + [1, 4, 1, 4, 1, 4], + [4, 7, 4, 7, 4, 7], + ] as Extent3D[])( + 'preserves the full-parent output for a tight mask at [%i, %i, %i, %i, %i, %i]', + (...extent) => { + const parentDimensions: Dims = [9, 9, 9]; + const dimensions: Dims = [ + extent[1] - extent[0] + 1, + extent[3] - extent[2] + 1, + extent[5] - extent[4] + 1, + ]; + const data = new Uint8Array( + dimensions[0] * dimensions[1] * dimensions[2] + ).fill(LABEL); + const croppedResult = smooth(data, dimensions, extent, parentDimensions); + const parentResult = smooth( + boxInParent(extent, parentDimensions), + parentDimensions + ); + expect(croppedResult).toEqual(parentResult); + } + ); + + it('retains all 62 voxels when mirroring grows a box toward parent faces', () => { + const smoothed = smooth( + new Uint8Array(64).fill(LABEL), + [4, 4, 4], + [1, 4, 1, 4, 1, 4], + [9, 9, 9] + ); + expect(smoothed.filter((value) => value === LABEL)).toHaveLength(62); + expect(smoothed[0 + 2 * 9 + 2 * 81]).toBe(LABEL); + }); + + it('erodes a cropped cube at its corners and keeps its centre', () => { + const smoothed = smooth( + new Uint8Array(27).fill(LABEL), + [3, 3, 3], + [1, 3, 1, 3, 1, 3], + [5, 5, 5] + ); + expect(smoothed[1 + 1 * 5 + 1 * 25]).toBe(0); + expect(smoothed[2 + 2 * 5 + 2 * 25]).toBe(LABEL); + }); + + it('mirrors only at parent faces, not at the mask allocation', () => { + const data = new Uint8Array(27).fill(LABEL); + const againstFace = smooth(data, [3, 3, 3], [0, 2, 2, 4, 2, 4], [9, 9, 9]); + const awayFromFace = smooth(data, [3, 3, 3], [1, 3, 2, 4, 2, 4], [9, 9, 9]); + expect(againstFace[0 + 3 * 9 + 3 * 81]).toBe(LABEL); + expect(awayFromFace[0 + 3 * 9 + 3 * 81]).toBe(0); + }); + + it('leaves a buffer with none of the label alone', () => { + expect(Array.from(smooth(new Uint8Array([0, 1, 0, 1]), [4, 1, 1]))).toEqual( + [0, 1, 0, 1] + ); + }); +}); diff --git a/src/segmentation/editing/algorithms/__tests__/gaussianSmoothGolden.spec.ts b/src/segmentation/editing/algorithms/__tests__/gaussianSmoothGolden.spec.ts new file mode 100644 index 000000000..be65c6a3d --- /dev/null +++ b/src/segmentation/editing/algorithms/__tests__/gaussianSmoothGolden.spec.ts @@ -0,0 +1,237 @@ +import { describe, it, expect } from 'vitest'; +import { gaussianSmoothLabelMapWorker } from '@/src/segmentation/editing/algorithms/gaussianSmooth.worker'; + +// A byte-exact record of the filter's output. The rest of the smoothing suite +// asserts on shape properties, which a change in the arithmetic can satisfy +// while every voxel moves; this catches that. +// +// Each entry is one row of the volume, '1' where the label survives, so a +// failure names the rows that moved. To refresh one after a deliberate change +// in behaviour, print `asRows(smooth(...))`. + +const LABEL = 3; +const DIMENSIONS: [number, number, number] = [9, 8, 7]; + +const SIGMA_0_6 = [ + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '011100100', + '011101100', + '011011100', + '000111100', + '001111100', + '000000000', + '000000000', + '000000000', + '011101100', + '011111100', + '011111100', + '001111100', + '011111100', + '000000000', + '000000000', + '000000000', + '011011100', + '011111100', + '001111100', + '011111100', + '011111000', + '000000000', + '000000000', + '000000000', + '000111100', + '001111100', + '011111100', + '011111000', + '011110000', + '000000000', + '000000000', + '000000000', + '001111100', + '011111100', + '011111000', + '011110000', + '011100100', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', +]; + +const SIGMA_1_0 = [ + '000000000', + '000000000', + '000000000', + '000001000', + '000011000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '001111000', + '001111100', + '000111000', + '000010000', + '000000000', + '000000000', + '000000000', + '001111000', + '011111100', + '011111100', + '001111100', + '001111000', + '000000000', + '000000000', + '000001000', + '001111100', + '011111100', + '011111100', + '011111100', + '001111000', + '000000000', + '000000000', + '000011000', + '000111000', + '001111100', + '011111100', + '011111000', + '001110000', + '000000000', + '000000000', + '000000000', + '000111000', + '001111000', + '011111000', + '011110000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000110000', + '001110000', + '001100000', + '000000000', + '000000000', + '000000000', +]; + +const ANISOTROPIC = [ + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '001100000', + '011110000', + '111111100', + '111111100', + '001111100', + '000111000', + '000000000', + '000000000', + '001000000', + '011111100', + '111111100', + '011111100', + '001111100', + '001111000', + '000000000', + '000000000', + '000001000', + '000111100', + '001111100', + '011111100', + '011111100', + '011111000', + '000000000', + '000000000', + '000011000', + '000111100', + '001111100', + '011111100', + '111111000', + '011110000', + '000000000', + '000000000', + '000111000', + '001111100', + '011111100', + '111111000', + '111111000', + '011100000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', + '000000000', +]; +/** An asymmetric blob with a pitted interior and two lone corner voxels. */ +const blob = () => { + const [di, dj, dk] = DIMENSIONS; + const data = new Uint8Array(di * dj * dk); + const mark = (i: number, j: number, k: number) => { + data[i + j * di + k * di * dj] = LABEL; + }; + for (let k = 1; k <= 5; k += 1) { + for (let j = 1; j <= 5; j += 1) { + for (let i = 1; i <= 6; i += 1) { + if ((i + j + k) % 7 !== 0) mark(i, j, k); + } + } + } + mark(8, 0, 0); + mark(0, 7, 6); + return data; +}; + +const smooth = (sigma: number, spacing: [number, number, number]) => + gaussianSmoothLabelMapWorker({ + data: blob(), + dimensions: DIMENSIONS, + spacing, + maskExtent: [0, 8, 0, 7, 0, 6], + parentDimensions: DIMENSIONS, + params: { sigma, label: LABEL }, + }).scalars; + +const asRows = (output: ArrayLike) => { + const bits = Array.from(output, (value) => (value ? '1' : '0')).join(''); + return bits.match(new RegExp(`.{${DIMENSIONS[0]}}`, 'g')); +}; + +describe('gaussian smooth golden output', () => { + it('matches byte for byte at a sigma below one voxel', () => { + expect(asRows(smooth(0.6, [1, 1, 1]))).toEqual(SIGMA_0_6); + }); + + it('matches byte for byte at a one-voxel sigma', () => { + expect(asRows(smooth(1.0, [1, 1, 1]))).toEqual(SIGMA_1_0); + }); + + it('matches byte for byte through anisotropic spacing', () => { + expect(asRows(smooth(1.0, [1, 1, 3]))).toEqual(ANISOTROPIC); + }); +}); diff --git a/src/segmentation/editing/algorithms/fillHoles.ts b/src/segmentation/editing/algorithms/fillHoles.ts new file mode 100644 index 000000000..f8f24deba --- /dev/null +++ b/src/segmentation/editing/algorithms/fillHoles.ts @@ -0,0 +1,165 @@ +import { TypedArray } from '@kitware/vtk.js/types'; + +// 4-connected neighbor offsets, shared so the flood-fill loops never allocate +// a neighbor array per visited voxel. +const NEIGHBOR_DU = [-1, 1, 0, 0]; +const NEIGHBOR_DV = [0, 0, -1, 1]; + +type MaskData = TypedArray | number[]; + +export type FillHolesOptions = { + // Flat label-map scalar array, indexed as i + j*dimI + k*dimI*dimJ. + data: MaskData; + // Label-map IJK dimensions [dimI, dimJ, dimK]. + dimensions: [number, number, number]; + // IJK axis perpendicular to the fill plane (the slice axis). + axis: 0 | 1 | 2; + // When set, only this slice index along `axis` is processed. + // When omitted, every slice along `axis` is processed. + sliceIndex?: number; + // The one label treated as foreground, and the one enclosed background is + // filled with. Storage is a mask per segment, so a fill is always one + // segment's own; running over several segments is one call each. + label: number; +}; + +type Plane = { + uDim: number; + vDim: number; + offset: (u: number, v: number) => number; +}; + +// `visited` and `stack` are reused across slices, so the whole run allocates +// them once. +// `visited`: 0 = unvisited, 1 = outside (border-connected), 2 = hole. +type Flood = { + plane: Plane; + out: MaskData; + label: number; + visited: Uint8Array; + stack: number[]; +}; + +const inPlane = (plane: Plane, u: number, v: number) => + u >= 0 && u < plane.uDim && v >= 0 && v < plane.vDim; + +// Drain `stack`, expanding the region into unvisited non-foreground neighbors +// (each marked with `mark`). `collect`, when given, receives the flat offset of +// every region cell. +function drain(flood: Flood, mark: number, collect?: number[]) { + const { plane, out, label, visited, stack } = flood; + while (stack.length) { + const p = stack.pop()!; + const u = p % plane.uDim; + const v = (p - u) / plane.uDim; + if (collect) collect.push(plane.offset(u, v)); + for (let n = 0; n < 4; n++) { + const nu = u + NEIGHBOR_DU[n]; + const nv = v + NEIGHBOR_DV[n]; + if (!inPlane(plane, nu, nv)) continue; + const np = nu + nv * plane.uDim; + if (out[plane.offset(nu, nv)] !== label && visited[np] === 0) { + visited[np] = mark; + stack.push(np); + } + } + } +} + +// Flood the non-foreground cells reachable from the slice border, marking them +// "outside". Whatever it does not reach is enclosed. +function markOutside(flood: Flood) { + const { plane, out, label, visited, stack } = flood; + const seed = (u: number, v: number) => { + const p = u + v * plane.uDim; + if (visited[p] === 0 && out[plane.offset(u, v)] !== label) { + visited[p] = 1; + stack.push(p); + } + }; + for (let u = 0; u < plane.uDim; u++) { + seed(u, 0); + seed(u, plane.vDim - 1); + } + for (let v = 0; v < plane.vDim; v++) { + seed(0, v); + seed(plane.uDim - 1, v); + } + drain(flood, 1); +} + +// Only fill background. A voxel another segment holds is not this segment's to +// take here; the write path decides that on confirm. +function fillBackground(out: MaskData, cells: number[], label: number) { + for (let c = 0; c < cells.length; c++) { + if (out[cells[c]] === 0) { + out[cells[c]] = label; + } + } +} + +// Any non-foreground cell not marked "outside" is part of a hole. Group each +// hole into a connected component and fill it. +function fillEnclosed(flood: Flood) { + const { plane, out, label, visited, stack } = flood; + for (let v = 0; v < plane.vDim; v++) { + for (let u = 0; u < plane.uDim; u++) { + const p = u + v * plane.uDim; + if (visited[p] !== 0 || out[plane.offset(u, v)] === label) continue; + + const holeCells: number[] = []; + visited[p] = 2; + stack.push(p); + drain(flood, 2, holeCells); + fillBackground(out, holeCells, label); + } + } +} + +// Fills enclosed background regions ("holes") on 2D slices of a label map. +// A hole is background that does not connect to the slice border. Only +// background (0) voxels are filled. Returns a copy of `data`; the input is left +// untouched. +export function fillHoles(opts: FillHolesOptions) { + const { data, dimensions, axis, sliceIndex, label } = opts; + const out = data.slice(); + + const strides = [1, dimensions[0], dimensions[0] * dimensions[1]]; + const sliceStride = strides[axis]; + const sliceCount = dimensions[axis]; + + // The two in-plane axes (everything that isn't the slice axis). + const [uAxis, vAxis] = [0, 1, 2].filter((a) => a !== axis); + const uDim = dimensions[uAxis]; + const vDim = dimensions[vAxis]; + const uStride = strides[uAxis]; + const vStride = strides[vAxis]; + + const visited = new Uint8Array(uDim * vDim); + const stack: number[] = []; + + const firstSlice = sliceIndex ?? 0; + const lastSlice = sliceIndex ?? sliceCount - 1; + + for (let slice = firstSlice; slice <= lastSlice; slice++) { + const base = slice * sliceStride; + visited.fill(0); + + const flood: Flood = { + plane: { + uDim, + vDim, + offset: (u, v) => base + u * uStride + v * vStride, + }, + out, + label, + visited, + stack, + }; + + markOutside(flood); + fillEnclosed(flood); + } + + return out; +} diff --git a/src/segmentation/editing/algorithms/fillHoles.worker.ts b/src/segmentation/editing/algorithms/fillHoles.worker.ts new file mode 100644 index 000000000..bf548bb4b --- /dev/null +++ b/src/segmentation/editing/algorithms/fillHoles.worker.ts @@ -0,0 +1,17 @@ +import * as Comlink from 'comlink'; +import { + fillHoles, + FillHolesOptions, +} from '@/src/segmentation/editing/algorithms/fillHoles'; + +// Runs the pure flood-fill off the main thread so whole-volume fills do not +// freeze the UI, mirroring gaussianSmooth.worker.ts. +export function fillHolesWorker(input: FillHolesOptions) { + return fillHoles(input); +} + +const workerApi = { + fillHolesWorker, +}; + +Comlink.expose(workerApi); diff --git a/src/segmentation/editing/algorithms/gaussianSmooth.worker.ts b/src/segmentation/editing/algorithms/gaussianSmooth.worker.ts new file mode 100644 index 000000000..54a130925 --- /dev/null +++ b/src/segmentation/editing/algorithms/gaussianSmooth.worker.ts @@ -0,0 +1,384 @@ +import * as Comlink from 'comlink'; +import { TypedArray } from '@kitware/vtk.js/types'; +import { createTypedArrayLike } from '@/src/utils'; +import { + extentSize, + extentUnion, + type Extent3D, +} from '@/src/segmentation/geometry'; + +export interface GaussianSmoothParams { + sigma: number; + label: number; +} + +export interface GaussianSmoothInput { + data: TypedArray | number[]; + dimensions: number[]; + spacing: [number, number, number]; + maskExtent: [number, number, number, number, number, number]; + parentDimensions: [number, number, number]; + params: GaussianSmoothParams; +} + +function generateGaussianKernel(sigma: number, radiusFactor = 1.5) { + const radius = Math.ceil(sigma * radiusFactor); + const size = 2 * radius + 1; + const kernel = new Float32Array(size); + const center = radius; + let sum = 0; + + for (let i = 0; i < size; i++) { + const x = i - center; + // VTK formula: exp(-(x * x) / (std * std * 2.0)) + const value = Math.exp(-(x * x) / (sigma * sigma * 2.0)); + kernel[i] = value; + sum += value; + } + + // Normalize kernel + for (let i = 0; i < size; i++) { + kernel[i] /= sum; + } + + return kernel; +} + +// Helper for robust boundary handling (mirroring) +function mirrorCoord(sampleCoord: number, axisDim: number) { + let finalCoord = sampleCoord; + if (sampleCoord < 0) { + finalCoord = -sampleCoord; // Reflect + } else if (sampleCoord >= axisDim) { + finalCoord = 2 * axisDim - sampleCoord - 2; // Reflect + } + // Clamp to ensure it's within bounds, useful if kernel is very large + return Math.max(0, Math.min(axisDim - 1, finalCoord)); +} + +/** + * What stays fixed for one axis pass: the two buffers, the kernel, and how a + * line along the convolved axis is addressed. Built once per pass, so walking + * a line allocates nothing. + */ +interface AxisPass { + inputData: TypedArray | number[]; + outputData: TypedArray | number[]; + kernel: Float32Array; + kernelCenter: number; + // Voxels along the convolved axis, and the index step between them. + count: number; + stride: number; +} + +// Convolves the one line of voxels that starts at `lineStart` and runs along +// the pass's axis. Every axis reduces to this, because a line differs only in +// where it starts, how long it is, and how far apart its voxels sit. +function convolveLine(pass: AxisPass, lineStart: number) { + const { inputData, outputData, kernel, kernelCenter, count, stride } = pass; + const kernelSize = kernel.length; + + for (let i = 0; i < count; i++) { + let sum = 0; + for (let k = 0; k < kernelSize; k++) { + const sample = mirrorCoord(i + k - kernelCenter, count); + sum += inputData[sample * stride + lineStart] * kernel[k]; + } + + outputData[i * stride + lineStart] = sum; + } +} + +function convolve1D( + inputData: TypedArray | number[], + outputData: TypedArray | number[], + kernel: Float32Array, + volume: { dimensions: number[]; axis: 0 | 1 | 2 } +) { + const { dimensions, axis } = volume; + const [dimX, dimY] = dimensions; + const strides = [1, dimX, dimX * dimY]; + const pass: AxisPass = { + inputData, + outputData, + kernel, + kernelCenter: Math.floor(kernel.length / 2), + count: dimensions[axis], + stride: strides[axis], + }; + + // The two axes the pass does not walk, the widest-striding one outermost: + // z, then y, then x. That is the loop order each axis wants for cache + // efficiency, so convolving along X visits z, y, x, along Y visits z, x, y, + // and along Z visits y, x, z. + const outer = axis === 2 ? 1 : 2; + const inner = axis === 0 ? 1 : 0; + + for (let o = 0; o < dimensions[outer]; o++) { + const outerOffset = o * strides[outer]; + for (let i = 0; i < dimensions[inner]; i++) { + convolveLine(pass, outerOffset + i * strides[inner]); + } + } +} + +function gaussianFilter3D( + inputData: TypedArray | number[], + dimensions: number[], + sigmaPixels: [number, number, number], + radiusFactor = 1.5 +) { + const totalSize = dimensions[0] * dimensions[1] * dimensions[2]; + const kernelX = generateGaussianKernel(sigmaPixels[0], radiusFactor); + const kernelY = generateGaussianKernel(sigmaPixels[1], radiusFactor); + const kernelZ = generateGaussianKernel(sigmaPixels[2], radiusFactor); + const temp = new Float32Array(totalSize); + const output = new Float32Array(totalSize); + + convolve1D(inputData, output, kernelX, { dimensions, axis: 0 }); + convolve1D(output, temp, kernelY, { dimensions, axis: 1 }); + convolve1D(temp, output, kernelZ, { dimensions, axis: 2 }); + + return output; +} + +// What a bounding-box scan holds still: the voxels being read, the bounds +// being widened, and the row addressing. Built once, so scanning a row +// allocates nothing. +interface RowScan { + data: TypedArray | number[]; + bounds: number[]; + dimX: number; + sliceSize: number; + label: number; +} + +// Its own function so that the per-voxel test sits two blocks deep rather than +// four. +function growBoundsOverRow(scan: RowScan, y: number, z: number) { + const { data, bounds, dimX, sliceSize, label } = scan; + const rowStart = y * dimX + z * sliceSize; + + for (let x = 0; x < dimX; x++) { + if (data[rowStart + x] !== label) continue; + bounds[0] = Math.min(bounds[0], x); + bounds[1] = Math.max(bounds[1], x); + bounds[2] = Math.min(bounds[2], y); + bounds[3] = Math.max(bounds[3], y); + bounds[4] = Math.min(bounds[4], z); + bounds[5] = Math.max(bounds[5], z); + } +} + +function calculateBoundingBox( + data: TypedArray | number[], + dimensions: number[], + label: number +) { + const [dimX, dimY, dimZ] = dimensions; + const bounds = [dimX, -1, dimY, -1, dimZ, -1]; + const scan: RowScan = { + data, + bounds, + dimX, + sliceSize: dimX * dimY, + label, + }; + + for (let z = 0; z < dimZ; z++) { + for (let y = 0; y < dimY; y++) { + growBoundsOverRow(scan, y, z); + } + } + + if (bounds[1] === -1) return null; + + return bounds; +} + +function expandBoundingBox({ + bounds, + maskExtent, + parentDimensions, + sigmaPixels, + radiusFactor = 1.5, +}: { + bounds: number[]; + maskExtent: GaussianSmoothInput['maskExtent']; + parentDimensions: GaussianSmoothInput['parentDimensions']; + sigmaPixels: [number, number, number]; + radiusFactor?: number; +}) { + return sigmaPixels.flatMap((sigma, axis) => { + const padding = Math.ceil(sigma * radiusFactor); + // The parent-image faces, stated in mask coordinates. Ending the + // convolution volume there keeps the established mirrored boundary + // wherever the mask sits, so the result does not depend on how much of + // the parent the mask happens to be allocated over. Crop faces are not + // clamped: outside the buffer reads as background. + const parentLow = -maskExtent[axis * 2]; + const parentHigh = parentDimensions[axis] - 1 - maskExtent[axis * 2]; + return [ + Math.max(parentLow, bounds[axis * 2] - padding), + Math.min(parentHigh, bounds[axis * 2 + 1] + padding), + ]; + }); +} + +/** + * Visits every voxel of `bounds` that the volume actually holds, giving each + * its offset in the volume and its offset in the padded sub-volume. The + * padding ring outside the volume is skipped by clipping the loops, not tested + * per voxel. + */ +function forEachClippedVoxel( + dimensions: number[], + bounds: number[], + visit: (origIndex: number, subIndex: number) => void +) { + const [dimX, dimY, dimZ] = dimensions; + const [minX, maxX, minY, maxY, minZ, maxZ] = bounds; + const subDimX = maxX - minX + 1; + const subDimY = maxY - minY + 1; + const lastX = Math.min(maxX, dimX - 1); + + for (let z = Math.max(minZ, 0); z <= Math.min(maxZ, dimZ - 1); z += 1) { + for (let y = Math.max(minY, 0); y <= Math.min(maxY, dimY - 1); y += 1) { + const rowOrig = y * dimX + z * dimX * dimY; + const rowSub = (y - minY) * subDimX + (z - minZ) * subDimX * subDimY; + for (let x = Math.max(minX, 0); x <= lastX; x += 1) { + visit(x + rowOrig, x - minX + rowSub); + } + } + } +} + +/** + * The label's own binary mask over `bounds`, which is the only thing the + * filter reads. Built in one pass rather than copying the labels out and + * thresholding them afterwards: the copy is a second volume-sized Float32 + * array, live at the same time as this one. + */ +function extractSubMask( + data: TypedArray | number[], + dimensions: number[], + bounds: number[], + label: number +) { + const [minX, maxX, minY, maxY, minZ, maxZ] = bounds; + const subDims = [maxX - minX + 1, maxY - minY + 1, maxZ - minZ + 1]; + // Zero filled, so everything outside the buffer stays background. + const subMask = new Float32Array(subDims[0] * subDims[1] * subDims[2]); + + forEachClippedVoxel(dimensions, bounds, (origIndex, subIndex) => { + subMask[subIndex] = data[origIndex] === label ? 255.0 : 0.0; + }); + + return { subMask, subDims }; +} + +// Output storage includes the padding ring: mirroring at a parent face can +// turn on voxels beyond the input mask's allocation. +function copySubVolumeBack( + subData: Float32Array, + originalData: TypedArray | number[], + region: { dimensions: number[]; bounds: number[] }, + label: number +) { + forEachClippedVoxel( + region.dimensions, + region.bounds, + (origIndex, subIndex) => { + const origLabel = originalData[origIndex]; + if (origLabel === label || origLabel === 0) { + originalData[origIndex] = subData[subIndex] > 127.5 ? label : 0; + } + } + ); +} + +export function gaussianSmoothLabelMapWorker(input: GaussianSmoothInput) { + const { + data: originalData, + dimensions, + spacing, + maskExtent, + parentDimensions, + params, + } = input; + const { sigma, label } = params; + + if (sigma <= 0) { + throw new Error('Sigma must be positive'); + } + + const sigmaPixels: [number, number, number] = [ + sigma / spacing[0], + sigma / spacing[1], + sigma / spacing[2], + ]; + + // Absent when the label is nowhere in the mask, which is also the whole + // answer for a mask with nothing to smooth: it comes back as it went in. + const bounds = calculateBoundingBox(originalData, dimensions, label); + if (!bounds) { + const outputData = createTypedArrayLike(originalData, originalData.length); + for (let i = 0; i < originalData.length; i++) { + outputData[i] = originalData[i]; + } + return { scalars: outputData, extent: maskExtent }; + } + + const expandedBounds = expandBoundingBox({ + bounds, + maskExtent, + parentDimensions, + sigmaPixels, + }); + const { subMask, subDims } = extractSubMask( + originalData, + dimensions, + expandedBounds, + label + ); + + const smoothedSubMask = gaussianFilter3D(subMask, subDims, sigmaPixels, 1.5); + + const expandedExtent = expandedBounds.map( + (value, axis) => value + maskExtent[axis - (axis % 2)] + ) as Extent3D; + const extent = extentUnion(maskExtent, expandedExtent); + const outputDimensions = extentSize(extent); + const outputData = createTypedArrayLike( + originalData, + outputDimensions[0] * outputDimensions[1] * outputDimensions[2] + ); + const outputBoundsInInput = extent.map( + (value, axis) => value - maskExtent[axis - (axis % 2)] + ); + forEachClippedVoxel( + dimensions, + outputBoundsInInput, + (origIndex, outIndex) => { + outputData[outIndex] = originalData[origIndex]; + } + ); + const smoothedBoundsInOutput = expandedExtent.map( + (value, axis) => value - extent[axis - (axis % 2)] + ); + + copySubVolumeBack( + smoothedSubMask, + outputData, + { dimensions: outputDimensions, bounds: smoothedBoundsInOutput }, + label + ); + + return { scalars: outputData, extent }; +} + +const workerApi = { + gaussianSmoothLabelMapWorker, +}; + +Comlink.expose(workerApi); diff --git a/src/segmentation/editing/coordinator.ts b/src/segmentation/editing/coordinator.ts new file mode 100644 index 000000000..0bb466539 --- /dev/null +++ b/src/segmentation/editing/coordinator.ts @@ -0,0 +1,29 @@ +import { defineStore } from 'pinia'; + +/** Coordinates temporary mask previews with competing edits and durable reads. */ +export const useSegmentationEditsStore = defineStore( + 'segmentationEdits', + () => { + let cancelPreview: (() => void) | undefined; + + function beforeEdit() { + const cancel = cancelPreview; + cancelPreview = undefined; + cancel?.(); + } + + function hold(cancel: () => void) { + beforeEdit(); + cancelPreview = cancel; + } + + function release(cancel: () => void) { + if (cancelPreview === cancel) cancelPreview = undefined; + } + + // Save and export read committed voxels, resolving an unconfirmed preview. + const beforeRead = beforeEdit; + + return { beforeEdit, beforeRead, hold, release }; + } +); diff --git a/src/segmentation/geometry.ts b/src/segmentation/geometry.ts new file mode 100644 index 000000000..7a2f4fe65 --- /dev/null +++ b/src/segmentation/geometry.ts @@ -0,0 +1,174 @@ +/** vtk.js index-space extent order: [iMin, iMax, jMin, jMax, kMin, kMax]. */ +export type Extent3D = [number, number, number, number, number, number]; + +/** Widens `box` in place to take in one more index. */ +export const growExtent = (box: Extent3D, i: number, j: number, k: number) => { + box[0] = Math.min(box[0], i); + box[1] = Math.max(box[1], i); + box[2] = Math.min(box[2], j); + box[3] = Math.max(box[3], j); + box[4] = Math.min(box[4], k); + box[5] = Math.max(box[5], k); +}; + +export function emptyExtent(): Extent3D { + return [0, -1, 0, -1, 0, -1]; +} + +/** vtk.js extents are inclusive, so an axis is empty only when max < min. */ +export function isEmptyExtent(extent: Extent3D) { + return ( + extent[1] < extent[0] || extent[3] < extent[2] || extent[5] < extent[4] + ); +} + +/** Voxel counts along i, j, k. Meaningless for an empty extent. */ +export function extentSize(extent: Extent3D) { + return [ + extent[1] - extent[0] + 1, + extent[3] - extent[2] + 1, + extent[5] - extent[4] + 1, + ] as [number, number, number]; +} + +export function extentContains(outer: Extent3D, inner: Extent3D) { + return ( + inner[0] >= outer[0] && + inner[1] <= outer[1] && + inner[2] >= outer[2] && + inner[3] <= outer[3] && + inner[4] >= outer[4] && + inner[5] <= outer[5] + ); +} + +export function extentContainsIndex( + extent: Extent3D, + i: number, + j: number, + k: number +) { + return ( + i >= extent[0] && + i <= extent[1] && + j >= extent[2] && + j <= extent[3] && + k >= extent[4] && + k <= extent[5] + ); +} + +/** A mask's extent with the row and plane strides that extent implies. */ +export type MaskBounds = { + extent: Extent3D; + mi: number; + mj: number; +}; + +/** Where a parent-index voxel sits in the buffer of a mask bounded that way. */ +export const maskOffset = ( + bounds: MaskBounds, + i: number, + j: number, + k: number +) => + i - + bounds.extent[0] + + (j - bounds.extent[2]) * bounds.mi + + (k - bounds.extent[4]) * bounds.mi * bounds.mj; + +/** + * The box `labelValue` actually occupies inside a mask bounded by `extent`, + * empty when it occupies nothing. A binding's extent is the allocation, padded + * and never shrunk by an erase, so it is not the segment's bounds. + */ +export function markedExtent( + scalars: ArrayLike, + extent: Extent3D, + labelValue: number +): Extent3D { + const ni = extent[1] - extent[0] + 1; + const nj = extent[3] - extent[2] + 1; + const nk = extent[5] - extent[4] + 1; + let bounds: Extent3D | undefined; + + const scanRow = (rowStart: number, j: number, k: number) => { + for (let index = 0; index < ni; index += 1) { + if (scalars[rowStart + index] !== labelValue) continue; + const i = extent[0] + index; + if (bounds) growExtent(bounds, i, j, k); + else bounds = [i, i, j, j, k, k]; + } + }; + + for (let row = 0; row < nj * nk; row += 1) { + scanRow(row * ni, extent[2] + (row % nj), extent[4] + Math.floor(row / nj)); + } + + return bounds ?? emptyExtent(); +} + +/** The parent-image slice indices containing `labelValue`, for i, j and k. */ +export function markedSlices( + scalars: ArrayLike, + extent: Extent3D, + labelValue: number +): [number[], number[], number[]] { + const ni = extent[1] - extent[0] + 1; + const nj = extent[3] - extent[2] + 1; + const occupied = [new Set(), new Set(), new Set()]; + + for (let offset = 0; offset < scalars.length; offset += 1) { + if (scalars[offset] !== labelValue) continue; + const i = extent[0] + (offset % ni); + const row = Math.floor(offset / ni); + const j = extent[2] + (row % nj); + const k = extent[4] + Math.floor(row / nj); + occupied[0].add(i); + occupied[1].add(j); + occupied[2].add(k); + } + + return occupied.map((slices) => [...slices]) as [ + number[], + number[], + number[], + ]; +} + +export function extentUnion(a: Extent3D, b: Extent3D): Extent3D { + return [ + Math.min(a[0], b[0]), + Math.max(a[1], b[1]), + Math.min(a[2], b[2]), + Math.max(a[3], b[3]), + Math.min(a[4], b[4]), + Math.max(a[5], b[5]), + ]; +} + +export function padExtent(extent: Extent3D, padding: number): Extent3D { + return [ + extent[0] - padding, + extent[1] + padding, + extent[2] - padding, + extent[3] + padding, + extent[4] - padding, + extent[5] + padding, + ]; +} + +export function clipExtent(extent: Extent3D, bounds: Extent3D): Extent3D { + return [ + Math.max(extent[0], bounds[0]), + Math.min(extent[1], bounds[1]), + Math.max(extent[2], bounds[2]), + Math.min(extent[3], bounds[3]), + Math.max(extent[4], bounds[4]), + Math.min(extent[5], bounds[5]), + ]; +} + +export function fullExtent(dimensions: number[] | Int32Array): Extent3D { + return [0, dimensions[0] - 1, 0, dimensions[1] - 1, 0, dimensions[2] - 1]; +} diff --git a/src/segmentation/io/__tests__/export.spec.ts b/src/segmentation/io/__tests__/export.spec.ts new file mode 100644 index 000000000..fbbb7a658 --- /dev/null +++ b/src/segmentation/io/__tests__/export.spec.ts @@ -0,0 +1,58 @@ +import { describe, expect, it } from 'vitest'; +import JSZip from 'jszip'; + +import { bundleExportFiles, layerFileName } from '@/src/segmentation/io/export'; + +// A labelmap file carries one label per voxel, so a segmentation with overlap +// leaves as several files. One file is the common case and stays the download +// it has always been: same name, same bytes, no archive around it. + +const bytes = (...values: number[]) => new Uint8Array(values); + +const readBlob = async (blob: Blob) => + Array.from(new Uint8Array(await blob.arrayBuffer())); + +describe('naming the file each group of segments writes', () => { + it('gives the first group the plain name', () => { + expect(layerFileName('Prostate', 'seg.nrrd', 0)).toBe('Prostate.seg.nrrd'); + }); + + it('numbers every later group after it', () => { + expect(layerFileName('Prostate', 'seg.nrrd', 1)).toBe( + 'Prostate_layer1.seg.nrrd' + ); + expect(layerFileName('Prostate', 'nii.gz', 2)).toBe( + 'Prostate_layer2.nii.gz' + ); + }); +}); + +describe('handing the written files to the browser', () => { + it('downloads a single file as itself', async () => { + const bundle = await bundleExportFiles('Prostate', [ + { name: 'Prostate.seg.nrrd', data: bytes(1, 2, 3) }, + ]); + + expect(bundle.name).toBe('Prostate.seg.nrrd'); + expect(await readBlob(bundle.blob)).toEqual([1, 2, 3]); + }); + + it('downloads several files as one archive of them', async () => { + const bundle = await bundleExportFiles('Prostate', [ + { name: 'Prostate.seg.nrrd', data: bytes(1, 2, 3) }, + { name: 'Prostate_layer1.seg.nrrd', data: bytes(4, 5) }, + ]); + + expect(bundle.name).toBe('Prostate.zip'); + const zip = await JSZip.loadAsync(bundle.blob); + expect(Object.keys(zip.files)).toEqual([ + 'Prostate.seg.nrrd', + 'Prostate_layer1.seg.nrrd', + ]); + expect( + Array.from( + await zip.files['Prostate_layer1.seg.nrrd'].async('uint8array') + ) + ).toEqual([4, 5]); + }); +}); diff --git a/src/segmentation/io/export.ts b/src/segmentation/io/export.ts new file mode 100644 index 000000000..08c1a7ed9 --- /dev/null +++ b/src/segmentation/io/export.ts @@ -0,0 +1,34 @@ +import JSZip from 'jszip'; + +export type ExportFile = { + name: string; + data: string | Uint8Array; +}; + +/** + * The name a group's file carries. The first keeps the plain stem, so a + * segmentation with no overlap saves as the one file it always did. + */ +export const layerFileName = (stem: string, format: string, layer: number) => + layer === 0 ? `${stem}.${format}` : `${stem}_layer${layer}.${format}`; + +export const archiveNameFor = (stem: string) => `${stem}.zip`; + +/** + * What a save hands to the browser: the single file itself, or every file in + * one archive. A labelmap file carries one label per voxel, so segments that + * overlap cannot share one and the save turns into several. + */ +export async function bundleExportFiles(stem: string, files: ExportFile[]) { + const [first] = files; + if (files.length === 1) { + return { name: first.name, blob: new Blob([first.data]) }; + } + + const zip = new JSZip(); + files.forEach((file) => zip.file(file.name, file.data)); + return { + name: archiveNameFor(stem), + blob: await zip.generateAsync({ type: 'blob' }), + }; +} diff --git a/src/segmentation/io/maskFileNaming.ts b/src/segmentation/io/maskFileNaming.ts new file mode 100644 index 000000000..873b4b76a --- /dev/null +++ b/src/segmentation/io/maskFileNaming.ts @@ -0,0 +1,27 @@ +const defaultName = (baseName: string, index: number) => + `Segment Group ${index} for ${baseName}`; + +/** + * Default names for mask files, counted per parent image. The + * count keeps rising so a deleted mask's name is not immediately handed to + * the next one, and `taken` skips a name something already holds. + */ +export function createMaskFileNamer(taken: () => Set) { + const nextIndex: Record = Object.create(null); + return { + pick(parentImageId: string, baseName: string) { + const held = taken(); + let name = ''; + do { + const index = nextIndex[parentImageId] ?? 1; + nextIndex[parentImageId] = index + 1; + name = defaultName(baseName, index); + } while (held.has(name)); + return name; + }, + /** Called by the deletion cascade, so a removed image stops counting. */ + forget(parentImageId: string) { + delete nextIndex[parentImageId]; + }, + }; +} diff --git a/src/segmentation/masks/labelValue.ts b/src/segmentation/masks/labelValue.ts new file mode 100644 index 000000000..5a72c3557 --- /dev/null +++ b/src/segmentation/masks/labelValue.ts @@ -0,0 +1,33 @@ +/** A composed labelmap is written as Uint8, so a label value has to fit in one byte. */ +export const LABELMAP_MAX_VALUE = 255; + +/** + * The value every mask marks its own voxels with. A mask holds one segment and + * nothing else, so the byte says only claimed or not; which segment it belongs + * to is the mask's identity, not its content. Export assigns its own values at + * write time, where one file does have to tell segments apart per voxel. + */ +export const SEGMENT_VALUE = 1; + +import { LABELMAP_BACKGROUND_VALUE } from '@/src/segmentation/model'; + +export function nextUnusedLabelValue( + used: Set, + maximum: number, + preferred?: number +) { + if ( + preferred !== undefined && + preferred > LABELMAP_BACKGROUND_VALUE && + preferred <= maximum && + !used.has(preferred) + ) { + return preferred; + } + let labelValue = LABELMAP_BACKGROUND_VALUE + 1; + while (used.has(labelValue)) labelValue += 1; + if (labelValue > maximum) { + throw new Error(`An image holds at most ${maximum} segments in a labelmap`); + } + return labelValue; +} diff --git a/src/segmentation/masks/overlap.ts b/src/segmentation/masks/overlap.ts new file mode 100644 index 000000000..4a87514c5 --- /dev/null +++ b/src/segmentation/masks/overlap.ts @@ -0,0 +1,211 @@ +import type vtkLabelMap from '@/src/vtk/LabelMap'; +import { + LABELMAP_BACKGROUND_VALUE, + maskScalars, +} from '@/src/segmentation/model'; +import { + clipExtent, + extentContainsIndex, + extentSize, + isEmptyExtent, + maskOffset, + type Extent3D, + type MaskBounds, +} from '@/src/segmentation/geometry'; + +// Bounded masks read as one parent-shaped picture: how a mask is written into +// that picture, and which masks can share one without losing a voxel. + +export type BoundedScalars = MaskBounds & { + mask: vtkLabelMap; + scalars: Uint8Array; +}; + +/** + * The extent is copied because the callers read it per voxel and a segment's + * own copy lives in the reactive tree. + */ +export function boundScalars( + mask: vtkLabelMap | undefined, + bounds: Extent3D +): BoundedScalars | undefined { + if (!mask || isEmptyExtent(bounds)) return undefined; + const extent = [...bounds] as Extent3D; + const [mi, mj] = extentSize(extent); + return { mask, scalars: maskScalars(mask), extent, mi, mj }; +} + +/** + * Clipping once here keeps a mask that misses the box out of the per-voxel + * containment test, and lets the caller skip the walk when none is left. + */ +const masksReaching = (masks: BoundedScalars[], within: Extent3D) => + masks.filter((bounded) => !isEmptyExtent(clipExtent(bounded.extent, within))); + +/** + * Whether any of these masks holds the voxel at PARENT indices i, j, k, over + * the box the caller is about to walk. Absent when no mask reaches that box. + * + * The answer is asked once per voxel the caller walks, so the sweep over the + * reaching masks is a plain indexed loop: a callback taking i, j, k would be a + * fresh closure per voxel. + */ +export function masksHolding(masks: BoundedScalars[], within: Extent3D) { + const reaching = masksReaching(masks, within); + if (reaching.length === 0) return undefined; + return (i: number, j: number, k: number) => { + for (let index = 0; index < reaching.length; index += 1) { + const bounded = reaching[index]; + if ( + extentContainsIndex(bounded.extent, i, j, k) && + bounded.scalars[maskOffset(bounded, i, j, k)] !== + LABELMAP_BACKGROUND_VALUE + ) + return true; + } + return false; + }; +} + +/** + * Clears the voxel at PARENT indices i, j, k from every one of these masks and + * answers true: nothing left here can refuse the write, which is the same + * per-voxel answer the occupancy test gives. Absent when no mask reaches + * `within`, the box the caller is about to walk. A mask that does not reach the + * voxel has nothing there to clear, so nothing grows. Finish the operation in + * a finally block to publish each changed mask once, including partial writes. + */ +export function masksClearing(masks: BoundedScalars[], within: Extent3D) { + const reaching = masksReaching(masks, within); + if (reaching.length === 0) return undefined; + const changed = new Set(); + // Indexed loop, as in masksHolding: claim runs once per voxel the caller + // walks, and a callback over the reaching masks would allocate per voxel. + const claim = (i: number, j: number, k: number) => { + for (let index = 0; index < reaching.length; index += 1) { + const bounded = reaching[index]; + if (extentContainsIndex(bounded.extent, i, j, k)) { + const offset = maskOffset(bounded, i, j, k); + if (bounded.scalars[offset] !== LABELMAP_BACKGROUND_VALUE) { + bounded.scalars[offset] = LABELMAP_BACKGROUND_VALUE; + changed.add(bounded.mask); + } + } + } + return true; + }; + return { + claim, + finish: () => { + changed.forEach((mask) => mask.modified()); + changed.clear(); + }, + }; +} + +/** + * A mask's buffer, positioned where one row of the shared box starts in it. + */ +type MaskRow = { scalars: Uint8Array; from: number }; + +const rowAt = ( + bounded: BoundedScalars, + i: number, + j: number, + k: number +): MaskRow => ({ + scalars: bounded.scalars, + from: maskOffset(bounded, i, j, k), +}); + +/** + * Whether both rows hold a voxel at the same step along i. Background is 0, so + * a claimed voxel is a truthy one. + */ +function rowsIntersect(a: MaskRow, b: MaskRow, count: number) { + const { scalars: av, from: ai } = a; + const { scalars: bv, from: bi } = b; + for (let n = 0; n < count; n += 1) { + if (av[ai + n] && bv[bi + n]) return true; + } + return false; +} + +/** + * Whether two masks claim one voxel in common. A mask is bounded to the voxels + * it holds, so boxes that miss cannot share a voxel and neither buffer is read. + * Boxes that meet are swept a row at a time: two segments that touch nowhere + * usually still share a box, so the sweep is the common case, and each row + * costs one offset per mask with a plain step along i from there. + */ +export function masksIntersect(a: BoundedScalars, b: BoundedScalars) { + const shared = clipExtent(a.extent, b.extent); + if (isEmptyExtent(shared)) return false; + + const [ni] = extentSize(shared); + const i = shared[0]; + for (let k = shared[4]; k <= shared[5]; k += 1) { + for (let j = shared[2]; j <= shared[3]; j += 1) { + if (rowsIntersect(rowAt(a, i, j, k), rowAt(b, i, j, k), ni)) return true; + } + } + return false; +} + +/** + * Items grouped so no group holds two masks that claim a voxel in common. + * Greedy first fit: an item takes the lowest group it does not intersect, so a + * segmentation with no overlap stays one group, in order. An item with no mask + * claims nothing and joins the first group. + */ +export function groupByLayer( + items: T[], + maskOf: (item: T) => BoundedScalars | undefined +) { + const layers: Array<{ items: T[]; masks: BoundedScalars[] }> = []; + const fits = ( + layer: { masks: BoundedScalars[] }, + mask: BoundedScalars | undefined + ) => !mask || layer.masks.every((other) => !masksIntersect(mask, other)); + + items.forEach((item) => { + const mask = maskOf(item); + const found = layers.find((layer) => fits(layer, mask)); + const layer = found ?? { items: [], masks: [] }; + if (!found) layers.push(layer); + layer.items.push(item); + if (mask) layer.masks.push(mask); + }); + + return layers.map((layer) => layer.items); +} + +/** + * Marks a bounded mask's voxels in a parent-shaped buffer as `labelValue`. The + * mask's own bytes only say claimed or not, so the value is the caller's. + * Masks are written in `order`, so a later one takes a voxel an earlier one + * also claims. + */ +export function writeMaskInto( + values: Uint8Array, + dimensions: readonly number[], + bounded: BoundedScalars, + labelValue: number +) { + const { extent, scalars } = bounded; + const [dx, dy] = dimensions; + const [ni, nj, nk] = extentSize(extent); + // Rows flat in one loop: a j loop inside a k loop would nest deeper than the + // style allows once the row test is in it. + for (let row = 0; row < nj * nk; row += 1) { + const j = extent[2] + (row % nj); + const k = extent[4] + Math.floor(row / nj); + const to = extent[0] + j * dx + k * dx * dy; + const from = maskOffset(bounded, extent[0], j, k); + for (let n = 0; n < ni; n += 1) { + // Background is 0, so a voxel this mask leaves unclaimed keeps whatever + // the buffer already holds there. + if (scalars[from + n]) values[to + n] = labelValue; + } + } +} diff --git a/src/segmentation/masks/storage.ts b/src/segmentation/masks/storage.ts new file mode 100644 index 000000000..1e7846132 --- /dev/null +++ b/src/segmentation/masks/storage.ts @@ -0,0 +1,101 @@ +import vtkDataArray from '@kitware/vtk.js/Common/Core/DataArray'; +import type vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; +import type { TypedArray, Vector3 } from '@kitware/vtk.js/types'; + +import { maskScalars } from '@/src/segmentation/model'; +import { + clipExtent, + extentSize, + isEmptyExtent, + maskOffset, + type Extent3D, +} from '@/src/segmentation/geometry'; +import vtkLabelMap from '@/src/vtk/LabelMap'; + +export const setMaskScalars = (mask: vtkLabelMap, values: Uint8Array) => + mask + .getPointData() + .setScalars(vtkDataArray.newInstance({ numberOfComponents: 1, values })); + +/** Places a bounded mask on the same index grid as its parent image. */ +export function placeMask( + mask: vtkLabelMap, + parent: vtkImageData, + extent: Extent3D +) { + const dimensions = isEmptyExtent(extent) ? [0, 0, 0] : extentSize(extent); + const origin = isEmptyExtent(extent) + ? Array.from(parent.getOrigin()) + : Array.from( + parent.indexToWorld([extent[0], extent[2], extent[4]] as Vector3) + ); + mask.setOrigin(origin as Vector3); + mask.setDimensions(dimensions as Vector3); + mask.computeTransforms(); + return dimensions; +} + +export function allocateMask(parent: vtkImageData, extent: Extent3D) { + const mask = vtkLabelMap.newInstance( + parent.get('spacing', 'origin', 'direction') + ); + const dimensions = placeMask(mask, parent, extent); + setMaskScalars( + mask, + new Uint8Array(dimensions[0] * dimensions[1] * dimensions[2]) + ); + return mask; +} + +/** + * Copies a mask onto another extent of its parent grid, padding with zero. + * An output buffer must match the destination size and not alias the source. + */ +export function reframeMaskScalars( + scalars: TypedArray | number[], + from: Extent3D, + to: Extent3D, + output?: Uint8Array +) { + const [mi, mj, mk] = extentSize(to); + const size = isEmptyExtent(to) ? 0 : mi * mj * mk; + const values = output ?? new Uint8Array(size); + if (values.length !== size) throw new Error('Mask output size mismatch'); + if (output) values.fill(0); + const shared = clipExtent(from, to); + const [si, sj] = extentSize(from); + // Bounds can be reactive; read them once before the per-row copy loop. + const source = { extent: [...from] as Extent3D, mi: si, mj: sj }; + const destination = { extent: [...to] as Extent3D, mi, mj }; + if (isEmptyExtent(shared)) return values; + const count = shared[1] - shared[0] + 1; + for (let k = shared[4]; k <= shared[5]; k += 1) { + for (let j = shared[2]; j <= shared[3]; j += 1) { + const start = maskOffset(source, shared[0], j, k); + const end = maskOffset(destination, shared[0], j, k); + if (count === 1) { + values[end] = scalars[start]; + } else { + const row = Array.isArray(scalars) + ? scalars.slice(start, start + count) + : scalars.subarray(start, start + count); + values.set(row, end); + } + } + } + return values; +} + +/** Preserves the vtk image instance while replacing its scalar storage. */ +export function regrowMask( + mask: vtkLabelMap, + parent: vtkImageData, + from: Extent3D, + to: Extent3D +) { + const previous = maskScalars(mask); + const values = reframeMaskScalars(previous, from, to); + placeMask(mask, parent, to); + setMaskScalars(mask, values); + mask.modified(); +} diff --git a/src/segmentation/masks/voxelAccess.ts b/src/segmentation/masks/voxelAccess.ts new file mode 100644 index 000000000..f29e642cf --- /dev/null +++ b/src/segmentation/masks/voxelAccess.ts @@ -0,0 +1,216 @@ +import type { TypedArray } from '@kitware/vtk.js/types'; + +import type { Maybe } from '@/src/types'; +import type { VoxelGesture } from '@/src/segmentation/model'; +import type { useImageCacheStore } from '@/src/store/image-cache'; +import { regrowMask } from '@/src/segmentation/masks/storage'; +import { + boundScalars, + masksClearing, + masksHolding, +} from '@/src/segmentation/masks/overlap'; +import { + listMasks, + maskScalars, + type LabelmapBinding, + type MaskVoxelAccessor, + type SegmentMask, + type Segmentation, + type VoxelStorage, +} from '@/src/segmentation/model'; +import { + clipExtent, + extentContains, + extentUnion, + fullExtent, + isEmptyExtent, + padExtent, + type Extent3D, +} from '@/src/segmentation/geometry'; + +export type VoxelAccessDeps = { + imageCacheStore: ReturnType; + findMask: (maskId: string) => SegmentMask | undefined; + getMask: (maskId: string) => SegmentMask; + segmentationOfMask: (maskId: string) => Segmentation | undefined; + ensureLabelmapBinding: (maskId: string) => LabelmapBinding; + maskLocked: (mask: SegmentMask) => boolean; +}; + +/** + * Reading and growing the voxels behind a mask. Split out so the store holds + * the records; every accessor re-resolves its binding rather than capturing a + * buffer, so none of them outlive a mask they were made for. + */ +export function createVoxelAccess(deps: VoxelAccessDeps) { + const { + imageCacheStore, + findMask, + getMask, + segmentationOfMask, + ensureLabelmapBinding, + maskLocked, + } = deps; + + function requireParentImage(maskId: string) { + const segmentation = segmentationOfMask(maskId); + if (!segmentation) throw new Error('No such segment'); + const parent = imageCacheStore.getVtkImageData(segmentation.parentImageId); + if (!parent) throw new Error('No such parent image'); + return parent; + } + + /** + * Grows one mask, in place, to cover `extent` in parent index space, with + * `padding` voxels of room beyond it when it has to grow at all. + */ + function ensureMaskContains(maskId: string, extent: Extent3D, padding = 0) { + if (isEmptyExtent(extent)) return false; + + const binding = getMask(maskId).representations.labelmap; + if (!binding) throw new Error('No storage: call materialize() first'); + const parent = requireParentImage(maskId); + // Refused before anything is touched, so a rejected growth leaves the mask + // exactly as it was. + const parentExtent = fullExtent(parent.getDimensions()); + if (!extentContains(parentExtent, extent)) + throw new Error('Extent leaves the parent image'); + + const current = binding.extent; + if (!isEmptyExtent(current) && extentContains(current, extent)) + return false; + + const requested = clipExtent(padExtent(extent, padding), parentExtent); + const grown = isEmptyExtent(current) + ? requested + : extentUnion(current, requested); + regrowMask(binding.image, parent, current, grown); + binding.extent = grown; + return true; + } + + /** + * The voxel half of the accessor seam, over whichever mask `findBinding` + * resolves. Resolution is deferred to every call so a stale accessor sees + * deletion or growth done through another one. `onMissing` names why storage + * is unreachable, so `exists()` can answer without throwing. + */ + function voxelStorage( + maskId: string, + findBinding: () => Maybe, + onMissing: () => never + ): VoxelStorage { + const findImage = () => findBinding()?.image; + const requireImage = () => findImage() ?? onMissing(); + const requireScalars = () => maskScalars(requireImage()); + + return { + exists: () => !!findImage(), + image: requireImage, + scalars: requireScalars, + snapshot: () => requireScalars().slice(), + apply: (scalars: TypedArray | number[]) => { + const image = requireImage(); + const data = maskScalars(image); + if (scalars.length !== data.length) { + throw new Error('Scalar length does not match storage'); + } + data.set(scalars); + image.modified(); + }, + ensureContains: (extent: Extent3D, padding = 0) => { + requireImage(); + return ensureMaskContains(maskId, extent, padding); + }, + }; + } + + function maskVoxels(maskId: string): MaskVoxelAccessor { + // Validates eagerly: an accessor for a nonexistent segment is refused up + // front, not just on first use. + getMask(maskId); + + const binding = () => getMask(maskId).representations.labelmap; + + // Deliberately tolerant where binding() is not: the segment itself can be + // deleted out from under an accessor, and that is an absent storage, not a + // lookup error. + const findBinding = () => findMask(maskId)?.representations.labelmap; + + const onMissing = (): never => { + throw new Error('No storage: call materialize() first'); + }; + + return { + binding, + materialize: () => ensureLabelmapBinding(maskId), + ...voxelStorage(maskId, findBinding, onMissing), + }; + } + + /** + * The accessor for consumers holding an id a segment may already have left: + * the renderer and the paint widget are computeds keyed on one that can + * vanish a tick before they do, so this stays constructible either way. + */ + const findMaskVoxels = (maskId: string) => + voxelStorage( + maskId, + () => findMask(maskId)?.representations.labelmap, + () => { + throw new Error('No such segment'); + } + ); + + /** A bound segment's buffer, absent when it has none or holds nothing. */ + const boundedMask = (binding?: LabelmapBinding) => + binding && boundScalars(binding.image, binding.extent); + + /** + * The masks of an image's other segments that `gesture` may take a voxel + * from, resolved once per run because the caller below runs per voxel. A + * locked segment is not editable, so an aimed gesture is not offered its mask + * at all. + */ + function siblingMasks(maskId: string, gesture: VoxelGesture) { + const segmentation = segmentationOfMask(maskId); + if (!segmentation) return []; + return listMasks(segmentation).flatMap((segment) => { + if (segment.id === maskId) return []; + if (gesture === 'aimed' && maskLocked(segment)) return []; + const bounded = boundedMask(segment.representations.labelmap); + return bounded ? [bounded] : []; + }); + } + + /** + * Whether the voxel at PARENT indices i, j, k is this segment's to write, + * taking it from the neighbours that have to yield it. Absent when no other + * segment reaches `within`, the box the caller is about to walk: every voxel + * in it is then uncontested and the question need not be asked per voxel. + * + * `gesture` is the whole of the policy, so see {@link VoxelGesture}. An aimed + * operation must call finish in a finally block after its last voxel write. + */ + function voxelClaim(maskId: string, gesture: VoxelGesture, within: Extent3D) { + const masks = siblingMasks(maskId, gesture); + if (gesture === 'aimed') return masksClearing(masks, within); + const held = masksHolding(masks, within); + return ( + held && { + claim: (i: number, j: number, k: number) => !held(i, j, k), + finish: () => undefined, + } + ); + } + + return { + requireParentImage, + ensureMaskContains, + maskVoxels, + findMaskVoxels, + boundedMask, + siblingMasks, + voxelClaim, + }; +} diff --git a/src/segmentation/model.ts b/src/segmentation/model.ts new file mode 100644 index 000000000..76f62d905 --- /dev/null +++ b/src/segmentation/model.ts @@ -0,0 +1,127 @@ +import type { Extent3D } from '@/src/segmentation/geometry'; +import type { ProcessingResultSource } from '@/src/types'; +import type { RGBAColor, TypedArray } from '@kitware/vtk.js/types'; + +import type vtkLabelMap from '@/src/vtk/LabelMap'; + +/** A fresh segmentation tints the anatomy under it rather than hiding it. */ +export const DEFAULT_SEGMENTATION_FILL_OPACITY = 0.3; + +export type LabelmapBinding = { + /** + * This mask's voxels, and no other segment's. Held raw: a vtk object must + * not be proxied, so every writer of a binding marks it. + */ + image: vtkLabelMap; + extent: Extent3D; // the mask's own bounds, in parent image index space + /** Reaches the saved archive's entry path, so a round trip keeps it. */ + name: string; + source?: ProcessingResultSource; +}; + +/** + * One image's mask for one segment type. Its id is its own, distinct from the + * type id: everything the user sees or sets, visibility and lock included, + * lives on the type, so this record is storage and nothing else. + */ +export type SegmentMask = { + id: string; + segmentId: string; + representations: { + // absent until voxels are allocated + labelmap?: LabelmapBinding; + }; +}; + +export const LABELMAP_BACKGROUND_VALUE = 0; + +export const makeDefaultSegmentName = (value: number) => `Segment ${value}`; + +/** + * One mask's label descriptor, derived from the segment it delineates. + * Identity lives on `Segment`; this is the value-keyed view the labelmap + * renderer and the .seg.nrrd writer consume. + */ +export type LabelmapSegment = { + value: number; + name: string; + color: RGBAColor; + visible: boolean; + locked?: boolean; + // Absent on descriptors that come off a file rather than off a segment. + fillOpacity?: number; + outlineOpacity?: number; +}; + +/** vtk declares getData() as number[] | TypedArray; mask storage is typed. */ +export const maskScalars = (mask: vtkLabelMap) => + mask.getPointData().getScalars().getData() as Uint8Array; + +export type Segmentation = { + id: string; + name: string; + parentImageId: string; + masks: Record; + order: string[]; + fillOpacity: number; + outlineOpacity: number; + outlineThickness: number; +}; + +/** The display multipliers every segment of a segmentation is scaled by. */ +export type SegmentationDisplayPatch = Partial< + Pick +>; + +/** + * The voxel operations every labelmap consumer routes through. Storage is one + * bounded mask per segment, sized to the region that segment covers. + * + * `ensureContains` may replace the scalar array, dimensions and strides: + * anything that cached those from `image()` or `scalars()` must re-fetch after + * calling it. + */ +export type VoxelStorage = { + /** + * Whether the storage is still reachable. An accessor outlives what it + * points at, so callers holding one across a deletion check this before a + * read or a write; every other method throws while it is false. + */ + exists(): boolean; + image(): vtkLabelMap; + /** Live mask buffer. Writers publish changes through apply() or image().modified(). */ + scalars(): TypedArray; + snapshot(): TypedArray; + /** Bulk copy-in; keeps image() and scalars() identity, marks it modified. */ + apply(scalars: TypedArray | number[]): void; + /** + * Ensures storage covers `extent`, growing to the union of what it has and + * what it was asked for. Returns whether storage was invalidated + * (scalars/dimensions/strides changed). An empty extent is already covered. + * Throws when the extent leaves the parent image, so callers clip. When the + * extent is not already covered, the mask grows by `padding` voxels beyond + * it on every face (clipped to the parent), so nearby requests that follow + * grow nothing. + */ + ensureContains(extent: Extent3D, padding?: number): boolean; +}; + +/** + * Voxel access for one segment. Re-resolves the binding on every call rather + * than capturing it, so a caller that holds an accessor across a segment + * deletion or a growth sees the current state, not a stale one. `exists()` is + * false, and every storage method throws, before `materialize()`. + */ +export type MaskVoxelAccessor = VoxelStorage & { + binding(): LabelmapBinding | undefined; + /** Allocates storage if needed and returns the binding. Idempotent. */ + materialize(): LabelmapBinding; +}; + +/** Segments in display order. `order` is the authority, `segments` the store. */ +export function listMasks(segmentation: Segmentation) { + return segmentation.order.map((id) => segmentation.masks[id]); +} + +/** Aimed writes clear unlocked neighbors; sweeps only grow into unclaimed voxels. */ +export type VoxelGesture = 'aimed' | 'sweep'; diff --git a/src/segmentation/rendering/__tests__/segmentDisplay.spec.ts b/src/segmentation/rendering/__tests__/segmentDisplay.spec.ts new file mode 100644 index 000000000..bfb018c2a --- /dev/null +++ b/src/segmentation/rendering/__tests__/segmentDisplay.spec.ts @@ -0,0 +1,129 @@ +import { describe, expect, it } from 'vitest'; + +import { + SEGMENT_ACTOR_OPACITY, + segmentFillAlpha, + segmentOutlineTables, +} from '@/src/segmentation/rendering/display'; +import type { LabelmapSegment } from '@/src/segmentation/model'; + +const makeMask = ( + value: number, + overrides: Partial = {} +): LabelmapSegment => ({ + value, + name: `Segment ${value}`, + color: [255, 0, 0, 255], + visible: true, + ...overrides, +}); + +describe('segmentFillAlpha', () => { + it('is the segment alpha when the fill is fully opaque', () => { + expect(segmentFillAlpha(makeMask(1, { fillOpacity: 1 }))).toBe(1); + }); + + it('scales the segment alpha by the fill opacity', () => { + expect(segmentFillAlpha(makeMask(1, { fillOpacity: 0.5 }))).toBe(0.5); + }); + + it('hides a fill the user set to zero', () => { + expect(segmentFillAlpha(makeMask(1, { fillOpacity: 0 }))).toBe(0); + }); + + it('hides an invisible segment whatever its fill opacity', () => { + expect( + segmentFillAlpha(makeMask(1, { visible: false, fillOpacity: 1 })) + ).toBe(0); + }); + + it('treats a descriptor without a fill opacity as opaque', () => { + expect(segmentFillAlpha(makeMask(1))).toBe(1); + }); + + it('scales the segment alpha by the segmentation\u2019s fill opacity', () => { + expect(segmentFillAlpha(makeMask(1, { fillOpacity: 0.5 }), 0.5)).toBe(0.25); + }); + + it('hides every fill when the segmentation\u2019s fill opacity is zero', () => { + expect(segmentFillAlpha(makeMask(1, { fillOpacity: 1 }), 0)).toBe(0); + }); +}); + +describe('segmentOutlineTables', () => { + it('indexes both tables by label value minus one', () => { + const tables = segmentOutlineTables( + [ + makeMask(1, { outlineOpacity: 0.25 }), + makeMask(2, { outlineOpacity: 0.5 }), + ], + 2, + 1 + ); + + expect(tables.opacities).toEqual([0.25, 0.5]); + expect(tables.thicknesses).toEqual([2, 2]); + }); + + it('hides an outline the user set to zero', () => { + const tables = segmentOutlineTables( + [makeMask(1, { outlineOpacity: 0 })], + 2, + 1 + ); + + expect(tables.opacities).toEqual([0]); + }); + + it('scales every segment by the group outline opacity', () => { + const tables = segmentOutlineTables( + [makeMask(1, { outlineOpacity: 0.5 })], + 2, + 0.5 + ); + + expect(tables.opacities).toEqual([0.25]); + }); + + it('leaves values no segment claims at the group defaults', () => { + const tables = segmentOutlineTables( + [makeMask(3, { outlineOpacity: 0.5 })], + 2, + 1 + ); + + expect(tables.opacities).toEqual([1, 1, 0.5]); + expect(tables.thicknesses).toEqual([2, 2, 2]); + }); + + it('drops the thickness of an invisible segment', () => { + const tables = segmentOutlineTables( + [makeMask(1, { visible: false }), makeMask(2)], + 2, + 1 + ); + + expect(tables.thicknesses).toEqual([0, 2]); + }); + + it('has no entries for an artifact with no bound segments', () => { + expect(segmentOutlineTables([], 2, 1)).toEqual({ + thicknesses: [], + opacities: [], + }); + }); +}); + +describe('SEGMENT_ACTOR_OPACITY', () => { + it('leaves the fill to the transfer functions', () => { + // A fully opaque segment reaches the screen at its own alpha, so the actor + // must not scale it down. + expect(SEGMENT_ACTOR_OPACITY).toBeGreaterThan(0.999); + }); + + it('stays out of the opaque render pass', () => { + // vtk.js treats an image slice at opacity 1 as opaque and restacks it + // against the base image and the sibling segment actors. + expect(SEGMENT_ACTOR_OPACITY).toBeLessThan(1); + }); +}); diff --git a/src/segmentation/rendering/__tests__/segmentRenderMask.spec.ts b/src/segmentation/rendering/__tests__/segmentRenderMask.spec.ts new file mode 100644 index 000000000..efb2930fb --- /dev/null +++ b/src/segmentation/rendering/__tests__/segmentRenderMask.spec.ts @@ -0,0 +1,114 @@ +import { describe, expect, it } from 'vitest'; +import vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; +import { allocateMask, regrowMask } from '@/src/segmentation/masks/storage'; +import { maskScalars } from '@/src/segmentation/model'; +import { type Extent3D } from '@/src/segmentation/geometry'; +import { segmentRenderMask } from '@/src/segmentation/rendering/renderMask'; + +function scene(extent: Extent3D) { + const parent = vtkImageData.newInstance({ + spacing: [2, 3, 4], + origin: [10, 20, 30], + direction: [0, 1, 0, -1, 0, 0, 0, 0, 1], + }); + parent.setExtent(4, 9, 5, 10, 6, 11); + const source = allocateMask(parent, extent); + maskScalars(source).fill(1); + source.modified(); + return { parent, source }; +} + +describe('render-only segment slices', () => { + it.each([0, 1, 2])('pads only the displayed plane on axis %i', (axis) => { + const extent: Extent3D = [6, 6, 7, 7, 8, 8]; + const { parent, source } = scene(extent); + const original = maskScalars(source); + const origin = [...source.getOrigin()]; + const rendered = segmentRenderMask(source, parent, extent, { + axis: axis, + index: extent[axis * 2], + })!; + const dimensions = [3, 3, 3]; + dimensions[axis] = 1; + expect(rendered.getDimensions()).toEqual(dimensions); + const center: [number, number, number] = [1, 1, 1]; + center[axis] = 0; + expect(Array.from(rendered.indexToWorld(center))).toEqual(origin); + expect([...maskScalars(rendered)]).toEqual([0, 0, 0, 0, 1, 0, 0, 0, 0]); + expect(source.getDimensions()).toEqual([1, 1, 1]); + expect(source.getOrigin()).toEqual(origin); + expect(maskScalars(source)).toBe(original); + expect([...original]).toEqual([1]); + expect( + segmentRenderMask(source, parent, extent, { + axis: axis, + index: extent[axis * 2] + 1, + }) + ).toBeNull(); + expect( + segmentRenderMask(source, parent, [0, -1, 0, -1, 0, -1], { + axis: axis, + index: 0, + }) + ).toBeNull(); + }); + + it.each([0, 1, 2, 3, 4, 5])('clips padding at parent face %i', (face) => { + const extent: Extent3D = [6, 6, 7, 7, 8, 8]; + const bounds = [4, 9, 5, 10, 6, 11]; + const boundaryAxis = Math.floor(face / 2); + const axis = (boundaryAxis + 1) % 3; + extent[boundaryAxis * 2] = bounds[face]; + extent[boundaryAxis * 2 + 1] = bounds[face]; + const { parent, source } = scene(extent); + const rendered = segmentRenderMask(source, parent, extent, { + axis: axis, + index: extent[axis * 2], + })!; + expect(rendered.getDimensions()[boundaryAxis]).toBe(2); + const corner: [number, number, number] = [0, 0, 0]; + if (face % 2) corner[boundaryAxis] = 1; + expect( + parent.worldToIndex(rendered.indexToWorld(corner))[boundaryAxis] + ).toBeCloseTo(bounds[face]); + expect([...maskScalars(rendered)].reduce((a, b) => a + b, 0)).toBe(1); + }); + + it('refreshes after edits, slice changes and growth', () => { + const extent: Extent3D = [6, 6, 7, 7, 8, 9]; + const { parent, source } = scene(extent); + const first = segmentRenderMask(source, parent, extent, { + axis: 2, + index: 8, + })!; + expect( + segmentRenderMask(source, parent, extent, { axis: 2, index: 8 }) + ).toBe(first); + maskScalars(source)[0] = 0; + source.modified(); + expect( + segmentRenderMask(source, parent, extent, { axis: 2, index: 8 }) + ).toBe(first); + expect([...maskScalars(first)].every((v) => v === 0)).toBe(true); + expect(first.getPointData().getScalars().getRange()).toEqual([0, 0]); + const nextSlice = segmentRenderMask(source, parent, extent, { + axis: 2, + index: 9, + })!; + expect([...maskScalars(nextSlice)]).toEqual([0, 0, 0, 0, 1, 0, 0, 0, 0]); + expect( + Array.from(parent.worldToIndex(nextSlice.indexToWorld([1, 1, 0]))) + ).toEqual([6, 7, 9]); + const grown: Extent3D = [5, 6, 7, 7, 8, 9]; + regrowMask(source, parent, extent, grown); + maskScalars(source)[0] = 1; + source.modified(); + const next = segmentRenderMask(source, parent, grown, { + axis: 2, + index: 8, + })!; + expect(next.getDimensions()).toEqual([4, 3, 1]); + expect(next.getPointData().getScalars().getRange()).toEqual([0, 1]); + expect([...maskScalars(next)].reduce((a, b) => a + b, 0)).toBe(1); + }); +}); diff --git a/src/segmentation/rendering/__tests__/segmentSliceVisibility.spec.ts b/src/segmentation/rendering/__tests__/segmentSliceVisibility.spec.ts new file mode 100644 index 000000000..d9d3d39fe --- /dev/null +++ b/src/segmentation/rendering/__tests__/segmentSliceVisibility.spec.ts @@ -0,0 +1,87 @@ +import { describe, expect, it } from 'vitest'; + +import { + segmentCoincidentOffset, + sliceWithinExtent, +} from '@/src/segmentation/rendering/display'; +import { emptyExtent, type Extent3D } from '@/src/segmentation/geometry'; + +// --------------------------------------------------------------------------- +// The two view-layer rules a mask per segment needs. +// +// `sliceWithinExtent` answers whether a segment's actor has anything to draw on +// the slice being viewed. A bounded mask covers only part of the volume, and +// vtkImageMapper clamps a slice outside its input to the nearest one, so an +// actor left visible off its own extent paints a stale slice over the image. +// The slice and the extent are both in the PARENT image's index space, on the +// index axis the view's LPS axis maps to. +// +// `segmentCoincidentOffset` gives each segment its own coincident-topology +// polygon offset, by its back-to-front stack index. Overlap is +// representable, so segments sharing one offset would z-fight. +// Greater stack indices draw in front. Registry order is mapped in reverse. +// --------------------------------------------------------------------------- + +const EXTENT: Extent3D = [1, 2, 0, 3, 2, 5]; + +describe('sliceWithinExtent', () => { + it('is true for a slice inside the extent', () => { + expect(sliceWithinExtent(EXTENT, 2, 3)).toBe(true); + }); + + it('includes both ends of the extent', () => { + expect(sliceWithinExtent(EXTENT, 2, 2)).toBe(true); + expect(sliceWithinExtent(EXTENT, 2, 5)).toBe(true); + }); + + it('is false below the extent', () => { + expect(sliceWithinExtent(EXTENT, 2, 1)).toBe(false); + }); + + it('is false above the extent', () => { + expect(sliceWithinExtent(EXTENT, 2, 6)).toBe(false); + }); + + it('reads the axis it is given', () => { + expect(sliceWithinExtent(EXTENT, 0, 3)).toBe(false); + expect(sliceWithinExtent(EXTENT, 1, 3)).toBe(true); + }); + + it('is false for a mask that covers nothing, on every axis', () => { + expect(sliceWithinExtent(emptyExtent(), 0, 0)).toBe(false); + expect(sliceWithinExtent(emptyExtent(), 1, 0)).toBe(false); + expect(sliceWithinExtent(emptyExtent(), 2, 0)).toBe(false); + }); +}); + +describe('segmentCoincidentOffset', () => { + it('puts the first segment in front of the base image', () => { + const [factor, units] = segmentCoincidentOffset(0); + + expect(factor).toBeLessThan(0); + expect(units).toBeLessThan(0); + }); + + it('puts a greater stack index in front of a smaller one', () => { + const [earlierFactor, earlierUnits] = segmentCoincidentOffset(0); + const [laterFactor, laterUnits] = segmentCoincidentOffset(1); + + expect(laterUnits).toBeLessThan(earlierUnits); + expect(laterFactor).toBeLessThanOrEqual(earlierFactor); + }); + + it('keeps that order all the way down a long list', () => { + const offsets = Array.from({ length: 64 }, (_, index) => + segmentCoincidentOffset(index) + ); + + expect( + offsets.every( + ([factor, units]) => Number.isFinite(factor) && Number.isFinite(units) + ) + ).toBe(true); + offsets.slice(1).forEach(([, units], index) => { + expect(units).toBeLessThan(offsets[index][1]); + }); + }); +}); diff --git a/src/segmentation/rendering/display.ts b/src/segmentation/rendering/display.ts new file mode 100644 index 000000000..484cae0e3 --- /dev/null +++ b/src/segmentation/rendering/display.ts @@ -0,0 +1,87 @@ +import type { LabelmapSegment } from '@/src/segmentation/model'; +import type { Extent3D } from '@/src/segmentation/geometry'; +import { isEmptyExtent } from '@/src/segmentation/geometry'; + +/** + * Whether a segment's actor has anything to draw on the slice being viewed. + * Extent and slice are both in the parent image's index space, on the index + * axis the view's LPS axis maps to. vtkImageMapper clamps a slice outside its + * input to the nearest one, so an actor left visible off its own extent paints + * a stale slice over the image. + */ +export function sliceWithinExtent( + extent: Extent3D, + axisIndex: number, + slice: number +) { + if (isEmptyExtent(extent)) return false; + return slice >= extent[axisIndex * 2] && slice <= extent[axisIndex * 2 + 1]; +} + +const SEGMENT_OFFSET_FACTOR = -4; + +/** + * Actor opacity for a segment's slice representation. Per-segment and + * per-segmentation opacity live in the transfer functions, so the actor itself + * carries none. It must stay below 1: vtk.js puts an image slice in the opaque + * render pass at an opacity of 1, which restacks it against the base image and + * the sibling segment actors. + */ +export const SEGMENT_ACTOR_OPACITY = 0.9999; + +/** + * A mask's coincident-topology polygon offset, by back-to-front stack index. + * Overlapping segments need distinct offsets to avoid z-fighting. Greater + * stack indices sit closer to the viewer; the registry maps its first entry + * to the greatest index. + */ +export function segmentCoincidentOffset(stackIndex: number) { + return [SEGMENT_OFFSET_FACTOR, SEGMENT_OFFSET_FACTOR - stackIndex] as [ + number, + number, + ]; +} + +/** + * Fill alpha in 0..1 for the slice representation's piecewise function: the + * segment's own alpha scaled by its fill opacity and by the segmentation's, + * the same way the outline tables compose theirs. + */ +export const segmentFillAlpha = ( + segment: LabelmapSegment, + segmentationOpacity = 1 +) => + segment.visible + ? ((segment.color[3] || 0) / 255) * + (segment.fillOpacity ?? 1) * + segmentationOpacity + : 0; + +/** + * The label outline tables vtk.js indexes by label value minus one, so both run + * from value 1 to the largest value in use. A value no segment claims keeps the + * segmentation defaults. + */ +export const segmentOutlineTables = ( + segments: LabelmapSegment[], + segmentationThickness: number, + segmentationOpacity: number +) => { + const byValue = new Map(segments.map((segment) => [segment.value, segment])); + const largestValue = segments.reduce( + (largest, segment) => Math.max(largest, segment.value), + 0 + ); + const at = (index: number) => byValue.get(index + 1); + + return { + thicknesses: Array.from({ length: largestValue }, (_, index) => { + const segment = at(index); + return !segment || segment.visible ? segmentationThickness : 0; + }), + opacities: Array.from( + { length: largestValue }, + (_, index) => segmentationOpacity * (at(index)?.outlineOpacity ?? 1) + ), + }; +}; diff --git a/src/segmentation/rendering/projection.ts b/src/segmentation/rendering/projection.ts new file mode 100644 index 000000000..70c2d57ef --- /dev/null +++ b/src/segmentation/rendering/projection.ts @@ -0,0 +1,60 @@ +import { computed } from 'vue'; + +import type { SegmentRegistry } from '@/src/segmentation/segmentRegistry'; +import { + sameLabelmapSegment, + toLabelmapSegment, +} from '@/src/segmentation/segment'; +import { SEGMENT_VALUE } from '@/src/segmentation/masks/labelValue'; +import { + listMasks, + type LabelmapSegment, + type Segmentation, +} from '@/src/segmentation/model'; + +export type SegmentProjectionDeps = { + segmentations: Record; + segmentRegistry: SegmentRegistry; +}; + +/** One descriptor per bound mask; unchanged appearances retain object identity. */ +export function createSegmentProjection({ + segmentations, + segmentRegistry, +}: SegmentProjectionDeps) { + function project() { + const byMask: Record = {}; + Object.values(segmentations).forEach((segmentation) => { + listMasks(segmentation).forEach((segment) => { + if (!segment.representations.labelmap) return; + byMask[segment.id] = toLabelmapSegment( + segmentRegistry.getSegment(segment.segmentId), + SEGMENT_VALUE + ); + }); + }); + return byMask; + } + + let projected: Record = {}; + + return computed(() => { + const fresh = project(); + const stable = Object.fromEntries( + Object.entries(fresh).map(([maskId, list]) => { + const previous = projected[maskId]; + return [ + maskId, + previous && sameLabelmapSegment(previous, list) ? previous : list, + ]; + }) + ); + // The record's own identity is what a consumer of the whole projection + // watches, so it survives a change that left every mask alone. + const unchanged = + Object.keys(stable).length === Object.keys(projected).length && + Object.entries(stable).every(([id, list]) => projected[id] === list); + if (!unchanged) projected = stable; + return projected; + }); +} diff --git a/src/segmentation/rendering/renderMask.ts b/src/segmentation/rendering/renderMask.ts new file mode 100644 index 000000000..78a2343c4 --- /dev/null +++ b/src/segmentation/rendering/renderMask.ts @@ -0,0 +1,83 @@ +import type vtkImageData from '@kitware/vtk.js/Common/DataModel/ImageData'; +import type vtkLabelMap from '@/src/vtk/LabelMap'; +import { + allocateMask, + reframeMaskScalars, +} from '@/src/segmentation/masks/storage'; +import { SEGMENT_VALUE } from '@/src/segmentation/masks/labelValue'; +import { maskScalars } from '@/src/segmentation/model'; +import { + clipExtent, + isEmptyExtent, + padExtent, + type Extent3D, +} from '@/src/segmentation/geometry'; + +// Keep only the current slice per axis, outside segmentation storage and export. +const renderMasks = new WeakMap< + vtkLabelMap, + Map< + number, + { + image: vtkLabelMap; + key: string; + mtime: number; + } + > +>(); + +/** Add known background within the scan, without inventing data beyond it. */ +export function segmentRenderMask( + source: vtkLabelMap, + parent: vtkImageData, + extent: Extent3D, + { axis, index: slice }: { axis: number; index: number } +) { + if (isEmptyExtent(extent)) return null; + const padded = clipExtent( + padExtent(extent, 1), + parent.getExtent() as Extent3D + ); + const index = Math.round(slice); + if (index < extent[axis * 2] || index > extent[axis * 2 + 1]) return null; + padded[axis * 2] = index; + padded[axis * 2 + 1] = index; + const key = [ + ...extent, + ...padded, + ...parent.getOrigin(), + ...parent.getSpacing(), + ...parent.getDirection(), + ].join(','); + let slices = renderMasks.get(source); + if (!slices) { + slices = new Map(); + renderMasks.set(source, slices); + } + let cached = slices.get(axis); + if (!cached || cached.key !== key) { + cached = { image: allocateMask(parent, padded), key, mtime: -1 }; + slices.set(axis, cached); + } + if (cached.mtime !== source.getMTime()) { + const values = reframeMaskScalars( + maskScalars(source), + [...extent], + padded, + maskScalars(cached.image) + ); + const scalars = cached.image.getPointData().getScalars(); + scalars.dataChange(); + // Each stored mask is binary. Supply the range to avoid another scan. + scalars.setRange( + { + min: values.includes(0) ? 0 : SEGMENT_VALUE, + max: values.includes(SEGMENT_VALUE) ? SEGMENT_VALUE : 0, + }, + 0 + ); + cached.image.modified(); + cached.mtime = source.getMTime(); + } + return cached.image; +} diff --git a/src/segmentation/segment.ts b/src/segmentation/segment.ts new file mode 100644 index 000000000..68c089fc3 --- /dev/null +++ b/src/segmentation/segment.ts @@ -0,0 +1,97 @@ +import type { RGBAColor } from '@kitware/vtk.js/types'; + +import { + STROKE_WIDTH_ANNOTATION_TOOL_DEFAULT, + TOOL_COLORS, +} from '@/src/config'; +import type { Maybe } from '@/src/types'; +import { cleanUndefined } from '@/src/utils'; +import type { LabelmapSegment } from '@/src/segmentation/model'; +import { cssColorToRGBA, rgbaToCssColor } from '@/src/segmentation/color'; + +/** + * Identity and shared appearance for everything drawn as one thing: a paint + * mask on any image, a rectangle, a polygon. Appearance fields are absent + * until set and mean "app default" while they are, so a configured or imported + * type that states nothing follows the default and a ruler type carries no + * meaningless opacity. + */ +export type Segment = { + id: string; + name: string; + color: RGBAColor; + // State the user sets on the thing itself, so it holds on every image. + visible: boolean; + locked: boolean; + fillOpacity?: number; + outlineOpacity?: number; + strokeWidth?: number; +}; + +export type SegmentInit = Partial>; + +export const DEFAULT_SEGMENT_COLOR = cssColorToRGBA(TOOL_COLORS[0]); + +const APPEARANCE_DEFAULTS = { + name: '', + color: DEFAULT_SEGMENT_COLOR, + visible: true, + locked: false, + fillOpacity: 1, + outlineOpacity: 1, + strokeWidth: STROKE_WIDTH_ANNOTATION_TOOL_DEFAULT, +}; + +/** + * The one resolver. Every renderer, editor and encoder reads a type through + * it; nothing reads the optional fields directly, so an absent field means the + * app default in exactly one place. + */ +export const resolveSegmentAppearance = (type: Maybe) => { + const stated: Partial = type ?? {}; + const resolved = { + ...APPEARANCE_DEFAULTS, + ...cleanUndefined({ + name: stated.name, + color: stated.color, + visible: stated.visible, + locked: stated.locked, + fillOpacity: stated.fillOpacity, + outlineOpacity: stated.outlineOpacity, + strokeWidth: stated.strokeWidth, + }), + }; + return { ...resolved, cssColor: rgbaToCssColor(resolved.color) }; +}; + +export type SegmentAppearance = ReturnType; + +/** + * The descriptor a record projects onto the label value its mask holds, all of + * it resolved from the segment. The labelmap renderer and the .seg.nrrd writer + * consume it. + */ +export const toLabelmapSegment = ( + type: Maybe, + labelValue: number +): LabelmapSegment => { + const resolved = resolveSegmentAppearance(type); + return { + value: labelValue, + name: resolved.name, + color: [...resolved.color] as RGBAColor, + visible: resolved.visible, + locked: resolved.locked, + fillOpacity: resolved.fillOpacity, + outlineOpacity: resolved.outlineOpacity, + }; +}; + +export const sameLabelmapSegment = (a: LabelmapSegment, b: LabelmapSegment) => + a.value === b.value && + a.name === b.name && + a.visible === b.visible && + a.locked === b.locked && + a.fillOpacity === b.fillOpacity && + a.outlineOpacity === b.outlineOpacity && + a.color.every((channel, index) => channel === b.color[index]); diff --git a/src/segmentation/segmentReferences.ts b/src/segmentation/segmentReferences.ts new file mode 100644 index 000000000..62541be05 --- /dev/null +++ b/src/segmentation/segmentReferences.ts @@ -0,0 +1,40 @@ +import { getActivePinia, type Pinia } from 'pinia'; + +// Who points at a segment, declared by the stores that hold references. +// +// The registry knows nothing about masks or shapes: it asks these declarations +// whether a segment is still referenced and hands them the removal. Declarations +// are registered from store setup and scoped to the application instance that +// ran it, so two applications, or two tests, never answer for each other. +// Dependency-free on purpose: a store import here would close a cycle back +// through the registry. + +export type SegmentReferenceHolder = { + has: (segmentId: string) => boolean; + remove: (segmentId: string) => void; +}; + +const holdersByApp = new WeakMap>(); + +const holdersOf = () => { + const pinia = getActivePinia(); + if (!pinia) return undefined; + const existing = holdersByApp.get(pinia); + if (existing) return existing; + const holders = new Map(); + holdersByApp.set(pinia, holders); + return holders; +}; + +export function declareSegmentReferences( + name: string, + holder: SegmentReferenceHolder +) { + holdersOf()?.set(name, holder); +} + +export const segmentIsReferenced = (segmentId: string) => + [...(holdersOf()?.values() ?? [])].some((holder) => holder.has(segmentId)); + +export const removeSegmentReferences = (segmentId: string) => + holdersOf()?.forEach((holder) => holder.remove(segmentId)); diff --git a/src/segmentation/segmentRegistry.ts b/src/segmentation/segmentRegistry.ts new file mode 100644 index 000000000..14ef5b2c6 --- /dev/null +++ b/src/segmentation/segmentRegistry.ts @@ -0,0 +1,276 @@ +import { computed, ref, type Ref } from 'vue'; + +import { TOOL_COLORS } from '@/src/config'; +import { useIdStore } from '@/src/store/id'; +import type { Maybe } from '@/src/types'; +import { cssColorToRGBA } from '@/src/segmentation/color'; +import { + resolveSegmentAppearance, + type Segment, + type SegmentInit, +} from '@/src/segmentation/segment'; +import { omit } from '@/src/utils'; +import { cleanUndefined } from '@/src/utils'; + +export type ConfiguredSegment = { + color?: string; + fillOpacity?: number; + outlineOpacity?: number; + strokeWidth?: number; +}; + +export type ConfiguredSegments = Record; + +export type SegmentRegistryOptions = { + hasReferences?: (segmentId: string) => boolean; + removeReferences?: (segmentId: string) => void; +}; + +const fromConfigured = ( + name: string, + configured: ConfiguredSegment +): SegmentInit => + cleanUndefined({ + name, + color: configured.color ? cssColorToRGBA(configured.color) : undefined, + fillOpacity: configured.fillOpacity, + outlineOpacity: configured.outlineOpacity, + strokeWidth: configured.strokeWidth, + }); + +const configuredAppearance = ({ + color, + fillOpacity, + outlineOpacity, + strokeWidth, +}: Segment) => ({ color, fillOpacity, outlineOpacity, strokeWidth }); + +/** + * Identity and shared appearance for a family of segments: one instance backs + * paint, rectangles, polygons and rulers together. Explicit order drives the + * picker, shortcuts, serialization and labelmap stacking. + */ +export const createSegmentRegistry = ({ + hasReferences = () => false, + removeReferences = () => {}, +}: SegmentRegistryOptions = {}) => { + const segmentById = ref>({}) as Ref< + Record + >; + + const segmentOrder = ref([]); + const segmentList = computed(() => + segmentOrder.value.map((id) => segmentById.value[id]) + ); + + const selectedSegmentId = ref>(); + const selectionRevision = ref(0); + + // A type that is gone is not selected. + const selectedSegment = computed(() => + selectedSegmentId.value + ? segmentById.value[selectedSegmentId.value] + : undefined + ); + + const selectSegment = (id: Maybe) => { + selectedSegmentId.value = id && segmentById.value[id] ? id : undefined; + // Reselecting the same segment can still request that its row be revealed. + if (selectedSegmentId.value) selectionRevision.value += 1; + }; + + const getSegment = (id: Maybe) => + id ? segmentById.value[id] : undefined; + + const appearanceOf = (id: Maybe) => + resolveSegmentAppearance(getSegment(id)); + + // Cached: the renderer asks for one index per mask per re-render, and the + // export sort asks twice per comparison. + const orderIndex = computed( + () => new Map(segmentOrder.value.map((id, index) => [id, index])) + ); + + const orderIndexOf = (id: Maybe) => + (id === undefined || id === null ? undefined : orderIndex.value.get(id)) ?? + -1; + + const findSegmentByName = (name: Maybe) => + segmentList.value.find((type) => type.name === name); + + const uniqueName = (stem: string) => { + const taken = new Set(segmentList.value.map((type) => type.name.trim())); + if (!taken.has(stem)) return stem; + let index = 2; + while (taken.has(`${stem} (${index})`)) index += 1; + return `${stem} (${index})`; + }; + + const defaultName = () => { + const taken = new Set(segmentList.value.map((type) => type.name.trim())); + let index = 1; + while (taken.has(`Segment ${index}`)) index += 1; + return `Segment ${index}`; + }; + + let nextColorIndex = 0; + const nextColor = () => { + const color = cssColorToRGBA(TOOL_COLORS[nextColorIndex]); + nextColorIndex = (nextColorIndex + 1) % TOOL_COLORS.length; + return color; + }; + + /** Mints a segment without touching the selection. Allocates no voxels. */ + const mintSegment = (init: SegmentInit = {}) => { + const id = useIdStore().nextId(); + segmentById.value = { + ...segmentById.value, + [id]: { + name: defaultName(), + color: nextColor(), + visible: true, + locked: false, + ...cleanUndefined(init), + id, + }, + }; + segmentOrder.value = [...segmentOrder.value, id]; + return id; + }; + + const addSegment = (init: SegmentInit = {}) => { + const id = mintSegment(init); + selectSegment(id); + return id; + }; + + const updateSegment = (id: string, patch: SegmentInit) => { + const type = segmentById.value[id]; + if (!type) return; + segmentById.value = { + ...segmentById.value, + [id]: { ...type, ...patch, id }, + }; + }; + + // Deleting a referenced segment takes its masks and shapes with it; the + // caller owns the confirmation. + const deleteSegment = (id: string) => { + if (!segmentById.value[id]) return; + removeReferences(id); + segmentOrder.value = segmentOrder.value.filter((key) => key !== id); + segmentById.value = omit(segmentById.value, id); + if (selectedSegmentId.value === id) { + selectSegment(segmentOrder.value[0]); + } + }; + + const moveSegment = (id: string, target: string, after = false) => { + if (id === target || !getSegment(id) || !getSegment(target)) return; + const order = segmentOrder.value.filter((key) => key !== id); + order.splice(order.indexOf(target) + Number(after), 0, id); + segmentOrder.value = order; + }; + + const ensureSelectedSegment = () => { + if (selectedSegment.value) return selectedSegment.value.id; + const first = segmentList.value[0]; + if (!first) return addSegment(); + selectSegment(first.id); + return first.id; + }; + + /** + * Exact-name lookup, minting on a miss. Import binds descriptors this way, + * so a file's segment lands on the one already carrying that name and the + * registry's own color wins. + */ + const segmentNamed = (name: string, init: SegmentInit = {}) => { + const existing = findSegmentByName(name); + if (existing) return existing.id; + return mintSegment({ ...init, name }); + }; + + // --- config overlay --- // + + // Keep each key's identity and appearance beneath its config contribution. + // New segments begin with automatic color and default optional appearance; + // a restored segment begins with its session appearance. + const configEntries = new Map< + string, + { id: string; appearance: ReturnType } + >(); + + const replaceConfigSegments = (configured: Maybe) => { + const next = configured ?? {}; + + Object.entries(next).forEach(([name, props]) => { + let entry = configEntries.get(name); + if (!entry || !getSegment(entry.id)) { + const id = findSegmentByName(name)?.id ?? mintSegment({ name }); + entry = { id, appearance: configuredAppearance(getSegment(id)!) }; + } + updateSegment(entry.id, { + ...entry.appearance, + ...fromConfigured(name, props), + }); + configEntries.set(name, entry); + }); + + [...configEntries.entries()] + .filter(([name]) => !(name in next)) + .forEach(([name, { id }]) => { + configEntries.delete(name); + // Content keeps the last configured appearance as session state. + if (segmentById.value[id] && !hasReferences(id)) deleteSegment(id); + }); + + // A configured registry offers a selection from the start; selecting + // creates nothing, so the first edit lands in a configured type rather + // than minting one beside it. + if (!selectedSegment.value) selectSegment(segmentList.value[0]?.id); + }; + + // --- wire --- // + + const serialize = () => segmentList.value.map((type) => ({ ...type })); + + /** + * Seats restored segments beside the ones already here. Ids are minted fresh + * and every incoming reference is remapped through the returned map, so an + * import into a populated scene overwrites nothing. + */ + const adopt = (incoming: Maybe) => + Object.fromEntries( + (incoming ?? []).map(({ id, ...init }) => [ + id, + mintSegment(init as SegmentInit), // side effect in Array.map + ]) + ); + + return { + segmentById, + segmentList, + selectedSegmentId, + selectionRevision, + selectedSegment, + selectSegment, + getSegment, + appearanceOf, + orderIndexOf, + findSegmentByName, + uniqueName, + mintSegment, + addSegment, + updateSegment, + moveSegment, + deleteSegment, + ensureSelectedSegment, + segmentNamed, + replaceConfigSegments, + serialize, + adopt, + }; +}; + +export type SegmentRegistry = ReturnType; diff --git a/src/utils/color.ts b/src/utils/color.ts index e71f644ea..80ca63bba 100644 --- a/src/utils/color.ts +++ b/src/utils/color.ts @@ -1,5 +1,15 @@ import type { RGBAColor } from '@kitware/vtk.js/types'; +/** A cursor over a palette: each call hands out the next opaque colour. */ +export function cycleColors(palette: readonly (readonly number[])[]) { + let index = 0; + return () => { + const color = palette[index]; + index = (index + 1) % palette.length; + return [...color, 255] as RGBAColor; + }; +} + /** * Converts an RGBA tuple to a hex string with alpha. *