Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Draft pre-init ingression #102

Merged
merged 19 commits into from
Sep 23, 2024
Merged

Draft pre-init ingression #102

Show file tree
Hide file tree
Changes from 1 commit
Commits
Show all changes
19 commits
Select commit Hold shift + click to select a range
eb86d69
Draft pre-init ingression.
tsalo Aug 30, 2024
21702d1
Add input_dir.
tsalo Aug 30, 2024
fb6bacf
Update pyproject.toml
tsalo Aug 30, 2024
bca8832
Update pyproject.toml
tsalo Aug 30, 2024
60e4cea
Update parser.py
tsalo Aug 30, 2024
747579a
Update pyproject.toml
tsalo Aug 30, 2024
742ff2b
Update parser.py
tsalo Aug 30, 2024
02e2dfe
Update parser.py
tsalo Aug 30, 2024
a316f94
Update parser.py
tsalo Aug 30, 2024
61d4f07
Update parser.py
tsalo Aug 30, 2024
4467831
add external ingression
smeisler Sep 12, 2024
08478ea
linting
smeisler Sep 12, 2024
b6305b9
linting: the sequel
smeisler Sep 12, 2024
897a780
update qsirecon_build
smeisler Sep 12, 2024
0df2933
use acres for loading data
smeisler Sep 13, 2024
297c1eb
update i2q version
smeisler Sep 19, 2024
04491af
update i2q for acres dependency
smeisler Sep 19, 2024
6e096ca
repin nilearn
smeisler Sep 23, 2024
69f3fd7
Update parser.py
smeisler Sep 23, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion pyproject.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ dependencies = [
"filelock",
"fury",
"indexed_gzip <= 1.8.7",
"Ingress2QSIRecon == 0.1.0",
"Ingress2QSIRecon == 0.1.3",
"jinja2 < 3.1",
"matplotlib",
"networkx ~= 2.8.8",
Expand Down
32 changes: 26 additions & 6 deletions qsirecon/cli/parser.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -475,16 +475,36 @@ def parse_args(args=None, namespace=None):
work_dir.mkdir(exist_ok=True, parents=True)

# Run ingression if necessary
if config.workflow.input_type in ("ukbb", "hcpya"):
from ingress2qsirecon.utils.workflows import ingress2qsirecon
if config.workflow.input_type in ("hcpya", "ukb"):
import os.path as op
import shutil

import ingress2qsirecon
from ingress2qsirecon.utils.functions import create_layout
from ingress2qsirecon.utils.workflows import create_ingress2qsirecon_wf

# Fake BIDS directory to be created
config.execution.bids_dir = work_dir / "bids"

wf = ingress2qsirecon(
config.execution.input_dir,
config.execution.bids_dir,
config.workflow.input_type,
# Make fake BIDS files
ingress2recon_dir = op.dirname(ingress2qsirecon.__file__)
mattcieslak marked this conversation as resolved.
Show resolved Hide resolved
if not op.exists(op.join(config.execution.bids_dir, "dataset_description.json")):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

could you switch this to use pathlib instead of op? I think we're going to try to consistently use it to be in line with nipreps

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Think I changed this, untested so far:

if config.workflow.input_type in ("hcpya", "ukb"):
        import shutil

        from ingress2qsirecon.data import load_resource
        from ingress2qsirecon.utils.functions import create_layout
        from ingress2qsirecon.utils.workflows import create_ingress2qsirecon_wf

        # Fake BIDS directory to be created
        config.execution.bids_dir = work_dir / "bids"

        # Make fake BIDS files
        bids_scaffold = load_resource("bids_scaffold")
        if not (config.execution.bids_dir / "dataset_description.json").exists():
            shutil.copytree(
                bids_scaffold,
                config.execution.bids_dir,
                dirs_exist_ok=True,
            )

shutil.copytree(op.join(ingress2recon_dir, "data", "bids_scaffold"),
config.execution.bids_dir, dirs_exist_ok=True)

if config.execution.participant_label is None:
participants_ingression = []
else:
participants_ingression = list(config.execution.participant_label)
layouts = create_layout(config.execution.input_dir, config.execution.bids_dir,
config.workflow.input_type, participants_ingression)

# Create the ingression workflow
wf = create_ingress2qsirecon_wf(
layouts,
base_dir=work_dir,
)

# Configure the nipype workflow
wf.config["execution"]["crashdump_dir"] = str(log_dir)
wf.run()
Expand Down
32 changes: 0 additions & 32 deletions qsirecon/interfaces/anatomical.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -26,9 +26,6 @@
)
from nipype.utils.filemanip import fname_presuffix

from ..utils.ingress import ukb_dirname_to_bids
from .images import to_lps

LOGGER = logging.getLogger("nipype.interface")


Expand Down Expand Up @@ -151,35 +148,6 @@ def _get_if_exists(self, name, pattern, excludes=None):
self._results[name] = files[0]


class UKBAnatomicalIngressInputSpec(QSIPrepAnatomicalIngressInputSpec):
recon_input_dir = traits.Directory(
exists=True, mandatory=True, help="directory containing a single subject's results"
)


class UKBAnatomicalIngress(QSIPrepAnatomicalIngress):
input_spec = UKBAnatomicalIngressInputSpec

def _run_interface(self, runtime):
# Load the Bias-corrected brain and brain mask
input_path = Path(self.inputs.recon_input_dir)
bids_name = ukb_dirname_to_bids(self.inputs.recon_input_dir)

ukb_brain = input_path / "T1" / "T1_unbiased_brain.nii.gz"
ukb_brain_mask = input_path / "T1" / "T1_brain_mask.nii.gz"

conformed_t1w_file = str(Path(runtime.cwd) / (bids_name + "_desc-preproc_T1w.nii.gz"))
conformed_mask_file = str(Path(runtime.cwd) / (bids_name + "_desc-brain_mask.nii.gz"))

to_lps(nb.load(ukb_brain)).to_filename(conformed_t1w_file)
to_lps(nb.load(ukb_brain_mask)).to_filename(conformed_mask_file)

self._results["t1_preproc"] = conformed_t1w_file
self._results["t1_brain_mask"] = conformed_mask_file

return runtime


"""

The spherical harmonic coefficients are stored as follows. First, since the
Expand Down
59 changes: 1 addition & 58 deletions qsirecon/interfaces/ingress.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -119,61 +119,4 @@ def _get_qc_filename(self, out_root, params, desc, suffix):
used_keys = ["subject_id", "session_id", "acq_id", "dir_id", "run_id"]
fname = "_".join([params[key] for key in used_keys if params[key]])
return out_root + "/" + fname + "_desc-%s_dwi.%s" % (desc, suffix)


class _UKBioBankDWIIngressInputSpec(QSIPrepDWIIngressInputSpec):
dwi_file = File(exists=False, help="The name of what a BIDS dwi file may have been")
data_dir = traits.Directory(
exists=True, help="The UKB data directory for a subject. Must contain DTI/ and T1/"
)


class UKBioBankDWIIngress(SimpleInterface):
input_spec = _UKBioBankDWIIngressInputSpec
output_spec = QSIPrepDWIIngressOutputSpec

def _run_interface(self, runtime):
runpath = Path(runtime.cwd)

# The UKB input files
in_dir = Path(self.inputs.data_dir)
dwi_dir = in_dir / "DTI" / "dMRI" / "dMRI"
ukb_bval_file = dwi_dir / "bvals"
ukb_bvec_file = dwi_dir / "bvecs" # These are the same as eddy rotated
ukb_dwi_file = dwi_dir / "data_ud.nii.gz"
ukb_dwiref_file = dwi_dir / "dti_FA.nii.gz"

# The bids_name is what the images will be renamed to
bids_name = Path(self.inputs.dwi_file).name.replace(".nii.gz", "")
dwi_file = str(runpath / (bids_name + ".nii.gz"))
bval_file = str(runpath / (bids_name + ".bval"))
bvec_file = str(runpath / (bids_name + ".bvec"))
b_file = str(runpath / (bids_name + ".b"))
btable_file = str(runpath / (bids_name + "btable.txt"))
dwiref_file = str(runpath / (bids_name.replace("_dwi", "_dwiref") + ".nii.gz"))

dwi_conform = ConformDwi(
dwi_file=str(ukb_dwi_file), bval_file=str(ukb_bval_file), bvec_file=str(ukb_bvec_file)
)

result = dwi_conform.run()
Path(result.outputs.dwi_file).rename(dwi_file)
Path(result.outputs.bvec_file).rename(bvec_file)
shutil.copyfile(result.outputs.bval_file, bval_file)
# Reorient the dwi file to LPS+
self._results["dwi_file"] = dwi_file
self._results["bvec_file"] = bvec_file
self._results["bval_file"] = bval_file

# Create a btable_txt file for DSI Studio
btable_from_bvals_bvecs(bval_file, bvec_file, btable_file)
self._results["btable_file"] = btable_file

# Create a mrtrix .b file
_convert_fsl_to_mrtrix(bval_file, bvec_file, b_file)
self._results["b_file"] = b_file

# Create a dwi ref file
to_lps(nb.load(ukb_dwiref_file)).to_filename(dwiref_file)
self._results["dwi_ref"] = dwiref_file
return runtime

98 changes: 31 additions & 67 deletions qsirecon/workflows/base.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -32,23 +32,12 @@ def init_qsirecon_wf():
qsirecon_wf = Workflow(name=f"qsirecon_{ver.major}_{ver.minor}_wf")
qsirecon_wf.base_dir = config.execution.work_dir

if config.workflow.input_type not in ("qsiprep", "ukb"):
if config.workflow.input_type not in ("qsiprep", "hcpya", "ukb"):
raise NotImplementedError(
f"{config.workflow.input_type} is not supported as recon-input yet."
)

if config.workflow.input_type == "qsiprep":
# This should work for --recon-input as long as the same dataset is in bids_dir
# or if the call is doing preproc+recon
to_recon_list = config.execution.participant_label
elif config.workflow.input_type == "ukb":
from ..utils.ingress import collect_ukb_participants, create_ukb_layout

# The ukb input will always be specified as the bids input - we can't preproc it first
ukb_layout = create_ukb_layout(config.execution.bids_dir)
to_recon_list = collect_ukb_participants(
ukb_layout, participant_label=config.execution.participant_label
)
to_recon_list = config.execution.participant_label

for subject_id in to_recon_list:
single_subject_wf = init_single_subject_recon_wf(subject_id=subject_id)
Expand Down Expand Up @@ -86,7 +75,7 @@ def init_single_subject_recon_wf(subject_id):
Single subject label
"""
from ..interfaces.bids import DerivativesDataSink
from ..interfaces.ingress import QSIPrepDWIIngress, UKBioBankDWIIngress
from ..interfaces.ingress import QSIPrepDWIIngress
from ..interfaces.interchange import (
ReconWorkflowInputs,
anatomical_workflow_outputs,
Expand Down Expand Up @@ -138,17 +127,16 @@ def init_single_subject_recon_wf(subject_id):

# This is here because qsiprep currently only makes one anatomical result per subject
# regardless of sessions. So process it on its
if config.workflow.input_type == "qsiprep":
anat_ingress_node, available_anatomical_data = init_highres_recon_anatomical_wf(
subject_id=subject_id,
extras_to_make=spec.get("anatomical", []),
needs_t1w_transform=needs_t1w_transform,
)
anat_ingress_node, available_anatomical_data = init_highres_recon_anatomical_wf(
subject_id=subject_id,
extras_to_make=spec.get("anatomical", []),
needs_t1w_transform=needs_t1w_transform,
)

# Connect the anatomical-only inputs. NOTE this is not to the inputnode!
config.loggers.workflow.info(
"Anatomical (T1w) available for recon: %s", available_anatomical_data
)
# Connect the anatomical-only inputs. NOTE this is not to the inputnode!
config.loggers.workflow.info(
"Anatomical (T1w) available for recon: %s", available_anatomical_data
)

aggregate_anatomical_nodes = pe.Node(
niu.Merge(len(dwi_recon_inputs)),
Expand All @@ -168,28 +156,11 @@ def init_single_subject_recon_wf(subject_id):
wf_name = _get_wf_name(dwi_file)

# Get the preprocessed DWI and all the related preprocessed images
if config.workflow.input_type == "qsiprep":
dwi_ingress_nodes[dwi_file] = pe.Node(
QSIPrepDWIIngress(dwi_file=dwi_file),
name=f"{wf_name}_ingressed_dwi_data",
)
anat_ingress_nodes[dwi_file] = anat_ingress_node

elif config.workflow.input_type == "ukb":
dwi_ingress_nodes[dwi_file] = pe.Node(
UKBioBankDWIIngress(dwi_file=dwi_file, data_dir=str(dwi_input["path"].absolute())),
name=f"{wf_name}_ingressed_ukb_dwi_data",
)
anat_ingress_nodes[dwi_file], available_anatomical_data = (
init_highres_recon_anatomical_wf(
subject_id=subject_id,
recon_input_dir=dwi_input["path"],
extras_to_make=spec.get("anatomical", []),
pipeline_source="ukb",
needs_t1w_transform=needs_t1w_transform,
name=f"{wf_name}_ingressed_ukb_anat_data",
)
)
dwi_ingress_nodes[dwi_file] = pe.Node(
QSIPrepDWIIngress(dwi_file=dwi_file),
name=f"{wf_name}_ingressed_dwi_data",
)
anat_ingress_nodes[dwi_file] = anat_ingress_node
mattcieslak marked this conversation as resolved.
Show resolved Hide resolved

# Aggregate the anatomical data from all the dwi files
workflow.connect([
Expand Down Expand Up @@ -404,26 +375,19 @@ def _get_iterable_dwi_inputs(subject_id):
the other files needed.

"""
from ..utils.ingress import create_ukb_layout

dwi_dir = config.execution.bids_dir
if config.workflow.input_type == "qsiprep":
if not (dwi_dir / f"sub-{subject_id}").exists():
raise Exception(f"Unable to find subject directory in {config.execution.bids_dir}")

layout = BIDSLayout(dwi_dir, validate=False, absolute_paths=True)
# Get all the output files that are in this space
dwi_files = [
f.path
for f in layout.get(
suffix="dwi", subject=subject_id, absolute_paths=True, extension=["nii", "nii.gz"]
)
if "space-T1w" in f.filename
]
config.loggers.workflow.info("found %s in %s", dwi_files, dwi_dir)
return [{"bids_dwi_file": dwi_file} for dwi_file in dwi_files]

if config.workflow.input_type == "ukb":
return create_ukb_layout(ukb_dir=config.execution.bids_dir, participant_label=subject_id)

raise Exception("Unknown pipeline " + config.workflow.input_type)
if not (dwi_dir / f"sub-{subject_id}").exists():
raise Exception(f"Unable to find subject directory in {config.execution.bids_dir}")

layout = BIDSLayout(dwi_dir, validate=False, absolute_paths=True)
# Get all the output files that are in this space
dwi_files = [
f.path
for f in layout.get(
suffix="dwi", subject=subject_id, absolute_paths=True, extension=["nii", "nii.gz"]
)
if "space-T1w" in f.filename
]
config.loggers.workflow.info("found %s in %s", dwi_files, dwi_dir)
return [{"bids_dwi_file": dwi_file} for dwi_file in dwi_files]
Loading
Loading