ENH: Adds populate_intended_for for fmaps

heudiconv/bids.py

@@ -19,6 +19,7 @@
     save_json,
     create_file_if_missing,
     json_dumps_pretty,
+    add_field_to_json,
     set_readonly,
     is_readonly,
     get_datetime,
@@ -46,6 +47,9 @@ class BIDSError(Exception):
 
 BIDS_VERSION = "1.4.1"
 
+SHIM_KEY = 'ShimSetting'
+lgr.debug('shim_key: {}'.format(SHIM_KEY))
+
 
 def populate_bids_templates(path, defaults={}):
     """Premake BIDS text files with templates"""
@@ -457,3 +461,145 @@ def convert_sid_bids(subject_id):
     lgr.warning('{0} contained nonalphanumeric character(s), subject '
                 'ID was cleaned to be {1}'.format(subject_id, sid))
     return sid, subject_id
+
+
+def get_shim_setting(json_file):
+    """
+    Gets the "ShimSetting" field from a json_file.
+    If not present:
+    - for fmap files: use the <acq> entity from the file name
+    - for other files: use the folder name ('anat', 'func', 'dwi', ...)
+
+    Parameters:
+    ----------
+    json_file : str
+
+    Returns:
+    -------
+    str with "ShimSetting" value or <acq> entity
+    """
+    data = load_json(json_file)
+    if SHIM_KEY in data.keys():
+        shims = data[SHIM_KEY]
+    else:
+        modality = op.basename(op.dirname(json_file))
+        if modality == 'fmap':
+            # extract the <acq> entity:
+            shims = re.search('(?<=[/_]acq-)\w+',json_file).group(0).split('_')[0]
+            if shims.lower() in ['fmri', 'bold']:
+                shims = 'func'
+        else:
+            shims = modality
+    return shims
+
+
+def populate_intended_for(path_to_bids_session):
+    """
+    Adds the 'IntendedFor' field to the fmap .json files in a session folder.
+    It goes through the session folders and checks what runs have the same
+    'ShimSetting' as the fmaps. If there is no 'ShimSetting' field in the json
+    file, we'll use the folder name ('func', 'dwi', 'anat') and see which fmap
+    with a matching '_acq' entity.
+
+    If several fmap runs have the same 'ShimSetting' (or '_acq'), it will use
+    the first one. Because fmaps come in groups (with reversed PE polarity,
+    or magnitude/phase), it adds the same runs to the 'IntendedFor' of the
+    corresponding fmaps by checking the '_acq' and '_run' entities.
+
+    Note: the logic behind the way we decide how to populate the "IntendedFor"
+    is: we want all images in the session (except for the fmap images
+    themselves) to have AT MOST one fmap.  (That is, a pair of SE EPI with
+    reversed polarities, or a magnitude a phase field map). If there are more
+    than one fmap (more than a fmap pair) with the same acquisition parameters
+    as, say, a functional run, we will just assign that run to the FIRST pair,
+    while leaving the other fmap pairs without any assigned images.  If the
+    user's intentions were different, he/she will have to manually edit the
+    fmap json files.
+
+    Parameters:
+    ----------
+    path_to_bids_session : str or os.path
+        path to the session folder (or to the subject folder, if there are no
+        sessions).
+    """
+    lgr.info('')
+    lgr.info('Adding "IntendedFor" to the fieldmaps in {}.'.format(path_to_bids_session))
+
+    # Resolve path (eliminate '..')
+    path_to_bids_session = op.abspath(path_to_bids_session)
+
+    # get the BIDS folder (if "data_folder" includes the session, remove it):
+    if op.basename(path_to_bids_session).startswith('ses-'):
+        bids_folder = op.dirname(path_to_bids_session)
+    else:
+        bids_folder = path_to_bids_session
+
+    fmap_dir = op.join(path_to_bids_session, 'fmap')
+    if not op.exists(fmap_dir):
+        lgr.warning('Fmap folder not found in {}.'.format(path_to_bids_session))
+        lgr.warning('We cannot add the IntendedFor field')
+        return
+
+    # Get a list of all fmap json files in the session:
+    # (we will remove elements later on, so don't just iterate)
+    fmap_jsons = sorted([j for j in glob(op.join(path_to_bids_session, 'fmap/*.json'))])
+
+    # Get a set with all non-fmap json files in the session (set is easier):
+    # We also exclude the SBRef files.
+    session_jsons = set(
+        j for j in glob(op.join(path_to_bids_session, '*/*.json')) if not (
+            j in fmap_jsons
+            # j[:-5] removes the '.json' from the end
+            or j[:-5].endswith('_sbref')
+        )
+    )
+
+    # Loop through all the fmap json files and, for each one, find which other
+    # non-fmap images in the session have the same shim settings. Those that
+    # match are added to the intended_for list and removed from the list of
+    # non-fmap json files in the session (since they have already assigned to
+    # a fmap).
+    # After finishing with all the non-fmap images in the session, we go back
+    # to the fmap json file list, and find any other fmap json files of the
+    # same acquisition type and run number (because fmaps have several files:
+    # normal- and reversed-polarity, or magnitude and phase, etc.) We add the
+    # same IntendedFor list to those other corresponding fmap json files, and
+    # remove them from the list of available fmap json files.
+    # Once we have gone through all the fmap json files, we are done.
+    runs_accounted_for = set()
+    fmaps_accounted_for = set()
+    for fm_json in fmap_jsons:
+        if fm_json not in fmaps_accounted_for:
+            lgr.debug('Looking for runs for {}'.format(fm_json))
+            fm_shims = get_shim_setting(fm_json)
+
+            intended_for = []
+            for image_json in session_jsons:
+                image_shims = get_shim_setting(image_json)
+                if image_shims == fm_shims:
+                    # BIDS specifies that the intended for are:
+                    # - **image** files
+                    # - path relative to the **subject level**
+                    image_json_relative_path = op.relpath(image_json, start=bids_folder)
+                    # image_json_relative_path[:-5] removes the '.json' extension:
+                    intended_for.append(
+                        image_json_relative_path[:-5] + '.nii.gz'
+                    )
+                    runs_accounted_for.add(image_json)
+            if len(intended_for) > 0:
+                intended_for = sorted([str(f) for f in intended_for])
+                # find all fmap json files with the same <acq> and <run> entities:
+                fm_json_name = op.basename(fm_json)
+                acq_match = re.findall('([/_]acq-([a-zA-Z0-9]*))', fm_json_name)
+                acq_str = acq_match[0][0] if acq_match else ''
+                run_match = re.findall('([/_]run-([a-zA-Z0-9]*))', fm_json_name)
+                run_str = run_match[0][0] if run_match else ''
+                # Loop through all the files that have the same "acq-" and "run-"
+                # Note: the following loop will also include 'fm_json'
+                for linked_fm_json in glob(op.join(path_to_bids_session, 'fmap/*' + acq_str + '*' + run_str + '*.json')):
+                    # add the IntendedFor field to the json file:
+                    add_field_to_json(linked_fm_json, {"IntendedFor": intended_for})
+                    fmaps_accounted_for.update({linked_fm_json})
+            # Remove the runs accounted for from the session_jsons list, so that
+            # we don't assign another fmap to this image:
+            session_jsons -= runs_accounted_for

heudiconv/cli/run.py

@@ -136,6 +136,7 @@ def get_parser():
             'heuristics', 'heuristic-info',
             'ls', 'populate-templates',
             'sanitize-jsons', 'treat-jsons',
+            'populate-intended-for'
         ),
         help='Custom action to be performed on provided files instead of '
              'regular operation.')

heudiconv/convert.py

@@ -28,6 +28,7 @@
 from .bids import (
     convert_sid_bids,
     populate_bids_templates,
+    populate_intended_for,
     save_scans_key,
     tuneup_bids_json_files,
     add_participant_record,
@@ -521,6 +522,28 @@ def convert(items, converter, scaninfo_suffix, custom_callable, with_prov,
         if custom_callable is not None:
             custom_callable(*item)
 
+    # Populate "IntendedFor" for fmap files.
+    # Because fmap files can only be used to correct for distortions in images
+    # collected within the same scanning session, find unique subject/session
+    # combinations from the outname in each item:
+    outnames = [item[0] for item in items]
+    # - grab "sub-<sID>[/ses-<ses>]", and keep only unique ones:
+    try:
+        sessions = set(
+            re.search(
+                'sub-(?P<subj>[a-zA-Z0-9]*)([{0}_]ses-(?P<ses>[a-zA-Z0-9]*))?'.format(op.sep),
+                oname
+            ).group(0)
+            for oname in outnames
+        )
+    except AttributeError:
+        # "sub-<sID>[/ses-<ses>]" is not present, so this is not BIDS compliant
+        # and it doesn't make sense to add "IntendedFor":
+        sessions = set()
+
+    for ses in sessions:
+        populate_intended_for(ses)
+
 
 def convert_dicom(item_dicoms, bids_options, prefix,
                   outdir, tempdirs, symlink, overwrite):

heudiconv/main.py

@@ -3,7 +3,7 @@
 import sys
 
 from . import __version__, __packagename__
-from .bids import populate_bids_templates, tuneup_bids_json_files
+from .bids import populate_bids_templates, tuneup_bids_json_files, populate_intended_for
 from .convert import prep_conversion
 from .due import due, Doi
 from .parser import get_study_sessions
@@ -47,13 +47,13 @@ def process_extra_commands(outdir, command, files, dicom_dir_template,
                            heuristic, session, subjs, grouping):
     """
     Perform custom command instead of regular operations. Supported commands:
-    ['treat-json', 'ls', 'populate-templates']
+    ['treat-json', 'ls', 'populate-templates', 'populate-intended-for']
 
     Parameters
     ----------
     outdir : str
         Output directory
-    command : {'treat-json', 'ls', 'populate-templates'}
+    command : {'treat-json', 'ls', 'populate-templates', 'populate-intended-for'}
         Heudiconv command to run
     files : list of str
         List of files
@@ -108,6 +108,12 @@ def process_extra_commands(outdir, command, files, dicom_dir_template,
         ensure_heuristic_arg(heuristic)
         from .utils import get_heuristic_description
         print(get_heuristic_description(heuristic, full=True))
+    elif command == 'populate-intended-for':
+        for subj in subjs:
+            session_path = op.join(outdir, 'sub-' + subj)
+            if session:
+                session_path = op.join(session_path, 'ses-' + session)
+            populate_intended_for(session_path)
     else:
         raise ValueError("Unknown command %s" % command)
     return