diff --git a/changelog.d/20251124_080614_markiewicz_bad_qform.md b/changelog.d/20251124_080614_markiewicz_bad_qform.md new file mode 100644 index 00000000..2cb16873 --- /dev/null +++ b/changelog.d/20251124_080614_markiewicz_bad_qform.md @@ -0,0 +1,5 @@ +### Fixed + +- NIfTI files with bad qform matrices, resulting from non-normalized quaternions, + would previously raise a NIFTI_HEADER_UNREADABLE error. Now only the axis codes + are disabled, preventing orientation checks, but not raising errors. diff --git a/src/files/nifti.test.ts b/src/files/nifti.test.ts index a9c44eb0..dcefacdd 100644 --- a/src/files/nifti.test.ts +++ b/src/files/nifti.test.ts @@ -96,6 +96,10 @@ Deno.test('Test extracting axis codes', async (t) => { const affine = [[0, 0, 1, 0], [1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1]] assertEquals(axisCodes(affine), ['A', 'S', 'R']) }) + await t.step('Fail gracefully on NaNs', async () => { + const affine = [[Number.NaN, 0, 1, 0], [1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1]] + assertEquals(axisCodes(affine), null) + }) }) testAsyncFileAccess('Test file access errors for loadHeader', loadHeader) diff --git a/src/files/nifti.ts b/src/files/nifti.ts index ed734191..5c66e051 100644 --- a/src/files/nifti.ts +++ b/src/files/nifti.ts @@ -34,9 +34,9 @@ async function extract(buffer: Uint8Array, nbytes: number): Promise { const buf = await readBytes(file, 1024) + let header try { const data = isCompressed(buf.buffer) ? await extract(buf, 540) : buf.slice(0, 540) - let header if (isNIFTI1(data.buffer)) { header = new NIFTI1() // Truncate to 348 bytes to avoid attempting to parse extensions @@ -48,29 +48,30 @@ export async function loadHeader(file: BIDSFile): Promise { if (!header) { throw { code: 'NIFTI_HEADER_UNREADABLE' } } - const ndim = header.dims[0] - return { - dim: header.dims, - // Hack: round pixdim to 3 decimal places; schema should add rounding function - pixdim: header.pixDims.map((pixdim) => Math.round(pixdim * 1000) / 1000), - shape: header.dims.slice(1, ndim + 1), - voxel_sizes: header.pixDims.slice(1, ndim + 1), - dim_info: { - freq: header.dim_info & 0x03, - phase: (header.dim_info >> 2) & 0x03, - slice: (header.dim_info >> 4) & 0x03, - }, - xyzt_units: { - xyz: ['unknown', 'meter', 'mm', 'um'][header.xyzt_units & 0x03], - t: ['unknown', 'sec', 'msec', 'usec'][(header.xyzt_units >> 3) & 0x03], - }, - qform_code: header.qform_code, - sform_code: header.sform_code, - axis_codes: axisCodes(header.affine), - } as NiftiHeader } catch (err) { throw { code: 'NIFTI_HEADER_UNREADABLE' } } + + const ndim = header.dims[0] + return { + dim: header.dims, + // Hack: round pixdim to 3 decimal places; schema should add rounding function + pixdim: header.pixDims.map((pixdim) => Math.round(pixdim * 1000) / 1000), + shape: header.dims.slice(1, ndim + 1), + voxel_sizes: header.pixDims.slice(1, ndim + 1), + dim_info: { + freq: header.dim_info & 0x03, + phase: (header.dim_info >> 2) & 0x03, + slice: (header.dim_info >> 4) & 0x03, + }, + xyzt_units: { + xyz: ['unknown', 'meter', 'mm', 'um'][header.xyzt_units & 0x03], + t: ['unknown', 'sec', 'msec', 'usec'][(header.xyzt_units >> 3) & 0x03], + }, + qform_code: header.qform_code, + sform_code: header.sform_code, + axis_codes: axisCodes(header.affine), + } as NiftiHeader } /** Vector addition */ @@ -124,13 +125,18 @@ function argMax(arr: number[]): number { * * @returns character codes describing the orientation of an image affine. */ -export function axisCodes(affine: number[][]): string[] { +export function axisCodes(affine: number[][]): string[] | null { // This function is an extract of the Python function transforms3d.affines.decompose44 // (https://github.com/matthew-brett/transforms3d/blob/6a43a98/transforms3d/affines.py#L10-L153) // // As an optimization, this only orthogonalizes the basis, // and does not normalize to unit vectors. + // Bad qforms result in NaNs in the rotation matrix + if (affine.some((row) => row.some((val) => !Number.isFinite(val)))) { + return null + } + // Operate on columns, which are the cosines that project input coordinates onto output axes const [cosX, cosY, cosZ] = [0, 1, 2].map((j) => [0, 1, 2].map((i) => affine[i][j])) @@ -148,7 +154,7 @@ export function axisCodes(affine: number[][]): string[] { // Check that indices are 0, 1 and 2 in some order if (maxIndices.toSorted().some((idx, i) => idx !== i)) { - throw { key: 'AMBIGUOUS_AFFINE' } + throw { code: 'AMBIGUOUS_AFFINE' } } // Positive/negative codes for each world axis