From d439e1319eed24751167f0240b3d3ab8dc264a40 Mon Sep 17 00:00:00 2001 From: "Gabriel A. Devenyi" Date: Sat, 8 Aug 2026 16:33:05 -0400 Subject: [PATCH 1/2] Apply ruff format to the repository `ruff format --check` reported twelve files as unformatted, so any change touching one of them dragged unrelated reflow into its diff. Run the formatter once, on its own, so subsequent changes show only what they actually change. Formatting only: no behaviour changes, and `ruff check` is clean before and after. The suite is unchanged at 2140 passed, 12 skipped. Co-Authored-By: Claude Opus 5 (1M context) Claude-Session: https://claude.ai/code/session_017SNrm3jFhpTPShP8xGkePC --- brukerapi/dataset.py | 32 +++----- brukerapi/schemas.py | 72 +++++------------ examples/read_2dseq.ipynb | 4 +- test/conftest.py | 3 +- test/loadtest.py | 90 +++++++++++++++++++++ test/test_dataset.py | 16 +--- test/test_geometry.py | 4 +- test/test_jcampdx.py | 135 ++++---------------------------- test/test_latent_conformance.py | 18 +---- test/test_paths.py | 5 +- 10 files changed, 142 insertions(+), 237 deletions(-) create mode 100644 test/loadtest.py diff --git a/brukerapi/dataset.py b/brukerapi/dataset.py index 49d8fed..435c2e9 100644 --- a/brukerapi/dataset.py +++ b/brukerapi/dataset.py @@ -245,11 +245,7 @@ def __init__(self, path, **state): (name for name in content if re.fullmatch(r"rawdata\.job\d+", name)), key=lambda name: int(name.rsplit("job", 1)[1]), ) - named_rawdata_jobs = sorted( - name - for name in content - if name not in rawdata_jobs and self.is_supported_path(self.path / name, {"rawdata"}) - ) + named_rawdata_jobs = sorted(name for name in content if name not in rawdata_jobs and self.is_supported_path(self.path / name, {"rawdata"})) if "fid" in content: self.path = self.path / "fid" elif "2dseq" in content: @@ -273,10 +269,7 @@ def __init__(self, path, **state): if self.type not in DEFAULT_STATES: raise UnsupportedDatasetType(self.type) if self.type == "fid" and self.subtype in FID_COMPANION_SUBTYPES and not state.get("_auxiliary"): - raise UnsupportedDatasetType( - f"{self.path.name} is a fid companion (spec §3.5); load it as " - f"Dataset({with_suffix(self.path, '')!s}).fid_companions[{self.subtype!r}]" - ) + raise UnsupportedDatasetType(f"{self.path.name} is a fid companion (spec §3.5); load it as Dataset({with_suffix(self.path, '')!s}).fid_companions[{self.subtype!r}]") if not self._is_supported_subtype() and not (state.get("_auxiliary") and self.type == "fid" and self.subtype in FID_COMPANION_SUBTYPES): raise UnsupportedDatasetType(self.path.name) @@ -596,9 +589,9 @@ def get(self, name, default=None): A name no configuration declares still raises, so a typo does not quietly become `default`:: - dataset.get("TE") # 4.0, or None where there is no echo time - dataset.get("TE", "n/a") # 4.0, or "n/a" - dataset.get("TEE") # AttributeError + dataset.get("TE") # 4.0, or None where there is no echo time + dataset.get("TE", "n/a") # 4.0, or "n/a" + dataset.get("TEE") # AttributeError Parameters are a separate layer, read with ``dataset[key]`` or ``dataset.parameters``. @@ -727,8 +720,7 @@ def rawdata_stored_scans(self): settings_stored = int(settings[9]) if settings_stored != job_stored and not getattr(self, "_warned_rawdata_stored_scans", False): warnings.warn( - f"ACQ_ScanPipeJobSettings nStoredScans={settings_stored} disagrees with " - f"ACQ_jobs nStoredScans={job_stored} for {self.path}; using the settings value", + f"ACQ_ScanPipeJobSettings nStoredScans={settings_stored} disagrees with ACQ_jobs nStoredScans={job_stored} for {self.path}; using the settings value", RuntimeWarning, stacklevel=2, ) @@ -779,9 +771,7 @@ def _infer_scheme_id(self): if descriptions and descriptions[0] == "Spectroscopic": return "SPECTROSCOPY" if int(dim) == 1 else "CSI" - if n_projections is not None and int(n_projections) > 0 and ( - (self.path.parent / "traj").exists() or "RADIAL" in family or "UTE" in family - ): + if n_projections is not None and int(n_projections) > 0 and ((self.path.parent / "traj").exists() or "RADIAL" in family or "UTE" in family): return "RADIAL" acq_size = np.atleast_1d(self._parameter_value("ACQ_size", [])) @@ -940,9 +930,7 @@ def _read_binary_file(self, path, dtype, shape): f"expected {expected_bytes} bytes for shape {tuple(shape)} and dtype {dtype}" ) - raise InvalidDataset( - f"Invalid dataset size at {path}: expected {expected_bytes} bytes for shape {tuple(shape)} and dtype {dtype}, got {actual_bytes} bytes" - ) + raise InvalidDataset(f"Invalid dataset size at {path}: expected {expected_bytes} bytes for shape {tuple(shape)} and dtype {dtype}, got {actual_bytes} bytes") return read_array(path, dtype, shape) @@ -1343,9 +1331,7 @@ def affine(self): """ descriptors = np.atleast_1d(np.asarray(self._parameter_value("VisuCoreDimDesc", []))).astype(str) if descriptors.size and any(descriptor != "spatial" for descriptor in descriptors): - raise UnsupportedDatasetType( - f"an image affine for {self.path}, whose frames are {sorted(set(descriptors))} rather than purely spatial (spec 7.2)," - ) + raise UnsupportedDatasetType(f"an image affine for {self.path}, whose frames are {sorted(set(descriptors))} rather than purely spatial (spec 7.2),") if len(self.slice_packages_index()) > 1: warnings.warn( f"{self.path} has multiple slice packages; a single affine cannot describe them -- use get_slice_packages() / affine_of_package(i)", diff --git a/brukerapi/schemas.py b/brukerapi/schemas.py index d851b73..c3d0846 100644 --- a/brukerapi/schemas.py +++ b/brukerapi/schemas.py @@ -51,7 +51,7 @@ "k_space", "encoded_dim", "shape_storage", - "dim_type" + "dim_type", ], "2dseq": [ "pv_version", @@ -223,9 +223,7 @@ def layouts(self): for name in ("encoding_space", "k_space"): logical_samples = int(np.prod(layouts[name])) if stored_samples % logical_samples: - raise InvalidDataset( - f"real-only AQ_mod=qf sample count {stored_samples} is incompatible with {name} layout {layouts[name]}" - ) + raise InvalidDataset(f"real-only AQ_mod=qf sample count {stored_samples} is incompatible with {name} layout {layouts[name]}") ratio = stored_samples // logical_samples if ratio > 1: layouts[name] = (layouts[name][0] * ratio,) + tuple(layouts[name][1:]) @@ -282,10 +280,7 @@ def raw(self): data = self._decode_raw_stream(stored, self.layouts) receivers = int(self._dataset.channels) if data.shape[0] % receivers: - raise InvalidDataset( - f"decoded FID sample count {data.shape[0]} is not divisible by " - f"the receiver count {receivers} for {self._dataset.path}" - ) + raise InvalidDataset(f"decoded FID sample count {data.shape[0]} is not divisible by the receiver count {receivers} for {self._dataset.path}") samples = data.shape[0] // receivers return np.transpose( np.reshape(data, (samples, receivers, data.shape[1]), order="F"), @@ -302,9 +297,7 @@ def to_kspace(self, data=None, *, bart=False): try: axes = tuple(BART_DIM_BY_TYPE[label] for label in self._dataset.dim_type) except KeyError as error: - raise UnknownAcqSchemeException( - f"cannot map FID axis {error.args[0]!r} to BART for {self._dataset.path}" - ) from error + raise UnknownAcqSchemeException(f"cannot map FID axis {error.args[0]!r} to BART for {self._dataset.path}") from error return self._as_bart(data, axes) def _reorder_objects(self, data, dir="FW"): @@ -391,9 +384,7 @@ def _reorder_fid_lines(self, data, dir="FW"): PVM_EncSteps1_sorted = self.permutation_inverse(PVM_EncSteps1_sorted) if data.shape[1] != len(PVM_EncSteps1_sorted): - raise InvalidDataset( - f"phase-encode reorder length {len(PVM_EncSteps1_sorted)} does not match k-space axis length {data.shape[1]} for scheme {self._dataset.scheme_id}" - ) + raise InvalidDataset(f"phase-encode reorder length {len(PVM_EncSteps1_sorted)} does not match k-space axis length {data.shape[1]} for scheme {self._dataset.scheme_id}") if np.array_equal(PVM_EncSteps1_sorted, np.arange(len(PVM_EncSteps1_sorted))): return data @@ -712,17 +703,13 @@ def to_kspace(self, data=None, *, bart=False): scheme_id = self._dataset._infer_scheme_id() if scheme_id is not None: raise UnknownAcqSchemeException( - f"rawdata-to-k-space is currently supported only for Cartesian PV-360 jobs, " - f"but {self._dataset.path} is {scheme_id}; use its acquisition-specific reader" + f"rawdata-to-k-space is currently supported only for Cartesian PV-360 jobs, but {self._dataset.path} is {scheme_id}; use its acquisition-specific reader" ) encoded_dim = self._dataset._parameter_value("ACQ_dim") matrix = np.atleast_1d(self._dataset._parameter_value("PVM_EncMatrix", [])) if encoded_dim not in (2, 3) or matrix.size < encoded_dim: - raise UnknownAcqSchemeException( - f"cannot establish a Cartesian rawdata layout for {self._dataset.path}; " - "pass data through an acquisition-specific reader" - ) + raise UnknownAcqSchemeException(f"cannot establish a Cartesian rawdata layout for {self._dataset.path}; pass data through an acquisition-specific reader") matrix = tuple(int(value) for value in matrix[:encoded_dim]) receivers = int(self._dataset.channels) @@ -747,22 +734,15 @@ def to_kspace(self, data=None, *, bart=False): permute = (0, 2, 3, 4, 1) else: if objects != 1: - raise UnknownAcqSchemeException( - f"cannot establish a 3-D Cartesian rawdata layout with NI={objects} for " - f"{self._dataset.path}" - ) + raise UnknownAcqSchemeException(f"cannot establish a 3-D Cartesian rawdata layout with NI={objects} for {self._dataset.path}") encoding_space = (readout, receivers, phase, matrix[2], repetitions) permute = (0, 2, 3, 4, 1) if data.ndim != 3 or data.shape[0] != readout or data.shape[1] != receivers: - raise InvalidDataset( - f"rawdata sample layout {data.shape} does not match Cartesian metadata " - f"(readout={readout}, receivers={receivers}) for {self._dataset.path}" - ) + raise InvalidDataset(f"rawdata sample layout {data.shape} does not match Cartesian metadata (readout={readout}, receivers={receivers}) for {self._dataset.path}") if data.size != int(np.prod(encoding_space)): raise InvalidDataset( - f"rawdata contains {data.size} complex samples but Cartesian layout {encoding_space} " - f"requires {int(np.prod(encoding_space))} for {self._dataset.path}" + f"rawdata contains {data.size} complex samples but Cartesian layout {encoding_space} requires {int(np.prod(encoding_space))} for {self._dataset.path}" ) k_space = np.transpose(np.reshape(data, encoding_space, order="F"), permute) @@ -779,10 +759,7 @@ def _self_gated_k_space(self, data, readout, phase, *, receivers, objects, repet """Arrange retrospectively gated Cartesian data before cine binning.""" steps = np.atleast_1d(self._dataset._parameter_value("PVM_EncGenSteps1", [])) if steps.size == 0 or steps.size % phase: - raise InvalidDataset( - f"self-gated phase-encode sequence has {steps.size} steps, which is incompatible " - f"with phase size {phase} for {self._dataset.path}" - ) + raise InvalidDataset(f"self-gated phase-encode sequence has {steps.size} steps, which is incompatible with phase size {phase} for {self._dataset.path}") acquired_frames = steps.size // phase output_frames = int(self._dataset._parameter_value("PVM_NMovieFrames", objects)) if output_frames != objects or acquired_frames % (objects * repetitions): @@ -793,10 +770,7 @@ def _self_gated_k_space(self, data, readout, phase, *, receivers, objects, repet acquisition_cycles = acquired_frames // (objects * repetitions) layout = (readout, receivers, phase, acquired_frames) if data.ndim != 3 or data.shape[0] != readout or data.shape[1] != receivers or data.size != int(np.prod(layout)): - raise InvalidDataset( - f"rawdata sample layout {data.shape} does not match self-gated Cartesian layout " - f"{layout} for {self._dataset.path}" - ) + raise InvalidDataset(f"rawdata sample layout {data.shape} does not match self-gated Cartesian layout {layout} for {self._dataset.path}") k_space = np.transpose(np.reshape(data, layout, order="F"), (0, 2, 3, 1)) order = np.argsort(np.reshape(steps, (phase, acquired_frames), order="F"), axis=0) @@ -813,10 +787,7 @@ def _reorder_phase_lines(self, data): return data indices = np.argsort(np.atleast_1d(steps)) if indices.size != data.shape[1]: - raise InvalidDataset( - f"phase-encode reorder length {indices.size} does not match k-space axis length " - f"{data.shape[1]} for {self._dataset.path}" - ) + raise InvalidDataset(f"phase-encode reorder length {indices.size} does not match k-space axis length {data.shape[1]} for {self._dataset.path}") return np.take(data, indices, axis=1) def _reorder_objects(self, data): @@ -827,12 +798,10 @@ def _reorder_objects(self, data): return data indices = np.argsort(np.atleast_1d(order).astype(int)) if indices.size != data.shape[2]: - raise InvalidDataset( - f"object-order length {indices.size} does not match k-space axis length " - f"{data.shape[2]} for {self._dataset.path}" - ) + raise InvalidDataset(f"object-order length {indices.size} does not match k-space axis length {data.shape[2]} for {self._dataset.path}") return np.take(data, indices, axis=2) + class Schema2dseq(Schema): """ Schema2dseq class @@ -948,8 +917,7 @@ def _apply_disk_slice_order(self, data): axis = 2 if axis is None or axis >= data.ndim: warnings.warn( - "VisuCoreDiskSliceOrder requests reversed slices but no slice axis is identifiable " - f"for {getattr(self._dataset, 'path', '')}; leaving order unchanged", + f"VisuCoreDiskSliceOrder requests reversed slices but no slice axis is identifiable for {getattr(self._dataset, 'path', '')}; leaving order unchanged", RuntimeWarning, stacklevel=2, ) @@ -974,9 +942,7 @@ def _complex_frame_axis(self, data): # component: there is nothing to combine, keep the axis as it is return None if data.shape[axis] != 2: - raise InvalidDataset( - f"complex 2dseq requires a two-element real/imag frame-group axis, got shape {data.shape} on axis {axis}" - ) + raise InvalidDataset(f"complex 2dseq requires a two-element real/imag frame-group axis, got shape {data.shape} on axis {axis}") return axis def _combine_complex_frames(self, data): @@ -1100,9 +1066,7 @@ def ra(self, slice_): # Frame-group selection above does not alter encoded dimensions. # Apply their requested selection after deserializing the selected # frames, preserving singleton axes for the final squeeze below. - encoded_slice = tuple( - slice(index, index + 1) if isinstance(index, int) else index for index in slice_[: self._dataset.encoded_dim] - ) + encoded_slice = tuple(slice(index, index + 1) if isinstance(index, int) else index for index in slice_[: self._dataset.encoded_dim]) array_ra = array_ra[encoded_slice + (slice(None),) * (array_ra.ndim - self._dataset.encoded_dim)] singletons = tuple(i for i, v in enumerate(slice_) if isinstance(v, int)) diff --git a/examples/read_2dseq.ipynb b/examples/read_2dseq.ipynb index 50827bf..e42ab1a 100644 --- a/examples/read_2dseq.ipynb +++ b/examples/read_2dseq.ipynb @@ -259,8 +259,8 @@ "print(dataset.TE)\n", "print(dataset.TR)\n", "print(dataset.imaging_frequency)\n", - "print(dataset.frame_group_values.get('VisuAcqEchoTime'))\n", - "print(dataset.metadata['visu_acq'].get('sequence_name'))" + "print(dataset.frame_group_values.get(\"VisuAcqEchoTime\"))\n", + "print(dataset.metadata[\"visu_acq\"].get(\"sequence_name\"))" ] } ], diff --git a/test/conftest.py b/test/conftest.py index f2e881b..ac404dc 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -78,8 +78,7 @@ def pytest_sessionstart(session): if _available_test_data_root() is None: pytest.exit( - "No Bruker test-data corpus is available. Provide test/test_data or resources/testdata, " - "or run pytest with --download_test_data.", + "No Bruker test-data corpus is available. Provide test/test_data or resources/testdata, or run pytest with --download_test_data.", returncode=1, ) diff --git a/test/loadtest.py b/test/loadtest.py new file mode 100644 index 0000000..4ab4baa --- /dev/null +++ b/test/loadtest.py @@ -0,0 +1,90 @@ +#!/usr/bin/env python +"""Empirical load-test harness for brukerapi against resources/testdata. + +Walks the testdata tree, finds every supported binary file (fid, 2dseq, rawdata.jobN, +traj, ser), attempts to construct a brukerapi.Dataset for it, and records the outcome: +success (+shape/dtype) or failure (+exception type + short message). +""" + +import json +import os +import sys +import traceback +from pathlib import Path + +from brukerapi.dataset import Dataset + +ROOT = Path("resources/testdata").resolve() + +# Binary stems brukerapi claims to support (see DEFAULT_STATES) +SUPPORTED_STEMS = {"fid", "2dseq", "rawdata", "traj", "ser"} + + +def classify(p: Path): + stem = p.name.split(".")[0] + return stem + + +def find_targets(root: Path): + targets = [] + for dirpath, dirnames, filenames in os.walk(root): + # skip git internals and provenance + parts = set(Path(dirpath).parts) + if ".git" in parts or "_sources" in parts: + dirnames[:] = [d for d in dirnames if d != ".git"] + continue + for fn in filenames: + stem = fn.split(".")[0] + if stem in SUPPORTED_STEMS: + targets.append(Path(dirpath) / fn) + return sorted(targets) + + +def try_load(path: Path): + rec = {"path": str(path.relative_to(ROOT)), "stem": classify(path), "size": path.stat().st_size} + try: + ds = Dataset(str(path)) + rec["ok"] = True + try: + rec["shape"] = list(ds.data.shape) + rec["dtype"] = str(ds.data.dtype) + except Exception as e: + rec["ok"] = False + rec["stage"] = "data-access" + rec["err"] = f"{type(e).__name__}: {e}" + return rec + except Exception as e: + rec["ok"] = False + rec["err_type"] = type(e).__name__ + rec["err"] = str(e)[:300] + rec["tb"] = traceback.format_exc()[-1500:] + return rec + + +def main(): + targets = find_targets(ROOT) + results = [try_load(p) for p in targets] + + ok = [r for r in results if r.get("ok")] + bad = [r for r in results if not r.get("ok")] + + print(f"TOTAL targets: {len(results)} OK: {len(ok)} FAIL: {len(bad)}") + print("\n=== FAILURE COUNTS BY (stem, err_type) ===") + from collections import Counter + + c = Counter((r["stem"], r.get("err_type", r.get("stage", "?"))) for r in bad) + for (stem, et), n in c.most_common(): + print(f" {n:4d} {stem:8s} {et}") + + print("\n=== OK COUNTS BY stem ===") + c2 = Counter(r["stem"] for r in ok) + for stem, n in c2.most_common(): + print(f" {n:4d} {stem}") + + out = Path(sys.argv[1]) if len(sys.argv) > 1 else Path("/tmp/loadtest_results.json") + out.write_text(json.dumps(results, indent=2)) + print(f"\nFull results -> {out}") + + +if __name__ == "__main__": + main() diff --git a/test/test_dataset.py b/test/test_dataset.py index 5f85317..c011cf3 100644 --- a/test/test_dataset.py +++ b/test/test_dataset.py @@ -768,10 +768,7 @@ def test_recipe_substitution_preserves_overlapping_identifiers(): substituted = dataset._sub_parameters(recipe) - assert substituted == ( - "self.Foo + self.FooBar + self['X'].value + " - "self['XY'].value + self['Matrix'].tuple" - ) + assert substituted == ("self.Foo + self.FooBar + self['X'].value + self['XY'].value + self['Matrix'].tuple") @pytest.mark.parametrize( @@ -1033,10 +1030,7 @@ def test_frame_group_values_align_echo_times_and_diffusion_matrices_to_data_axes continue candidate = Dataset(path) try: - if ( - "VisuAcqEchoTime" in candidate.frame_group_values - and np.atleast_1d(candidate["VisuAcqEchoTime"].value).size > 1 - ): + if "VisuAcqEchoTime" in candidate.frame_group_values and np.atleast_1d(candidate["VisuAcqEchoTime"].value).size > 1: echo = candidate break except KeyError: @@ -1123,11 +1117,7 @@ def test_data_load(test_data): reference = np.transpose(reference, (0, 2, 1)) if dataset.type == "2dseq" and np.iscomplexobj(actual) and not np.iscomplexobj(reference): complex_axis = next( - ( - axis - for axis, dim_type in enumerate(dataset.dim_type) - if str(dim_type).upper() == "FG_COMPLEX" - ), + (axis for axis, dim_type in enumerate(dataset.dim_type) if str(dim_type).upper() == "FG_COMPLEX"), None, ) if complex_axis is None: diff --git a/test/test_geometry.py b/test/test_geometry.py index edb2db3..f86d4ae 100644 --- a/test/test_geometry.py +++ b/test/test_geometry.py @@ -19,9 +19,7 @@ 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")) -) +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): diff --git a/test/test_jcampdx.py b/test/test_jcampdx.py index 3717fca..f0fe95e 100644 --- a/test/test_jcampdx.py +++ b/test/test_jcampdx.py @@ -50,13 +50,7 @@ def test_jcampdx(test_jcampdx_data): def test_jcampdx_iteration_and_length_follow_its_mapping_interface(tmp_path): path = tmp_path / "visu_pars" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$VisuCoreDim=2\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$VisuCoreDim=2\n##END=\n") jcamp = JCAMPDX(path) assert list(jcamp) == list(jcamp.keys()) @@ -142,13 +136,7 @@ def test_parallel_lists_do_not_split_a_string_on_a_delimiter_it_contains(): def test_jcampdx_get_value_keeps_the_whole_enum_display_name(tmp_path): path = tmp_path / "configscan" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$CONFIG_SCAN_operation_mode=(operation, <[1H] TX Volume, RX Surface Array>)\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$CONFIG_SCAN_operation_mode=(operation, <[1H] TX Volume, RX Surface Array>)\n##END=\n") assert JCAMPDX(path).get_value("CONFIG_SCAN_operation_mode") == [ "operation", @@ -179,15 +167,7 @@ def test_run_length_expansion_handles_multiple_and_nested_runs(): def test_jcampdx_data_parameter_parses_multiline_xy_pairs(tmp_path): path = tmp_path / "data" - path.write_text( - "##TITLE=XY Data\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$POINTS=(XY..XY)\n" - "1.0, 2.0\n" - "3.0, 4.0\n" - "##END=\n" - ) + path.write_text("##TITLE=XY Data\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$POINTS=(XY..XY)\n1.0, 2.0\n3.0, 4.0\n##END=\n") assert np.array_equal( JCAMPDX(path).get_value("POINTS"), @@ -197,15 +177,7 @@ def test_jcampdx_data_parameter_parses_multiline_xy_pairs(tmp_path): def test_jcampdx_float_and_list_serialization_round_trip(tmp_path): source = tmp_path / "source" - source.write_text( - "##TITLE=Serialization Test\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$FLOAT=0.0\n" - "##$VALUES=( 2 )\n" - "0.0 0.0\n" - "##END=\n" - ) + source.write_text("##TITLE=Serialization Test\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$FLOAT=0.0\n##$VALUES=( 2 )\n0.0 0.0\n##END=\n") jcamp = JCAMPDX(source) jcamp.get_parameter("FLOAT").value = 1.25 @@ -222,15 +194,7 @@ def test_jcampdx_float_and_list_serialization_round_trip(tmp_path): def test_jcampdx_data_parameter_setter_round_trip(tmp_path): source = tmp_path / "source-data" - source.write_text( - "##TITLE=XY Data\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$POINTS=(XY..XY)\n" - "1.0, 2.0\n" - "3.0, 4.0\n" - "##END=\n" - ) + source.write_text("##TITLE=XY Data\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$POINTS=(XY..XY)\n1.0, 2.0\n3.0, 4.0\n##END=\n") jcamp = JCAMPDX(source) expected = np.array([[5.0, 6.0], [7.0, 8.0]]) @@ -337,18 +301,7 @@ def test_parse_value_does_not_treat_unclosed_parenthesis_as_list(): def test_jcampdx_size_parsing_accepts_compact_and_padded_brackets(tmp_path): path = tmp_path / "sizes" - path.write_text( - "##TITLE=Size Test\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$COMPACT=(2)\n" - "1 2\n" - "##$PADDED=( 2 )\n" - "3 4\n" - "##$MATRIX=(2, 3)\n" - "1 2 3 4 5 6\n" - "##END=\n" - ) + path.write_text("##TITLE=Size Test\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$COMPACT=(2)\n1 2\n##$PADDED=( 2 )\n3 4\n##$MATRIX=(2, 3)\n1 2 3 4 5 6\n##END=\n") jcamp = JCAMPDX(path) assert jcamp.get_parameter("COMPACT").size == (2,) @@ -385,12 +338,7 @@ def test_jcampdx_detects_whitespace_padded_supported_versions(tmp_path, header, def test_jcampdx_keeps_double_hash_inside_bracketed_value(tmp_path): path = tmp_path / "configscan" - path.write_text( - "##TITLE=Config Scan\n" - "##JCAMPDX= 5.0\n" - "##$PULPROG=\n" - "##END=\n" - ) + path.write_text("##TITLE=Config Scan\n##JCAMPDX= 5.0\n##$PULPROG=\n##END=\n") assert JCAMPDX(path).get_value("PULPROG") == "HpMode,On##$EndBis,04,FA#" @@ -402,13 +350,7 @@ def test_jcampdx_record_without_assignment_raises_typed_error(): def test_load_parameter_allows_hash_and_dollar_in_value(tmp_path): path = tmp_path / "special-value" - path.write_text( - "##TITLE=Special Value\n" - "##JCAMPDX=5.0\n" - "##$VALUE=\n" - "##$NEXT=2\n" - "##END=\n" - ) + path.write_text("##TITLE=Special Value\n##JCAMPDX=5.0\n##$VALUE=\n##$NEXT=2\n##END=\n") key, parameter = JCAMPDX.load_parameter(path, "VALUE") @@ -444,16 +386,7 @@ def test_jcampdx_round_trip_preserves_comments_and_end_marker(tmp_path): def test_jcampdx_version_detection_is_label_based_within_header(tmp_path): path = tmp_path / "reordered-header" - path.write_text( - "##TITLE=Reordered Header\n" - "##DATATYPE=Parameter Values\n" - "##ORIGIN=Test\n" - "$$ header comment\n" - "##OWNER=Tester\n" - "##JCAMPDX=4.24\n" - "##$VALUE=42\n" - "##END=\n" - ) + path.write_text("##TITLE=Reordered Header\n##DATATYPE=Parameter Values\n##ORIGIN=Test\n$$ header comment\n##OWNER=Tester\n##JCAMPDX=4.24\n##$VALUE=42\n##END=\n") jcamp = JCAMPDX(path) @@ -491,15 +424,7 @@ def test_a_wrap_inside_a_string_does_not_invent_a_space(tmp_path): def test_a_wrap_at_a_space_keeps_exactly_one_space(tmp_path): path = tmp_path / "acqp" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$ACQ_size=( 4 )\n" - "128 64 \n" - "32 16\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$ACQ_size=( 4 )\n128 64 \n32 16\n##END=\n") assert np.array_equal(JCAMPDX(path)["ACQ_size"].value, np.array([128, 64, 32, 16])) @@ -523,16 +448,7 @@ def test_write_does_not_wrap_comment_records(tmp_path): """ comment = "$$ /opt/PV6.0.1/data/imag/20200913_160003_In_situ_experiment_with_a_very_long_name/74/acqp" path = tmp_path / "acqp" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##OWNER=imag\n" - f"{comment}\n" - "##$ACQ_size=( 2 )\n" - "128 64\n" - "##END=\n" - ) + path.write_text(f"##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##OWNER=imag\n{comment}\n##$ACQ_size=( 2 )\n128 64\n##END=\n") original = JCAMPDX(path) original.write(tmp_path / "acqp.written") @@ -654,15 +570,7 @@ def test_a_trailing_backslash_in_a_string_is_content_not_an_escape(tmp_path): would drop the record. """ path = tmp_path / "visu_pars" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$VisuStudyDescription=( 2048 )\n" - "<\\\n" - ">\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$VisuStudyDescription=( 2048 )\n<\\\n>\n##END=\n") assert JCAMPDX(path)["VisuStudyDescription"].value == "\\" @@ -675,13 +583,7 @@ def test_the_last_parameter_survives_a_file_without_an_end_marker(tmp_path): exception and no warning. """ path = tmp_path / "visu_pars" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$VisuCoreDim=2\n" - "##$VisuRespSynchUsed=No\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$VisuCoreDim=2\n##$VisuRespSynchUsed=No\n") parameters = JCAMPDX(path) @@ -693,16 +595,7 @@ def test_a_malformed_record_raises_a_typed_error(tmp_path): """Spec 2.3: an element count that does not fill the declared size is a diagnosable condition, not a raw numpy ValueError.""" path = tmp_path / "acqp" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$SHORT=( 3, 3 )\n" - "1 2 3 4\n" - "##$BADSIZE=( a )\n" - "1 2\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$SHORT=( 3, 3 )\n1 2 3 4\n##$BADSIZE=( a )\n1 2\n##END=\n") parameters = JCAMPDX(path) with pytest.raises(InvalidJcampdxFile, match="do not fill the declared size"): diff --git a/test/test_latent_conformance.py b/test/test_latent_conformance.py index 16c326a..fbe0ad8 100644 --- a/test/test_latent_conformance.py +++ b/test/test_latent_conformance.py @@ -218,15 +218,7 @@ def test_rawdata_settings_control_stored_scans_receivers_and_discarded_jobs(tmp_ def test_a_dollar_comment_inside_a_string_is_data(tmp_path): """Spec 2.2: the text inside `<...>` is free-form, `$$` included.""" path = tmp_path / "method" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$PVM_Comment=( 64 )\n" - "\n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$PVM_Comment=( 64 )\n\n##END=\n") assert JCAMPDX(path)["PVM_Comment"].value == "a$$b" @@ -235,12 +227,6 @@ def test_a_scalar_struct_is_not_eaten_as_a_size_bracket(tmp_path): """`##$VisuCoreSlicePacksDef=(0, 1)` is a value, not a size -- including when the line has a trailing blank.""" path = tmp_path / "visu_pars" - path.write_text( - "##TITLE=Parameter List\n" - "##JCAMPDX=4.24\n" - "##DATATYPE=Parameter Values\n" - "##$VisuCoreSlicePacksDef=(0, 1) \n" - "##END=\n" - ) + path.write_text("##TITLE=Parameter List\n##JCAMPDX=4.24\n##DATATYPE=Parameter Values\n##$VisuCoreSlicePacksDef=(0, 1) \n##END=\n") assert JCAMPDX(path)["VisuCoreSlicePacksDef"].value == [0, 1] diff --git a/test/test_paths.py b/test/test_paths.py index febe765..479e103 100644 --- a/test/test_paths.py +++ b/test/test_paths.py @@ -84,8 +84,7 @@ def test_dataset_reads_from_archive_identically(study_dir, study_zip): def test_parameters_resolve_through_relative_paths(study_dir, study_zip): """``../../acqp`` resolves inside an archive, where ``..`` is not collapsed.""" - dataset = Dataset(study_zip / "1" / "pdata" / "1" / "2dseq", scale=False, - parameter_files=["acqp"]) + dataset = Dataset(study_zip / "1" / "pdata" / "1" / "2dseq", scale=False, parameter_files=["acqp"]) assert dataset["ACQ_scan_name"].value == "demo" @@ -111,7 +110,7 @@ def test_get_value_default_for_absent_key(study_dir): def test_path_helpers_work_for_both_kinds(study_dir, study_zip): assert isinstance(as_path(str(study_dir)), Path) - assert as_path(study_zip) is study_zip # passed through, not coerced + assert as_path(study_zip) is study_zip # passed through, not coerced for root in (study_dir, study_zip): proc = root / "1" / "pdata" / "1" From a158f8c79a851797359572feacd5caf497bbaa0d Mon Sep 17 00:00:00 2001 From: "Gabriel A. Devenyi" Date: Sat, 8 Aug 2026 16:24:51 -0400 Subject: [PATCH 2/2] Nest the object loop inside the phase loop in SchemaRawdata.to_kspace The 2-D branch reshaped the PV360 scan stream as (readout, receivers, phase, objects, repetitions) in Fortran order, which reads the scan index as phase + Nphase*object. Spec 5.2 "Acquisition loop nesting (default)" puts it the other way round: NS > ACQ_phase_factor > NSLICES > NI or NSLICES > NA > ACQ_size[1]/ACQ_phase_factor > ACQ_size[1] > ACQ_size[2] > NAE > NR The object level is INSIDE the phase-encode-group level, so consecutive scans step the object first. The library's own fid path already has this right -- properties_fid_core.json CART_2D encoding_space is (EncMatrix[0], EncNReceivers, ACQ_phase_factor, NI, EncMatrix[1]//ACQ_phase_factor, NR) -- so to_kspace disagreed with both the spec and the rest of the library, and ACQ_phase_factor was missing entirely, mis-nesting RARE-style echo trains too. Nothing raised: the array kept its shape and dtype, so every 2-D multi-slice PV360 job -- all of Bruker's standard FLASH/RARE/MSME/MGE protocols -- came back with its slice and phase-encode indices interleaved. Adopt the fid nesting. On the PV360 standard data the phase profile now peaks at exactly N/2 for every 2-D scan, where before it was flat and off-centre: T1_FLASH 160/288, peak/median 1.26 -> 144/288, 6.46 T1_RARE 170/256, 1.79 -> 128/256, 4.06 T2_TurboRARE 143/256, 3.44 -> 128/256, 10.63 T2map_MSME 187/192, 2.37 -> 96/192, 10.44 and a zero-filled sum-of-squares reconstruction correlates with the vendor 2dseq far better: T1_FLASH 0.07-0.36 -> 0.47-0.56 T1_RARE 0.03-0.17 -> 0.62-0.99 T2_TurboRARE 0.06-0.22 -> 0.61-0.77 The 3-D branch is unchanged; with NI = 1 there is no object level to interleave. Closes #192 Co-Authored-By: Claude Opus 5 (1M context) Claude-Session: https://claude.ai/code/session_017SNrm3jFhpTPShP8xGkePC --- brukerapi/schemas.py | 16 +++++++++++++-- test/test_rawdata.py | 48 ++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 62 insertions(+), 2 deletions(-) diff --git a/brukerapi/schemas.py b/brukerapi/schemas.py index c3d0846..25e0e5a 100644 --- a/brukerapi/schemas.py +++ b/brukerapi/schemas.py @@ -730,13 +730,24 @@ def to_kspace(self, data=None, *, bart=False): return self._as_bart(k_space, axes) if bart else k_space if encoded_dim == 2: - encoding_space = (readout, receivers, phase, objects, repetitions) - permute = (0, 2, 3, 4, 1) + # Spec 5.2, "Acquisition loop nesting (default)": + # NS > ACQ_phase_factor > NSLICES > NI or NSLICES > NA + # > ACQ_size[1]/ACQ_phase_factor > ACQ_size[1] > ACQ_size[2] > NAE > NR + # The object level sits INSIDE the phase-encode-group level, so the scan + # index runs object + NI*group, not phase + Nphase*object. This is the + # same nesting the fid path uses (properties_fid_core.json, CART_2D). + phase_factor = int(self._dataset._parameter_value("ACQ_phase_factor", 1)) or 1 + if phase % phase_factor: + raise UnknownAcqSchemeException(f"ACQ_phase_factor={phase_factor} does not divide the {phase} phase-encode steps of {self._dataset.path}") + encoding_space = (readout, receivers, phase_factor, objects, phase // phase_factor, repetitions) + permute = (0, 2, 4, 3, 5, 1) + k_space_shape = (readout, phase, objects, repetitions, receivers) else: if objects != 1: raise UnknownAcqSchemeException(f"cannot establish a 3-D Cartesian rawdata layout with NI={objects} for {self._dataset.path}") encoding_space = (readout, receivers, phase, matrix[2], repetitions) permute = (0, 2, 3, 4, 1) + k_space_shape = (readout, phase, matrix[2], repetitions, receivers) if data.ndim != 3 or data.shape[0] != readout or data.shape[1] != receivers: raise InvalidDataset(f"rawdata sample layout {data.shape} does not match Cartesian metadata (readout={readout}, receivers={receivers}) for {self._dataset.path}") @@ -746,6 +757,7 @@ def to_kspace(self, data=None, *, bart=False): ) k_space = np.transpose(np.reshape(data, encoding_space, order="F"), permute) + k_space = np.reshape(k_space, k_space_shape, order="F") k_space = self._reorder_phase_lines(k_space) k_space = self._reorder_objects(k_space) axes = ( diff --git a/test/test_rawdata.py b/test/test_rawdata.py index 03fa2ef..9a55fb3 100644 --- a/test/test_rawdata.py +++ b/test/test_rawdata.py @@ -5,6 +5,7 @@ from brukerapi.dataset import Dataset from brukerapi.exceptions import InvalidDataset, UnknownAcqSchemeException +from test.synthetic import Verbatim, write_binary, write_jcampdx RAWDATA_JOB_PATHS = sorted(Path("test/test_data").rglob("rawdata.job*")) @@ -100,3 +101,50 @@ def test_pv360_epi_rawdata_job_requires_an_acquisition_specific_reader(test_data with pytest.raises(UnknownAcqSchemeException, match="is EPI"): dataset.to_kspace() + + +def test_to_kspace_nests_the_object_loop_inside_the_phase_loop(tmp_path): + """Spec 5.2, "Acquisition loop nesting (default)": + + NS > ACQ_phase_factor > NSLICES > NI or NSLICES > NA + > ACQ_size[1]/ACQ_phase_factor > ACQ_size[1] > ACQ_size[2] > NAE > NR + + The object level sits INSIDE the phase-encode-group level, so consecutive + scans step the object first. Reading it the other way round interleaves the + slice and phase-encode indices of every 2-D multi-slice job, silently: the + array keeps its shape and dtype. + """ + readout, phase, objects = 2, 3, 2 + # scan s = object + NI*group, and each scan carries 1000*object + group so a + # misplaced sample is visible + scans = [1000 * obj + group for group in range(phase) for obj in range(objects)] + samples = np.zeros((2 * readout, 1, len(scans)), dtype=""], + "ACQ_word_size": "_32_BIT", + "BYTORDA": "little", + "ACQ_dim": 2, + "ACQ_dim_desc": Verbatim("( 2 )\nSpatial Spatial"), + # spec 13.1: ACQ_size[0] need not equal the job scan size + "ACQ_size": np.array([8, phase]), + "ACQ_phase_factor": 1, + "NI": objects, + "NR": 1, + "ACQ_jobs": Verbatim(f"( 1 )\n({2 * readout}, 1, 0, {len(scans)}, 101, 5000, {len(scans)}, 1, )"), + }, + ) + write_jcampdx(experiment / "method", {"PVM_EncNReceivers": 1, "PVM_EncMatrix": np.array([readout, phase])}) + write_binary(experiment / "rawdata.job0", samples, np.dtype("