From 56ea4d0d897831da04d80bcba026d0d91bc8291c Mon Sep 17 00:00:00 2001
From: Fidel Alfaro Almagro <falmagro@fmrib.ox.ac.uk>
Date: Fri, 24 Feb 2023 17:47:44 +0000
Subject: [PATCH] Adding last IDP scripts
---
bip/data/FileTree.tree | 4 +-
bip/data/IDPs/IDPs.json | 3 +
bip/pipelines/IDPs_gen/IDP_diff_autoPtx.py | 130 +++++++--------------
bip/pipelines/IDPs_gen/IDPs_gen.py | 24 +++-
bip/pipelines/dMRI_diff/diff_autoptx.py | 12 +-
5 files changed, 73 insertions(+), 100 deletions(-)
create mode 100755 bip/data/IDPs/IDPs.json
diff --git a/bip/data/FileTree.tree b/bip/data/FileTree.tree
index d7f4238..d87f4a5 100644
--- a/bip/data/FileTree.tree
+++ b/bip/data/FileTree.tree
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:30abf95478fe8611d0162a4bce17d4cabae1471ea0d322398d1071a8417e0009
-size 23361
+oid sha256:6f9005fcde92e75653faa54de4e45f5bb29367989e01dadd2344bf52a67bd74e
+size 23447
diff --git a/bip/data/IDPs/IDPs.json b/bip/data/IDPs/IDPs.json
new file mode 100755
index 0000000..0b5725d
--- /dev/null
+++ b/bip/data/IDPs/IDPs.json
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:8005353f07fe3d5cf09f29f9faf45e42182ab1317f45d1370cb0e3623a9883c9
+size 244014
diff --git a/bip/pipelines/IDPs_gen/IDP_diff_autoPtx.py b/bip/pipelines/IDPs_gen/IDP_diff_autoPtx.py
index 786285e..31b3d40 100755
--- a/bip/pipelines/IDPs_gen/IDP_diff_autoPtx.py
+++ b/bip/pipelines/IDPs_gen/IDP_diff_autoPtx.py
@@ -1,6 +1,6 @@
#!/usr/bin/env python
#
-# IDP_diff_autoPtx.py - Generating IDP file with autoPtx metrics
+# IDP_diff_autoPtx.py - Generating IDP file with autoPtx metrics
#
# Author: Fidel Alfaro Almagro <fidel.alfaroalmagro@ndcn.ox.ac.uk>
# Author: Paul McCarthy <pauldmccarthy@gmail.com>
@@ -10,111 +10,59 @@
#
import os
-import json
import logging
-from pipe_tree import In, Out, Ref
-from bip.utils.log_utils import redirect_logging
+import nibabel as nib
+from fsl import wrappers
+from pipe_tree import In, Out, Ref, Var
+from bip.utils.log_utils import redirect_logging, tempdir
log = logging.getLogger(__name__)
def run(ctx,
- T1_orig: In,
logs_dir: Ref,
- IDP_diff_autoPtx: Out):
-
- with redirect_logging('IDP_diff_autoPtx', outdir=logs_dir):
+ TBSS_prefix: Ref,
+ aptx_txt_prefix: Ref,
+ tmp_dir: Ref,
+ aptx_tract_tmp: In(no_iter=True),
+ IDP_diff_autoPtx: Out):
+ with redirect_logging('IDP_diff_autoPtx', outdir=logs_dir),\
+ tempdir(tmp_dir):
-set -x
+ autoPtx_all = tmp_dir + '/autoPtx_all.nii.gz'
+ autoPtx_tmp = tmp_dir + '/autoPtx_tmp.nii.gz'
-origDir=`pwd`
-scriptName=`basename "$0"`
-direc=$1
+ result_NaN = ("NaN " * 27).strip()
+ result = ""
-cd $direc
+ # TODO: See if it is possible to generate the aptx_tract_tmp files here
+ wrappers.fslmerge(autoPtx_all, *list(aptx_tract_tmp.data))
+ tractnorm = wrappers.fslstats(autoPtx_all, t=True).m.run()
-basedMRI="dMRI"
-#if [ -d $basedMRI/unusable ] ; then
-# basedMRI="$basedMRI/unusable"
-#elif [ -d $basedMRI/incompatible ] ; then
-# basedMRI="$basedMRI/incompatible"
-#fi
+ for suffix in ['FA','L1','L2','L3','MO','MD','ICVF','OD','ISOVF']:
+ all_file = TBSS_prefix + suffix + ".nii.gz"
+ APTX_file = aptx_txt_prefix + suffix + ".txt"
-#Setting the string of NaN in case there is a problem.
-numVars="27"
-nanResult="";
-for i in $(seq 1 $numVars) ; do
- nanResult="NaN $nanResult" ;
-done
+ if os.path.exists(all_file):
+ wrappers.fslmaths(autoPtx_all).mul(all_file).run(autoPtx_tmp)
+ tractmeas = wrappers.fslstats(autoPtx_tmp, t=True).m.run()
-result=""
+ img = nib.load(autoPtx_all)
+ jmax = int(round(img.header['dim'][4]))
-allgood=1
-for i in FA MD MO L1 L2 L3 ICVF OD ISOVF ; do
- if [ ! -f $basedMRI/autoptx_preproc/APTX_${i}.txt ] || [ `cat $basedMRI/autoptx_preproc/APTX_${i}.txt | wc -w` != 27 ] ; then
- allgood=0
- fi
-done
+ cad = ""
+ for i in range(jmax):
+ cad += " " + str(tractmeas[i] / tractnorm[i])
+ cad = cad.strip()
-if [ $allgood = 0 ] ; then
- if [ -f $basedMRI/autoptx_preproc/tracts/mcp/tracts/tractsNorm.nii.gz ] ; then
+ with open(APTX_file, 'wt', encoding="utf-8") as f:
+ f.write(f'{cad}\n')
- ntracts=`echo $basedMRI/autoptx_preproc/tracts/*/tracts/tractsNorm.nii.gz | wc -w`
- if [ $ntracts != 27 ] ; then
- >&2 echo only found $ntracts tracts
- fi
-
- # hack to ameliorate partial misalignment of mcp
- fslmaths $FSLDIR/data/standard/FMRIB58_FA_1mm -thr 1000 -bin -ero -ero -mul $basedMRI/autoptx_preproc/tracts/mcp/tracts/tractsNorm $basedMRI/autoptx_preproc/tracts/mcp/tracts/tractsNorm -odt float
-
- for p in ar_l ar_r atr_l atr_r cgc_l cgc_r cgh_l cgh_r cst_l cst_r fma fmi ifo_l ifo_r ilf_l ilf_r mcp ml_l ml_r ptr_l ptr_r slf_l slf_r str_l str_r unc_l unc_r ; do
- if [ -f $basedMRI/autoptx_preproc/tracts/$p/tracts/tractsNorm.nii.gz ] ; then
- t=`fslstats $basedMRI/autoptx_preproc/tracts/$p/tracts/tractsNorm -P 99`
- fslmaths $basedMRI/autoptx_preproc/tracts/$p/tracts/tractsNorm -min $t -mul 1e6 -div $t -thr 1e5 $basedMRI/autoptx_preproc/tmp_$p
- fi
- done
-
- cad=""
-
- fslmerge -t $basedMRI/autoptx_preproc/tmp_all \
- $basedMRI/autoptx_preproc/tmp_{ar_l,ar_r,atr_l,atr_r,cgc_l,cgc_r,cgh_l,cgh_r,cst_l,cst_r,fma,fmi,ifo_l,ifo_r,ilf_l,ilf_r,mcp,ml_l,ml_r,ptr_l,ptr_r,slf_l,slf_r,str_l,str_r,unc_l,unc_r}.nii.gz
- tractnorm=`fslstats -t $basedMRI/autoptx_preproc/tmp_all -m`
-
- if [ -f $basedMRI/autoptx_preproc/tmp_all.nii.gz ] ; then
- for i in FA MD MO L1 L2 L3 ICVF OD ISOVF ; do
- /bin/rm -f $basedMRI/autoptx_preproc/APTX_${i}.txt
- if [ -f $basedMRI/TBSS/stats/all_${i}.nii.gz ] ; then
- fslmaths $basedMRI/autoptx_preproc/tmp_all -mul $basedMRI/TBSS/stats/all_${i} $basedMRI/autoptx_preproc/tmp
- tractmeas=`fslstats -t $basedMRI/autoptx_preproc/tmp -m | sed 's/-/_/g'`
- jmax=`fslnvols $basedMRI/autoptx_preproc/tmp_all`
- for (( j=1; j <= $jmax; ++j )) ; do
- tm=`echo $tractmeas | cut -d " " -f$j`
- tn=`echo $tractnorm | cut -d " " -f$j`
- echo "10 k $tm $tn / p" | dc - >> $basedMRI/autoptx_preproc/APTX_${i}.txt
- done
- fi
- done
- rm $basedMRI/autoptx_preproc/tmp_*.nii.gz
- fi
- else
- >&2 echo "mcp tract missing"
- fi
-fi
-
-for i in FA MD MO L1 L2 L3 ICVF OD ISOVF ; do
- if [ -f $basedMRI/autoptx_preproc/APTX_${i}.txt ] && [ `cat $basedMRI/autoptx_preproc/APTX_${i}.txt | wc -w` = 27 ] ; then
- miniResult=`cat $basedMRI/autoptx_preproc/APTX_${i}.txt`
- else
- miniResult="$nanResult"
- fi
- result="$result $miniResult"
-done
-
-mkdir -p IDP_files
-
-echo $result > IDP_files/$scriptName.txt
-echo $result
-
-cd $origDir
+ result += cad
+ else:
+ result += " " + result_NaN
+ result = result.strip()
+ with open(IDP_diff_autoPtx, 'wt', encoding="utf-8") as f:
+ f.write(f'{result}\n')
diff --git a/bip/pipelines/IDPs_gen/IDPs_gen.py b/bip/pipelines/IDPs_gen/IDPs_gen.py
index 769243a..a1b41ce 100755
--- a/bip/pipelines/IDPs_gen/IDPs_gen.py
+++ b/bip/pipelines/IDPs_gen/IDPs_gen.py
@@ -11,17 +11,24 @@
import logging
from bip.utils.log_utils import redirect_logging
-from bip.pipelines.IDPs_gen import IDP_subject_ID, IDP_subject_centre
+from bip.pipelines.IDPs_gen import IDP_subject_ID
+from bip.pipelines.IDPs_gen import IDP_subject_centre
from bip.pipelines.IDPs_gen import IDP_subject_COG_table
from bip.pipelines.IDPs_gen import IDP_all_align_to_T1
-from bip.pipelines.IDPs_gen import IDP_T1_FIRST_vols, IDP_T1_SIENAX
-from bip.pipelines.IDPs_gen import IDP_T1_GM_parcellation, IDP_T1_SIENAX
-from bip.pipelines.IDPs_gen import IDP_T1_align_to_std, IDP_T1_noise_ratio
+from bip.pipelines.IDPs_gen import IDP_T1_FIRST_vols
+from bip.pipelines.IDPs_gen import IDP_T1_SIENAX
+from bip.pipelines.IDPs_gen import IDP_T1_GM_parcellation
+from bip.pipelines.IDPs_gen import IDP_T1_SIENAX
+from bip.pipelines.IDPs_gen import IDP_T1_align_to_std
+from bip.pipelines.IDPs_gen import IDP_T1_noise_ratio
from bip.pipelines.IDPs_gen import IDP_T2_FLAIR_WMH
from bip.pipelines.IDPs_gen import IDP_SWI_T2star
-from bip.pipelines.IDPs_gen import IDP_func_head_motion, IDP_func_TSNR
+from bip.pipelines.IDPs_gen import IDP_func_head_motion
+from bip.pipelines.IDPs_gen import IDP_func_TSNR
from bip.pipelines.IDPs_gen import IDP_func_task_activation
-from bip.pipelines.IDPs_gen import IDP_diff_eddy_outliers, IDP_diff_TBSS
+from bip.pipelines.IDPs_gen import IDP_diff_eddy_outliers
+from bip.pipelines.IDPs_gen import IDP_diff_TBSS
+from bip.pipelines.IDPs_gen import IDP_diff_autoPtx
log = logging.getLogger(__name__)
@@ -102,6 +109,11 @@ def add_to_pipeline(ctx, pipe, tree, targets):
kwargs={'ctx' : ctx})
targets.append('IDP_diff_TBSS')
+ pipe(IDP_diff_autoPtx.run,
+ submit=dict(jobtime=200, name="BIP_IDP_diff_autoPtx_"+ subj),
+ kwargs={'ctx' : ctx})
+ targets.append('IDP_diff_autoPtx')
+
pipe(IDP_subject_COG_table.run,
submit=dict(jobtime=200, name="BIP_IDP_subject_COG_table_" + subj),
kwargs={'ctx' : ctx})
diff --git a/bip/pipelines/dMRI_diff/diff_autoptx.py b/bip/pipelines/dMRI_diff/diff_autoptx.py
index 894dc2d..017bd35 100755
--- a/bip/pipelines/dMRI_diff/diff_autoptx.py
+++ b/bip/pipelines/dMRI_diff/diff_autoptx.py
@@ -39,7 +39,8 @@ def run(ctx,
aptx_exclude: Out,
density: Out,
waytotal: Out,
- tractsNorm: Out):
+ tractsNorm: Out,
+ aptx_tract_tmp: Out):
with redirect_logging('diff_autoptx_' + autoptx_tract.value,
outdir=logs_dir):
@@ -49,6 +50,7 @@ def run(ctx,
t_seeds = int(ctx.tract_struct[t_name]["num_seeds"] * 300)
orig_tract_dir = ctx.get_data('dMRI/autoptx/protocols/' + t_name + '/')
+ FMRIB58_FA_1mm = ctx.get_standard('FMRIB58_FA_1mm.nii.gz')
# Does the protocol defines a second
# run with inverted seed / target masks?
@@ -120,3 +122,11 @@ def run(ctx,
f.write(f'{way}\n')
wrappers.fslmaths(density).div(way).range().run(tractsNorm)
+
+ # Hack to ameliorate partial misalignment of mcp
+ if t_name == "mcp":
+ wrappers.fslmaths(FMRIB58_FA_1mm).thr(1000).bin().ero(2).mul(tractsNorm).run(tractsNorm)
+
+ # This will be needed for IDP generation
+ val = wrappers.fslstats(tractsNorm).P(99).run()
+ wrappers.fslmaths(tractsNorm).min(val).mul(1e6).div(val).thr(1e5).run(aptx_tract_tmp)
--
GitLab