{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Parallel reconstruction using Q-Ball\n\nWe show an example of parallel reconstruction using a Q-Ball Constant Solid\nAngle model (see Aganj et al. (MRM 2010)) and `peaks_from_model`.\n\nImport modules, fetch and read data, and compute the mask.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import time\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\nfrom dipy.reconst.shm import CsaOdfModel\nfrom dipy.direction import peaks_from_model\nfrom dipy.segment.mask import median_otsu\n\n\nhardi_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)\n\nmaskdata, 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": [
        "sphere = get_sphere('repulsion724')\n\nstart_time = time.time()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We will first run `peaks_from_model` using parallelism with 2 processes. If\n`num_processes` is None (default option) then this function will find the total\nnumber of processors from the operating system and use this number as\n`num_processes`. Sometimes it makes sense to use only a few of the processes in\norder to allow resources for other applications. However, most of the times\nusing the default option will be sufficient.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "csapeaks_parallel = peaks_from_model(model=csamodel,\n                                     data=maskdata,\n                                     sphere=sphere,\n                                     relative_peak_threshold=.5,\n                                     min_separation_angle=25,\n                                     mask=mask,\n                                     return_odf=False,\n                                     normalize_peaks=True,\n                                     npeaks=5,\n                                     parallel=True,\n                                     num_processes=2)\n\ntime_parallel = time.time() - start_time\nprint(\"peaks_from_model using 2 processes ran in : \" +\n      str(time_parallel) + \" seconds\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "peaks_from_model using 2 process ran in  : 114.333221912 seconds, using 2\nprocess\n\nIf we don't use parallelism then we need to set `parallel=False`:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "start_time = time.time()\ncsapeaks = peaks_from_model(model=csamodel,\n                            data=maskdata,\n                            sphere=sphere,\n                            relative_peak_threshold=.5,\n                            min_separation_angle=25,\n                            mask=mask,\n                            return_odf=False,\n                            normalize_peaks=True,\n                            npeaks=5,\n                            parallel=False,\n                            num_processes=None)\n\ntime_single = time.time() - start_time\nprint(\"peaks_from_model ran in : \" + str(time_single) + \" seconds\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "peaks_from_model ran in : 196.872478008 seconds\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(\"Speedup factor : \" + str(time_single / time_parallel))"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Speedup factor : 1.72191839533\n\nIn Windows if you get a runtime error about frozen executable please start\nyour script by adding your code above in a ``main`` function and use::\n\n    if __name__ == '__main__':\n        import multiprocessing\n        multiprocessing.freeze_support()\n        main()\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
}