Merge pull request #491 from cbinyu/bvecs_just_for_dwi

BF: Ensure bvec/bval files are only created for dwi output

Author: yarikoptic

heudiconv/convert.py

@@ -41,6 +41,7 @@
 )
 
 LOCKFILE = 'heudiconv.lock'
+DW_IMAGE_IN_FMAP_FOLDER_WARNING = 'Diffusion-weighted image saved in non dwi folder ({folder})'
 lgr = logging.getLogger(__name__)
 
 
@@ -674,9 +675,22 @@ def save_converted_files(res, item_dicoms, bids_options, outtype, prefix, outnam
         return
 
     if isdefined(res.outputs.bvecs) and isdefined(res.outputs.bvals):
-        outname_bvecs, outname_bvals = prefix + '.bvec', prefix + '.bval'
-        safe_movefile(res.outputs.bvecs, outname_bvecs, overwrite)
-        safe_movefile(res.outputs.bvals, outname_bvals, overwrite)
+        if prefix_dirname.endswith('dwi'):
+            outname_bvecs, outname_bvals = prefix + '.bvec', prefix + '.bval'
+            safe_movefile(res.outputs.bvecs, outname_bvecs, overwrite)
+            safe_movefile(res.outputs.bvals, outname_bvals, overwrite)
+        else:
+            if bvals_are_zero(res.outputs.bvals):
+                os.remove(res.outputs.bvecs)
+                os.remove(res.outputs.bvals)
+                lgr.debug("%s and %s were removed since not dwi", res.outputs.bvecs, res.outputs.bvals)
+            else:
+                lgr.warning(DW_IMAGE_IN_FMAP_FOLDER_WARNING.format(folder= prefix_dirname))
+                lgr.warning(".bvec and .bval files will be generated. This is NOT BIDS compliant")
+                outname_bvecs, outname_bvals = prefix + '.bvec', prefix + '.bval'
+                safe_movefile(res.outputs.bvecs, outname_bvecs, overwrite)
+                safe_movefile(res.outputs.bvals, outname_bvals, overwrite)
+
 
     if isinstance(res_files, list):
         res_files = sorted(res_files)
@@ -806,3 +820,23 @@ def  add_taskname_to_infofile(infofiles):
 
         # write to outfile
         save_json(infofile, meta_info)
+
+
+def bvals_are_zero(bval_file):
+    """Checks if all entries in a bvals file are zero (or 5, for Siemens files).
+    Returns True if that is the case, otherwise returns False
+
+    Parameters
+    ----------
+    bval_file : file with the bvals
+
+    Returns
+    -------
+    True if all are zero; False otherwise.
+    """
+
+    with open(bval_file) as f:
+        bvals = f.read().split()
+
+    bvals_unique = set(float(b) for b in bvals)
+    return bvals_unique == {0.} or bvals_unique == {5.}

heudiconv/heuristics/test_b0dwi_for_fmap.py

@@ -0,0 +1,33 @@
+"""Heuristic to extract a b-value=0 DWI image (basically, a SE-EPI)
+both as a fmap and as dwi
+
+It is used just to test that a 'DIFFUSION' image that the user
+chooses to extract as fmap (pepolar case) doesn't produce _bvecs/
+_bvals json files, while it does for dwi images
+"""
+
+
+def create_key(template, outtype=('nii.gz',), annotation_classes=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):
+    """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
+    """
+    fmap = create_key('sub-{subject}/fmap/sub-{subject}_acq-b0dwi_epi')
+    dwi = create_key('sub-{subject}/dwi/sub-{subject}_acq-b0dwi_dwi')
+
+    info = {fmap: [], dwi: []}
+    for s in seqinfo:
+        if 'DIFFUSION' in s.image_type:
+            info[fmap].append(s.series_id)
+            info[dwi].append(s.series_id)
+    return info

heudiconv/tests/data/b0dwiForFmap/b0dwi_for_fmap+00001.dcm

@@ -1,11 +1,18 @@
 """Test functions in heudiconv.convert module.
 """
+import os.path as op
+from glob import glob
+
 import pytest
+from .utils import TESTS_DATA_PATH
 
 from heudiconv.convert import (update_complex_name,
                                update_multiecho_name,
-                               update_uncombined_name)
+                               update_uncombined_name,
+                               DW_IMAGE_IN_FMAP_FOLDER_WARNING,
+                               )
 from heudiconv.bids import BIDSError
+from heudiconv.cli.run import main as runner
 
 
 def test_update_complex_name():
@@ -76,3 +83,31 @@ def test_update_uncombined_name():
     out_fn_true = 'sub-X_ses-Y_task-Z_run-01_ch-04_bold'
     out_fn_test = update_uncombined_name(metadata, fn, channel_names)
     assert out_fn_test == out_fn_true
+
+
+def test_b0dwi_for_fmap(tmpdir, capfd):
+    """Make sure we raise a warning when .bvec and .bval files
+    are present but the modality is not dwi.
+    We check it by extracting a few DICOMs from a series with
+    bvals: 5 5 1500
+    """
+    tmppath = tmpdir.strpath
+    subID = 'b0dwiForFmap'
+    args = (
+        "-c dcm2niix -o %s -b -f test_b0dwi_for_fmap --files %s -s %s"
+        % (tmpdir, op.join(TESTS_DATA_PATH, 'b0dwiForFmap'), subID)
+    ).split(' ')
+    runner(args)
+
+    # assert that it raised a warning that the fmap directory will contain
+    # bvec and bval files.
+    output = capfd.readouterr().err.split('\n')
+    expected_msg = DW_IMAGE_IN_FMAP_FOLDER_WARNING.format(folder=op.join(tmppath, 'sub-%s', 'fmap') % subID)
+    assert [o for o in output if expected_msg in o]
+
+    # check that both 'fmap' and 'dwi' directories have been extracted and they contain
+    # *.bvec and a *.bval files
+    for mod in ['fmap', 'dwi']:
+        assert op.isdir(op.join(tmppath, 'sub-%s', mod) % (subID))
+        for ext in ['bval', 'bvec']:
+            assert glob(op.join(tmppath, 'sub-%s', mod, 'sub-%s_*.%s') % (subID, subID, ext))