diff --git a/heudiconv/bids.py b/heudiconv/bids.py index ab32feb1..43871ee5 100644 --- a/heudiconv/bids.py +++ b/heudiconv/bids.py @@ -1074,12 +1074,25 @@ class BIDSFile: order matters """ + # The full list of entities, in the order mandated by the BIDS entity table + # (``rules/entities.yaml`` of the BIDS schema). Entities we do not produce + # ourselves are listed as well, so that a filename which carries them (e.g. + # coming from a heuristic) survives a parse/serialize round-trip in the + # right order. _known_entities = [ "sub", + "tpl", "ses", + "cohort", + "sample", "task", + "tracksys", "acq", + "nuc", + "voi", "ce", + "trc", + "stain", "rec", "dir", "run", @@ -1089,7 +1102,22 @@ class BIDSFile: "inv", "mt", "part", + # not a BIDS entity: heudiconv's own, for uncombined multi-channel data. + # Kept at the position update_uncombined_name() places it at. + "ch", + "proc", + "hemi", + "space", + "split", "recording", + "chunk", + "atlas", + "seg", + "scale", + "res", + "den", + "label", + "desc", ] def __init__( @@ -1113,15 +1141,32 @@ def __eq__(self, other: Any) -> bool: @classmethod def parse(cls, filename: str) -> BIDSFile: - """Parse the filename for BIDS entities, suffix and extension""" - # use re.findall to find all lower-case-letters + '-' + alphanumeric + '_' pairs: - entities_list = re.findall("([a-z]+)-([a-zA-Z0-9]+)[_]*", filename) - # keep only those in the _known_entities list: - entities = {k: v for k, v in entities_list if k in BIDSFile._known_entities} - # get whatever comes after the last key-value pair, and remove any '_' that + """Parse the filename for BIDS entities, suffix and extension + + Raises + ------ + ValueError + If no ``key-value`` pair could be found at all, i.e. the name is + not a BIDS one. + """ + # Entities are the leading run of lower-case-letters + '-' + alphanumeric + # pairs; everything after it is the suffix (+ extension). Matching the + # run as a whole, rather than every such pair anywhere in the name, keeps + # us from mistaking a part of the suffix for an entity: reproin marks + # duplicate series with a trailing '__dup-01', and a 'T1w-mod' suffix + # contains a 'w-mod' pair. + match = re.match( + r"((?:[a-z]+-[a-zA-Z0-9]+_)*[a-z]+-[a-zA-Z0-9]+)(?=_|\.|$)", filename + ) + if not match: + raise ValueError(f"No BIDS entities found in {filename!r}") + # keep all of them: dropping the ones we do not know about would silently + # lose information from the filename (see __str__, which puts the unknown + # ones back at the end). + entities = dict(re.findall("([a-z]+)-([a-zA-Z0-9]+)", match.group(1))) + # get whatever comes after the entities, and remove any '_' that # might come in front: - ending = filename.split("-".join(entities_list[-1]))[-1] - ending = remove_prefix(ending, "_") + ending = remove_prefix(filename[match.end() :], "_") # the first dot ('.') separates the suffix from the extension: if "." in ending: suffix, extension = ending.split(".", 1) @@ -1136,17 +1181,24 @@ def __str__(self) -> str: # reconstitute the ending for the filename: suffix = "_" + self.suffix if self.suffix else "" extension = "." + self.extension if self.extension else "" - return ( - "_".join( - [ - "-".join([e, self._entities[e]]) - for e in self._known_entities - if e in self._entities - ] - ) + ordered = [e for e in self._known_entities if e in self._entities] + # entities we do not know about cannot be placed within the entity table + # order, so keep them (in the order they were given) right before the + # suffix rather than dropping them on the floor + unknown = [e for e in self._entities if e not in self._known_entities] + out = ( + "_".join(["-".join([e, self._entities[e]]) for e in ordered + unknown]) + suffix + extension ) + if unknown: + lgr.warning( + "Unknown BIDS entities (%s) in %s: keeping them, but their " + "placement within the filename might not be BIDS-compliant.", + ", ".join(unknown), + out, + ) + return out def __getitem__(self, entity: str) -> Optional[str]: return self._entities[entity] if entity in self._entities else None diff --git a/heudiconv/convert.py b/heudiconv/convert.py index 03de785a..b6a4f3ae 100644 --- a/heudiconv/convert.py +++ b/heudiconv/convert.py @@ -20,6 +20,7 @@ from .bids import ( BIDS_VERSION, BIDSError, + BIDSFile, add_participant_record, populate_bids_templates, populate_intended_for, @@ -521,6 +522,79 @@ def update_uncombined_name( return filename +def update_multiorient_name( + metadata: dict[str, Any], + filename: str, + iops: list[str], +) -> str: + """ + Insert `_chunk-` entity into filename if data are from a sequence + that outputs multiple FoV (localizer, multi-FoV bold) + + The index is the position of this file's orientation within ``iops``, which + lists the distinct orientations of the series in the order in which the + files carrying them are converted, i.e. sorted by the name dcm2niix gave + them. Indexing the *orientations* rather than the files means that files + which share an orientation share a chunk, which matters when the series is + additionally split by echo or by magnitude/phase. Taking the order from + the files rather than from the orientation values themselves keeps the + index stable: sorting the orientations would give an order which is both + arbitrary and liable to change between subjects, since a slightly different + obliquity would reshuffle it. + + Parameters + ---------- + metadata : dict + Scan metadata dictionary from BIDS sidecar file. + filename : str + Incoming filename + iops : list of str + The distinct ``ImageOrientationPatientDICOM`` values of the series, as + strings, in the order in which they were encountered. + + Returns + ------- + filename : str + Updated filename with chunk entity added, if appropriate. + """ + iop = metadata.get("ImageOrientationPatientDICOM") + if iop is None: + lgr.warning( + "Not embedding multi-orientation information into %r: it has no " + "ImageOrientationPatientDICOM while other files of the series do.", + filename, + ) + return filename + if str(iop) not in iops: # should not happen: iops comes from this metadata + lgr.warning( + "Not embedding multi-orientation information into %r: its " + "orientation is not among the %d collected for the series.", + filename, + len(iops), + ) + return filename + chunk = str(iops.index(str(iop)) + 1) + try: + bids_file = BIDSFile.parse(filename) + if bids_file["chunk"]: + lgr.warning( + "Not embedding multi-orientation information as %r already uses " + "the chunk- entity; falling back to appending an index to the " + "suffix, which is not BIDS-compliant.", + filename, + ) + return filename + bids_file["chunk"] = chunk + return str(bids_file) + except ValueError as exc: + # not a name we can take apart (no entities at all, or no sub-); + # leave it to the caller's fallback rather than aborting the conversion + lgr.warning( + "Not embedding multi-orientation information into %r: %s", filename, exc + ) + return filename + + def convert( items: list[tuple[str, tuple[str, ...], list[str]]], converter: str, @@ -1029,6 +1103,9 @@ def rename_files() -> None: echo_times: set[float] = set() channel_names: set[str] = set() image_types: set[str] = set() + # the distinct image orientations, kept in the order dcm2niix emitted + # them rather than sorted -- see update_multiorient_name + iops: list[str] = [] for metadata in bids_metas: if not metadata: continue @@ -1044,6 +1121,13 @@ def rename_files() -> None: image_types.update(metadata["ImageType"]) except KeyError: pass + try: + iop = str(metadata["ImageOrientationPatientDICOM"]) + except KeyError: + pass + else: + if iop not in iops: + iops.append(iop) is_multiecho = ( len(set(filter(bool, echo_times))) > 1 @@ -1054,10 +1138,16 @@ def rename_files() -> None: is_complex = ( "M" in image_types and "P" in image_types ) # Determine if data are complex (magnitude + phase) + is_multiorient = len(iops) > 1 echo_times_lst = sorted(echo_times) # also converts to list channel_names_lst = sorted(channel_names) # also converts to list - ### Loop through the bids_files, set the output name and save files + ### Loop through the bids_files and set the output names. We do not + ### save anything yet: the renaming below is only of use if it gives + ### every file a distinct name, and we cannot tell until we have them + ### all (dcm2niix splits a series on more criteria than we handle, so + ### e.g. two of the images can well share an orientation). + renamed: list[tuple[str, str, Optional[str]]] = [] for fl, suffix, bids_file, bids_meta in zip( res_files, suffixes, bids_files, bids_metas ): @@ -1084,12 +1174,38 @@ def rename_files() -> None: bids_meta, this_prefix_basename, channel_names_lst ) + if is_multiorient: + this_prefix_basename = update_multiorient_name( + bids_meta, this_prefix_basename, iops + ) + # Fallback option: # If we have failed to modify this_prefix_basename, because it didn't fall # into any of the options above, just add the suffix at the end: if this_prefix_basename == prefix_basename: this_prefix_basename += suffix + renamed.append((fl, this_prefix_basename, bids_file)) + + # If the renaming did not manage to tell the files apart, we would + # overwrite (or, with overwrite=False, fail on) our own output, so use + # the plain numeric suffix for all of them instead. + basenames = [basename for _, basename, _ in renamed] + if len(set(basenames)) != len(basenames): + lgr.warning( + "Renaming of the %d files converted for %s did not give them " + "unique names; falling back to appending an index to the suffix, " + "which is not BIDS-compliant.", + len(basenames), + prefix_basename, + ) + renamed = [ + (fl, prefix_basename + suffix, bids_file) + for (fl, _, bids_file), suffix in zip(renamed, suffixes) + ] + + ### Now save the files under the names we settled on + for fl, this_prefix_basename, bids_file in renamed: # Finally, form the outname by stitching the directory and outtype: outname = op.join(prefix_dirname, this_prefix_basename) outfile = outname + "." + outtype diff --git a/heudiconv/heuristics/bids_localizer.py b/heudiconv/heuristics/bids_localizer.py new file mode 100644 index 00000000..d13dff49 --- /dev/null +++ b/heudiconv/heuristics/bids_localizer.py @@ -0,0 +1,52 @@ +"""Heuristic demonstrating conversion of a multi-orientation localizer. + +It only cares about converting sequences which have "localizer" in their +series_description and outputs to BIDS. + +Note that BIDS has no suffix for localizers/scouts, so the `_localizer` +name used below is *not* BIDS-compliant -- this heuristic exists to +exercise the `chunk-` naming of the multiple orientations dcm2niix +produces for such a series, not to serve as a model for what to do with +localizers (`reproin` converts them to DICOMs only). It also makes no +attempt to tell two localizer series apart, so it is only usable on a +session which has a single one. +""" + +from __future__ import annotations + +from typing import Optional + +from heudiconv.utils import SeqInfo + + +def create_key( + template: Optional[str], + outtype: tuple[str, ...] = ("nii.gz",), + annotation_classes: None = None, +) -> tuple[str, tuple[str, ...], None]: + if template is None or not template: + raise ValueError("Template must be a valid format string") + return (template, outtype, annotation_classes) + + +def infotodict( + seqinfo: list[SeqInfo], +) -> dict[tuple[str, tuple[str, ...], None], list[str]]: + """Heuristic evaluator for determining which runs belong where + + allowed template fields - follow python string module: + + item: index within category + subject: participant id + seqitem: run number during scanning + subindex: sub index within group + """ + localizer = create_key("sub-{subject}/anat/sub-{subject}_localizer") + + info: dict[tuple[str, tuple[str, ...], None], list[str]] = { + localizer: [], + } + for s in seqinfo: + if "localizer" in s.series_description: + info[localizer].append(s.series_id) + return info diff --git a/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015350928736278242 b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015350928736278242 new file mode 100755 index 00000000..7f944f28 Binary files /dev/null and b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015350928736278242 differ diff --git a/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.201811301535098526678240 b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.201811301535098526678240 new file mode 100755 index 00000000..989f6c68 Binary files /dev/null and b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.201811301535098526678240 differ diff --git a/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015351140807678244 b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015351140807678244 new file mode 100755 index 00000000..fa7cc8da Binary files /dev/null and b/heudiconv/tests/data/01-localizer_64ch/MR.1.3.12.2.1107.5.2.43.167006.2018113015351140807678244 differ diff --git a/heudiconv/tests/test_bids.py b/heudiconv/tests/test_bids.py index 38ddd215..ff5ae11f 100644 --- a/heudiconv/tests/test_bids.py +++ b/heudiconv/tests/test_bids.py @@ -8,11 +8,13 @@ from datetime import datetime, timedelta from glob import glob import itertools +import logging import os import os.path as op from pathlib import Path from random import choice, random, seed, shuffle import re +import shutil import string from typing import Any, Dict, List, Optional, Tuple @@ -1575,6 +1577,103 @@ def test_BIDSFile() -> None: assert my_bids_file["echo"] == "2" +@pytest.mark.parametrize( + "shuffled,expected", + [ + # entities which heudiconv itself adds, all of them at once + ( + "sub-1_chunk-2_part-mag_echo-1_acq-A_dir-AP_ses-B_run-3_rec-C_ce-D_task-E", + "sub-1_ses-B_task-E_acq-A_ce-D_rec-C_dir-AP_run-3_echo-1_part-mag_chunk-2", + ), + # entities from the non-MR modalities, which the entity table + # interleaves with the ones above + ( + "sub-1_stain-A_nuc-1H_trc-C_voi-D_tracksys-E_acq-F_task-G", + "sub-1_task-G_tracksys-E_acq-F_nuc-1H_voi-D_trc-C_stain-A", + ), + # entities we never produce ourselves still have to round-trip + ("sub-1_space-T1w_hemi-L_run-1", "sub-1_run-1_hemi-L_space-T1w"), + ], +) +@pytest.mark.ai_generated +def test_BIDSFile_entity_order(shuffled: str, expected: str) -> None: + """__str__ must order entities as the BIDS entity table mandates""" + assert str(BIDSFile.parse(shuffled + "_T1w.nii.gz")) == expected + "_T1w.nii.gz" + + +@pytest.mark.ai_generated +def test_BIDSFile_unknown_entities(caplog: pytest.LogCaptureFixture) -> None: + """Entities we do not know about must be kept, not silently dropped""" + with caplog.at_level(logging.WARNING): + out = str(BIDSFile.parse("sub-1_run-1_madeup-X_T1w.nii.gz")) + assert out == "sub-1_run-1_madeup-X_T1w.nii.gz" + assert "madeup" in caplog.text + + +@pytest.mark.ai_generated +def test_convert_multiorient( + tmp_path: Path, + heuristic: str = "bids_localizer.py", + subID: str = "loc", +) -> None: + """Test conversion of a series which dcm2niix splits into several images + because they were acquired with different orientations (a 3-plane + localizer): each of them should get its own `chunk-` index. + """ + datadir = op.join(TESTS_DATA_PATH, "01-localizer_64ch") + outdir = tmp_path / "out" + outdir.mkdir() + args = gen_heudiconv_args(datadir, str(outdir), subID, heuristic) + runner(args) + + anatdir = outdir / f"sub-{subID}" / "anat" + expected = { + f"sub-{subID}_chunk-{chunk}_localizer.{ext}" + for chunk in (1, 2, 3) + for ext in ("nii.gz", "json") + } + # nothing more, nothing less -- in particular no file left with the + # dcm2niix `_i0000` postfix or with a mangled suffix + assert {f.name for f in anatdir.iterdir()} == expected + + +@pytest.mark.ai_generated +def test_convert_multiorient_nonunique( + tmp_path: Path, + heuristic: str = "bids_localizer.py", + subID: str = "loc", +) -> None: + """dcm2niix splits a series on more than the orientation, so it can produce + two images which share one. `chunk-` cannot tell those apart, and we must + not end up overwriting (or failing to move onto) our own output. + """ + import pydicom + + datadir = tmp_path / "dicoms" + shutil.copytree(op.join(TESTS_DATA_PATH, "01-localizer_64ch"), datadir) + # add a 4th image repeating the 1st one's orientation, at a different + # matrix size so that dcm2niix writes it out separately + ds = pydicom.dcmread(sorted(datadir.iterdir())[0]) + arr = ds.pixel_array[::2, ::2] + ds.PixelData = arr.tobytes() + ds.Rows, ds.Columns = arr.shape + ds.PixelSpacing = [float(v) * 2 for v in ds.PixelSpacing] + ds.SOPInstanceUID = pydicom.uid.generate_uid() + ds.InstanceNumber = 99 + ds.save_as(datadir / "extra.dcm") + + outdir = tmp_path / "out" + outdir.mkdir() + runner(gen_heudiconv_args(str(datadir), str(outdir), subID, heuristic)) + + niftis = sorted((outdir / f"sub-{subID}" / "anat").glob("*.nii.gz")) + # every converted image is still there, under a name of its own + assert len(niftis) == 4 + assert len({f.name for f in niftis}) == 4 + # and since chunk- could not do it, none of them claims to be a chunk + assert not any("chunk-" in f.name for f in niftis) + + @pytest.mark.skipif(not have_datalad, reason="no datalad") def test_ME_mag_phase_conversion( monkeypatch: pytest.MonkeyPatch, diff --git a/heudiconv/tests/test_convert.py b/heudiconv/tests/test_convert.py index 7e888678..fa259273 100644 --- a/heudiconv/tests/test_convert.py +++ b/heudiconv/tests/test_convert.py @@ -3,6 +3,7 @@ from __future__ import annotations from glob import glob +import json import os.path as op from pathlib import Path from typing import Optional @@ -17,6 +18,7 @@ bvals_are_zero, update_complex_name, update_multiecho_name, + update_multiorient_name, update_uncombined_name, ) from heudiconv.utils import load_heuristic @@ -143,6 +145,53 @@ def test_update_uncombined_name() -> None: update_uncombined_name(metadata, base_fn, set(channel_names)) # type: ignore[arg-type] +SAG = "[0, 1, 0, 0, 0, -1]" +COR = "[1, 0, 0, 0, 0, -1]" +TRA = "[1, 0, 0, 0, 1, 0]" +# as encountered in file order, i.e. *not* sorted -- COR sorts before TRA +THREE_PLANES = [SAG, TRA, COR] + + +@pytest.mark.parametrize( + "iop,base_fn,expected", + [ + # index is the position of the orientation within the series, so the + # last-listed orientation is chunk-3 even though it sorts second + (SAG, "sub-X_task-Z_run-01_bold", "sub-X_task-Z_run-01_chunk-1_bold"), + (TRA, "sub-X_task-Z_run-01_bold", "sub-X_task-Z_run-01_chunk-2_bold"), + (COR, "sub-X_task-Z_run-01_bold", "sub-X_task-Z_run-01_chunk-3_bold"), + # chunk- goes after part- and heudiconv's own ch-, per the entity table + (SAG, "sub-X_acq-A_part-mag_T2starw", "sub-X_acq-A_part-mag_chunk-1_T2starw"), + (SAG, "sub-X_ch-03_T1w", "sub-X_ch-03_chunk-1_T1w"), + # names we cannot take apart are handed back untouched, so that the + # caller's fallback applies instead of the conversion blowing up + (SAG, "sub-X_chunk-7_T1w", "sub-X_chunk-7_T1w"), + (SAG, "localizer", "localizer"), + (SAG, "run-01_T1w", "run-01_T1w"), + # reproin's __dup-NN marker must not be mistaken for an entity and eat + # the suffix along with it + (SAG, "sub-X_task-Z_bold__dup-01", "sub-X_task-Z_chunk-1_bold__dup-01"), + ], +) +@pytest.mark.ai_generated +def test_update_multiorient_name(iop: str, base_fn: str, expected: str) -> None: + """Unit testing for heudiconv.convert.update_multiorient_name(), which updates + filenames with the chunk field if appropriate. + """ + metadata = {"ImageOrientationPatientDICOM": json.loads(iop)} + assert update_multiorient_name(metadata, base_fn, THREE_PLANES) == expected + + +@pytest.mark.ai_generated +def test_update_multiorient_name_no_orientation() -> None: + """A sidecar without an orientation while its siblings have one used to + raise a KeyError; it should just be left alone.""" + assert ( + update_multiorient_name({}, "sub-X_task-Z_bold", THREE_PLANES) + == "sub-X_task-Z_bold" + ) + + def test_b0dwi_for_fmap(tmp_path: Path, caplog: pytest.LogCaptureFixture) -> None: """Make sure we raise a warning when .bvec and .bval files are present but the modality is not dwi.