{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Continuous and analytical diffusion signal modelling with 3D-SHORE\n\nWe show how to model the diffusion signal as a linear combination\nof continuous functions from the SHORE basis [Merlet2013]_, [\u00d6zarslan2008]_,\n[\u00d6zarslan2009]_. We also compute the analytical Orientation Distribution\nFunction (ODF).\n\nFirst import the necessary modules:\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "from dipy.reconst.shore import ShoreModel\nfrom dipy.viz import window, actor\nfrom dipy.core.gradients import gradient_table\nfrom dipy.data import get_fnames, get_sphere\nfrom dipy.io.gradients import read_bvals_bvecs\nfrom dipy.io.image import load_nifti"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Download and read the data for this tutorial.\n\n``fetch_isbi2013_2shell()`` provides data from the [ISBI HARDI contest 2013](http://hardi.epfl.ch/static/events/2013_ISBI/) acquired for two shells at\nb-values 1500 $s/mm^2$ and 2500 $s/mm^2$.\n\nThe six parameters of these two functions define the ROI where to reconstruct\nthe data. They respectively correspond to ``(xmin,xmax,ymin,ymax,zmin,zmax)``\nwith x, y, z and the three axis defining the spatial positions of the voxels.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fraw, fbval, fbvec = get_fnames('isbi2013_2shell')\n\ndata, affine = load_nifti(fraw)\nbvals, bvecs = read_bvals_bvecs(fbval, fbvec)\ngtab = gradient_table(bvals, bvecs)\ndata_small = data[10:40, 22, 10:40]\n\nprint('data.shape (%d, %d, %d, %d)' % data.shape)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "``data`` contains the voxel data and ``gtab`` contains a ``GradientTable``\nobject (gradient information e.g. b-values). For example, to show the b-values\nit is possible to write::\n\n    print(gtab.bvals)\n\nInstantiate the SHORE Model.\n\n``radial_order`` is the radial order of the SHORE basis.\n\n``zeta`` is the scale factor of the SHORE basis.\n\n``lambdaN`` and ``lambdaL`` are the radial and angular regularization\nconstants, respectively.\n\nFor details regarding these four parameters see [Cheng2011]_ and [Merlet2013]_.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "radial_order = 6\nzeta = 700\nlambdaN = 1e-8\nlambdaL = 1e-8\nasm = ShoreModel(gtab, radial_order=radial_order,\n                 zeta=zeta, lambdaN=lambdaN, lambdaL=lambdaL)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Fit the SHORE model to the data\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "asmfit = asm.fit(data_small)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Load an odf reconstruction sphere\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "sphere = get_sphere('repulsion724')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Compute the ODFs\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "odf = asmfit.odf(sphere)\nprint('odf.shape (%d, %d, %d)' % odf.shape)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Display the ODFs\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "# Enables/disables interactive visualization\ninteractive = False\n\nscene = window.Scene()\nsfu = actor.odf_slicer(odf[:, None, :], sphere=sphere, colormap='plasma',\n                       scale=0.5)\nsfu.RotateX(-90)\nsfu.display(y=0)\nscene.add(sfu)\nwindow.record(scene, out_path='odfs.png', size=(600, 600))\nif interactive:\n    window.show(scene)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        ".. figure:: odfs.png\n   :align: center\n\n   Orientation distribution functions.\n\n## References\n.. [\u00d6zarslan2008] \u00d6zarslan E. et al., \"Simple harmonic oscillator based\n   estimation and reconstruction for one-dimensional q-space MR,\u201d in Proc Intl\n   Soc Mag Reson Med, 16, Toronto, Canada, 2008\n\n.. [\u00d6zarslan2009] \u00d6zarslan E. et al., \"Simple harmonic oscillator based\n   reconstruction and estimation for three-dimensional q-space MRI,\" in Proc\n   Intl Soc Mag Reson Med, 17, Honolulu, HI, 2009\n\n.. [Merlet2013] Merlet S. et al., \"Continuous diffusion signal, EAP and ODF\n   estimation via Compressive Sensing in diffusion MRI\", Medical Image\n   Analysis, 2013.\n\n.. [Cheng2011] Cheng J. et al., \"Theoretical Analysis and Practical Insights on\n   EAP Estimation via Unified HARDI Framework\", MICCAI workshop on\n   Computational Diffusion MRI, 2011.\n\n.. include:: ../links_names.inc\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
}