{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulating and fitting data using a simple Stick model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As in the previous tutorial we will use the wu-minn acquisition scheme to do our experiments. Instead of loading it from scratch, we load it from dmipy.data.saved_acquisition_schemes, which contains some saved dmipy acquisition schemes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# load the necessary modules\n",
    "from dmipy.signal_models import cylinder_models\n",
    "from dmipy.core import modeling_framework\n",
    "from dmipy.data import saved_acquisition_schemes\n",
    "import numpy as np\n",
    "\n",
    "acq_scheme = saved_acquisition_schemes.wu_minn_hcp_acquisition_scheme()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using this acquisition scheme, we will simulate data using a simple Stick model, and then use the same Stick model to fit the signal again. First, we instantiate the model(s) we need.\n",
    "\n",
    "NOTE: this example the same for any other dmipy model. One only needs to change the model and appropriate input parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "stick = cylinder_models.C1Stick()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In dmipy, all the simulation and fitting functionality is contained in the MultiCompartmentModel module. To simulate some data, we therefore make a MultiCompartmentModel that just contains the Stick model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from dmipy.core.modeling_framework import MultiCompartmentModel\n",
    "stick_model = MultiCompartmentModel(models=[stick])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can visualize the flow diagram of the model using the `visualize_model_setup` command:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAR4AAAD7CAYAAAC8N/tTAAAABmJLR0QA/wD/AP+gvaeTAAAgAElE\nQVR4nOzdd1RU1/o38O8whV6lizQL3QKCiigW8NoTU+2JiZqrN2o09phfmrHEGEvUXHOvLmONMRoT\nNdFYCVLUICpNRBGkSe8wMOV5/+DlXEdAUZk5M7A/a82COW0/Z86ZZ07ZZ28BEREYhmE0SI/vABiG\n6XhY4mEYRuNY4mEYRuNEfAeg7SoqKpCTk4OysjJUV1dDJpOhqqoKenp6MDc3h0AggIWFBezs7GBn\nZwexWMx3yMwLkslkqKysRFlZGSorKyGXy1FeXg6lUslN0zi8kUQigbGxMfdeJBLB1NQUEokEJiYm\nsLCwgKmpKUQi9pUDWOLh3L17F3FxcUhISEBSUhJSUlKQlZWFmpqaVi9DIBDAzs4Ozs7O8PPzg4+P\nD/z8/NC/f3+YmJioMXrmSerr65GZmYmsrCzk5eWhsLAQBQUFePjwIQoLC7lXeXk5qqqqUFdXp7ZY\nDA0NYWJiAnNzc9jY2MDGxga2trawt7fn3js6OsLZ2RlOTk7t9odM0FHvahUUFOC3337DxYsXcenS\nJeTm5kIsFqN79+7w9fWFl5cXXF1dYWdnBycnJ1hZWcHIyAhisRgmJiZQKpUoLy8HAJSUlCA/Px8P\nHz5ETk4O7t69i6SkJCQmJiI/Px8ikQj+/v4YPHgwRo0ahdDQUAiFQp4/gfalqqoKKSkpSE5Oxt27\nd5GRkYGMjAzcv38feXl53NGKWCzmvuAODg7c/zY2NjA3N4eJiQlMTU1hYmICS0tLmJiYQCwWw9jY\nGBKJhCvPyMgI+vr63Pva2lpIpVLufV1dHWpqalBXV4eqqiqUl5ejoqICVVVV3PuCgoJmk6BCoQAA\nCIVCODk5wdXVlXv16NEDXl5e8PT0hKGhoYY+3bbXoRJPcXEx9u/fj6NHjyI6Ohr6+voYNGgQBg8e\njNDQUAQGBqrsXG0hLy8PkZGR+Ouvv3Dp0iUkJSWhU6dOGD9+PCZNmoSwsDAIBII2LbM9UygUSElJ\nwbVr15CYmIikpCTcvn0bmZmZAAB9fX1069YNbm5ucHV1hYuLC/eldXZ2hq2tLc9r8GREhPz8fGRm\nZnLJ89HXvXv3IJPJoKenB1dXV3h7e8Pb2xs9e/ZEYGAgunfvrhP7U4dIPLGxsdixYweOHDkCiUSC\n8ePHY8KECRg5ciSMjIw0Gsu9e/dw7Ngx/PLLL4iJiUG3bt0we/ZsvPPOO+jUqZNGY9EF2dnZiIqK\nwrVr13Dt2jVcv34dVVVVMDQ0hI+PD7y9veHl5QUvLy/4+PjAzc2tXR9NymQy3Lt3j0u4iYmJuH37\nNpKSkiCTyWBhYYG+ffsiMDAQgYGBGDhwoFYm23adeC5fvoz169fj5MmTCAgIwOzZszF58mStud5y\n584d7N69G//9738hlUrx7rvvYsWKFbC3t+c7NN4UFBQgIiICly9fRlRUFOLi4iAUCuHh4YGAgADu\nFRgYqHKq09HJ5XKkpqYiLi4OcXFxiIqKwo0bN6BQKODu7o6wsDCEhYVh2LBh2vEDR+1QSkoKhYWF\nEQAaOnQoRURE8B3SE1VWVtL69evJxsaGTExM6Msvv6S6ujq+w9IIhUJBsbGxtHz5cvLx8SEAJBaL\naeDAgfTRRx/RuXPnqLq6mu8wdVJFRQWdOnWKFi9eTAEBASQUCklPT48CAgLos88+o/j4eN5ia1eJ\np6amhlauXEkSiYT69Omj9QnncVVVVbR69WoyMjIiDw8POn/+PN8hqYVcLqfTp0/Te++9Rw4ODgSA\n3N3dadGiRfTHH39QZWUl3yG2S6WlpfTrr7/SnDlzqHPnzgSAXF1dacGCBRQREUFKpVJjsbSbxJOU\nlER+fn5kbm5OW7duJblczndIzy0jI4NeeuklEggEtHjx4nZz9JOSkkKffPIJubi4EADy9vamZcuW\nUWRkpEZ3eqZBYmIiffLJJxQQEEAAyMnJiZYtW0ZpaWlqL7tdJJ49e/aQkZERDRgwgO7fv893OG1m\nz549ZGJiQoGBgZSVlcV3OM9FKpXSrl27KDAwkACQs7MzffzxxxrZuZnWS0hIoA8//JDs7e1JIBBQ\naGgoHT58WG0/4DqfeD7//HMSCAS0dOlSkslkfIfT5m7fvk0+Pj7UpUsXSkxM5DucVistLaW1a9eS\ng4MDSSQSmjp1Kp0/f54UCgXfoTFPIJPJ6MSJEzRhwgQSCoXk7u5O3377LVVVVbVpOTqdeN5//30S\nCoW0c+dOvkNRq5KSEho8eDBZWlrS1atX+Q7nicrLy2nFihVkampK5ubmtHTpUsrJyeE7LOY53L17\nl+bOnUtGRkbUqVMnWrduHdXW1rbJsnU28Xz88cckEono6NGjfIeiEbW1tTR69Giytram27dv8x1O\nEzKZjL777juytbUlKysr+uqrr6i8vJzvsJg2UFhYSB9//DEZGxuTi4sLHTx48IWvyelk4vn+++9J\nIBDQ7t27+Q5Fo6qrq2nAgAHk4uJChYWFfIfDuXbtGvn4+JBEIqEPPviAiouL+Q6JUYOcnByaMWMG\n6enp0YABA17oB1DnEk9iYiIZGhrSqlWr+A6FF8XFxeTi4kLjxo3j/U6QQqGgtWvXklgspuHDh9Od\nO3d4jYfRjPj4eOrbty8ZGxs/92UOnUo89fX15OfnRwMHDmyXF5Jb6/LlyyQSiWjHjh28xVBcXExD\nhgwhfX192rBhA7to3MHU19fTihUrSCgU0oQJE565kqdOJZ6tW7eSgYEBpaen8x0K75YvX05WVlZU\nUlKi8bJzcnLIx8eHXF1dea39yvAvIiKCbGxsaODAgVRaWtrq+XQm8ZSVlZG1tTUtXbqU71C0Qnl5\nOdnY2NDixYs1Wm56ejq5ubmRt7c3ZWdnq62chw8f0uHDh2n16tVqK+NxZWVlzzXfs3zh2qOUlBTq\n0qUL9erViwoKClo1j84kns2bN5OpqWmH38iP2rhxI5mYmGjsEYPa2lrq3bs39enTh4qKitRWTnJy\nMs2dO5cAkIeHh9rKIWq4G7d27VoaOHAgCYXCVs9XW1tLq1evpv79+5Oenp4aI9QNmZmZ5ObmRuHh\n4a067daZxNOrVy+aNWsW32FoleLiYjIwMNDY3b1//vOfZGFhQffu3VN7WbW1tRpJPEQNz/hZWlrS\nsz4z/bzztVc3b94kQ0ND+vzzz586rU58YsnJyQSALl++zHcoWueNN96g8PBwtZdz7tw5AkDHjh1T\ne1mNNJV4iIg8PDyeK4E873zt1fbt20koFFJCQsITp9OJXiYiIyNhYmKCfv368R2K1hk+fDiio6NV\nGh5Xh08//RSjR4/GhAkT1FoOo9vmzJkDPz8/rF69+onT6UTiiY2NRVBQkNpa6K+ursb+/fsxadIk\nBAcHIyYmBn369IGLiwsuX76M1NRUvPzyy7C2toanpyf+/vtvbt7vv/8eAoGAa26yoqICGzduVBmm\nTgMHDkR1dTUSEhLUVsbdu3cRFRWFDz74QG1ltMadO3fw2muvYdmyZZg2bRoGDRqEW7duAXixbfio\ntLQ0jBs3DpaWlggMDMTFixe5cTU1NVi0aBFmz56NVatWYcWKFaiurm51jE9DRIiJicGHH34IV1dX\nPHz4EK+++iqsrKzg6+uLo0ePPrUchUKBS5cu4YMPPoCrqytycnIQGhoKZ2dnlJaWPutH/swEAgE+\n+OADHD9+/MnlaeYA7MUMHTqU5syZo7blKxQKSktLIwBkZmZGJ0+epKSkJAJALi4u9NVXX1FZWRld\nv36dAFBoaKjK/O7u7k0Ot5sbpg719fUkEAjol19+UVsZ33//PZmYmGi8qRE8dqrVrVs3cnd3J6KG\n9TY3NycfHx8ievFt2HjKtGDBAvrzzz/p3//+NxkZGZGenh7dvHmTZDIZBQUF0cyZM7mKm3fv3iWh\nUKiynZ8U49PI5XI6ceIEGRgYEAB6//33KSIigg4cOEAmJiYqlxtaKkcqlVJUVBQZGhoSAFqzZg2d\nPXuW3n33XY3dhCguLiaBQEAnTpxocRqdSDw9e/ZUe01lpVLZZEd3dHRU2amUSiXZ2NiQubm5yrzN\nnedr8tzfzMyM/vOf/6ht+YsXL6bAwEC1Lb8lj2+PjRs30sGDB4moIdG4u7uTSCTixrfFNnz0+bLN\nmzcTAJo+fTp9++23BICSk5NV5uvevbvK8p8WY2s0LvPRJ8I3bdpEAOjNN99sVTk9evQgALw9vuLk\n5EQbN25scbxO9KtVU1Oj9q48mjstMjU1bTKNlZUVUlNT1RrLszI2Nm5yyN+WamtrNd4ofnMWLVqE\nqqoqbN++HSUlJairq1O5ttUW29DMzIz7/+WXX8YHH3yA5ORk7rTB1dVVZXo9PdWrFU+LsTUal/lo\nB4Hjx4/HwoULkZaW1qpyGj8LKyurZyq7rTxtn9SJazyWlpYaOT/VVSUlJWrdwaysrFBYWKi25bfW\n1atX4efnB3d3d3z88cdqb7Tfzs4OAODs7IycnBwADV0k8RGjo6MjAKBLly5qLactEBEKCgqe2Ki8\nTiQea2trrdjxW9L469LYA+Wjnf2RmjvxqKioQF1dHaytrdVWRp8+fZCamsp78p8+fTpkMhlGjRoF\nAFwnfer6jLOysgAAY8eOhaenJwDg1KlTvMTYmPDCwsLUWk5buHPnDkpLS+Hv79/iNDqReHx8fBAf\nH6/WMmprawGobjiZTAagoZ/sRo29RTb29giA2ylXr16NtLQ0bNmyhUtCZ86cUZm2rV2/fh0A4Ovr\nq7Yyhg8fDn19fRw+fFhtZTyusevoR3vnzMvLQ05ODv78808cOHAAZWVlABp+/bOysl5oGzb+eJSU\nlHDL2LRpE8aPH4+3334bS5YsgVAoxMqVK3H69GnU1NTgwoULyM3NBQDcv3+/VTE+i0dPnc6dOwd/\nf3+89957rSqncR2rqqqeqcy2cODAATg4OKBv374tT6SB60wv7NixYyQUCtXWsNTDhw9p4cKFBIAk\nEgmdPXuWTp8+zd2xmDdvHhUVFdHWrVsJAAGg9evXc23ipKamUlBQEBkZGVF4eDilpqZSSEgITZ06\nlQ4dOkRSqVQtcRMRrVmzhpycnNS2/EazZs0iNzc3ta5Lo3v37tG8efO4z3rTpk1UUlJC27ZtIzMz\nMwoMDKSYmBjavHkzWVhY0Pjx4ykpKemFtuGff/5JY8eOpdDQUJo1axbNmzePtm3bpnInLyIigoKD\ng8nExITc3Nxo7dq1NGjQIHrvvffo3LlzJJfLnxhjax8zabzQvWHDBiosLKT8/Hxau3atyl2plsoZ\nPnw4zZ8/n1vHWbNm0fXr19t2Az1BSUkJWVpa0qeffvrE6XQi8RQWFpJYLKb9+/fzHYrWCQoKorfe\nekvt5WRlZZGZmRm9//77ai+ro9Pl2tCvvvoqOTo6PvWBW51Zu5dffpmGDRvGdxha5datWwSALl26\npJHyfvrpJwJABw4c0Eh57U3jUciTXikpKTqbeDZv3kxCoZDOnTv31Gl1Zu1OnjxJAoGAtf/yiBkz\nZlD37t012hLh+++/T0ZGRnT69GmNldnRODk5EQCd6thw3759JBKJaM2aNa2aXmcSj1KppODgYBo6\ndCjfoWiFW7dukVAo1PjRh1wup3fffZckEgkdPnxYo2W3d5WVlbRixQru6GfGjBkUHR3Nd1hPtWPH\nDtLT06Nly5a1+kdQZxIPUcPFPWj4CWltpFAoaNiwYeTv789Lk6NKpZI++OADEgqF9Mknn3ToZmg7\nspqaGpo3bx4JBAJat27dM82rU4mHiOjdd98lS0tLyszM5DsU3jQ2sB4bG8trHNu3bydDQ0Pq378/\n3b17l9dYGM2Kj48nb29vsrCwoEOHDj3z/DqXeKqqqsjT05NCQkLarHMxXXLx4kUSi8X09ddf8x0K\nETX0Wd+7d28yMTGhtWvXdsht0pGUl5fT8uXLSSKRUGho6HMfAOhc4iFq6OLG0tKSJkyYoPEnpvkU\nHx9P5ubm9Oabb/Letc2jpFIpffrpp2RsbEyurq506NAhrYqPeXFyuZx27txJdnZ2ZGVlRZs3b36h\n03ydTDxEDV28GBkZ0dSpU6m+vp7vcNQuPj6e7OzsaPjw4VRXV8d3OM16tMO3gIAAOnz4cIf6YWiP\npFIp7dq1i7y8vEgsFtOCBQva5Il3nU08RA21TU1NTWnEiBFUUVHBdzhqc/bsWTIzM6Phw4frxHrG\nx8fT66+/TkKhkNzc3Gjr1q0qTTww2q+kpITWrFlDDg4OJJFIaMaMGW3aYaNOJx4iori4OLK3tycf\nH5+ntvOqa5RKJX399dckFotp8uTJWnuk05L09HSaP38+GRsbk5mZGU2bNo3Onj3LTsO0lEKhoMjI\nSJo9ezYZGxuTqakpzZ8/n7Kystq8LJ1PPEQNXWsMHDiQDA0N6bvvvmsXO/bDhw9p1KhRJBaLad26\ndTq9ToWFhfTNN99Qz549CQD16NGDVq9eTWlpaXyHxhBRQkICrVq1irp06UIAqH///vTdd9+p9ei6\nXSQeoob+kT766CMSCoUUEhJCN2/e5Duk5yKXy+nbb78lCwsLcnNzo5iYGL5DalNxcXE0f/58sra2\nJgDk7e1NK1asoCtXruh0ctUlcrmcIiIiaNGiRdS1a1cCQI6OjrR06dImLSyqS7tJPI2uX79O/fv3\nJ5FIRHPmzKEHDx7wHVKrKJVK+vXXX6lXr14kkUhoxYoV7fq6iEwmowsXLtCCBQvIzc2N2/nfeust\n2rNnT4eup6UOd+7coZ07d9LEiRO5pO/p6UnLli2j6OhojVdEbXeJh6jhXHX37t3k6upK+vr6NGfO\nHEpNTeU7rGbV19fTkSNHqE+fPiQQCGjChAl0+/ZtvsPSuBs3btDq1atp6NChXGPnXbt2pZkzZ9Le\nvXspJSWFl1raukgmk9HNmzfpv//9L02bNo179svExIRGjhxJGzZs4P37ICDSgibL1EQmk2HPnj1Y\nt24d7t+/j6FDh2L27NkYP3682ttwfpr79+9j9+7d2LVrF/Lz8/HSSy/hk08+Qa9evXiNSxtIpVLE\nxMTg4sWLuHDhAv7++2/U1dXBzMwMAQEB6Nu3LwIDA9G7d2+4u7tDKBTyHTJvZDIZ0tLSEB8fj7//\n/hvXrl1DfHw8ampqYGxsjH79+mHo0KEYOnQogoKCIBaL+Q4ZANCuE08jpVKJ06dPY+fOnTh16hQM\nDAwwevRovPrqqwgPD9dIg9hEhKSkJPz22284evQorl+/DgcHB7zzzjuYOXNmk0bEmf+RyWS4efMm\n98W6du0akpOToVAooK+vDy8vL3h6esLX1xeenp7w8PCAm5ubSmPpuq68vBwZGRlITU1FUlISUlJS\nkJSUhLS0NMhkMojFYvTs2ROBgYEIDAxE37594e3trba+6F5Uh0g8j8rLy8Px48dx7NgxXLp0CUql\nEr6+vggNDUX//v3h4+MDLy8vSCSSFyqnsLAQCQkJuHnzJiIjIxEZGYmioiJIJBKMHTsW//znPzF0\n6FCt3TG0XXV1Nffla/ybnJyMjIwMrv1ha2truLq6wsXFhftrb28POzs72NjYwM7OjrdeGB5VWFiI\nwsJCFBQUIC8vDwUFBbh//z4yMzORkZGBzMxMrr1rkUgEd3d3Lsn6+vrCy8sLXl5e0NfX53lNWq/D\nJZ5HlZaWIjIyEhEREYiMjMSNGzcgk8kgEonQrVs3ODs7w8HBAZ07d4aFhQXMzc2hp6cHc3NzyGQy\nVFVVQSaToby8HPn5+cjNzUVeXh7S0tJQUFAAALCxsUFwcDBCQ0MREhKC//u//0NaWhpu3LihVT0D\ntBe1tbW4d+8eMjIyuFfjFzgrKwsFBQUqbTKLRCLY2trC1tYWZmZmMDExgYmJCSwsLGBqagoTExMY\nGhpy272RWCxW2X6VlZUqbSSXlZWBiFBdXY2qqipUVVWhrKwMlZWVqKqqQmVlJQoKClBYWKgyn1Ao\nhI2NTZOE6erqCldXV3Tt2vWFfxS1QYdOPI+TyWQqh7JZWVl4+PAhsrOzUVFRgfLycigUClRUVEAk\nEsHU1BQSiQSmpqaws7ODg4MDHB0d4ebmBl9fX/j6+sLe3l6ljJycHPTs2RNvvvkmduzYwdOadlwK\nhYI7wjh8+DDWrFmDf/7zn7CwsEBlZSWXGMrLy1FRUYHKykrU1dWhvr5epZ8oqVTKNS4PNPQj9WhC\nMDU1hUgkgqGhIUxMTGBqagpLS0vuf1NTU9jY2MDGxgb29vawtbXl3j/eV1d7xBIPD44dO4bXXnsN\nv/32G8aOHct3OB1S4w/AxIkTsX37dr7D6XBY4uHJ1KlTcfbsWdy6dYvrOI7RDKVSiREjRiAzMxPx\n8fHslJcH7f+YTktt374dhoaGXD9JjOZs3rwZERER2L9/P0s6PGGJhyfm5ubYv38/Tp48iV27dvEd\nToeRnJyMVatW4dNPP0W/fv34DqfDYqdaPFu6dCl27NiB+Ph4dO/ene9w2rW6ujr069cPxsbG+Ouv\nvzp0xUO+scTDs8Yvg5GRESIjI9mXQY2WLFmC7777jiV5LcBOtXimr6+PgwcP4saNG1i7di3f4bRb\nkZGR2LRpE7Zu3cqSjhZgRzxaYtOmTVi6dCkuX77Mrj20sfLycvTq1Qu9e/fG8ePH+Q6HAUs8WoOI\nMGbMGNy5c4fVam5jU6ZMwblz51jVBS3CTrW0hEAgwH/+8x+UlpZi6dKlfIfTbhw9ehSHDh3C7t27\nWdLRIuyIR8uwWs1th9VO1l4s8WghVqv5xbHaydqNnWppIVar+cWx2snajSUeLWRubo59+/axWs3P\nidVO1n7sVEuLsVrNz47VTtYNLPFoMVar+dmx2sm6gZ1qaTFWq/nZsNrJuoMd8egAVqv56VjtZN3C\nEo8OICKMHj0ad+/eZbeGW8BqJ+sWdqqlAwQCAf773/+ipKQEy5Yt4zscrcNqJ+sedsSjQ44ePYrX\nX3+d1Wp+BKudrJtY4tExrFbz/7DaybqLnWrpGFar+X9Y7WTdxRKPjnm0VvPu3bv5Doc3rHaybmOn\nWjqqI1eUY7WTdR9LPDqqI3/5OnLSbS/YqZaOaqzVHB8f36FqNbPaye0DO+LRcR2pVjOrndx+sMSj\n4zpSrWZWO7n9YKdaOq6j1GpmtZPbF3bE00401mo+ceIExowZw3c4bYrVTm5/WOJpR9rjqQirndw+\nsVOtdmTHjh3trlYzq53cPrHE0448rVZzYWEhKisreYjsyZRKJe7fv99kOKud3I4R0+4sXryYjI2N\n6c6dO9ywH374gYyMjOjtt9/mMbLm7dmzhwQCAS1evJikUikREUmlUurVqxcFBweTXC7nOUKmrbHE\n0w49+qXNz8+nV155hQAQALKysiKFQsF3iComTJhAenp6JBQKycvLixISEppNnkz7wS4ut1MJCQkI\nCAiAsbExqqurIZPJuHFxcXHw9/fnMbr/kcvlsLS0RFVVFQBAJBJBIBBAIBBg+/btmDlzJs8RMurA\nrvG0Q1KpFAcOHIBcLkdFRYVK0pFIJDh9+jSP0amKiYnhkg7QkIhkMhnkcjn27duHrKwsHqNj1IUl\nnnYmMTER/v7+2LhxI4gISqVSZbxMJsPJkyd5iq6pM2fOQCKRNBmuVCoRExMDLy8v7Nu3j4fIGHVi\np1rtSFRUFIYMGQIigkKhaHE6oVCI4uJimJubazC65vXq1Qu3bt1qcbxAIAARYcOGDVi8eLEGI2PU\niR3xtCPdu3dH3759nzqdUqnE+fPnNRDRkxUUFCAhIeGJ0wiFQjg4OCA0NFRDUTGawBJPO2Jra4vL\nly/jyy+/hFAobLGNHpFIhDNnzmg4uqbOnDkDgUDQ4niBQIAxY8YgMTERgYGBGoyMUTeWeNoZoVCI\nZcuWITo6Gp07d4ZIJGoyjUwmw6+//spDdKpOnz7dbHIUi8WQSCTYtGkTjh8/DisrKx6iY9SJXeNp\nxyoqKjBnzhwcPHiQu1byqKSkJHh7e/MSm1KphLW1NUpLS1WGC4VC+Pj44MiRI+jRowcvsTHqx454\n2jEzMzMcOHAAP/30E4yNjVWOfsRiMa+31ePi4lSSjp6eHgQCAebOnYtr166xpNPOscTTAbz++utI\nSEhAnz59uFMbuVyOU6dO8RbT6dOnIRaLATQkQUtLS5w6dQpbt25t9vY6076wU60ORCaT4f/+7/+w\nfv16AA0XmcvKymBkZNTs9BUVFSgpKUFJSQmkUilqamoANFRQrK2tBdCQNBqfGpdIJDA2NoaZmRms\nrKxgZWXV4sXj/v3748qVKxAIBBg1ahT27NkDGxubtl5lRkuxxNMBXbhwARMnTkRhYSE+/PBDmJub\n48GDB8jKysKDBw9QXFyMkpISyOXyFy7L0tIS1tbWcHR0hLOzM1xcXGBtbY1FixZBKBRi06ZNmDt3\n7hPvbjHtD0s87RwRISUlBVeuXMHNmzeRmJiIpKQkPHz4EEDDUYqLiwucnZ25l42NDXfE0vgyMjKC\ngYEBgIYeLhqPkurr61FdXQ2g4fStsrKSO1JqTGCFhYXIzc3FgwcPkJGRgczMTO4xCTMzM3h7e8PP\nzw++vr7o168f+vTpw0632jmWeNoZuVyO6OhoXLhwAbGxsbhy5Qp3OuXr6ws/Pz/ui+7l5QUnJyde\n4iwpKUFycjKSkpKQkJCA5ORk3LhxA6WlpTAwMEBAQAD69euHIUOGYNiwYTA2NuYlTkY9WOJpB3Jz\nc3Hq1CmcPn0a58+fR3l5OVxdXTFw4ED0798fAwYMQK9evZqt06NNiAh37txBbGwsYmNjERMTg4SE\nBIhEIoSEhGDkyJEYM2YMb1UAmLbDEo+OKi0txYkTJ3DkyBHuDtHAgQMRFhaGsLAwBAQE8B1imygq\nKsLFixdx7tw5nDx5Erm5ufD29sbrr7+OyZMns9vuOoolHh1CRDh37hy2b9+OU6dOQSKRYOzYsZg4\ncSJGjRrFXYNpr5RKJaKiovDjjz/i559/RkFBAQYMGIA5c+bgjTfegL6+PsDF2gkAACAASURBVN8h\nMq3EEo8OqK6uxu7du7F9+3akpqZi8ODBmD17Nl566aUO2wC6QqHA+fPnsXv3bhw7dgwWFhaYOXMm\n/vWvf6Fz5858h8c8BUs8WqympgbfffcdvvrqK1RXV2PKlCn417/+hZ49e/IdmlbJy8vDzp078f33\n36O0tBSzZ8/G8uXL4eDgwHdoTEs00Lwq84wUCgXt2LGD7O3tydjYmJYuXUqFhYV8h6X1pFIpbdu2\njTp37kyGhoa0ZMkSqqio4Dssphks8WiZuLg4CgwMJLFYTIsWLaL8/Hy+Q9I5tbW1tGXLFurUqRM5\nOTnR0aNH+Q6JeQx7VktLyGQyLFu2DEFBQTAwMEB8fDw2btwIW1tbvkPTOQYGBpg/fz5u376NsLAw\nvPbaa3j55ZdRXFzMd2jM/8eu8WiB7OxsTJw4ETdu3MCmTZswc+ZM9ghBG4qIiMD06dMBAIcPH0b/\n/v15johhRzw8i4mJgb+/P0pKSnDlyhXMmjWLJZ02FhoaiuvXr8PPzw+DBw/Grl27+A6pw2OJh0eX\nLl3CiBEjEBwcjKtXr8LHx4fvkNqtTp064cSJE1ixYgVmzZqFLVu28B1Sh6bddejbsUuXLmHMmDEY\nN24c9u3bx7VNw6iPQCDAZ599BnNzcyxcuBByuRwffvgh32F1SOwaDw8yMjIQGBiIYcOG4eDBgy02\nys6oz6ZNm7B48WKcPHkSo0aN4jucDoclHg2TSqUIDg6GUqlEdHR0i41wdVT19fW4evUqQkJC1F7W\n22+/jd9++w1///033N3d1V4e8z8s8WjY2rVrsXbtWty8eRNubm5tvnwiwpEjR7B3717k5OTAxsYG\nBgYG6NKlC7p06YKioiJ8/fXX3PQ5OTk4c+YMTp8+jaysLMTExDRZ3u7du/HNN99AJBKhqKgIubm5\nABoaFBs6dOgT4+nXrx8GDx6MDRs2PHG6kpISbNiwAVu3bkVNTU2ThunVoba2Fv369YOLiwtOnDih\n9vKYR/BUf6hDKigoIHNzc/riiy/UtvwhQ4ZQ165dKTY2lpRKJRE11ITet28fWVlZ0TvvvNNkvszM\nTAJAHh4eTcbt2rWLANChQ4e4YceOHSMzMzPau3evyrQPHjxoMv+bb75Jq1atalX8SqWSbGxsSJO7\n5YULFwgAnT17VmNlMqzmskatWrWKHBwcqKqqqs2XrVAoKDg4mCwtLamoqKjZaS5evEhvvvlms+Na\nSjyhoaEEgMrKylSG//jjj7RmzRrufXp6OoWEhLzAGjTw8PDQaOIhIho9ejQNGTJEo2V2dOx2ugb9\n/PPPmDJlilpa0zt27Biio6OxfPlydOrUqdlphgwZgtdff/2ZlqtUKgE0XIylR05/Xn31VXh6egJo\nqAA5duxYFBYWPmf0/Jo1axb++usv7hSS0QC+M19HkZCQQAAoNjZWLcufNGkSAaC///77ueZHC0c8\nP/30EwEgADRu3DjKy8trMs0XX3xBAMjc3Jzee+89IiKSy+V0+PBhmj59Og0aNIibtrKykj7//HOa\nMmUKzZs3jwYPHkybNm3iTgsfP+LZsGEDSSQSWrRoEUVGRj7Xuj1NTU0NmZqa0nfffaeW5TNNscSj\nIYcOHSKRSERyuVwty+/bt2+zp0St1VLiISLau3cvmZubEwCytLSk7777rsl6NDf/49eO6uvrKTQ0\nlKZOnUoKhYKIiHbv3k0A6LfffiMi1cRTXFxMU6dOpZs3bz7XOj2LwYMH0/vvv6/2cpgG7FRLQ7Kz\ns+Ho6Ki2OjuNy23s+6otTZs2DXfv3sWcOXNQXl6OOXPm4KWXXuJ6imhJly5dVN5/++23iIiIwKpV\nq6Cnp8cte/fu3U1un6enp2PJkiX45ptvNNL+kJOTE7KystReDtOAJR4NqaiogJmZmdqW39gAekpK\nilqWb21tjR07diAuLg5dunTBqVOnsHTp0ifO8/gzZ5cuXQIAlZ4tRCIRZsyYAUtLS5Vpx4wZg+rq\nalhbW7fNCjyFhYUFysvLNVIWwxKPxjg4OCAvL09tyw8NDQUAxMbGttkyIyIicP36dZVhvXv35hLI\njz/++EzLy8/PBwCkpaU9ddqvv/4ahw8f5no9VbecnBw4OjpqpCyGJR6NcXJyQklJyVNPT57X1KlT\n4e/vjy1btrR4d0YqleKHH35o9TJNTU2xaNEiKBQKleHu7u6ws7Nr0lbQ03oe7dWrFwDgyy+/5O6W\nAQ2PkPz+++8q044ZMwYrV67EypUrm4xTh8zMTN76GOuQ+L7I1FEUFhaSSCSiI0eOqK2M5ORkcnZ2\nJjc3Nzp69CjJZDIiIqqurqbz58/TsGHDKCYmpsl81dXVBIC6deumMryiooIA0FtvvaXShOhvv/1G\nAGjXrl3csK5du5KRkRFlZmY2md/BwYGIiO7du0dGRkYEgIYOHUrbtm2jVatW0ezZs7mLza6urgSA\nFAoFyWQyGjp0KJmbm9P169fb7oN6TG5uLunp6dGJEyfUVgajiiUeDQoLC6M33nhDrWVUVFTQunXr\naPTo0eTq6ko+Pj7Uq1cvWrlyZbMVCy9cuECzZs0iACQSiWj9+vUUHx/Pjbe3tycAZGVlRWFhYRQW\nFkYDBgygY8eOqSxn+fLlZG9vTz///DMREVVVVdHy5cu5W/EbN26k8vJyunXrFo0YMYIsLCzI0dGR\nFixYQGVlZVRcXEyff/45N/2XX35J2dnZ9MMPPxAAMjU1pTVr1lBpaWmbf2Zbtmwhc3Nzkkqlbb5s\npnnsWS0N2rt3L2bOnInExETWEZ2WqKurg5eXF8LDw7Fz506+w+kwWOLRIKVSib59+8LV1RXHjh3j\nOxwGwFdffYVPP/0UqampTW7/M+rDLi5rkJ6eHtauXYtffvkFP/30E9/hdHgpKSlYvXo1li5dypKO\nhrEjHh4sWLAAu3btQkxMDPz8/PgOp0OqrKxE//79YW5ujkuXLkEikfAdUofCEg8PZDIZwsLCkJWV\nhQsXLsDV1ZXvkDqUmpoaTJgwAQkJCfj7779Z/R0esFMtHojFYhw9ehQWFhYYPHgw7ty5w3dIHUZl\nZSVGjx6N69ev4+TJkyzp8IQlHp5YW1vjwoUL6Ny5MwYPHoyIiAi+Q2r3MjIyMHToUNy5cweXLl2C\nv78/3yF1WCzx8MjCwgJ//vkngoODERYWhrVr12qkyc+O6MSJEwgICIBcLkdkZCTrSohnLPHwzNTU\nFEePHsVXX32FTz75BCNGjMDdu3f5DqvdKCsrw9y5c/HSSy/h5ZdfRkxMDLp27cp3WB0eSzxaQCAQ\nYOHChYiMjERBQQH8/Pzw+eefo66uju/QdNrBgwfh5eWFY8eO4cCBA9i1axcMDQ35DosB2LNa2kYm\nk9HmzZvJ1NSUnJ2dafPmzawq/zOKjIyk0NBQEggENG3atBbboGb4w454tIxIJMKCBQuQnJyMkSNH\nYsmSJfDx8cGePXvYEdATEBHOnTuHwYMHY9CgQTAyMsLVq1exd+/eFtugZvjDEo+WcnJyws6dO5Ga\nmorQ0FDMnj0bzs7O+Oijj1hLeY+oqKjAtm3b4OPjg/DwcOjr6yM6Ohq///47+vbty3d4TAtYBUId\nkZubi507d+I///kPCgsLMWrUKEycOBHjx4+HiYkJ3+FplEKhwPnz5/Hjjz/i559/hlKpxJQpUzB3\n7lyuzR9Gu7HEo2NkMhmOHTuGvXv34uzZsxCLxRg7dixee+01hIWFNWlCtL2QSqWIjIzE8ePH8fPP\nP6OgoABBQUGYMmUKpk+fDgsLC75DZJ4BSzw6rLi4GD///DN+/PFHREZGAgCCgoIwcuRIhIeHIyAg\nQGefQSIipKSk4Ny5czhz5gwuXbqEmpoa9OzZE2+88QYmTpzIbovrMJZ42onS0lKcP38ep0+fxpkz\nZ5CdnQ0DAwMEBARgwIABCA4ORp8+feDi4tKkEXZtUFBQgBs3biA2NhYxMTGIjY1FWVkZLCwsEBYW\nhn/84x/4xz/+wZ4ibydY4mmnUlNTuS9xTEwMkpKSoFAoYGJiAi8vL/j5+cHLywtdu3aFs7MznJ2d\nYWNjo9aYKioq8ODBA2RkZCAzMxPJyclITk5GYmIiioqKAACurq4IDg7GgAEDMGDAAPTu3VttXQIx\n/GGJp4OorKxEYmIiEhMTuS97SkoKcnNzucc0DA0N4eLiAmtra1hZWam8JBIJzM3NATTc8jc1NQUA\n1NbWQiqVAgCqq6tRX1+PkpISlVdxcTGysrJQVlbGxWNlZQVPT0/4+vrC29sbPj4+8PPzg52dnYY/\nGYYPLPF0cHV1dcjKysKDBw/w4MEDZGZmori4GCUlJSgtLeWSh0wm4xKHXC5HZWUlAMDAwICrDWxq\nagqRSMQlK0tLS1hZWaFTp07o3LkznJ2d4eLiAhcXlw53J45RxRIP81yWLFmCv/76C1euXOE7FEYH\nsQqEDMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8\nDMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNo\nHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8DMNoHEs8\nDMNoHEs8DMNoHEs8DMNoHEs8DMNonICIiO8gGO328OFDvPTSS6ipqeGGFRYWora2Fs7OztwwoVCI\njRs3Yvjw4XyEyegQEd8BMNrP0NAQ8fHxkMlkTcYlJiaqvFcqlZoKi9Fh7FSLeSpzc3OMHj0aItGT\nf6c6deqEYcOGaSgqRpexxMO0ytSpU6FQKFocL5FIMHXqVAiFQg1Gxegqdo2HaRWpVApra2tUV1e3\nOE1sbCz69eunwagYXcWOeJhWMTAwwKuvvgqxWNzseCcnJwQFBWk4KkZXscTDtNrkyZObvcAskUjw\n9ttvQyAQ8BAVo4vYqRbTagqFAra2tigpKWkyLjExET4+PjxExegidsTDtJpQKMTkyZMhkUhUhnt5\nebGkwzwTlniYZzJp0iTU19dz78ViMd566y0eI2J0ETvVYp4JEcHZ2RnZ2dkAAIFAgPT0dLi6uvIb\nGKNT2BEP80wEAgGmTZsGsVgMgUCAoKAglnSYZ8YSD/PMJk2axN3dmjZtGs/RMLqIPavVQZSVlYGI\nUFdXxz3s2TgMAORyOSorK1ucX6FQoKKignvv6OiIhw8fQl9fH0eOHAEAGBkZQV9fv8VlGBoawsDA\ngHtvZmYGoVAIiUQCY2NjlWFM+8au8WiJuro6lJWVoaysDKWlpaiqqkJpaSlqa2shlUq5/2tra1FW\nVoaamhrU1taivLxc5X/gfwlFKpWitraW5zV7PmKxGCYmJgD+l4yMjY1haGgIMzMzmJiYwMDAgPvf\n0NAQpqamMDU1haGhIUxMTLhxFhYWsLCwgKWlJYyMjHheMwZgiafNSaVSFBUVoaCgAPn5+SgqKkJR\nURFKSkpUEsujfxsTSXMajxIsLS1b9T/wvy9qc1/e5oY1Mjc3h55ey2ffTxtfXl7+xKfTHx9fWlrK\nfWaNCbK5pNk4rLKyElKpFJWVla36vzkSiYRLQo8mpMb/raysYG1tDWtra9jY2MDW1ha2trbcZ8a0\nDZZ4WqGurg55eXnIzs5GdnY28vLymiSWxvdVVVUq8+rr68PGxgaWlpbcDt7cTv/4MBMTEy6RMM+n\nsrKSO3J8PNE3939paSlKS0tRWFjY5Jk0Q0NDWFtbc4nIxsYG1tbWsLOzg4ODAzp37gxHR0d06dKF\nO21kWtbhE09dXR0yMjKQkZGBvLw8PHjwoEmSyc/P56YXCoWws7ODra0t7O3tVX4dG9/b2NjAxsYG\ndnZ2MDU15XHtmOdVW1uLwsJCFBQUoKCgAEVFRSgsLER+fj4KCwtRWFiIoqIiPHz4EA8fPkRdXR03\nr7m5OTp37gwnJyc4ODjA2dkZDg4OcHJyQpcuXeDm5gZzc3Me145/HSLxlJaWIj09vdlXZmYm19yD\ngYEBHB0d4eDgoPLX3d2d+9/Z2fmp7dIwHU9paSlyc3ORl5fH/U1PT1cZ9vDhQ+5ivqWlpcp+5e7u\nzr28vLza/bWodpN4ysvLcfv2baSkpKj8zcjI4GraGhgYwM3NDe7u7nBzc1N5ubq6wsLCgue1YNqz\nmpoaZGZm4v79+0hPT8f9+/dVXo03B/T09ODk5IQePXrAy8sLXl5e8PT0hJeXF+zt7Xlei7ahc4mn\nrKwM8fHxSElJQVJSElJTU5GSkoLc3FwADcnF09MTHh4e8PLyQteuXblk4+DgwHP0DNOykpISLgml\np6cjNTUVycnJuH37NsrKygA0HCl5eHjA29sbnp6e8Pb2Ru/evdG5c2eeo382Wp14cnNzkZycjKSk\nJMTFxSEuLg4pKSkgIlhYWKBr167w9vaGj48P3N3duY3B6oEw7U1paSmSkpKQnJyM9PR07v+MjAwo\nlUpYWFjAx8cHAQEB3MvLy+uJdyH5pDWJp7y8HNHR0YiKisK1a9dw48YNFBQUAADc3d3Rp08f9O7d\nG3369EGfPn3g6OjIc8QMw7/y8nLcuHED8fHx3CslJQVyuRwmJibo1asXAgICEBwcjJCQEK05MuIt\n8WRlZSEyMhLR0dGIjIxEYmIilEolPDw8EBQUxCWY3r17s2svDPMMpFIpEhMTcf36ddy4cYP7IZfL\n5XBzc0NISAgGDhyIkJAQeHt789KAm8YST2lpKf7880/88ccfuHDhArKysiAWi7lsPGjQIAQHB8PW\n1lYT4TBMh1JVVYUrV67g8uXLiIqKQmxsLCorK2FlZYVBgwZh5MiRGDVqFFxcXDQSj9oSDxHhxo0b\n+OOPP/D7778jNjYWAoEAAwcORFhYGEJCQhAUFNTubxsyjDZSKBS4efMmLl++jAsXLuD8+fOoqqqC\nt7c3Ro8ejVGjRiEkJKRJo29tpc0TT0xMDPbv349ffvkFeXl5cHBwwKhRozBq1CiEh4d3+IpTDKON\n6uvrERkZiT/++AN//PEHkpOTYWpqipEjR2Lq1KkYOXJkmyahNkk89+7dw/79+7F//37cvXsXvr6+\nmDRpEkaPHo1evXqxRsAZRsdkZGTg9OnT+OmnnxAREQFLS0tMnDgRU6ZMwYABA154+c+deOrr63H4\n8GH8+9//RkxMDOzs7DBp0iRMnz4dvXv3fuHAGIbRDllZWThw4AD279+PpKQkdOvWDTNnzsTs2bOf\n+3nCZ048VVVV2L59O7Zs2YKioiK88sorePvttxEeHs7qzzBMO3f9+nXs3bsXe/bsgVwux4wZM7Bs\n2TI4OTk924Kolerr62nTpk1ka2tLpqamtGzZMsrKymrt7O3Cw4cP6fDhw7R69Wq+Q2HagdLS0lZP\nq237XkVFBW3ZsoVcXFzIwMCAFixYQMXFxa2ev1WJJyoqivz8/MjAwICWLFlChYWFzx1wWzt69Ci9\n9tprBIAA0MWLF1uc9vLly9x0r7zyCl24cKHV5SQnJ9PcuXMJAHl4eDx1+qCgIFq8eLHKMKVSSYcP\nH6YxY8ZQ7969KTw8nMaNG0dz586ltWvX0ocfftjqeJgXo6n95nG1tbW0evVq6t+/P+np6bVqnmfd\n9zSprq6Otm3bRvb29mRjY0N79+5t1XxPTDxKpZLWrVtHIpGI/vGPf1BaWlqbBNvWqquruR1j3Lhx\nLU43ceJEMjQ0JACUl5f3zOXU1tY2u/EfPHjQZNo333yTVq1axb0vKCigIUOGUNeuXSk2NpaUSiUR\nESkUCtq3bx9ZWVnRO++888wxaavmPhO+PR6Tpvabx9XU1JClpSU9wwlHi/uetigrK6P333+f9PT0\naPLkyVRZWfnE6Vtcc4VCQe+88w6JRCJat24d90XRVgAoODiYBAIB3blzp8n43NxcGjFiBHl4eDzT\nBm+unEc3fnp6OoWEhDxxHoVCQcHBwWRpaUlFRUXNTnPx4kV68803nzsubdKaz0TTWopJU/vN455n\nedqceBpduHCB7O3tyd/fv8V9negJiWfBggVkYGBAp06dUkuAbQ0A/fTTTwSA/vWvfzUZ/+mnn9Iv\nv/zSpoknKyuLvL29n7ozHDlyhADQ+vXrnzjdzz///NxxaYvWfiaa9KSYNLXfPK69Jh4iort375Kz\nszMNGDCA6uvrm52m2TU/e/YsCQQC+vHHH9UaYFsCQDKZjJydncnIyIhKSkq4cXV1dTRgwACSy+XN\nbvCdO3dyh9xEROXl5fT111+rDHu0nMaN/8UXXxAAMjc3p/fee4+IiORyOR0+fJimT59OgwYNIiKi\nSZMmEQD6+++/W70+ZWVltGTJElq2bBktXLiQwsPDaeHChdx6VVVV0b59+2jixIk0YMAAio6Opt69\ne5OzszNFRkbS7du36aWXXqJOnTqRh4cHXbt2jYgaTp+jo6Np0aJF5OLiQnl5efTKK6+QpaUl+fj4\nqCS/1NRUevXVV2np0qU0depUCgkJoZs3b5JcLqeLFy/SggULyMXFhbKzs2nw4MHUpUsXWrRoUZPP\n5HljbVRTU0Pr1q2jd955hwICAmj48OF069YtUiqVdPz4cZo1axZ17tyZSkpKaPr06WRlZUU+Pj7c\ncprbTm2x37RmOxE1nNItXLiQZs2aRR999BEtX76cnJycVJbX0jo+HqsuJB4iopSUFDI2NqZPP/20\n2fHNJp5hw4bRmDFj1BpYW2vciBs2bGhydHHo0CHasGEDEbX8S+Pu7t5keHPDHt/4ze0MmZmZKsP7\n9u1LAKisrKxV61JRUUHdu3enTz75hBuWn59P3bt3Jzc3NyotLSWFQkFpaWkEgMzMzOjkyZOUlJRE\nAMjFxYW++uorKisro+vXrxMACg0NJaKGxHjixAkyMDAgAPT+++9TREQEHThwgExMTAgAXb58mYiI\nunXrRu7u7kTUcFfT3NycfHx8SCqVUlRUFHfdY82aNXT27Fl69913qbKyssln8ryxNpo5cyalpKRw\n78PDw8nW1pbKysooKyuLjI2NCQCtXr2aMjIyaN++fQSAgoKCnridGocTPd9+05rtJJPJKCgoiGbO\nnMldrrh79y4JhUKV5bW0juXl5U9dB221Zs0asrS0pKqqqibjmnwD6+vrSSQS6dTRDtH/dqDS0lIy\nNjYmJycn7jAvPDycu9XXUuJpbnhzw1qTeJRKpcrwfv36EQDKzc1t1bqsXLmy2el/+OEHAkBLlixp\nthwiIkdHR5WYlUol2djYkLm5ucqyunfvTgBUdopNmzYRAO5a08aNG+ngwYNE1JA83N3dSSQScdP3\n6NGDADS5jdqaz6S1scbGxnJHno+/Tpw4oRLHo8uxtbUliUTyxJgahxM9337Tmu307bffEgBKTk5W\nmabx82/tOj5pHbRVdnY2AaC//vqrybgmrQSVlZVBLpfr7FPiFhYWmDFjBrKzs3H06FHcuHED7u7u\nsLKy0lgMjz8i4u3tDQBISUlp1fxRUVEA0KSh+MGDBwMAoqOjmy2nuXkEAgGsrKy4ZjUbNTYQ9WiP\nCOPHjwcApKWlAQAWLVqEcePGYfv27fjyyy9RV1cHuVyusmwArfpsnzfWa9euwcfHB9TwI6nyGjt2\nbLPLFggEsLS05Jq8bY3n2W9as53+/PNPAGjSzfOjDXS1Zh11ka2tLQQCAQoLC5uMa5J4rK2tYW5u\njlu3bmkkOHWYP38+BAIBNm3ahG3btmHevHm8xhMaGgoAiI2NbdX0jTtlRkaGynA7OzsAUNuDto2N\nq3Xp0gUAcPXqVfj5+cHd3R0ff/wxL31LFRcXIz09vUl3MwC4RvrbyrPuN63ZTjk5OQAa1qMlmlxH\nTbp58yaICN26dWsyrkniEQgEmDx5MrZt26YzvVA2bpzGv927d8fYsWNx9epV5OTkwMfHh5uWWnhC\npPFXs7GbEqVSyf3ytjRPo0ePApozdepU+Pv7Y8uWLVzb0I+TSqX44YcfAPzvF/PUqVMq02RlZQEA\nwsLCnlje82r8cjQuf/r06ZDJZBg1ahQAcJ3xPe3zAJ7+mbSWp6cnamtrsX79epXhycnJ2LZt2zMt\n6/GYXnS/ac128vT0bHaaR7XlOmqTDRs2wM/PDz179mw6srlzs6ysLLKysqK3335b6+vvEDXUtQBA\nOTk53LALFy4QAPrtt99Upu3cuTMBoJqaGpXhL7/8MgGgVatW0Z07d+ibb77hKnn98ccfJJfLuQpn\nLi4u3Hxdu3YlIyMjyszM5IZVVFQQAHJwcOCGJScnk7OzM7m5udHRo0dJJpMRUcMdj/Pnz9OwYcMo\nJiaGG+bj40OdO3dWuX4wf/58Cg4O5q5B1NTUEADq0aMHN03jBfGKigpumIuLCwEguVzODWu8ZtEY\nBxHRnj17yN/fn1u+mZkZAaAzZ87Q/v37ycbGhgBQbGwsPXjwgFvu45XFmvtMnjfW2tpacnNzIwA0\nY8YM2r9/P3300UcUHh7OXXhtnOfRfbXx+lHjujQX04vuN63ZTvHx8SQUCsnKyor++OMPbnubmpoS\nAEpPT2/VOja372mz77//ngQCAZ0+fbrZ8S1WJDh16hRJJBKaOXNmi/fitcHx48dp7NixBIDGjBlD\n586dI6KGC4yvvPIKtwMnJSVxFwMB0Ouvv65S9T01NZWCgoLIyMiIwsPDKTU1lUJCQmjq1Kl06NAh\nSk5Opnnz5nHzb9q0iUpKSmj58uVkb2/P3Yauqqqi5cuXc9Nt3LiR23kqKipo3bp1NHr0aHJ1dSUf\nHx/q1asXrVy5skllq4qKClqyZAmFh4fTokWLaMmSJfTZZ5+RVColooZndxYuXEgASCKR0NmzZ+n0\n6dPc3ZJ58+ZRUVERbd26lYtl/fr13OMujYlnw4YNVFhYSPn5+bR27VqVJLJt2zYyMzOjwMBAiomJ\noc2bN5OFhQUNHz6c5s+fzy131qxZdP36dW6+xz+TF431/v37NG7codevWAAAB3JJREFUOLK0tCQ7\nOzuaNWsWFRQUcDE2zvPFF19QWVkZd5EcAC1btoxqamqaxNRW+83TthMRUUREBAUHB5OJiQm5ubnR\n2rVradCgQfTee+/RuXPnSC6XP3Ed79271+y+p6127txJenp6Knf7HvfEGkwnT54kY2NjCg4OpoyM\njLaOj+FRW1eIY5iqqip65513SCAQ0Oeff/7EaZ/Y98WYMWNw9epVlJeXw8fHBxs2bIBMJnv6yR3D\nMB3KL7/8Am9vbxw/fhy//vorPv744ydO/9ROd7y9vREfH48vv/wSn332Gbp164bvv/++zS4eMvxo\nvINSVVXFcySMLrt8+TKGDh2KV155Bf369UNKSgrGjRv31Pla1duXWCzGggULkJycjBEjRuBf//oX\nevTogc2bN6OysvKFg2c0p6qqCitXrkR2djaAhlvIMTExPEfF6BKFQoEjR45gwIABGDRoECQSCWJj\nY/HTTz+1uv7fczV9eu/ePWzevBl79uyBQCDAq6++imnTpmHIkCFa23MhwzAvJiEhAfv378fBgweR\nm5uLl19+GR9++CGCg4OfeVkv1Nh7aWkp9u3bh/379+PatWtwcnLC5MmTMW3aNPj6+j7vYhmG0RK5\nubk4dOgQ9u3bh5s3b8LV1RVTpkzBjBkz0LVr1+debpt1b3P79m2uQeiMjAz4+vpy/fMMHDgQYrG4\nLYphGEbNbt68yXVzExUVBTMzM7z++uuYOnUqQkJC2qTXmDbvV4uIcPnyZfzyyy/4/fffkZqaCjMz\nM4SHh3P9a7F+zxlGe1RUVODcuXNcssnJyYGdnR1GjRqF8ePHY/To0dDX12/TMtXehXF6ejrXm+jF\nixchlUrh7e2t0n+zm5ubOkNgGOYRRUVFiI6ORmRkJKKjo3Ht2jUolUoEBQVxZyn+/v5q7Q9PY32n\nA0BtbS0iIiJw8eJFXL58GXFxcairq4Ojo6NKIurZsydEIpGmwmKYdi0tLQ1RUVFcv+mpqakQCATc\nAcDgwYMxYsQIdOrUSWMxaTTxPE4ul3P9N0dFReHixYsoKiqCWCxG9+7dERAQwL38/f1ZP+sM8wQK\nhQKZmZlISkpCXFwc4uLicOXKFRQWFkIsFqNnz57cj/uwYcM0mmgex2vieZxSqURKSgquXbuG+Ph4\nxMfH48aNG6isrIRIJIKnpyf69OmDPn36wM/PDx4eHlwTDgzTkZSUlCAlJQXJycncd+XWrVuoqamB\nRCKBr68v913x9/dHQEBAm/Z9/qK0KvE0h4hw9+5d7sNtfBUUFABoaITJw8MDXl5e8PLygoeHB7y9\nvdG1a1d2J43RaUSEzMxMpKamIjk5Gbdv3+b+b2xcy9TUFD179uSSTJ8+feDr66v1+77WJ56WFBUV\nISUlBbdv38bt27eRnJyM1NRUZGZmQqlUQiwWw93dHd26dYObmxvc3Nzg7u7O/W9mZsb3KjAMpFIp\nMjIycP/+faSnp+P+/fvc/3fu3EFNTQ2AhsbFvLy84OnpCU9PT+5H1tnZWa0XgdVFZxNPS2pra7lf\nhtu3byM9PZ3boI82wmVtba2SkFxcXNC5c2c4OTnBwcEBdnZ2OrlBGe1SWlqKnJwcZGdnIy8vj0sy\njcklLy+Pa2TMyspK5cexR48e8Pb2hqenJywtLXlek7bV7hLPk0ilUm6jP/pKT0/HgwcPUFJSwk0r\nkUhgb2+PLl26wNHREY6OjlxScnZ2hrW1NWxsbGBtbc3jGjF8qaysRH5+PgoK/l9757qrKAxF4aXi\npQhSjYD390/mrYxGESIV0KJGnV90BuV4MsnxMnF/SZOyEqMlu6tFYO8VFotFzlym0ykWiwVms1ku\ni6eu65hMJje776z/qJS278hHGc93SCkLg2c+nyvN87xcahBN05QBOY4D13WVKbmuC8dx0O120el0\nwDkH5xyMsReOkrjmeDxCCAEhBNbrNYIggO/78H0fy+VSHa9WK3iehyAIkKap+nypVILruhgMBrld\n83g8Rr/fx2g0wnA4BOf8haN8L8h4/pHz+ayCLwtE3/cRBIFaAbNA9TzvproDANTrdbTbbWVEWbvW\nDMMAYwymaeb6pmmi0WjcVDf4NKSUSNMUYRhCSgkpJYQQ2O12SNMUQggkSQIhBMIwVOaStUwrSrKu\n6zps20av11MLyfXCkh07jvNWd4z+B8h4Hsx+v0cQBDeBf28iZG273X5bosU0TTDGYBgGTNOEpmmo\n1+vqmSfOOUqlUk6zLAvlcjmnZbRaLVQqlcLvqtVquXI41+PM/ggtQgiRS5h+OBzUhN9sNjifzzkt\niiKcTqecFoYh0jSFlBJhGN49L+VyGZZlwTCMQlO/p9m2Tc+MPRgynjfndDohiiIkSQIpJeI4RhzH\nkFIiSRJEUQQpJbbbrZrAu90O+/0el8sFQggAf3YHRVrG35U1irhnhNlE/4pms5nbFWiapnZshmGg\nWq0WapVKRd2BtCwLjDHouq4uWRljaLfbqs85h67rP/5uEfGzkPEQBPF0KGsXQRBPh4yHIIinQ8ZD\nEMTT0QD8evWPIAjis/gNcb/bEoTA1KIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.display import Image\n",
    "stick_model.visualize_model_setup(view=False, cleanup=False, with_parameters=True)\n",
    "Image('Model Setup.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To find more information on the model parameters check the function documentation of the Stick model. It is also possible to print the parameter cardinality to figure out the parameter names, and their input format."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "OrderedDict([('C1Stick_1_mu', 2), ('C1Stick_1_lambda_par', 1)])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stick_model.parameter_cardinality"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First thing to notice is name of the model \"C1Stick\" and the parameter \"mu\" is separated by a number \"\\_1\\_\". If multiple sticks are given to the MultiCompartmentModel, then this number will increase as more Sticks are added.\n",
    "\n",
    "\n",
    "The number after the parameter name indicates the cardinality of the parameter, meaning that the orientation of the stick \"mu\" takes two angles on the sphere as [theta, phi], and one value as parallel diffusivity lambda_par.\n",
    "\n",
    "For the example we align the Stick with some angle and give it a diffusivity of 1.7e-9 m^2/s. We obtain the right ordering for the input of the function by using the model's parameters_to_parameter_vector() function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mu = (np.pi / 2., np.pi / 2.)  # in radians\n",
    "lambda_par = 1.7e-9  # in m^2/s\n",
    "parameter_vector = stick_model.parameters_to_parameter_vector(\n",
    "    C1Stick_1_lambda_par=lambda_par, C1Stick_1_mu=mu)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, this produces a parameter vector with the 'correct' order for the model to understand it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1.57079633e+00, 1.57079633e+00, 1.70000000e-09])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "parameter_vector"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can generate the diffusion-weighted signals for these model parameters and the wu-minn acquisition scheme as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "288"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "E = stick_model.simulate_signal(acq_scheme, parameter_vector)\n",
    "len(E)  # See that this produces the signal attenuation for the entire acquisition scheme"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's assume this signal is now unknown, and we want to fit the Stick model to this signal to find best fitting model parameters. the model.fit(scheme, data) is the easiest way to fit some data. As a default, dmipy uses a global optimizer that we call brute2fine, which does exactly what the name implies: first to a global brute-force optimization and then refine solution to a local minimum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using parallel processing with 8 workers.\n",
      "Setup brute2fine optimizer in 0.0109260082245 seconds\n",
      "Fitting of 1 voxels complete in 0.267587184906 seconds.\n",
      "Average of 0.267587184906 seconds per voxel.\n"
     ]
    }
   ],
   "source": [
    "res = stick_model.fit(acq_scheme, E)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the data and the acquisition scheme we fit the stick_model using the following one-liner. We can see the correct model parameters are obtained."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimized result: [[ 1.57080726e+00 -1.57079593e+00  1.70010076e-09]]\n",
      "Ground truth:     [1.57079633e+00 1.57079633e+00 1.70000000e-09]\n"
     ]
    }
   ],
   "source": [
    "print 'Optimized result:', res.fitted_parameters_vector\n",
    "print 'Ground truth:    ', parameter_vector"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}