MRtrix pipeline_2

s='5770';

setenv('s', s)

workingdir=['/Volumes/Data2/Analysis/Monk/analysis/dti/' s];

cd(workingdir)

%% Generation of tissue-segmented image for "anatomically constrained tractography"

img_t1w=[s '_nu_brain_mni_warped.nii.gz'];

cmd=['5ttgen fsl ' img_t1w ' mr_5TT.mif -premasked']

unix(cmd)

%% labelconvert

if ~exist ('mr_parcels.mif'),

lut=fullfile('/Volumes/Data2/Analysis/Monk/analysis/dti','mr_roi_lut.txt');

img_parcels=fullfile(workingdir,'dmri','roi','infant-neo-aal_warped_diff.nii.gz');

cmd=['labelconvert ' img_parcels ' ' lut ' ' lut ' mr_parcels.mif'];

unix(cmd)

end

%% setup

img_dwi=fullfile(workingdir,'dmri','dwi.nii.gz');

% img_dwi_denoised=fullfile(workingdir,'dwi_denoised.nii.gz');

img_bvecs=fullfile(workingdir,'dmri','bvecs');

img_bvals=fullfile(workingdir,'dmri','bvals');

%% DWI processing2-converting nifti to mif

cmd = (['mrconvert ' img_dwi ' -force mr_dwi.mif -fslgrad 'img_bvecs ' ' img_bvals ' -datatype float32 -stride 0,0,0,1'])

unix(cmd)

%% denoising

cmd = (['dwidenoise mr_dwi.mif -force mr_dwi_denoised.mif -nthreads 4']);

unix(cmd)

%% mask and bias field correction

% img_dwi_biasCorr=fullfile(workingdir,'dwi_denoised_biasCorr.nii.gz');

unix(['dwi2mask mr_dwi_denoised.mif - | maskfilter erode -npass 5 -force mr_eroded_mask.mif'])

cmd = (['dwibiascorrect mr_dwi_denoised.mif -force mr_dwi_denoised_biasCorr.mif -fsl -mask mr_eroded_mask.mif -fslgrad ' img_bvecs ' ' img_bvals ' -nthreads 4']);

unix(cmd)

 

%% DWI processing3-bad vol exclusion

%%%% if stripe noise is found in a volume, make sure you delete %%%%

cmd=['mrconvert -coord 3 1:5,7:end mr_dwi_denoised_biasCorr.mif mr_dwi_denoised_biasCorr_reduced.mif'];

unix(cmd)

%index=1:14; badindex=str2num(Index_of_bad_volumes1{1}); index(badindex)=[];

%cmd_test=['echo mrconvert -coord 3 ' num2str(index) ' mr_DWI.mif mr_DWI_reduced.mif'];

 

 

%%%% generating b0

unix('dwiextract mr_dwi_denoised_biasCorr_reduced.mif - -bzero | mrmath - mean mr_meanb0.mif -axis 3')

 

%% dwi2response

% dwi2response tournier <Input DWI> <Output response text file>

% shview <Output response text file>

%%%% eroded mask

% cmd=['dwi2mask mr_dwi_denoised_reduced.mif - | maskfilter - erode -npass 5 mr_eroded_mask.mif']

% unix(cmd)

 

%if ~exist('mr_dwi_response_tournier.txt'),

cmd=['dwi2response tournier -mask mr_eroded_mask.mif -voxels mr_voxels_reduced_eroded.mif mr_dwi_denoised_biasCorr_reduced.mif mr_dwi_response_tournier_reduced_eroded.txt'];

unix(cmd)

%end

 

%% upsampling

filelist={'mr_dwi_denoised_biasCorr_reduced','mr_eroded_mask','mr_meanb0'};

parfor i=1:length(filelist)

    im=[filelist{i} '.mif'];

    newim=[filelist{i} '_upsample.mif'];

    cmd=['mrresize ' im ' -scale 2.0 ' newim];

    unix(cmd)

end

%% FOD

% make sure to use "DILATED MASK" for FOD generation

unix(['dwi2mask mr_dwi_denoised_biasCorr_reduced.mif - | maskfilter - dilate -npass 5 -force mr_dilate_mask.mif'])

unix(['mrresize mr_dilate_mask.mif -scale 2.0 mr_dilate_mask_upsample.mif'])

cmd1=['dwi2fod csd mr_dwi_denoised_biasCorr_reduced_upsample.mif mr_dwi_response_tournier_reduced_eroded.txt mr_fod_upsample.mif -mask mr_eroded_mask_upsample.mif -nthreads 4'];

unix(cmd1)

%% streamline tractography

cmd=['tckgen mr_fod_upsample.mif mr_fod_upsample_tckgen_100M.tck -seed_image mr_parcels.mif -mask mr_eroded_mask_upsample.mif -number 100M -maxlength 250 -nthreads 4'];

unix(cmd)

zip('mr_fod_upsample_tckgen.zip', 'mr_fod_upsample_tckgen.tck')

%% sift

% example: tcksift 100M.tck WM_FODs.mif 10M_SIFT.tck -act 5TT.mif -term_number 10M