diff --git a/brukerapi/dataset.py b/brukerapi/dataset.py index 2d60898..c6da217 100644 --- a/brukerapi/dataset.py +++ b/brukerapi/dataset.py @@ -4,6 +4,7 @@ import os.path import re import warnings +from contextlib import suppress from copy import deepcopy from pathlib import Path @@ -445,7 +446,10 @@ def add_parameter_file(self, file): dataset["PVM_DwDir"].value """ - path = traverse(self.path.parent, RELATIVE_PATHS[self.type][file]) + if file in RELATIVE_PATHS[self.type]: + path = traverse(self.path.parent, RELATIVE_PATHS[self.type][file]) + else: + path = as_path(file) if not hasattr(self, "_parameters") or self._parameters is None: self._parameters = {path.name: JCAMPDX(path)} @@ -474,6 +478,16 @@ def _read_parameters(self): # if jcampdx file is not found, but not required, pass pass + # The vendor's first reconstruction is the conventional default. A + # fid can have several reconstructions, however, so callers can select + # the one whose input-order declaration applies with ``reco_path=``. + reco_path = self._state.get("reco_path") + if reco_path is not None: + self.add_parameter_file(reco_path) + elif self.type in {"fid", "fid_proc"}: + with suppress(FileNotFoundError): + self.add_parameter_file("reco") + def _write_parameters(self, parent): for type_, jcampdx in self._parameters.items(): jcampdx.write(parent / type_) diff --git a/brukerapi/schemas.py b/brukerapi/schemas.py index b831b5a..18f30d0 100644 --- a/brukerapi/schemas.py +++ b/brukerapi/schemas.py @@ -162,6 +162,45 @@ def acquisition_factor(self): """Number of stored scalar samples per logical sample.""" return 1 if str(self._dataset._parameter_value("AQ_mod", "")).lower() == "qf" else 2 + def _declared_or_inferred(self, parameter, declared_value, inferred, warning_flag): + """Prefer a loaded format declaration, retaining scheme inference.""" + if declared_value is None: + return inferred + if declared_value != inferred and not getattr(self._dataset, warning_flag, False): + warnings.warn( + f"{parameter}={self._dataset._parameter_value(parameter)!r} disagrees with " + f"scheme_id={self._dataset.scheme_id!r} for {self._dataset.path}; " + f"using the declared value", + RuntimeWarning, + stacklevel=3, + ) + setattr(self._dataset, warning_flag, True) + return declared_value + + @property + def continuous_train(self): + """Whether §5.3 declares a phase-factor continuous acquisition train.""" + value = self._dataset._parameter_value("ACQ_scan_size") + declared = None if value is None else str(value).strip("<>").casefold() == "acq_phase_factor_scans" + return self._declared_or_inferred( + "ACQ_scan_size", + declared, + "EPI" in self._dataset.scheme_id, + "_warned_scan_size_disagreement", + ) + + @property + def mirror_odd_lines(self): + """Whether a selected reconstruction requests §6.2 odd-line reversal.""" + value = self._dataset._parameter_value("RECO_inp_order") + declared = None if value is None else str(value).strip("<>").upper() == "REV_ALT_ROWS" + return self._declared_or_inferred( + "RECO_inp_order", + declared, + "EPI" in self._dataset.scheme_id, + "_warned_reco_input_order_disagreement", + ) + @property def layouts(self): """Dictionary of possible logical layouts of data @@ -193,7 +232,7 @@ def layouts(self): layouts["encoding_permuted"] = tuple(np.array(layouts["encoding_space"])[np.array(layouts["permute"])]) - if "EPI" in self._dataset.scheme_id: + if self.continuous_train: discarded = self._dataset.block_size - self._dataset.acq_length block_format = self._dataset._parameter_value("GO_block_size") scan_shift = self._dataset._parameter_value("ACQ_scan_shift", 0) @@ -221,7 +260,7 @@ def deserialize(self, data, layouts): # Typically for RARE, or EPI data = self._reorder_fid_lines(data, dir="FW") - if "EPI" in self._dataset.scheme_id: + if self.mirror_odd_lines: data = self._mirror_odd_lines(data) data = self._reorder_objects(data, dir="FW") @@ -385,7 +424,7 @@ def _mirror_odd_lines(self, data): def serialize(self, data, layouts): data = self._reorder_objects(data, dir="BW") - if "EPI" in self._dataset.scheme_id: + if self.mirror_odd_lines: data = self._mirror_odd_lines(data) data = self._reorder_fid_lines(data, dir="BW") diff --git a/test/test_dataset.py b/test/test_dataset.py index aa62dde..afbf232 100644 --- a/test/test_dataset.py +++ b/test/test_dataset.py @@ -691,7 +691,8 @@ def test_monotonic_phase_encode_steps_skip_reordering(): @pytest.mark.skipif(not PV51_STUDY_PATH.is_dir(), reason="PV51 test data is not available") def test_epi_standard_kblock_trims_trailing_padding(): - dataset = Dataset(PV51_STUDY_PATH / "13" / "fid") + with pytest.warns(RuntimeWarning, match="RECO_inp_order=.*disagrees with scheme_id"): + dataset = Dataset(PV51_STUDY_PATH / "13" / "fid") dataset["GO_block_size"].val_str = "Standard_KBlock_Format" dataset.block_size = 256 dataset.acq_length = 200 @@ -706,7 +707,8 @@ def test_epi_standard_kblock_trims_trailing_padding(): @pytest.mark.skipif(not PV51_STUDY_PATH.is_dir(), reason="PV51 test data is not available") def test_epi_standard_kblock_warns_on_nonzero_discarded_samples(): - dataset = Dataset(PV51_STUDY_PATH / "13" / "fid") + with pytest.warns(RuntimeWarning, match="RECO_inp_order=.*disagrees with scheme_id"): + dataset = Dataset(PV51_STUDY_PATH / "13" / "fid") dataset["GO_block_size"].val_str = "Standard_KBlock_Format" dataset.block_size = 8 dataset.acq_length = 6 diff --git a/test/test_raw_layouts.py b/test/test_raw_layouts.py index 89a9fbb..56c78b2 100644 --- a/test/test_raw_layouts.py +++ b/test/test_raw_layouts.py @@ -5,9 +5,10 @@ """ import numpy as np +import pytest from brukerapi.dataset import LOAD_STAGES, Dataset -from test.synthetic import Verbatim, write_binary, write_fid +from test.synthetic import Verbatim, write_binary, write_fid, write_jcampdx # pv6 EPSI, shrunk: 6 read points, 4 phase steps, 2 segments, 4 spectral points # per read line -- the real scan is 96 x 64 x 4 x 64. @@ -84,6 +85,48 @@ def test_epsi_keeps_every_stored_complex_sample(tmp_path): assert np.array_equal(np.sort_complex(dataset.data.reshape(-1)), np.sort_complex(expected)) +def test_fid_uses_pdata_one_reco_unless_a_different_one_is_selected(tmp_path): + """`pdata/1` is conventional; `reco_path` selects another reconstruction.""" + path = write_fid( + tmp_path / "30", + acqp={ + "GO_block_size": "continuous", + "GO_raw_data_format": "GO_32BIT_SGN_INT", + "BYTORDA": "little", + "ACQ_dim": 2, + "ACQ_dim_desc": Verbatim("( 2 )\nSpatial Spatial"), + "ACQ_size": np.array([8, 2]), + "ACQ_phase_factor": 1, + "ACQ_scan_size": "ACQ_phase_factor_scans", + "PULPROG": "", + "NI": 1, + "NR": 1, + }, + method={"Method": "", "PVM_DigNp": 4, "PVM_EncMatrix": np.array([4, 2]), "PVM_EncNReceivers": 1}, + blocks=2, + ) + write_jcampdx(path.parent / "pdata" / "1" / "reco", {"RECO_inp_order": "NO_REORDERING"}) + reco = write_jcampdx(path.parent / "pdata" / "2" / "reco", {"RECO_inp_order": "REV_ALT_ROWS"}) + + default = Dataset(path, load=LOAD_STAGES["properties"]) + default.load_schema() + assert "reco" in default.parameters + assert default["RECO_inp_order"].value == "NO_REORDERING" + assert not default._schema.mirror_odd_lines + + dataset = Dataset(path, reco_path=reco, load=LOAD_STAGES["properties"]) + dataset.load_schema() + assert "reco" in dataset.parameters + with pytest.warns(RuntimeWarning, match="ACQ_scan_size=.*disagrees with scheme_id"): + assert dataset._schema.continuous_train + with pytest.warns(RuntimeWarning, match="RECO_inp_order=.*disagrees with scheme_id"): + assert dataset._schema.mirror_odd_lines + + lines = np.arange(8).reshape((4, 2), order="F") + expected = np.column_stack((lines[:, 0], lines[::-1, 1])) + assert np.array_equal(dataset._schema._mirror_odd_lines(lines), expected) + + def kblock_fid(tmp_path, acqp, method, *, blocks, samples_per_block): """A Standard_KBlock_Format fid: real samples first, then zero padding.""" block_size = int(np.ceil(int(np.atleast_1d(acqp["ACQ_size"])[0]) * int(method["PVM_EncNReceivers"]) * 4 / 1024.0) * 1024 // 4)