Skip to content
Merged
81 changes: 52 additions & 29 deletions fmriprep/workflows/bold/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,8 +31,6 @@
"""
import typing as ty

import nibabel as nb
import numpy as np
from nipype.interfaces import utility as niu
from nipype.pipeline import engine as pe
from niworkflows.utils.connections import listify
Expand Down Expand Up @@ -294,44 +292,20 @@ def init_bold_wf(
fieldmap_id=fieldmap_id,
omp_nthreads=omp_nthreads,
)
bold_anat_wf = init_bold_volumetric_resample_wf(
metadata=all_metadata[0],
fieldmap_id=fieldmap_id if not multiecho else None,
omp_nthreads=omp_nthreads,
mem_gb=mem_gb,
name='bold_anat_wf',
)
bold_anat_wf.inputs.inputnode.resolution = "native"

workflow.connect([
(inputnode, bold_native_wf, [
("fmap_ref", "inputnode.fmap_ref"),
("fmap_coeff", "inputnode.fmap_coeff"),
("fmap_id", "inputnode.fmap_id"),
]),
(inputnode, bold_anat_wf, [
("t1w_preproc", "inputnode.target_ref_file"),
("t1w_mask", "inputnode.target_mask"),
("fmap_ref", "inputnode.fmap_ref"),
("fmap_coeff", "inputnode.fmap_coeff"),
("fmap_id", "inputnode.fmap_id"),
]),
(bold_fit_wf, bold_native_wf, [
("outputnode.coreg_boldref", "inputnode.boldref"),
("outputnode.bold_mask", "inputnode.bold_mask"),
("outputnode.motion_xfm", "inputnode.motion_xfm"),
("outputnode.boldref2fmap_xfm", "inputnode.boldref2fmap_xfm"),
("outputnode.dummy_scans", "inputnode.dummy_scans"),
]),
(bold_fit_wf, bold_anat_wf, [
("outputnode.coreg_boldref", "inputnode.bold_ref_file"),
("outputnode.boldref2fmap_xfm", "inputnode.boldref2fmap_xfm"),
("outputnode.boldref2anat_xfm", "inputnode.boldref2anat_xfm"),
]),
(bold_native_wf, bold_anat_wf, [
("outputnode.bold_minimal", "inputnode.bold_file"),
("outputnode.motion_xfm", "inputnode.motion_xfm"),
]),
]) # fmt:skip

boldref_out = bool(nonstd_spaces.intersection(('func', 'run', 'bold', 'boldref', 'sbref')))
Expand Down Expand Up @@ -399,6 +373,35 @@ def init_bold_wf(
if config.workflow.level == "resampling":
return workflow

# Resample to anatomical space
bold_anat_wf = init_bold_volumetric_resample_wf(
metadata=all_metadata[0],
fieldmap_id=fieldmap_id if not multiecho else None,
omp_nthreads=omp_nthreads,
mem_gb=mem_gb,
name='bold_anat_wf',
)
bold_anat_wf.inputs.inputnode.resolution = "native"

workflow.connect([
(inputnode, bold_anat_wf, [
("t1w_preproc", "inputnode.target_ref_file"),
("t1w_mask", "inputnode.target_mask"),
("fmap_ref", "inputnode.fmap_ref"),
("fmap_coeff", "inputnode.fmap_coeff"),
("fmap_id", "inputnode.fmap_id"),
]),
(bold_fit_wf, bold_anat_wf, [
("outputnode.coreg_boldref", "inputnode.bold_ref_file"),
("outputnode.boldref2fmap_xfm", "inputnode.boldref2fmap_xfm"),
("outputnode.boldref2anat_xfm", "inputnode.boldref2anat_xfm"),
]),
(bold_native_wf, bold_anat_wf, [
("outputnode.bold_minimal", "inputnode.bold_file"),
("outputnode.motion_xfm", "inputnode.motion_xfm"),
]),
]) # fmt:skip

# Full derivatives, including resampled BOLD series
if nonstd_spaces.intersection(('anat', 'T1w')):
ds_bold_t1_wf = init_ds_volumes_wf(
Expand Down Expand Up @@ -506,7 +509,11 @@ def init_bold_wf(
]) # fmt:skip

if config.workflow.cifti_output:
from .resampling import init_bold_fsLR_resampling_wf, init_bold_grayords_wf
from .resampling import (
init_bold_fsLR_resampling_wf,
init_bold_grayords_wf,
init_goodvoxels_bold_mask_wf,
)

bold_MNI6_wf = init_bold_volumetric_resample_wf(
metadata=all_metadata[0],
Expand All @@ -517,12 +524,29 @@ def init_bold_wf(
)

bold_fsLR_resampling_wf = init_bold_fsLR_resampling_wf(
estimate_goodvoxels=config.workflow.project_goodvoxels,
grayord_density=config.workflow.cifti_output,
omp_nthreads=omp_nthreads,
mem_gb=mem_gb["resampled"],
)

if config.workflow.project_goodvoxels:
goodvoxels_bold_mask_wf = init_goodvoxels_bold_mask_wf(mem_gb["resampled"])

workflow.connect([
(inputnode, goodvoxels_bold_mask_wf, [("anat_ribbon", "inputnode.anat_ribbon")]),
(bold_anat_wf, goodvoxels_bold_mask_wf, [
("outputnode.bold_file", "inputnode.bold_file"),
]),
(goodvoxels_bold_mask_wf, bold_fsLR_resampling_wf, [
("outputnode.goodvoxels_mask", "inputnode.volume_roi"),
]),
]) # fmt:skip

bold_fsLR_resampling_wf.__desc__ += """\
A "goodvoxels" mask was applied during volume-to-surface sampling in fsLR space,
excluding voxels whose time-series have a locally high coefficient of variation.
"""

bold_grayords_wf = init_bold_grayords_wf(
grayord_density=config.workflow.cifti_output,
mem_gb=1,
Expand Down Expand Up @@ -571,7 +595,6 @@ def init_bold_wf(
("midthickness_fsLR", "inputnode.midthickness_fsLR"),
("sphere_reg_fsLR", "inputnode.sphere_reg_fsLR"),
("cortex_mask", "inputnode.cortex_mask"),
("anat_ribbon", "inputnode.anat_ribbon"),
]),
(bold_anat_wf, bold_fsLR_resampling_wf, [
("outputnode.bold_file", "inputnode.bold_file"),
Expand Down
59 changes: 16 additions & 43 deletions fmriprep/workflows/bold/resampling.py
Original file line number Diff line number Diff line change
Expand Up @@ -501,7 +501,6 @@ def _calc_lower_thr(in_stats):

def init_bold_fsLR_resampling_wf(
grayord_density: ty.Literal['91k', '170k'],
estimate_goodvoxels: bool,
omp_nthreads: int,
mem_gb: float,
name: str = "bold_fsLR_resampling_wf",
Expand All @@ -520,7 +519,6 @@ def init_bold_fsLR_resampling_wf(

from fmriprep.workflows.bold.resampling import init_bold_fsLR_resampling_wf
wf = init_bold_fsLR_resampling_wf(
estimate_goodvoxels=True,
grayord_density='92k',
omp_nthreads=1,
mem_gb=1,
Expand All @@ -530,9 +528,6 @@ def init_bold_fsLR_resampling_wf(
----------
grayord_density : :class:`str`
Either ``"91k"`` or ``"170k"``, representing the total *grayordinates*.
estimate_goodvoxels : :class:`bool`
Calculate mask excluding voxels with a locally high coefficient of variation to
exclude from surface resampling
omp_nthreads : :class:`int`
Maximum number of threads an individual process may use
mem_gb : :class:`float`
Expand All @@ -544,35 +539,33 @@ def init_bold_fsLR_resampling_wf(
------
bold_file : :class:`str`
Path to BOLD file resampled into T1 space
surfaces : :class:`list` of :class:`str`
Path to left and right hemisphere white, pial and midthickness GIFTI surfaces
morphometrics : :class:`list` of :class:`str`
Path to left and right hemisphere morphometric GIFTI surfaces, which must include thickness
white : :class:`list` of :class:`str`
Path to left and right hemisphere white matter GIFTI surfaces.
pial : :class:`list` of :class:`str`
Path to left and right hemisphere pial GIFTI surfaces.
midthickness : :class:`list` of :class:`str`
Path to left and right hemisphere midthickness GIFTI surfaces.
midthickness_fsLR : :class:`list` of :class:`str`
Path to left and right hemisphere midthickness GIFTI surfaces in fsLR space.
sphere_reg_fsLR : :class:`list` of :class:`str`
Path to left and right hemisphere sphere.reg GIFTI surfaces, mapping from subject to fsLR
anat_ribbon : :class:`str`
Path to mask of cortical ribbon in T1w space, for calculating goodvoxels
cortex_mask : :class:`list` of :class:`str`
Path to left and right hemisphere cortical masks.
volume_roi : :class:`str` or Undefined
Pre-calculated goodvoxels mask. Not required.

Outputs
-------
bold_fsLR : :class:`list` of :class:`str`
Path to BOLD series resampled as functional GIFTI files in fsLR space
goodvoxels_mask : :class:`str`
Path to mask of voxels, excluding those with locally high coefficients of variation

"""
import templateflow.api as tf
from niworkflows.engine.workflows import LiterateWorkflow as Workflow
from niworkflows.interfaces.utility import KeySelect
from smriprep import data as smriprep_data
from smriprep.interfaces.workbench import SurfaceResample

from fmriprep.interfaces.gifti import CreateROI
from fmriprep.interfaces.workbench import (
MetricFillHoles,
MetricRemoveIslands,
VolumeToSurfaceMapping,
)
from fmriprep.interfaces.workbench import VolumeToSurfaceMapping

fslr_density = "32k" if grayord_density == "91k" else "59k"

Expand All @@ -587,13 +580,13 @@ def init_bold_fsLR_resampling_wf(
niu.IdentityInterface(
fields=[
'bold_file',
'anat_ribbon',
'white',
'pial',
'midthickness',
'midthickness_fsLR',
'sphere_reg_fsLR',
'cortex_mask',
'volume_roi',
]
),
name='inputnode',
Expand All @@ -612,7 +605,7 @@ def init_bold_fsLR_resampling_wf(
)

outputnode = pe.Node(
niu.IdentityInterface(fields=['bold_fsLR', 'goodvoxels_mask']),
niu.IdentityInterface(fields=['bold_fsLR']),
name='outputnode',
)

Expand Down Expand Up @@ -687,6 +680,7 @@ def init_bold_fsLR_resampling_wf(
# Resample BOLD to native surface, dilate and mask
(inputnode, volume_to_surface, [
('bold_file', 'volume_file'),
('volume_roi', 'volume_roi'),
]),
(select_surfaces, volume_to_surface, [
('midthickness', 'surface_file'),
Expand All @@ -713,27 +707,6 @@ def init_bold_fsLR_resampling_wf(
(joinnode, outputnode, [('bold_fsLR', 'bold_fsLR')]),
]) # fmt:skip

if estimate_goodvoxels:
workflow.__desc__ += """\
A "goodvoxels" mask was applied during volume-to-surface sampling in fsLR space,
excluding voxels whose time-series have a locally high coefficient of variation.
"""

goodvoxels_bold_mask_wf = init_goodvoxels_bold_mask_wf(mem_gb)

workflow.connect([
(inputnode, goodvoxels_bold_mask_wf, [
("bold_file", "inputnode.bold_file"),
("anat_ribbon", "inputnode.anat_ribbon"),
]),
(goodvoxels_bold_mask_wf, volume_to_surface, [
("outputnode.goodvoxels_mask", "volume_roi"),
]),
(goodvoxels_bold_mask_wf, outputnode, [
("outputnode.goodvoxels_mask", "goodvoxels_mask"),
]),
]) # fmt:skip

return workflow


Expand Down