Skip to content

Commit

Permalink
feat(MultiscaleSpatialImage): add getImageInImageSpace
Browse files Browse the repository at this point in the history
  • Loading branch information
@PaulHax
PaulHax committed Jun 4, 2024
1 parent f176658 commit d59023e
Show file tree
Hide file tree
Showing 3 changed files with 123 additions and 40 deletions.
156 changes: 120 additions & 36 deletions packages/io/src/MultiscaleSpatialImage.ts
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@ import {
} from './dimensionUtils.js';
import { transformBounds } from './transformBounds.js';
import {
ChunkParameter,
Dimension,
Extent,
ReadOnlyDimensionBounds,
Expand All @@ -42,6 +43,11 @@ function setMatrixElement(
matrixData[column + row * columns] = value;
}

type ReadonlyFloatArray =
| Readonly<Float32Array>
| Readonly<Float64Array>
| ReadonlyMat3;

type ImageDataFromChunksWorkerArgs = {
scaleInfo: {
chunkSize: Map<Dimension, number>;
Expand Down Expand Up @@ -80,7 +86,8 @@ function inflate(bounds: Bounds, delta: number) {
}

// code modified from vtk.js/ImageData
const extentToBounds = (ex: Extent, indexToWorld: ReadonlyMat4) => {
const extentToBounds = (extent: Extent, indexToWorld: ReadonlyMat4) => {
const ex = extent;
// prettier-ignore
const corners = [
ex[0], ex[2], ex[4],
Expand Down Expand Up @@ -130,7 +137,10 @@ const extentToBounds = (ex: Extent, indexToWorld: ReadonlyMat4) => {
};

// returns a copy
export const ensure3dDirection = (d: Float64Array): mat3 => {
export const ensure3dDirection = (
maybe2dDirection: ReadonlyFloatArray,
): mat3 => {
const d = maybe2dDirection;
if (d.length >= 9) {
return mat3.fromValues(
d[0],
Expand Down Expand Up @@ -173,21 +183,28 @@ const makeMat4 = ({
return mat4.scale(mat, mat, spacing);
};

const makeIndexToWorld = ({
type Maybe2dVec =
| [number, number]
| [number, number, number]
| vec3
| Array<number>;

const maybe2dSpatialToMat4 = ({
direction: inDirection,
origin,
spacing,
}: {
direction: ReadonlyMat3;
origin: Array<number>;
spacing: Array<number>;
direction: ReadonlyFloatArray;
origin: Maybe2dVec;
spacing: Maybe2dVec;
}) => {
const inDirection3d = ensure3dDirection(inDirection);
// ITK (and VTKMath) uses row-major index axis, but gl-matrix uses column-major. Transpose.
const DIMENSIONS = 3;
const direction = Array(inDirection.length) as mat3;
const direction = Array(inDirection3d.length) as mat3;
for (let idx = 0; idx < DIMENSIONS; ++idx) {
for (let col = 0; col < DIMENSIONS; ++col) {
direction[col + idx * 3] = inDirection[idx + col * DIMENSIONS];
direction[col + idx * 3] = inDirection3d[idx + col * DIMENSIONS];
}
}

Expand Down Expand Up @@ -240,6 +257,53 @@ export const worldBoundsToIndexBounds = ({
return new Map([...spaceBounds, ...ctBounds]);
};

// Ensures CXYZT dimensions are present
const ensureBoundsCXYZT = ({
indexBounds,
fullIndexBounds,
}: {
indexBounds: ReadonlyBounds;
fullIndexBounds: ReadOnlyDimensionBounds;
}) => {
const fullIndexBoundsWithZCT = ensuredDims(
[0, 1] as [number, number],
CXYZT,
fullIndexBounds,
);
// clamp to existing integer indexes
const imageBoundsByDim = chunk(2, indexBounds);
const spaceBounds = (['x', 'y', 'z'] as const).map((dim, idx) => {
// eslint-disable-next-line @typescript-eslint/no-non-null-assertion
const [min, max] = fullIndexBoundsWithZCT.get(dim)!;
const [bmin, bmax] = imageBoundsByDim[idx];
return [
dim,
[
Math.floor(Math.min(max, Math.max(min, bmin))),
Math.ceil(Math.min(max, Math.max(min, bmax))),
],
] as const;
});
const ctBounds = (['c', 't'] as const).map(
(dim) => [dim, fullIndexBoundsWithZCT.get(dim)!] as const,
);
return new Map([...spaceBounds, ...ctBounds]);
};

const normalizedImageBoundsToIndexBounds = (
arrayShape: ChunkParameter,
normalizedImageBounds: ReadonlyBounds,
) => {
const toIndexScale = XYZ.map((axis) => arrayShape.get(axis) ?? 1);
const normalizedToIndex = maybe2dSpatialToMat4({
direction: mat3.identity(mat3.create()),
origin: vec3.create(),
spacing: toIndexScale,
}) as ReadonlyMat4;
const indexBounds = transformBounds(normalizedToIndex, normalizedImageBounds);
return indexBounds;
};

function isContained(
benchmarkBounds: ReadOnlyDimensionBounds,
testedBounds: ReadOnlyDimensionBounds,
Expand Down Expand Up @@ -328,7 +392,7 @@ export class MultiscaleSpatialImage {
const info = this.scaleInfos[scale];
if (info.origin) return info.origin;

const origin = new Array(this.spatialDims.length);
const origin = new Array<number>(this.spatialDims.length);
for (let index = 0; index < this.spatialDims.length; index++) {
const dim = this.spatialDims[index];
if (info.coords.has(dim)) {
Expand All @@ -347,7 +411,7 @@ export class MultiscaleSpatialImage {
const info = this.scaleInfos[scale];
if (info.spacing) return info.spacing;

const spacing = new Array(this.spatialDims.length);
const spacing = new Array<number>(this.spatialDims.length);
for (let index = 0; index < this.spatialDims.length; index++) {
const dim = this.spatialDims[index];
const dimCoords = await info.coords.get(dim);
Expand Down Expand Up @@ -391,12 +455,6 @@ export class MultiscaleSpatialImage {
}
}

// const rotation = quat.create();
// quat.fromEuler(rotation, 0, 55, 30);
// quat.fromEuler(rotation, 0, 45, 0);
// quat.fromEuler(rotation, 35, 0, 0);
// mat3.fromQuat(direction, rotation);

return direction;
}

Expand Down Expand Up @@ -513,17 +571,39 @@ export class MultiscaleSpatialImage {
this.scaleSpacing(scale),
]);

const direction = ensure3dDirection(this.direction);

const indexToWorld = makeIndexToWorld({
direction,
const indexToWorld = maybe2dSpatialToMat4({
direction: this.direction,
origin,
spacing,
}) as ReadonlyMat4;
this.scaleInfos[scale].indexToWorld = indexToWorld;
return indexToWorld;
}

async buildAndCacheImage(scale: number, indexBounds: Bounds) {
const indexBoundsCXYZT = ensureBoundsCXYZT({
indexBounds,
fullIndexBounds: this.getIndexBounds(scale),
});
const { dims } = this.scaleInfos[scale];
const indexBoundsByDimension = orderBy(dims)(indexBoundsCXYZT);
const cachedImage = findImageInBounds({
cache: this.cachedImages,
scale,
bounds: indexBoundsByDimension,
});
if (cachedImage) return cachedImage;

const image = await this.buildImage(scale, indexBoundsByDimension);
storeImage({
cache: this.cachedImages,
scale,
bounds: indexBoundsByDimension,
image,
});
return image;
}

/* Retrieve bounded image at scale. */
async getImage(
requestedScale: number,
Expand All @@ -532,25 +612,29 @@ export class MultiscaleSpatialImage {
const scale = Math.min(requestedScale, this.scaleInfos.length - 1);
const indexToWorld = await this.scaleIndexToWorld(scale);

const { dims } = this.scaleInfos[scale];
const indexBounds = orderBy(dims)(
worldBoundsToIndexBounds({
bounds: worldBounds,
fullIndexBounds: this.getIndexBounds(scale),
worldToIndex: mat4.invert(mat4.create(), indexToWorld),
}),
const fullIndexBounds = ensuredDims(
[0, 0],
XYZ,
this.getIndexBounds(scale),
);
let indexBounds = XYZ.flatMap((dim) => fullIndexBounds.get(dim)) as Bounds;
if (worldBounds) {
const worldToIndex = mat4.invert(mat4.create(), indexToWorld);
indexBounds = transformBounds(worldToIndex, worldBounds);
}

const cachedImage = findImageInBounds({
cache: this.cachedImages,
scale,
bounds: indexBounds,
});
if (cachedImage) return cachedImage;
return this.buildAndCacheImage(scale, indexBounds);
}

const image = await this.buildImage(scale, indexBounds);
storeImage({ cache: this.cachedImages, scale, bounds: indexBounds, image });
return image;
async getImageInImageSpace(
scale: number,
normalizedImageBounds: ReadonlyBounds = [0, 1, 0, 1, 0, 1],
) {
const indexBounds = normalizedImageBoundsToIndexBounds(
this.scaleInfos[scale].arrayShape,
normalizedImageBounds,
);
return this.buildAndCacheImage(scale, indexBounds);
}

getIndexBounds(scale: number) {
Expand Down
5 changes: 2 additions & 3 deletions packages/io/src/dimensionUtils.ts
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
import { Dimension, SpatialDimensions } from './types.js';

export const CXYZT = Object.freeze(['c', 'x', 'y', 'z', 't'] as const); // viewer indexing
export const XYZ = Object.freeze(['x', 'y', 'z']) as SpatialDimensions;

export const ensuredDims = <T>(
defaultValue: T,
Expand Down Expand Up @@ -30,7 +31,7 @@ export const orderBy =
}),
);

export const chunk = <T>(chunkSize: number, array: Array<T>) => {
export const chunk = <T>(chunkSize: number, array: ReadonlyArray<T>) => {
const chunks = [];
for (let i = 0; i < array.length; i += chunkSize) {
chunks.push(array.slice(i, i + chunkSize));
Expand All @@ -41,5 +42,3 @@ export const chunk = <T>(chunkSize: number, array: Array<T>) => {
export const nonNullable = <T>(value: T): value is NonNullable<T> => {
return value != null;
};

export const XYZ = Object.freeze(['x', 'y', 'z']) as SpatialDimensions;
2 changes: 1 addition & 1 deletion packages/io/src/types.ts
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ export type TypedArrayConstructor = Constructor<TypedArray>;

export type SpatialDimensions = ['x'] | ['x', 'y'] | ['x', 'y', 'z'];

type ChunkParameter = Map<Dimension, number>;
export type ChunkParameter = Map<Dimension, number>;
type LazyCoords = {
get: (dim: Dimension) => Awaitable<Float32Array | undefined>;
has: (dim: Dimension) => boolean;
Expand Down

0 comments on commit d59023e

Please sign in to comment.