Matrix | Python Fiddle

import os
import subprocess

rootd='/home/Benson/Documents/Data_Analysis/NKI/Data/1898228/'
os.chdir(rootd)
%% Prepare freesurfer for Tractography
os.mkdir(anat)
subprocess.call(['mri_convert', 'mri/rawavg.mgz',  'anat/str.nii.gz'])
subprocess.call(['mri_convert', 'mri/orig.mgz',  'anat/fs.nii.gz'])
subprocess.call(['mri_binarize', '--i', 'mri/aparc+aseg.mgz', '--ventricles', '--o', 'anat/ventricles.nii.gz'])
subprocess.call(['mri_binarize', '--i', 'mri/aparc+aseg.mgz', '--match', '2', '--o', 'anat/wm.lh.nii.gz'])
subprocess.call(['mri_binarize', '--i', 'mri/aparc+aseg.mgz', '--match', '41', '--o', 'anat/wm.rh.nii.gz'])
%% Binarized wm filenames are put into txt files
subprocess.call(['ls', '-1', 'anat/wm* > waypoint_files.txt']);
%% Copy over label files & white surfaces
(['rsync label/*.label label/']); unix(['rsync surf/{l,r}h.white surf/']);

%% File transformation (Generation of transformation matrices + file transformation
unix(['fs2str=dmri.bedpostX/xfms/fs2str.mat']);unix(['str2fs=dmri.bedpostX/xfms/str2fs.mat']);
unix(['fa2fs=dmri.bedpostX/xfms/fa2fs.mat']);unix(['fs2fa=dmri.bedpostX/xfms/fs2fa.mat']);
unix(['fa2str=dmri.bedpostX/xfms/fa2str.mat']);unix(['str2fa=dmri.bedpostX/xfms/str2fa.mat']);
%% 1 - Register structurual to Freesurfer
%% $- all defined above!
% unix(['tkregister2 --mov anat/fs.nii.gz --targ anat/str.nii.gz --regheader --reg /tmp/junk --fslregout $fs2str  --noedit']);
%% a- Create freesurfer to anatomical space
unix([' tkregister2 --mov mri/orig.mgz --targ mri/rawavg.mgz --regheader --reg junk --fslregout dmri.bedpostX/xfms/fs2str.mat --noedit ']);
%% b - Create Inverse matrices.
% # invert to create str2fs
unix(['convert_xfm -omat dmri.bedpostX/xfms/str2fs.mat -inverse dmri.bedpostX/xfms/fs2str.mat ']);
%% c- Now transform FA to structural:
unix(['flirt -in dmri/dtifit_FA.nii.gz -ref anat/str.nii.gz -omat dmri.bedpostX/xfms/fa2str.mat -dof 6']);
%% d invert to create str2fa
unix(['convert_xfm -omat dmri.bedpostX/xfms/str2fa.mat -inverse dmri.bedpostX/xfms/fa2str.mat']);
%% e  Concatenation and inverse
% # Concatenate and inverse
unix(['convert_xfm -omat dmri.bedpostX/xfms/fa2fs.mat -concat dmri.bedpostX/xfms/str2fs.mat dmri.bedpostX/xfms/fa2str.mat']);
unix(['convert_xfm -omat dmri.bedpostX/xfms/fs2fa.mat -inverse dmri.bedpostX/xfms/fa2fs.mat']);
%% Generate label files.
unix(['mri_annotation2label --subject 1898228 --hemi lh --labelbase lhtracto']);
unix(['mri_annotation2label --subject 1898228 --hemi rh --labelbase rhtracto']);
nodes1=dir('label/lhtracto*');
%% Generate  annotation files
for i=1:length(nodes1)
unix(['mris_label2annot --s 1898228 --h lh --l label/', nodes1(i).name , '  --ctab label/aparc.annot.ctab  --a ' , nodes1(i).name(1:12)]);
clc
end
nodes2=dir('label/rhtracto*');
for i=1:length(nodes2)
unix(['mris_label2annot --s 1898228 --h rh --l label/', nodes2(i).name , '  --ctab label/aparc.annot.ctab  --a ' , nodes2(i).name(1:12)]);
clc
end

%% Annotation to volume
%Left hemisphere
for i=1:length(nodes1)
k1=deblank([nodes1(i).name(1:12),'.annot']); k2=deblank([nodes1(i).name(1:12),'.nii.gz']);
unix(['mri_label2vol --annot label/lh.', k1 , ' --temp anat/fs.nii.gz --subject 1898228 --hemi lh --identity --fillthresh .5 --o anat/',k2 , ' --proj frac 0 1 .1']);
end
%Right hemisphere
for i=1:length(nodes2)-1
k1=deblank([nodes2(i).name(1:12),'.annot']); k2=deblank([nodes2(i).name(1:12),'.nii.gz']);
unix(['mri_label2vol --annot label/rh.', k1 , ' --temp anat/fs.nii.gz --subject 1898228 --hemi rh --identity --fillthresh .5 --o anat/',k2 , ' --proj frac 0 1 .1']);
end

%%
gz=deblank([nodes1(1).name(1:12),'.nii.gz']);

probtrackx2 --network -x anat/masks.txt -V 1 --onewaycondition -c 0.2 -S 2000 --steplength=0.5 -P 5000 --fibthresh=0.01 --distthresh=0.0 --sampvox=0.0 --forcedir --opd -s merged -m dmri.bedpostX/nodif_brain_mask --xfm=dmri.bedpostX/xfms/fs2fa.mat --dir=trackresults --targetmasks=seeds.txt --avoid=anat/ventricles.nii.gz --seedref=anat/fs.nii.gz -s dmri.bedpostX/merged -m dmri.bedpostX/nodif_brain_mask -l --usef --s2tastext  --os2t

% Normalize with waytotal
load networkmatrix.mat
load waytotal.mat
for i=1:length(networkmatrix)
row(i,:)=networkmatrix(i,:)./waytotal(i);
end

Python Fiddle

Python Cloud IDE