/NII_RES/maps.
Package: mdm
Multidimensional data management (MDM). This package contains functions
for building the local data structure.
We work with a structure often referred to as 's'. This structure has the
following fields
s.nii_fn - full path to a nifti file
s.xps - the eXperimental Parameter Structure
s.mask_fn - (optional) path to a nifti file with a 3D mask
that determines which parts of the data that contains
actual data and not just background
EXPERIMENTAL PARAMETER STRUCTURE (xps)
The xps has the following fields. All are in SI units. Not all fields must
be present. The model code checks that the necessary fields are present.
All parameters in the xps should be are stored as variables of size n x m,
where n is the number of image volumes and m is the number of parameters.
For example, a 30 direction DTI set with 6 b0 images would render an xps
with many fields, for example, a 'bt' field of size 36 x 6, where 36
is the total number of images and 6 are the number of parameters required
to describe the b-tensor.
A) General parameters
- n: number of images/signal values (fourth dimension). All
parameters need to have n entries.
- t_ex: cumulative time from start of the experiment to the time of
the excitation pulse
- t_acq: time from excitation to readout of echo, i.e., k = 0
(typically equals te in a SE sequence)
- intention: a string that describes the intention of the experiment,
for example
- PGSE/DTI
- DDE/uFA
- Angular-DDE
- FEXI
- slice_order: not yet defined
- a_ind: Averaging index. After averaging (arithmetic or geometric),
there will be max(a_ind) number of images left.
- s_ind: Series index. Refers to data acquired in different series,
for example, with different prescans (e.g. gain adjustment).
Can also index acquisitions with different echo times et c.
Potentially but not necessarily automatically calculated:
- b_ind: Indexes measurements according to total b-values
- bd_ind: Indexes measurements according to b-anisotropy (b_delta)
- be_ind: Indexes measurements according to b-asymmetry (b_eta)
- br_ind: Indexes measurements according to b-tensor rotations
B) Parameters describing the total diffusion encoding effects
- b: total b-value
- bt: b-matrix, n x 6 (see tm_* for format)
- alpha: flow compensation factor (see Ahlgren 16)
- alpha2: flow attenuation factor (see Ahlgren 16)
i) Parameters derived automatically:
- bt2: outer product of b-matrix (fourth order tensor)
- u: symmetry axis of bt, n x 3
- b_delta: anisotropy of b-tensor
- b_eta: asymmetry of b-tensor
- b_s: spherical component of b-tensor (Martins et al, PRL 16)
- b_p: planar component of b-tensor (Martins et al, PRL 16)
- b_l: linear component of b-tensor (Martins et al, PRL 16)
- b_paszz: b-eigenvalue furthest from the mean ("symmetry axis").
- b_pasyy: b-eigenvalue closest from the mean ("symmetry axis")
- b_pasxx: b-eigenvalue that is not zz or yy (Eriksson et al, JCP 15)
- b_alpha: Euler angle of the b-tensor, from LAB TO PAS
- b_beta: Euler angle of the b-tensor
- b_gamma: Euler angle of the b-tensor
ii) Fields needed where we have not yet decided on a format
- gradient waveform
- nex? number of averages per acquisition, or perhaps also noise sigma?
C) Diffusion-related parameters in a multiple diffusion encoding setting
i) The following parameters are valid for SDE, DDE et c
- mde_delta1, mde_delta2: Duration of diffusion encoding
gradients.
Note: For SDE, we have only mde_delta1
- mde_capital_delta1, delta2, et c: Time between leading edged of encoding
gradients. For SDE, only
mde_capital_delta1 is defined
- mde_ramp_time Assumed to be the same throughout
- mde_g1, mde_g2: n x 3 vector, gradient amplitude
Derived parameters:
- mde_q1, mde_q2, ...: q-vectors
- mde_td1, mde_td2, ...: Diffusion times
ii) Parameters valid only for DDE, TDE or more
- mde_tm12, mde_tm23, ...: diffusion-related mixing times between
diffusion-encoding block 1 and 2, block 2 and 3,
in a DDE, TDE, et c setting. In DDE-based FEXI,
we have only mde_tm1.
- mde_b1, mde_b2, ...: b-value per block in a DDE, TDE et c sequence
- mde_bt1, mde_bt2, ...: b-tensor per block in a DDE, TDE et c sequence
Optional, may be calculated automatically:
- mde_tm12_ind: Index according to mixing times
- mde_b1_ind, mde_b2_ind: Index according to b1-values, b2-values, ...
D) Relaxation parameters dealing with total relaxation weighting
- te: echo time, total time with T2 weighting from excitation to readout
- tm: mixing time, total time with T1 relaxation from excitation to
readout
- ts: saturation recovery time, time between magnetization was zero and
excitation
- ti: inversion time, time between inversion of magnetization and
excitation
- tr: repetition time, time between excitations
In addition, we see the need for managing more complex sequence with
multiple RF pulses. We suggest to store such timings in fields called
- ste_tm: a vector of mixing times in a STE experiment
- ste_te: a vector of echo times in a STE experiment
Functions
function res = mdm_bruker_acqus2mat(data_path)
% Read Bruker acquisition parameters in directory data_path
% Includes acqus, acqu2s, and (if defined) vd, vc, vp, and fq1 lists.
% Save experimental parameter structures NMRacqus and NMRacqu2s in
% data_path.
function xps = mdm_bruker_dt_axderare2d_acqus2xps(data_path, xps_fn)
% Calculation of b-tensors for the Bruker pulse program DT_axderare2d
% as used in Topgaard, Phys. Chem. Chem. Phys., 2016.
% http://dx.doi.org/10.1039/c5cp07251d
%
% data_path: directory for gradient text files
% xps_fn (optional) filename of xps
function out_path = mdm_bruker_dt_rare2d_recon(data_path, out_path, rps, opt)
% Recon images acquired with Bruker pulse sequence DT_axderare2d.
%
% Works in progress
function mdm_bruker_dt_rare2d_ser2nii(data_path, nii_fn, rps)
% image reconstruction from Bruker DT_**rare2d pulse programs
% saves complex image as nifti file
% image resolution in field n.pixdim in nifti header
%
% data_path: folder where the Bruker ser file is located
% nii_fn: nifti file name (including complete path and extension)
% rps: image recon parameters structure
% rps.smooth : Gaussian smooting [m]
% rps.npix.read : image size in read dimension
% rps.npix.phase : image size in phase dimension
function gamma = mdm_bruker_gamma(NMRacqus)
% NMRacqus: Bruker acquistion parameter structure
% NMRacqus.nuc1: nucleus in acquisiton channel 1
function g = mdm_bruker_grad_read(fn)
% Read Bruker gradient shape file.
% fn: file name
% g: column vector with gradient values
function Gmax = mdm_bruker_maxgradient(NMRacqus)
% NMRacqus: Bruker acquistion parameter structure
% NMRacqus.probhd: probehead name
% Max gradients for the Bruker imaging probes in Lund
function [o_fn, tpm_fn] = mdm_coreg(i_fn, r_fn, p_fn, o_path, opt)
% Coregisters input file 'i_fn' to reference 'r_fn' using elastix parameters
% 'p_fn'. Saves the files to 'o_path' as 'o_fn'
%
% For details on motion correction, see
%
% Nilsson M, Szczepankiewicz F, van Westen D, Hansson O (2015)
% Extrapolation-Based References Improve Motion and Eddy-Current
% Correction of High B-Value DWI Data: Application in Parkinson's
% Disease Dementia. PLoS One 10(11):e0141825.
function mfs_fn = mdm_data2fit(fun_4d_data2fit, s, mfs_fn, opt)
% fun_4d_data2fit - model fit function
% s - input data structure
% mfs_fn - the model fit structure is written to mfs_fn
% opt - options structure
function dps = mdm_dps_load(dps_fn)
% load the derived parameter structure and check a few things
function dps_fn = mdm_dps_save(dps, s, dps_fn, opt)
% Saves derived parameter structure
function status = mdm_fit(varargin)
% initiate
function dps_fn = mdm_fit2param(fun_4d_fit2param, mfs_fn, dps_fn, opt)
% fun_4d_fit2param - parameter derivation function
% mfs_fn - model fit structure read from dps_fn
% dps_fn - derived parameter structure is written to dps_fn
% opt - options structure
function gdir_fn = mdm_fn_nii2gdir(nii_fn)
% convert a filename ending in .nii or .nii.gz to one ending in _gdir.txt
function g = mdm_gwf_read(gwf_fn)
% Read a gradient waveform
function [bt, q] = mdm_gwf_to_btensor(gwf, dt, g_max) (mdm_gwf_to_btensor.m)
% assume max(abs(gwf)) = 1
function o = mdm_iter_lund(input_path, f_handle, f_args, opt)
function s = mdm_mask(s, mask_fun, path, opt)
% Mask the data in s.nii_fn using mask_fun (e.g. mio_mask_simple)
%
% Save the mask as a nifti to s.mask_fn. This field is created if it does
% not exist already
% init
function M = mdm_mask_load(s, opt)
% Load the mask from s.mask_fn
%
% Optionally, mask using i_range, j_range and k_range in opt
function mdm_mat2nii(I, path, sdim)
% Makes nifti from I using path and header h
function s = mdm_mec_b0(s, p_fn, o_path, opt)
% Perform motion and eddy currect correction by registering to b0
%
% s - input structure OR a nii_fn with xps computed from it
% p_fn - parameter filename, to elastix registration scheme
% o_path - output path for the new files
% opt - options (optional)
%
% Output
% s - s.nii_fn will be updated to refer to the corrected volume
function s = mdm_mec_eb(s_target, s_source, p_fn, o_path, opt)
% Perform motion and eddy currect correction by registering to extrapolated
% references
%
% s_target - input structure for target data OR a nifti filename
% s_source - input structure for source data (should be low b-values)
% OR a nifti filename
% p_fn - parameter filename, to elastix registration scheme
% o_path - output path for the new files
% opt - options (optional)
%
% Output
% s - s.nii_fn will be updated to refer to the corrected volume
% init
function mfs = mdm_mfs_load(mfs_fn)
% Load model fit structure
function mfs_fn = mdm_mfs_save(mfs, s, mfs_fn, opt)
% Saves the model fit structure
function mdm_nii2pdf(nii_fn, pdf_fn, opt)
% Converts one or several nii files in nii_fn to pdf files
function datatype = mdm_nii_datatype(value, t)
% Returns the datatype of a nifti based on the header info in a
% matlab-friendly format
%
% set t = 1 to change float64 --> double and float32 --> single
function [nii_fn, tmp_path, tmp_fn] = mdm_nii_gunzip(nii_fn, h_only)
% Extract the file to the system tempdir (faster if using an external HD)
function h = mdm_nii_h_empty()
% create a header where all apppropriate fields for the nifti header exists
function out_nii_fn = mdm_nii_merge(nii_fn_cell, out_nii_fn, opt)
function o = mdm_nii_oricode(h)
% Find the principal orientation of the nifti header
function [I,h] = mdm_nii_read(nii_fn)
function [I,h] = mdm_nii_read_and_rescale(nii_fn)
function h = mdm_nii_read_header(nii_fn)
% Possibly unzip the file
function out_fn = mdm_nii_subsample(in_fn, ind, out_fn) (mdm_nii_subsample.m)
function s = mdm_nii_to_s(nii_fn)
% converts a nii_fn to an input structure with two fields
%
% s.nii_fn
% s.xps
%
% assumes the xps filename can be constructred from the nii_fn
function mdm_nii_write(I, nii_fn, h, is_colour)
function opt = mdm_opt(opt)
% Specifies default options
function [a_ind, c_list, id_ind] = mdm_pa_ind_from_xps(xps)
% Volumes with identical rotations is defined from xps.a_ind
%
% If this does not exist, we try to compute it, but beware of errors in
% this step
function [nii_fn, xps] = mdm_par2nii_and_xps(par_fn, nii_fn)
function fn = mdm_param2nii(dps_fn, o_path, fig_opt, opt)
function paths = mdm_paths(tmp, prefix, suffix)
% paths will get the required fields:
% paths.(mfs_fn/dps_fn/nii_path)
%
% if tmp is a string, these will be created assuming tmp is a path to the
% folder where the mfs, dps and nii:s are going to be stored
%
% if tmp is a struct, we verify that these fields are present
function xps = mdm_philips_fexi_prot2xps(prot_fn, xps)
% Parse a Philips protocol text file in prot_fn, and store it in an xps
%
% The protocol file must have been generated with the FEXI patch
%
% An xps must be provided and some fields must have been filled in
% (see this file for details)
% associate tm so different s_ind's
function s_protocol = mdm_philips_parse_txt_protocol(prot_fn)
function s = mdm_powder_average(s, o_path, opt)
% Average over rotations. Image volumes with identical rotations is defined
% from s.xps.a_ind
%
% To do: find a way of keeping track of number of averages per step
% Init
function s = mdm_s_from_lund_nii_plus_gdir(nii_fn)
function s_new = mdm_s_subsample(s, ind, path, opt)
function s = mdm_smooth(s, filter_sigma, o_path, opt)
function txt = mdm_txt_read(txt_fn, opt) (mdm_txt_read.m)
% This function returns each line in txt_fn as a cell, skipping empty
% lines and those starting with #
%
function txt = mdm_txt_write(txt, txt_fn, opt)
% This function write each line in txt, a cell, as separate lines in txt_fn
function res = mdm_unix2txt(in_fn, out_fn)
% Converts the unix text file 'in_fn' to the text file
% 'out_fn'. If 'out_fn' is omitted, the original file is
% overwritten. res = 1 in case of success.
% Open the files
function xps = mdm_xps_calc_btpars(xps)
function mdm_xps_check(xps)
% Check that necessary fields are present in the xps, and that all other
% fields are well formatted
function xps_fn = mdm_xps_fn_from_nii_fn(nii_fn)
function xps = mdm_xps_from_bt(bt)
function xps = mdm_xps_from_btens(btens_fn)
% assume b-tensors are stored in a file where each row corresponds to a
% b-tensor in voigt format, i.e.
%
% b_xx b_yy b_zz sqrt(2) * b_xy sqrt(2) * b_xz sqrt(2) * b_yz
%
% see tm_1x3_to_1x6.m and tm_1x6_to_3x3.m for conversion between tensors
% on 3x3 or 1x6 (voigt) format.
function xps = mdm_xps_from_bval_bvec(bval_fn, bvec_fn, b_delta) (mdm_xps_from_bval_bvec.m)
function xps = mdm_xps_from_gdir(gdir_fn, delimeter, b_delta)
% read a gradient textfile in the Lund format, which is defined as follows
% n_x, n_y, n_z, b
%
% units of b are in s/mm2, i.e., b = 1000 s/mm2 for a standard DTI
%
% the xps is always in SI units
%
% optional arguments
% delimeter, asumed to be ',' unless specified
% b_delta, assumed to be 1 unless specified (e.g. PGSE / LTE / SDE)
function xps = mdm_xps_load(xps_fn)
function xps = mdm_xps_merge(xps_cell, opt)
% xps_cell could be a cell of xps structres or filenames
function xps_pa = mdm_xps_pa(xps, opt) (mdm_xps_pa.m)
function mdm_xps_save(xps, xps_fn, opt) (mdm_xps_save.m)
function xps = mdm_xps_subsample(xps, ind)
% help the user
function fn = ut_mdm(c_ut)
% Run unit tests on the files in this package
Package: mio
Multidimensional Image Operations (MIO)
Functions for manipulating image volumes. These functions should take an
image volume (3D or 4D) as first input and yield another image volume as
the first output
Functions
function [I_res,tp,h_res,elastix_t] = mio_coreg(I_mov, I_ref, p, opt, h_mov, h_ref) (mio_coreg.m)
function mfs_fn = mio_fit_model(fun, s, o_fn, opt)
% Fits the model specified in 'fun' to the data referred to in 's',
% and saves the output in 'o_fn', which should have '.mat' as its extension
%
% Typically, 'fun' is a local function defined within model_4d_data2fit
function M = mio_mask_cv(I, opt)
% all images should have the same contrast
%
function M = mio_mask_expand(I, n, opt) (mio_mask_expand.m)
% expand the mask by n voxels
function M = mio_mask_fill(M,d)
% fill holes in a mask on a slice-by-slice basis
function M = mio_mask_keep_largest(M)
function M = mio_mask_mic(I, opt)
function M = mio_mask_pca(I, opt)
% pca-based filtering, inspired by an MRM paper (to do: find ref)
function M = mio_mask_simple(I, opt)
% creates a mask in a very simple way by looking a signal variation along
% the fourth dimension
% define the mask from the variation of the signal
function M = mio_mask_thresh(I, opt)
function I = mio_min_max_cut(I, min_max)
% cut values below min or above max, where [min max] = min_max
function opt = mio_opt(opt)
function A = mio_pa(I, xps, opt)
function I = mio_pad(I, pad_xyz)
% Adds zero around 'I' according to the three numbers in 'pad_xyz',
% which may be negative (this leads to trimming instead)
function I_ref = mio_ref_extrapolate(I, xps_source, xps_target, M, ind)
% Extrapolate references according to Nilsson et al, 2015, Plos One
% doi:10.1371/journal.pone.0141825
function [I_B,sc] = mio_rescale(I_B, I_A, M)
% rescale B to A
function I = mio_smooth_4d(I, filter_sigma, opt)
% Gaussian smoothing with a width controlled by 'filter_sigma'
% init
function I_out = mio_transform(I_in, t, h, opt)
% I: Image volume (x,y,z,c) to be transformed (looping over 4th dim)
% t: Elastix Transform structure
% h: Nifti header
%
% Note that the contents of h and t must be matched
function p = mio_volume_loop(fun, I, M, opt, S)
% Applies the function fun to all voxels where M is positive
%
% In order to execute the loop in parallel, initiate parpool before
% executing this function. If number of workers > 1, volume loop will be
% executed in parallel
%
% fun - accepts one or two input arguments (signal vector + optionally
% supplementary data)
%
% I - 4D matrix
%
% M - 3D mask
%
% opt - option structure
%
% S - 4D supplementary data
function fn = ut_mio(c_ut)
% Run unit tests on the files in this package
Package: mio/elastix
These functions allows for interaction with the ElastiX program for
image registration.
To get ElastiX running, follow these instructions
1) Download ElastiX from http://elastix.isi.uu.nl
2) Unzip and rename the downloaded folder as 'elastix'
MAC
3) Move the program folder into /usr/local/elastix, e.g. by
'sudo mv ~/Downloads/elastix /usr/local'
4) Typing 'ls /usr/local/elasmtix' should now show to folder, bin and lib
5) To run, you need to be able to start elastix from your terminal by
typing 'elastix'. Do this by editing ~/.bash_profile, and add the following
lines
export PATH=/usr/local/elastix/bin:/usr/local/elastix/parameters:$PATH
export DYLD_LIBRARY_PATH=/usr/local/elastix/lib:$DYLD_LIBRARY_PATH
The editing can be done from within Matlab by
typing 'edit ~./bash_profile'
WINDOWS
3) Move the program folder into Functions
function p = elastix_p_3dof(n_iter)
% Use an affine transform, but with all parameters except
% translation x, y and rotation around z set to zero
function p = elastix_p_6dof(n_iter)
% Use an affine transform, but with all parameters except
% translation x, y and rotation around z set to zero
function p = elastix_p_affine(n_iter)
function p = elastix_p_read(fn)
% Reads an elastix parameter structure from file
function fn = elastix_p_write(p, fn)
% Writes the elastix parameter structure to a file
function [res_fn, tp_fn] = elastix_run_elastix(i_fn, ref_fn, p_fn, o_path)
% Runs elastix. For help, see readme.txt in the elastix folder.
function res_fn = elastix_run_transformix(i_fn, t_fn, o_path)
% Run transformix. For help, see readme.txt in the elastix folder.
function [p,tmat,str] = elastix_tp_read(fn)
% Reads a transform parameter file
% Read the whole file
function tpm = elastix_tpm_read(tpm_fn)
% Read a set of transform parameters from the Affine DTI Transform
%
%
function elastix_tpm_write(tpm, fn)
% Writes the 'transform parameter matrix' to the file named 'fn'
Package: msf
Multidimensional support functions (msf). Functions that help with everyday
things...Functions
function gamma = msf_const_gamma(nucleus)
function msf_delete(fn)
% Deletes the file fn, or multiple files if fn is a cell array
function s = msf_ensure_field(s, f, v)
% Sets the field 'f' in the structure 's' to the value 'v', unless there is
% already a field 'f' in 's'.
function [path, name, ext] = msf_fileparts(fn)
% Just as the built-in 'fileparts', except that this run fileparts twice in
% order to count '.nii.gz' as one extension
function [t,ss] = msf_fit(f, s, l, u, n, opt)
% Perform lsqcurvefit with function 'f' on data 's', repeating it 'n' times
%
% f - function handle
% s - signal
% l - lower bound
% u - upper bound
% n - number of repetitions
% opt - options for lsqcurvefit
function msf_fprintf(opt, str, varargin)
% prints to terminal if opt.verbose = 1
function msf_imagesc(I,d,k,c)
% Displays a 2D slice through a 3D or 4D image volume I
%
% d - dimension
% k - slice
% c - volume
function val = msf_isfield(s, field_name)
function msf_log(str, opt)
% Displays 'str' if opt.verbose is true
function msf_mkdir(folder_path)
% Recursively makes folder_path
function a = msf_nanmean(a,d)
function s = msf_rmfield(s, f)
% Removes field 'f' from 's' is possible. If 'f' is a cell array, all
% fields in 'f' are removed
function sz = msf_size(m, d)
% a size function that makes sure to return data with the correct
% dimensionality
%
% m - matrix
% d - dimensions, i.e. length of sz
function o = msf_str2double(str, delim)
% Splits 'str' according to 'delim' and converts to double
function tmp_path = msf_tmp_path(do_mkdir)
% Obtains a temporary path with a random name
Package: methods
Collection of models for multidimensional diffusion MRI Package: methods/dtd_pa
Diffusion tensor distribution from powder averaged data by Topgaard.Functions
function m = dtd_pa_1d_data2fit(signal, xps, opt, ind)
% Size-shape diffusion tensor distribution
% Assumes powder averaging
% de Almeida Martins and Topgaard, Phys. Rev. Lett. 116, 087601 (2016).
% http://dx.doi.org/10.1103/PhysRevLett.116.087601
% Modified for general b-tensor shapes
function s = dtd_pa_1d_fit2data(m, xps) (dtd_pa_1d_fit2data.m)
function res = dtd_4d_data2fit(s, mfs_fn, opt) (dtd_pa_4d_data2fit.m)
function dps = dtd_pa_4d_fit2param(mfs_fn, dps_fn, opt)
function dtd = dtd_pa_data2dtd(stemp,b,b_delta,dtd_nodes) (dtd_pa_data2dtd.m)
function dtd_nodes = dtd_pa_dist2nodes(dtd) (dtd_pa_dist2nodes.m)
function [n,par,perp,w] = dtd_pa_dist2par(dtd_pa) (dtd_pa_dist2par.m)
function m = dtd_pa_dtd2m(dtd,opt) (dtd_pa_dtd2m.m)
function dtd = dtd_pa_extinction(stemp, b, b_delta, dtd, opt)
function dtd_pa = dtd_pa_m2dtd(m) (dtd_pa_m2dtd.m)
function dtd_pa_mkpdf(dps_fn, pdf_path, opt)
function [n,par,perp] = dtd_pa_nodes2par(dtd_nodes) (dtd_pa_nodes2par.m)
function dtd_nodes = dtd_pa_nodes_merge(dtd_nodes1,dtd_nodes2) (dtd_pa_nodes_merge.m)
function dtd_nodes_out = dtd_pa_mutate(dtd_nodes,opt) (dtd_pa_nodes_mutate.m)
function dtd_nodes_out = dtd_pa_nodes_select(dtd_nodes,ind) (dtd_pa_nodes_select.m)
function dtd = dtd_pa_nodesw2dist(dtd_nodes,w) (dtd_pa_nodesw2dist.m)
function opt = dtd_pa_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dtd_pa_pipe(s, paths, opt)
% s - input structure
% paths - either a pathname or a path structure (see mdm_paths)
% opt - (optional) options that drive the pipeline
%
% fn - a cell arary with filenames to generated nii files
function dtd_pa_plot(m, opt) (dtd_pa_plot.m)
function dtd = dtd_pa_proliferation(stemp, b, b_delta, opt) (dtd_pa_proliferation.m)
function dtd_nodes = dtd_pa_rand(n,dmin,dmax) (dtd_pa_rand.m)
function dtd_out = dtd_pa_sort(dtd_in) (dtd_pa_sort.m)
Package: methods/dti_euler
Diffusion tensor imaging with the diffusion tensor represented by
euler angles. Fitting is done in a weighted mannor to reduce influence of
high b-values.
By Topgaard.Functions
function m = dti_euler_1d_data2fit(signal, xps, opt, ind)
function s = dti_euler_1d_fit2data(m, xps)
% convert to readable parameters
function mfs_fn = dti_euler_4d_data2fit(s, mfs_fn, opt) (dti_euler_4d_data2fit.m)
function dps = dti_euler_4d_fit2param(mfs_fn, dps_fn, opt)
function dti_euler_mic_check_xps(xps) (dti_euler_check_xps.m)
% checks that all required fields are found in the xps
function opt = dti_euler_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dti_euler_pipe(s, paths, opt)
% s - input structure
% paths - either a pathname or a path structure (see mdm_paths)
% opt - (optional) options that drive the pipeline
%
% fn - a cell arary with filenames to generated nii files
Package: methods/dti_nls
A simple DTI model with NLS fitting. Utilizes Cholesky decomposition to
ensure positive eigenvalues of the diffusion tensorFunctions
function m = dti_nls_1d_data2fit(signal, xps, opt)
% Yields a 1x7 vector 'm' with the fit parameters
% m(1) - s0
% m(2:7) - diffusion tensor
%
% This file features a number of functions important for this framework.
% As input to the lsqcurvefit, we use a local vector t. Locally and during
% the fit, we represent the diffusion tensor by its cholesky factorization
% to ensure that all eigenvalues are positive. We also make sure the units
% of 't' are in the same range as 's0' in order to improve the fitting's
% stop constrains. Conversion between local 't' and true model variables
% 'm' are done by the function 't2m'.
function s = dti_nls_1d_fit2data(m, xps)
% Predict the signal 's' as a 1 x xps.n vector from the model parameters
% in 'm'.
function mfs_fn = dti_nls_4d_data2fit(s, o, opt)
% Loops over a 4D volume to produce fit parameters with the DTI nls model
%
% Input:
%
% s - data structure
% s.nii_fn - full path to nifti filename with data
% s.xps - experimental parameter structure
%
% o - output folder
%
% opt - options structure, optional
%
% Output:
%
% mfs_fn - path to .mat file with the ModelFitStructure (mfs)
function fn = dti_nls_4d_fit2param(mfs_fn, o_path, opt)
% Creates meaningful parameters from the model fit structure
%
% Input:
%
% mfs_fn - Path to .mat file with model fit structure
%
% o_path - Output path, where parameters maps are stored
%
% opt - Options, optional argument
%
%
% Output:
%
% fn - A cell array with paths to parameter maps that were written
function dti_nls_check_xps(xps)
% checks that all required fields are found in the xps
function opt = dti_nls_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dti_nls_pipe_example(s, o_path, opt)
% s - input structure
% o_path - output path
function fn = dti_nls_pipe_mic(s, o_path, opt)
% s - input structure
% o_path - output path
%
% Pipeline for running dti_nls in a microimaging setting
function fn = ut_dti_nls(c_ut)
% if c_ut is not supplied, the function returns the number of unit tests
% n_ut = number of unit tests
Package: methods/dtd_pake
Estimate microscopic anisotropy by utilizing the dtd_pake model.
Assumes powder averaging and axisymmetric b-tensors
Eriksson et al., J. Chem. Phys. 142, 104201 (2015).
http://dx.doi.org/10.1063/1.4913502
Functions
function m = dtd_pake_1d_data2fit(signal, xps, opt, ind)
function s = dtd_pake_1d_fit2data(m, xps)
% m = [s0 d_iso d_delta]
% convert to readable parameters
function msf_fn = dtd_pake_4d_data2fit(s, mfs_fn, opt) (dtd_pake_4d_data2fit.m)
% Error function fit
% Assumes powder averaging and axisymmetric b-tensors
% Eriksson et al., J. Chem. Phys. 142, 104201 (2015).
% http://dx.doi.org/10.1063/1.4913502
function dps = dtd_pake_4d_fit2param(mfs_fn, dps_fn, opt) (dtd_pake_4d_fit2param.m)
function dtd_pake_check_xps( xps )
function opt = dtd_pake_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dtd_pake_pipe(s, paths, opt)
% s - input structure
% paths - either a pathname or a path structure (see mdm_paths)
% opt - (optional) options that drive the pipeline
%
% fn - a cell arary with filenames to generated nii files
function dtd_pake_plot(S, xps, h, h2)
Package: methods/fexi11
Filter exchange imaging model. Requires DDE with variable mixing times.
Works best with powder averaged data.
Lasi? S, Nilsson M, L?tt J, St?hlberg F, Topgaard D (2011) Apparent
exchange rate mapping with diffusion MRI. Magn Reson Med 66(2):356?365.
Nilsson M, et al. (2013) Noninvasive mapping of water diffusional
exchange in the human brain using filter-exchange imaging.
Magn Reson Med 69(6):1573?1581.
Lampinen B, et al. (2016) Optimal experimental design for filter exchange
imaging: Apparent exchange rate measurements in the healthy brain and in
intracranial tumors. Magn Reson Med. doi:10.1002/mrm.26195.Functions
function m = fexi11_1d_data2fit(signal, xps, opt, ind)
% Yields a 1xN vector 'm' with the fit parameters
% m(1) - ADC0
% m(2) - sigma
% m(3) - AXR
% m(4:end) - vector of S0, length(m) = 3+max(s.xps.s_ind)
function s = fexi11_1d_fit2data(m, xps)
function mfs_fn = fexi11_4d_data2fit(s, o, opt)
% Loops over a 4D volume to produce fit parameters with the FEXI nls model
%
% Input:
%
% s - data structure
% s.nii_fn - full path to nifti filename with data
% s.xps - experimental parameter structure
%
% o - output folder
%
% opt - options structure, optional
%
% Output:
%
% mfs_fn - path to .mat file with the ModelFitStructure (mfs)
function fn = fexi11_nls_4d_fit2param(mfs_fn, o_path, opt) (fexi11_4d_fit2param.m)
% Creates meaningful parameters from the model fit structure
%
% Input:
%
% mfs_fn - Path to .mat file with model fit structure
%
% o_path - Output path, where parameters maps are stored
%
% opt - Options, optional argument
%
%
% Output:
%
% fn - A cell array with paths to parameter maps that were written
%
% 2do: Update this function to agree with structure of dti_euler
function fexi11_check_xps(xps)
% checks that all required fields are found in the xps
function fn = fexi11_invivo_pipeline(s, o_path, opt) (fexi11_invivo_pipeline.m)
% s - input structure
% o_path - output path
% init stuff
function opt = fexi11_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = ut_fexi11(c_ut)
% if c_ut is not supplied, the function returns the number of unit tests
% n_ut = number of unit tests
Package: methods/dtd_gamma
Gamma model fit.
Does powder averaging
Lasic et al, Front. Phys. 2, 11 (2014).
http://dx.doi.org/10.3389/fphy.2014.00011
Modified for general b-tensor shapesFunctions
function m = dtd_gamma_1d_data2fit(signal, xps, opt, ind)
% m = [s0 d_iso mu2_iso mu2_aniso]
function s = dtd_gamma_1d_fit2data(m, xps)
% s0 = m(1);
% d_iso = m(2);
% mu2_iso = m(3);
% mu2_aniso = m(4);
%
% convert to readable parameters
function mfs_fn = dtd_gamma_4d_data2fit(s, mfs_fn, opt)
function dps = dtd_gamma_4d_fit2param(mfs_fn, dps_fn, opt)
function dtd_gamma_check_xps( xps )
function opt = dtd_gamma_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dtd_gamma_pipe(s, paths, opt)
% s - input structure
% paths - either a pathname or a path structure (see mdm_paths)
% opt - (optional) options that drive the pipeline
%
% fn - a cell arary with filenames to generated nii files
%
% Gamma fit
% Does powder averaging
% Lasic et al, Front. Phys. 2, 11 (2014).
% http://dx.doi.org/10.3389/fphy.2014.00011
% Modified for general b-tensor shapes
function dtd_gamma_plot(S, xps, h, h2)
Package: methods/quick_dti
Quick DTI fitting using least squares fitting without adjusting for
heteroscedasticity.Functions
function m = dti_1d_data2fit(signal, xp, opt) (qdti_1d_data2fit.m)
function p = dti_1d_fit2param(m, xp, opt, signal) (qdti_1d_fit2param.m)
function mfs_fn = qdti_4d_data2fit(s, o, opt) (qdti_4d_data2fit.m)
function fn = qdti_4d_fit2param(mfs_fn, o_path, opt)
% 2do: Update this function to agree with structure of dti_euler
function opt = qdti_opt(opt)
function fn = qdti_pipe_example(s, o_path, opt) (qdti_pipe_example.m)
% s - input structure
% o_path - output path
function status = ut_qdti(c_ut) (ut_qdti.m)
% n_ut = number of unit tests
Package: methods/dtd_saupe
Functions
function dps_fn = dtd_saupe_4d_fit2param(mfs_fns, dps_fn, opt) (dtd_saupe_4d_fit2param.m)
% Saupe order tensor imaging
% Combination of DTI, gamma, and erf
% Topgaard, Phys. Chem. Chem. Phys. (2016).
% http://dx.doi.org/10.1039/c5cp07251d
function opt = dtd_saupe_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = dtd_saupe_pipe(s, paths, opt)
Package: methods/vasco16
VascoPule model. Aim is to extract the blood fraction by the use of two
datasets obtained with and without flow compensation (stored in xps.alpha2)
This model is described in Ahlgren et al. (2016) NMR Biomed.
WARNING: This code is note yet verified. It is not the same code as
was used in the Ahlgren et al paper. Functions
function m = vasco16_1d_data2fit(signal, xps, opt, ind)
function s = vasco16_1d_fit2data(m, xps)
% convert to readable parameters
function mfs_fn = vasco16_4d_data2fit(s, o, opt)
function fn = vasco16_4d_fit2param(mfs_fn, o_path, opt)
% 2do: Update this function to agree with structure of dti_euler
function vasco16_check_xps(xps)
function opt = vasco16_opt(opt)
% Makes sure that all needed fields in the options structure are present
function fn = vasco16_pipeline(s, o_path, opt)
% unit of v2 seems wrong
% calculation of alpha must be checked
function vasco16_plot(S, xps, h, h2)
Package: tools
This folder contains various tools that are useful in the context of
designing or analysing multidimensional dMRI data, e.g., math-related
functions.
Package: tools/tensor_maths
Functions for transforming first, second, and fourth-order tensors.
The functions were used extensively in the following papers:
Westin CF, et al (2016). "Q-space trajectory imaging for multidimensional
diffusion MRI of the human brain". NeuroImage.
Westin CF, et al. (2014) "Measurement Tensors in Diffusion MRI:
Generalizing the Concept of Diffusion Encoding."
Med Image Comput Comput Assist Interv 8675:209?216.Functions
function t = tm_1x15_to_6x6(N)
% Converts a fourth-order tensor on the 1x15 format to the 6x6 format
function [E_bulk, E_shear, E2_iso] = tm_1x21_iso()
% Obtain the isotropic fourth-order tensors in the 1x21 Voigt-like format
function t = tm_1x21_to_6x6(t)
% Convert fourth-order tensor in 1x21 format to 6x6 format
function t = tm_1x3_to_1x15(N)
% Convert a vector in 1x3 format to a fourth-order tensor in 1x15 format
function t = tm_1x3_to_1x6(ad, rd, n)
% Convert a vector (1x3) to a tensor (1x6), but enable some shaping of the
% tensor, so that its longest eigenvalue is 'ad' and its shortest to
% eigenvalues are 'rd'
function [L, v_1] = tm_1x6_eigvals(v_1x6, c_method)
%
% c_method
% 1 - Cardano's fast method
% 2 - eigs in Matlab
function t = tm_1x6_to_1x21(n2)
% Convert a second-order tensor on Voigt notation (1x6) to a
% fourth-order tensor in Voight-like notation (1x21). This correspons
% to taking the outer product of 'n2'
function t = tm_1x6_to_1x56(B2)
% Take the third outer product of B2, making the outcome a sixth order
% tensor
function t = tm_1x6_to_3x3(t)
% convert a second order tensor from Voigt format to 3x3 matrix format
function t = tm_1x6_to_tpars(t_1x6)
function t = tm_2d_1x2_to_1x3(ad, rd, n)
% Two-dimensional version of tm_1x3_to_1x6
function t = tm_2d_1x3_to_1x6(t)
% Two-dimensional version of tm_1x6_to_1x21
function t = tm_2d_1x3_to_2x2(t)
% Two-dimensional version of tm_1x6_to_3x3
function t = tm_2d_1x6_to_3x3(t)
% Two-dimensional version of tm_1x21_to_6x6
function E_iso = tm_3x3_iso()
% Get the isotropic 3x3 tensor
function t = tm_3x3_to_1x6(t)
% Convert a second-order 3x3 tensor to Voigtm-format 1x6
function t = tm_3x3_to_tpars(t3x3)
% calculate derived tensor parameters
function [E_bulk, E_shear, E2_iso] = tm_6x6_iso()
% Get the three different isotropic fourth-order tensors
function t = tm_6x6_to_1x15(N)
% Convert a 6x6 tensor to a Voigt-like 1x15 representation
% Warning: degrees of freedom are lost in this conversion.
function t = tm_6x6_to_1x21(t)
% Convert a fourth-order tensor on the 6x6 format to a Voigt-like format
function t_nx6 = tm_ax2nx6(par,perp,theta,phi)
% Produce a set of 1x6 tensors from parallel and perpendicular
% eigenvalues, plus theta and phi angles of the parallel part
function [rotmat,rotmatinv] = tm_euler_angles2rotmat(alpha,beta,gamma)
% Active right-handed rotation
% First gamma around z, then beta around y, and finally alpha around z
% Input must be scalars.
function fa = tm_fa(a,b)
% With one input argumnet: calculate the FA of the second order tensor 'a'
% With two arguments: assume 'a' to be MD and 'b' to be the variance in
% eigenvectors, and calculate FA from that
function x = tm_inner(x,y)
% Take the inner product of 'x' and 'y'.
% However, rewrite required
function md = tm_md(d2)
% Calculate the mean diffusivity of the second-order tensor 'd2'
function x = tm_outer(x,y)
% Take the outer product of 'x' and 'y'
function t = tm_t2tpars(t3x3)
% calculate derived tensor parameters
function v_lambda = tm_v_lambda(d2)
% Calculate the eigenvalue variance of the second-order tensor in 'd2'
function fn = ut_tensor_maths(c_ut)
% if c_ut is not supplied, the function returns the number of unit tests
%
% n_ut = number of unit tests
Package: tools/uvec
Functions to generate/retrieve set of vectorsFunctions
function fn = ut_uvec(c_ut)
% if c_ut is not supplied, the function returns the number of unit tests
% n_ut = number of unit tests
function u = uvec_dodeca()
% Optimal dodeca direction set
function [u, phi,theta] = uvec_elstat(n)
% Loads electrostatically optimized angles for n = 10-1000
function u = uvec_icosa()
% Returns an optimal icosahedral direction set
function u = uvec_tricosa()
% Truncated icosa direction set