Source code for thesis.workflows.mrtrix3.operations.density_map
"""Density-map post-processing for MRtrix3 tractography outputs.
Bundles the steps that do not have Nipype interfaces:
* ``run_tcksift2_task`` — wraps ``tcksift2`` (per-streamline weighting)
and captures the proportionality coefficient ``mu`` via ``-out_mu``.
``mu`` is required to make SIFT2 weight sums comparable across
subjects and methods (units: fibre cross-sectional area).
* ``rename_density_map_task`` — moves the ``ComputeTDI`` output to the
``fdt_paths.nii.gz`` filename consumed by the existing tract_similarity
workflow, preserving compatibility with the ProbTrackX2 layout.
* ``write_waytotal_task`` — writes a single-line ``waytotal`` file. When a
SIFT2 weights file is supplied, the value is the sum of those weights
(matches the units of the SIFT2-weighted TDI). Otherwise the raw
streamline count from ``tckinfo`` is written.
* ``read_sift2_mu`` — module-level helper that reads ``sift2_mu.txt``
produced by ``tcksift2 -out_mu``.
"""
# NOTE: Function-node bodies run in subprocesses and cannot pickle loguru;
# diagnostics use ``print``.
[docs]
def run_tcksift2_task(
tracks_file: str,
fod_file: str,
fivett_file: str,
output_dir: str,
) -> tuple[str, str]:
"""Run ``tcksift2`` to produce per-streamline weights and mu.
Args:
tracks_file: Streamline file (``tracks.tck``).
fod_file: WM FOD image used by SIFT2 (typically the mtnormalise output).
fivett_file: 5TT image for ACT-aware proportionality.
output_dir: Directory where ``sift2_weights.txt`` and
``sift2_mu.txt`` are written.
Returns:
``(weights_file, mu_file)`` — absolute paths to the per-streamline
weights file and the single-float ``mu`` proportionality
coefficient file.
"""
import shutil
import subprocess
from pathlib import Path
out_dir = Path(output_dir)
out_dir.mkdir(parents=True, exist_ok=True)
weights = out_dir / "sift2_weights.txt"
mu_path = out_dir / "sift2_mu.txt"
if shutil.which("tcksift2") is None:
raise RuntimeError("tcksift2 not found on PATH; ensure MRtrix3 is installed.")
cmd = [
"tcksift2",
str(tracks_file),
str(fod_file),
str(weights),
"-act",
str(fivett_file),
"-out_mu",
str(mu_path),
"-force",
]
print(f"[mrtrix3] {' '.join(cmd)}")
# capture_output so the actual MRtrix error (e.g. "no streamlines to
# weight" when tracks.tck is empty) reaches the user instead of a bare
# CalledProcessError. Trial 2026-05-05_trials_30 spent a day finding
# the cause precisely because stderr was swallowed here.
result = subprocess.run(cmd, capture_output=True, text=True)
if result.returncode != 0:
raise RuntimeError(
f"tcksift2 exited with code {result.returncode}.\n"
f"Command: {' '.join(cmd)}\n"
f"stderr:\n{result.stderr}\n"
f"stdout (tail):\n{(result.stdout or '')[-2000:]}"
)
return str(weights), str(mu_path)
[docs]
def read_sift2_mu(mu_file):
"""Read the single-float ``mu`` proportionality coefficient.
Args:
mu_file: Path to ``sift2_mu.txt`` written by ``tcksift2 -out_mu``.
Returns:
``mu`` as a float, or ``None`` if the file is missing or empty.
"""
from pathlib import Path
if mu_file is None:
return None
p = Path(mu_file)
if not p.is_file():
return None
try:
text = p.read_text().strip()
except OSError:
return None
if not text:
return None
for token in text.split():
try:
return float(token)
except ValueError:
continue
return None
[docs]
def rename_density_map_task(in_file: str, output_dir: str) -> str:
"""Copy the TDI output to ``fdt_paths.nii.gz`` for downstream compatibility.
The MRtrix tckmap interface writes to a temporary path inside the
Nipype run-dir; downstream consumers (tract_similarity) expect the
file to live next to a ``waytotal`` under the run output directory.
Args:
in_file: Path to the ComputeTDI output NIfTI.
output_dir: Directory where ``fdt_paths.nii.gz`` is written.
Returns:
Absolute path to ``fdt_paths.nii.gz``.
"""
import shutil
from pathlib import Path
out_dir = Path(output_dir)
out_dir.mkdir(parents=True, exist_ok=True)
target = out_dir / "fdt_paths.nii.gz"
shutil.copyfile(str(in_file), str(target))
return str(target)
[docs]
def write_waytotal_task(
tracks_file: str,
output_dir: str,
weights_file: str | None = None,
) -> str:
"""Write a single-line ``waytotal`` file matching the TDI denominator.
The downstream atlas / tract_similarity readers divide the TDI
(``fdt_paths.nii.gz``) by this value to get a streamline fraction. The
written value therefore must match the units of the TDI numerator:
* When ``weights_file`` points at a SIFT2 weights text file
(``sift2_weights.txt``), the TDI was generated with
``-tck_weights_in`` and each voxel holds a sum of weights. The
correct denominator is the sum of all per-streamline weights; we
compute it here.
* When ``weights_file`` is ``None`` or missing, the TDI is a raw
streamline count, so we fall back to ``tckinfo -count``.
Args:
tracks_file: Path to ``tracks.tck``.
output_dir: Directory where ``waytotal`` is written.
weights_file: Optional path to a SIFT2 weights file (one float per
line). When supplied and readable, the sum is written instead
of the raw streamline count.
Returns:
Absolute path to the written ``waytotal`` file.
"""
import shutil
import subprocess
from pathlib import Path
out_dir = Path(output_dir)
out_dir.mkdir(parents=True, exist_ok=True)
target = out_dir / "waytotal"
if weights_file:
weights_path = Path(weights_file)
if weights_path.is_file():
total = 0.0
with weights_path.open() as fh:
for line in fh:
for token in line.split():
try:
total += float(token)
except ValueError:
continue
target.write_text(f"{total:.6f}\n")
return str(target)
print(f"[mrtrix3] weights_file {weights_file!r} missing; using tckinfo count")
count = 0
if shutil.which("tckinfo") is not None:
try:
result = subprocess.run(
["tckinfo", str(tracks_file), "-count"],
check=True,
capture_output=True,
text=True,
)
for line in (result.stdout or "").splitlines():
stripped = line.strip()
if not stripped:
continue
last_token = stripped.split()[-1]
if last_token.isdigit():
count = int(last_token)
except (subprocess.CalledProcessError, ValueError) as exc:
print(f"[mrtrix3] tckinfo failed for {tracks_file}: {exc}")
target.write_text(f"{count}\n")
return str(target)
