In [1]:
import numpy as np
import nibabel as nib
import dsiadapt as dsi
from dipy.viz import fvtk
from dipy.data import get_sphere
import dsitool as dst

%pylab inline
np.set_printoptions(threshold=numpy.nan)


Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.kernel.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.

In [2]:
# Load ODF sphere directions
sphere = get_sphere('symmetric724')
sphere = sphere.subdivide(3)

In [3]:
files = np.array(['DSI11_invivo_xfib', 'DSI11_exvivo_xfib'])

In [4]:
# Load fov and mdd needed for ODF reconstruction
# fov and mdd data are generated by script s-figure2-pdf
# The four elements are in order to in vivo DSI11, ex vivo DSI11, ex vivo DSI15, ex vivo DSI17
fovarr = np.load('results/fov.npy');
mddarr = np.load('results/mdd.npy');

In [ ]:
qgridsz = 201;
for ii in np.arange(files.shape[0]):
    filename = files[ii]
    data = dst.loaddata(filename);
    gtab = dst.loadgtab(filename);
    fov = fovarr[ii];
    mdd = mddarr[ii];
    
    hannwidth = np.inf;
    integdistance = mdd;
    
    wgtarr = np.array([0, 1, 2, 4, 8])
    for pdfwgtorder in wgtarr: 
        rend = np.round(integdistance / fov * qgridsz);
        dsimodel = dsi.DiffusionSpectrumModel(gtab, qgrid_size=qgridsz, filter_width=hannwidth, r_start=0, r_end=rend, r_step=0.2);
        dsifit = dsimodel.fit(data);
        odf = dsifit.odf(sphere, pdfwgtorder);
        odfvox = odf[2, 0, 0]
        
        # save odf
        r = fvtk.ren();
        fvtk.add(r, fvtk.sphere_funcs(odfvox, sphere));
        odfname = 'figure/figure5/' + filename + '_odf_order_' + str(pdfwgtorder) + '.png';
        fvtk.record(r, n_frames=1, out_path=odfname, size=(600,600))
        #fvtk.show(r)