{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Reconstruct with Constant Solid Angle (Q-Ball)\n\nWe show how to apply a Constant Solid Angle ODF (Q-Ball) model from Aganj et\nal. [Aganj2010]_ to your datasets.\n\nFirst import the necessary modules:\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\nfrom dipy.core.gradients import gradient_table\nfrom dipy.data import get_fnames, default_sphere\nfrom dipy.io.gradients import read_bvals_bvecs\nfrom dipy.io.image import load_nifti\nfrom dipy.reconst.shm import CsaOdfModel\nfrom dipy.direction import peaks_from_model\nfrom dipy.segment.mask import median_otsu" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Download and read the data for this tutorial and load the raw diffusion data\nand the affine.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "hardi_fname, hardi_bval_fname, hardi_bvec_fname = get_fnames('stanford_hardi')\n\ndata, affine = load_nifti(hardi_fname)\n\nbvals, bvecs = read_bvals_bvecs(hardi_bval_fname, hardi_bvec_fname)\ngtab = gradient_table(bvals, bvecs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "img contains a nibabel Nifti1Image object (data) and gtab contains a\nGradientTable object (gradient information e.g. b-values). For example to\nread the b-values it is possible to write print(gtab.bvals).\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print('data.shape (%d, %d, %d, %d)' % data.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "data.shape ``(81, 106, 76, 160)``\n\nRemove most of the background using DIPY's mask module.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "maskdata, mask = median_otsu(data, vol_idx=range(10, 50), median_radius=3,\n numpass=1, autocrop=True, dilate=2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We instantiate our CSA model with spherical harmonic order of 4\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "csamodel = CsaOdfModel(gtab, 4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`Peaks_from_model` is used to calculate properties of the ODFs (Orientation\nDistribution Function) and return for\nexample the peaks and their indices, or GFA which is similar to FA but for ODF\nbased models. This function mainly needs a reconstruction model, the data and a\nsphere as input. The sphere is an object that represents the spherical discrete\ngrid where the ODF values will be evaluated.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "csapeaks = peaks_from_model(model=csamodel,\n data=maskdata,\n sphere=default_sphere,\n relative_peak_threshold=.5,\n min_separation_angle=25,\n mask=mask,\n return_odf=False,\n normalize_peaks=True)\n\nGFA = csapeaks.gfa\n\nprint('GFA.shape (%d, %d, %d)' % GFA.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "GFA.shape ``(81, 106, 76)``\n\nApart from GFA, csapeaks also has the attributes peak_values, peak_indices and\nODF. peak_values shows the maxima values of the ODF and peak_indices gives us\ntheir position on the discrete sphere that was used to do the reconstruction of\nthe ODF. In order to obtain the full ODF, return_odf should be True. Before\nenabling this option, make sure that you have enough memory.\n\nLet's visualize the ODFs of a small rectangular area in an axial slice of the\nsplenium of the corpus callosum (CC).\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "data_small = maskdata[13:43, 44:74, 28:29]\n\nfrom dipy.viz import window, actor\n\n# Enables/disables interactive visualization\ninteractive = False\n\nscene = window.Scene()\n\ncsaodfs = csamodel.fit(data_small).odf(default_sphere)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is common with CSA ODFs to produce negative values, we can remove those\nusing ``np.clip``\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "csaodfs = np.clip(csaodfs, 0, np.max(csaodfs, -1)[..., None])\ncsa_odfs_actor = actor.odf_slicer(csaodfs, sphere=default_sphere,\n colormap='plasma', scale=0.4)\ncsa_odfs_actor.display(z=0)\n\nscene.add(csa_odfs_actor)\nprint('Saving illustration as csa_odfs.png')\nwindow.record(scene, n_frames=1, out_path='csa_odfs.png', size=(600, 600))\nif interactive:\n window.show(scene)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: csa_odfs.png\n :align: center\n\n Constant Solid Angle ODFs.\n\n.. include:: ../links_names.inc\n\n## References\n\n.. [Aganj2010] Aganj I, Lenglet C, Sapiro G, Yacoub E, Ugurbil K, Harel N.\n \"Reconstruction of the orientation distribution function in single- and\n multiple-shell q-ball imaging within constant solid angle\", Magnetic\n Resonance in Medicine. 2010 Aug;64(2):554-66. doi: 10.1002/mrm.22365\n\n\n\n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.13" } }, "nbformat": 4, "nbformat_minor": 0 }