diff --git a/.github/workflows/python-ci.yml b/.github/workflows/python-ci.yml index 2111a8b..4f50b24 100644 --- a/.github/workflows/python-ci.yml +++ b/.github/workflows/python-ci.yml @@ -78,7 +78,10 @@ jobs: uses: actions/cache@v5 with: path: test/test_data/PV360_StdData - key: pv360-stddata-6f1b67e5dbc3d7b3646a6315959ccf6d4bd02237 + # The suffix tracks which LFS media conftest smudges. Bumping it when + # that set grows rebuilds the cache once, instead of re-downloading the + # newly required files on every run, because a cache hit is not re-saved. + key: pv360-stddata-6f1b67e5dbc3d7b3646a6315959ccf6d4bd02237-nifti - name: Run all dataset tests run: | diff --git a/test/conftest.py b/test/conftest.py index 67f4114..f2e881b 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -98,7 +98,10 @@ def _resolve_requested_datasets(opt: str | None): def _is_required_github_data_file(path: Path): - return path.name in {"2dseq", "traj"} or path.name.startswith("rawdata.job") + # The NIfTI exports are the independent geometry oracle of test_geometry.py. + # Leaving them as LFS pointers made that check skip itself in CI, so it was + # green whatever it would have found. + return path.name in {"2dseq", "traj"} or path.name.startswith("rawdata.job") or path.suffix == ".nii" def _is_git_lfs_pointer(path: Path): diff --git a/test/test_geometry.py b/test/test_geometry.py index 97c978b..bc00ef8 100644 --- a/test/test_geometry.py +++ b/test/test_geometry.py @@ -17,10 +17,20 @@ SAGITTAL_ORIENTATION = np.array([0.0, 1.0, 0.0, 0.0, 0.0, 1.0, -1.0, 0.0, 0.0]) PV360_NIFTI_ROOT = Path("test/test_data/PV360_StdData") +# Every reconstruction that ships a NIfTI export, rather than a fixed list, so +# a dataset added to the corpus is checked without editing this file. +PV360_NIFTI_EXPORTS = sorted( + directory.parent.relative_to(PV360_NIFTI_ROOT).as_posix() for directory in PV360_NIFTI_ROOT.glob("*/pdata/*/nifti") if any(directory.glob("*.nii")) +) + +def nifti_affine(path): + """Read the NIfTI-1 image affine without making nibabel a test dependency. -def nifti_sform(path): - """Read the NIfTI-1 sform without making nibabel a test dependency.""" + The standard stores it two independent ways, and a file need only carry one. + ParaVision writes the quaternion (`qform_code > 0`) and leaves the sform + zeroed, so a reader that consults the sform alone finds nothing at all. + """ header = path.read_bytes()[:348] if header.startswith(b"version https://git-lfs.github.com/spec/"): pytest.skip("ParaVision NIfTI exports are Git LFS pointers; fetch the LFS assets to validate their affines") @@ -34,10 +44,37 @@ def nifti_sform(path): else: pytest.fail(f"{path} is not a NIfTI-1 image") - sform_code = struct.unpack_from(f"{byte_order}h", header, 254)[0] - assert sform_code > 0, f"{path} has no sform" - sform = np.asarray(struct.unpack_from(f"{byte_order}12f", header, 280)).reshape(3, 4) - return np.vstack((sform, (0.0, 0.0, 0.0, 1.0))) + qform_code, sform_code = struct.unpack_from(f"{byte_order}2h", header, 252) + if sform_code > 0: + sform = np.asarray(struct.unpack_from(f"{byte_order}12f", header, 280)).reshape(3, 4) + return np.vstack((sform, (0.0, 0.0, 0.0, 1.0))) + if qform_code > 0: + return nifti_qform_affine(header, byte_order) + pytest.fail(f"{path} declares neither an sform nor a qform") + + +def nifti_qform_affine(header, byte_order): + """NIfTI-1 "method 2": the affine as a rotation quaternion and a scaling. + + Only the vector part of the unit quaternion is stored; `pixdim[0]` is the + handedness factor applied to the third column. + """ + pixdim = struct.unpack_from(f"{byte_order}8f", header, 76) + handedness = pixdim[0] if pixdim[0] else 1.0 + b, c, d, *origin = struct.unpack_from(f"{byte_order}6f", header, 256) + a = np.sqrt(max(0.0, 1.0 - (b * b + c * c + d * d))) + rotation = np.array( + [ + [a * a + b * b - c * c - d * d, 2 * (b * c - a * d), 2 * (b * d + a * c)], + [2 * (b * c + a * d), a * a + c * c - b * b - d * d, 2 * (c * d - a * b)], + [2 * (b * d - a * c), 2 * (c * d + a * b), a * a + d * d - b * b - c * c], + ] + ) + + affine = np.eye(4) + affine[:3, :3] = rotation @ np.diag((pixdim[1], pixdim[2], handedness * pixdim[3])) + affine[:3, 3] = origin + return affine def load(tmp_path, **kwargs): @@ -135,29 +172,24 @@ def test_report_carries_the_affine_and_omits_it_for_spectroscopy(tmp_path): assert "affine" not in spectroscopy.to_dict() -@pytest.mark.parametrize( - "relative_path", - [ - "T1_FLASH/pdata/1", - "T1_FLASH_3D_iso/pdata/1", - "T1_RARE/pdata/1", - "T2_TurboRARE/pdata/1", - "UTE3D/pdata/1", - ], -) -def test_affine_agrees_with_the_paravision_nifti_sform(relative_path): +@pytest.mark.parametrize("relative_path", PV360_NIFTI_EXPORTS) +def test_affine_agrees_with_the_paravision_nifti_export(relative_path): """PV's NIfTI export is an independent geometry oracle. - VisuCore geometry is in the DICOM patient frame, whereas NIfTI sforms are + VisuCore geometry is in the DICOM patient frame, whereas a NIfTI affine is in RAS. Therefore the equivalent NIfTI affine negates the patient x/y axes (spec 12). + + Every volume of an export shares one geometry -- a diffusion series writes + one file per direction -- so all of them are checked, which also pins the + export as internally consistent. """ reconstruction = PV360_NIFTI_ROOT / relative_path - nifti_path = next((reconstruction / "nifti").glob("*.nii")) dataset = Dataset(reconstruction / "2dseq", load=LOAD_STAGES["properties"]) - expected_sform = np.diag((-1.0, -1.0, 1.0, 1.0)) @ dataset.affine - assert np.allclose(nifti_sform(nifti_path), expected_sform, atol=1e-5) + expected = np.diag((-1.0, -1.0, 1.0, 1.0)) @ dataset.affine + for nifti_path in sorted((reconstruction / "nifti").glob("*.nii")): + assert np.allclose(nifti_affine(nifti_path), expected, atol=1e-5), nifti_path.name def test_slice_spacing_is_centre_to_centre_not_thickness_plus_distance(tmp_path):