{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Finding faces, with processed data"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab inline\n",
      "import numpy as np\n",
      "import matplotlib.pyplot as plt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Populating the interactive namespace from numpy and matplotlib\n"
       ]
      }
     ],
     "prompt_number": 2
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "You'll need the processed data from the website.  Use your web browser to download http://nipy.bic.berkeley.edu/practical_neuroimaging/srabold.nii.gz into this current directory, or do:\n",
      "\n",
      "`curl -O http://nipy.bic.berkeley.edu/practical_neuroimaging/srabold.nii.gz`\n",
      "\n",
      "in the terminal / shell (for OSX and Linux)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import nibabel as nib"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "img = nib.load('srabold.nii.gz')"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print(img.header)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "<class 'nibabel.nifti1.Nifti1Header'> object, endian='<'\n",
        "sizeof_hdr      : 348\n",
        "data_type       : \n",
        "db_name         : \n",
        "extents         : 0\n",
        "session_error   : 0\n",
        "regular         : \n",
        "dim_info        : 0\n",
        "dim             : [  4  40  64  64 121   1   1   1]\n",
        "intent_p1       : 0.0\n",
        "intent_p2       : 0.0\n",
        "intent_p3       : 0.0\n",
        "intent_code     : none\n",
        "datatype        : float64\n",
        "bitpix          : 64\n",
        "slice_start     : 0\n",
        "pixdim          : [-1.    3.5   3.75  3.75  1.    1.    1.    1.  ]\n",
        "vox_offset      : 352.0\n",
        "scl_slope       : 1.0\n",
        "scl_inter       : 0.0\n",
        "slice_end       : 0\n",
        "slice_code      : unknown\n",
        "xyzt_units      : 10\n",
        "cal_max         : 0.0\n",
        "cal_min         : 0.0\n",
        "slice_duration  : 0.0\n",
        "toffset         : 0.0\n",
        "glmax           : 0\n",
        "glmin           : 0\n",
        "descrip         : \n",
        "aux_file        : \n",
        "qform_code      : scanner\n",
        "sform_code      : scanner\n",
        "quatern_b       : 0.0\n",
        "quatern_c       : 1.0\n",
        "quatern_d       : 0.0\n",
        "qoffset_x       : 68.25\n",
        "qoffset_y       : -118.125\n",
        "qoffset_z       : -118.125\n",
        "srow_x          : [ -3.5    0.     0.    68.25]\n",
        "srow_y          : [   0.       3.75     0.    -118.125]\n",
        "srow_z          : [   0.       0.       3.75  -118.125]\n",
        "intent_name     : \n",
        "magic           : n+1\n"
       ]
      }
     ],
     "prompt_number": 5
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Defaults for plotting"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.rcParams['image.cmap'] = 'gray'\n",
      "plt.rcParams['image.interpolation'] = 'nearest'"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Loading the image may take a little while:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "data = img.get_data()\n",
      "data.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 7,
       "text": [
        "(40, 64, 64, 121)"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.imshow(data.mean(axis=-1)[:, :, 30])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 8,
       "text": [
        "<matplotlib.image.AxesImage at 0x10542fad0>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAWwAAADwCAYAAAAkTF41AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnW9sVeUdx7+nUGH+AQHpbaFo6+RfS2kv1mEy3dRSnHMw\nCMZoNkYA94JlyRy+0XduL7RuWQw6Xxjjkm4mTLOEjSyzEWJBwWmn9IYoaIm2WmpbKVCUv4XesxeE\nm9LzfW6fp6f34gPfT9JEfj3nOc/znHN/3p7v708QhmEIIYQQ33oKLvUEhBBC2CGHLYQQniCHLYQQ\nniCHLYQQniCHLYQQniCHLYQQnjBqh93U1IR58+Zh9uzZeOaZZ8ZyTkIIIQjBaOKwBwcHMXfuXGzf\nvh0zZ87Ebbfdhs2bN2P+/Pm5mKMQQgiM8ht2S0sLbrnlFpSVlaGwsBAPPfQQ/vWvf4313IQQQgxh\n/GhO6urqwqxZszL/Li0txXvvvXfRMUEQxJuZEEJcoZhefIzKYY/WGbPzxo0bZ33s4OBgxMYWFvd/\nFmxMXzL4TWtn9qG2dDqNgoICp3WyMdPptPX5tmMy20j3/cJ6TMfaPjcXxrC1D8dlP7LNc+h6GC7P\nvO1+ujwLtsf68jkykW2fwzDM/D5X6xzVK5GZM2eis7Mz8+/Ozk6UlpaO2aSEEEJEGZXDrq2txYED\nB9DR0YGBgQG8+uqrWL58+VjPTQghxBBG9Upk/Pjx+POf/4x7770Xg4ODWL9+PY0QGfrnQRAEuPrq\nqyPHXHPNNfQa7FXJiRMnIrazZ89GbC5/xrI/WUd69TL0z1N2ftw/h2zPZ3s0fjy/pYWFhRHbVVdd\nlfnvgYEBXHXVVXTtpj/r2ZinT5+2Oo7ZAOA73/kOtQ/nzJkzWcc8ffo0Jk6cCAA4d+5c5Fi2TvZ8\nmJ7PC2MPhd039sy6PB8Xnu8zZ85gwoQJAPj9MD3zpn0eDttPdi9N17f9HA21jfQ5sn2dlIvXQaPZ\nz8HBwcxn0vaZY3Ni515gVA4bAO677z7cd999VsdebgKk7ftLnxjqvC8HmEP1mQvO+nLicvscmfS4\nseTy2jEhhLiMkcMWQghPkMMWQghPkMMWQghPGLXoaMNwsdFFCGIKta2aa1LhWUSJbbKESThl14qj\neJuuNVKUxwVMwge7PlOjXSIQ2LWmTZsWsd1yyy0R2/Tp0+mYxcXFERu7bwcPHozYTBEmbJ8GBgYi\ntpMnT0Zsx44do2P29PREbGw/WZSJKdKAHcvm7hKZxJ4l9tnq7++3HpNFlDDY82H6HLAIimzREkNx\nCWqwFTpNnyN2P5iNPbOm9Qzfk2+++cY4L33DFkIIT5DDFkIIT5DDFkIIT5DDFkIITxhVAwOrgYMg\nIgaw1HSTCGBbPcxF4LNdqm3FuLi4CBtMqGIihks6MUtjnzRpUsRWVlZGx2TC4Y033hixVVRURGwm\ngdA2W4wJM0zoAbjYza7DhK+PP/6Yjrlv376Ira+vL2IziZYMJrgmEomILW5WKntuvv7664itq6uL\nnv/VV19FbEyIjPs5sq3Oafq8M3vcIAP2mWGfTdvAATZmR0eH0VfpG7YQQniCHLYQQniCHLYQQniC\nHLYQQnhCTkXH4bASkS4lCZk44JIhFUd0dDmW2ZhQZMr8ZPO0zZwy1dxle8+6BFVWVkZsyWSSjjll\nypSI7brrrovY2DptawMD/L67ZM1ee+21VtdnczdlnR0+fDhia29vj9j++9//Rmz79++nY7L7WVRU\nFLFNnjw5Yrv++uvpmOy5s81wPXr0KB2TiY5McDUJ4Lawz4GtEGmys7UzMdAlc5Sdz8RJk18Zvvfd\n3d0SHYUQwnfksIUQwhNiFX8qKyvDpEmTMG7cOBQWFqKlpWWs5iWEEGIYsRx2EATYsWMHpk6dOlbz\nEUIIYSB2eVUXzdIlS8m2DKJLNlUcfdVFiGQihItIxgRGJs6yMU29/1i24ve///2Irby8PGIz/Q/Z\nNuPOpbGvrVjDMiVNAjbbJzYmy/JkAh/AMxBvuummiM223CwAfPjhhxFbb29vxMbEPFYaFuDZxez5\ncvls2n7mmMDnIubFyXYG7DMQbYVu07G2ZWRNQruLOBvrHXYQBFiyZAlqa2vx0ksvxRlKCCHECMT6\nhr17926UlJTg0KFDqK+vx7x583DnnXeO1dyEEEIMIdY37JKSEgDniwCtXLlSoqMQQjhy7tw5nDlz\nJvOTjVE77JMnT2aSCk6cOIE33ngDVVVVox1OCCGuSMaPH48JEyZkfrIeO9qL9Pb2YuXKlQDO/x/i\nZz/7GZYuXTra4YQQQozAqB12eXk5UqmU0zkuURq2da7jppHHTVdnUQ2mWs/DMdVvZs1Y586dG7HN\nnDkzYrvhhhvomDNmzIjYWOSIbUQFwNfJIgNM62Sw/WQqPFPxXWqr26YYmyJP2DrZ9VnkiCnqhq29\no6MjYmPRBqYoJPaZYWnkrKnx8ePH6ZisgTGLdHBJI2ew8+N8hk3nM0z33TaNnT0fbN9Mx5pQpqMQ\nQniCHLYQQniCHLYQQniCHLYQQnhC7NT0uNjWiDUdG1eItJ2TSdBiYg87lsVXmlK7b7755oitpqYm\nYmNNW1l6NcBTrFkIkW2zX4Dvva0YaLpHTOyJez9txUS2Ry4pyux+slrcLmIemzuzsesAwKlTpyK2\n7u5uqzmZ5mn7mWP7YUojZ58PUyq3LbbCMNtPU+AAExjZfWNCe9zm4IC+YQshhDfIYQshhCfIYQsh\nhCfIYQshhCfkVHQcLkTYZqwB8URHl4w3diybp23tZ4ALPex8ln0I8Ea4LKuRCU2mjDcmrLC1s2w7\n034yUcj2fNN+xqmrbFo7axbMBEZm+/rrr+mYLBuVNfFl6+np6aFjsueBiXH9/f0Rm0kkY58Pdo/Y\n2k11mpmdiWzsONOYbJ1s7rZZqy4w8d0kBDKBkdniZnma0DdsIYTwBDlsIYTwBDlsIYTwBDlsIYTw\nhJyKjsNfsrs00GTYHuvSPJTBREcm2gH2pRHnzZsXsS1ZsoQee6GTz1CYSOeSQWhblpIJUi4ZZy4i\nMINdn+0xE4pMDXNZ9idrTsvuu0nMs22UzOZuyka1LSnMrn2hmchw2L1j82QiqkkYZtdnwptLI9o4\nuDT2tc16NZUEZuuMW6LZBX3DFkIIT5DDFkIIT5DDFkIITxjRYa9btw6JROKiBrtHjhxBfX095syZ\ng6VLl9JAfiGEEGPLiA577dq1aGpqusjW0NCA+vp6tLW1oa6uDg0NDTmboBBCiPOMGCVy5513RhqB\nbt26FTt37gQArFmzBnfdddeYO20WWWCrxsZtymlqwMmwbYg6e/bsiI1FgwA8xdk20sJFhWfHujTR\nZefHjQJg0R+2a2eRDqbz2TpPnjxpfW3b6CB2HIuEAXipAbafLBrGdI9Yw90TJ05EbOyvZFNavm0a\nuu3zBYxN2rYN7PPO1mN6jl0a5saZk4lRvcPu7e1FIpEAACQSCfT29o5mGCGEEA7EjsMOgmBM4guF\nEOJKZHBw0DrHZFTfsBOJRKbiWHd3N4qKikYzjBBCXPGMGzcOhYWFmZ9sjMphL1++HI2NjQCAxsZG\nrFixYjTDCCGEcGDEVyIPP/wwdu7cib6+PsyaNQu///3v8fjjj+PBBx/Eyy+/jLKyMrz22mujnoDp\ndYqt6GibzmvCVtBi9XpN12LpzExMMwkbLP2VpQm7vIqyFQhZirNLujsTbF0EJSbqsDFZzXGTSHbk\nyBGrObF7ZEpNZ3vHniXbms4ux7J5mgRXlgbP5j5t2rSIzVS7mu0na9hrEkIZ7B6zebrsJ8O23rpL\n42nbz6Fpni7zH9Fhb968mdq3b99ufREhhBDxUaajEEJ4ghy2EEJ4ghy2EEJ4Qk7rYQ8nbu1r2wxE\n00t8W4HRRchkdZVZBiMTf5h4BJgFzjgwAYcJmddff72VDeBiC8sWZMcxYdU0T9smvIcOHaJjMoGS\nCWrsXpoa+zI7y2C0FbkAviY2TxeRjH1m2DyLi4sjthtuuIGO+dlnn1ld/+jRoxGbSaBjc2L74ZKV\naCuAs7mbxmTnx83SdDlf37CFEMIT5LCFEMIT5LCFEMIT5LCFEMITcio6DhcYTGUlbclXCUZ2HZP4\nVF5eHrEtWLAgYps+fbrVdQAujjKRjokqpsapLBOOlelk8zSJjqxMJxNrmOBqyqI7duyY1ZhMsDWJ\nuKx0KLsOy9Yz1XZgpVDZftpmPwJcHGUCNJsns5nOt20ka5one8bYPWYCtGlMW2GZYQomYP7GNqDA\nJUCCYdv0OpudoW/YQgjhCXLYQgjhCXLYQgjhCXLYQgjhCXLYQgjhCXlNTY8b5WFbC9el3jBTmFm0\nwY033kjHvP322yO2GTNmRGyssa5JibatM+2SSs2iGth1WFSBqa4xi6Bg62Q1pU0RQ6zWsm16tSla\ngM2JRYmwqAZTOQTb585Up5rBom5Yje+4DXMZLpESttFJU6dOjdhM82R12G0bOpvuO5unbTSKyYfE\nqX2tKBEhhLiCkMMWQghPGNFhr1u3DolEAlVVVRnbk08+idLSUiSTSSSTSTQ1NeV0kkIIISwc9tq1\nayMOOQgCbNy4Ea2trWhtbcWPfvSjnE1QCCHEeUYUHe+88050dHRE7DYvyocf49LA0vZFPDvOpbEv\nEyaKiooiturqajomO5aJgbZNhQEu4LAxmfBmEh2ZkMruB0uPZmIYwEVHW6HIZe1MKGJzN9URZ3Z2\nP1zSiVlqPUv5ZmOaBFd2vm1deNOYbO22ademGtu2deXZ+aZ52jaZtr1vpmNtn0+XBtcMFyExL6Lj\n888/j+rqaqxfv56q1kIIIcaWUTnsDRs2oL29HalUCiUlJXjsscfGel5CCCGGMSqHXVRUhCAIEAQB\nHnnkEbS0tIz1vIQQ4oogDEOk0+nMTzZG5bC7u7sz/71ly5aLIkiEEELYEwQBCgoKMj/ZGFF0fPjh\nh7Fz50709fVh1qxZ+N3vfocdO3YglUohCAKUl5fjxRdftJoYEyHivtxn57sImexY1nzUlOnIxBJb\nscMkwLDsOCY0MZtLbWBT7ew4sOu7NI1l4ioTN9l1WKYiYK4VbXMdl6xZllXIhETTvjPB2LYpsam+\nuO19N90Phq0Qymxx62Ez8dz0ebddk0vta1sfkqva/SM67M2bN0ds69aty8lkhBBCmFGmoxBCeIIc\nthBCeIIcthBCeEJey6sy8lnGkMHEK9ZQlJUnBbiow+bOBBBWdtQ0J1uB0WXtbO7sOqYxWRYdE+56\ne3ut58SEN9uMTpO4yIQ7dj+ZGOdSXpVlhLIGxqYxmfjFREu2x6ZmwewZY3NnZU9ZCVrAvvkyu7Yp\nG9U2e5Od7yIMxxUD45w/FkKkvmELIYQnyGELIYQnyGELIYQnyGELIYQnyGELIYQn5DVKxCWN/FLC\n1GlT+ipTfm1r7pqiBZidpeS6pKazaANW+5pFX/T19dExWWQBO58dx+qIA0BJSUnExvbTNpIGACZP\nnhyxsQiGKVOmRGym8gFs71hUBXtuTGnkDNvUcpbSb7o+a3TM7jE7DuD3g6Xls8+26XPExmQ2tnem\nZ97Wt7j4INvILJcx1YRXCCEuQ+SwhRDCE+SwhRDCE+SwhRDCE/IqOro00IzzIt8kQthe/9ChQ1Y2\ngAta1113XcTGxDjTPJnYEldAOXz4cMTW1tZmdb5JdPz8888jNiZ+MYFv0aJFdEwmqDHBlO3ntGnT\n6Jhsn1nNcZbebRLJmJ0Jw998803ExgRLgAt3tgIfuw7ARTqWwn/06FGr4wAuhNoK7exeAvzzbtt8\nOVcCHyNOkIStr8u2l/qGLYQQniCHLYQQnpDVYXd2duLuu+9GZWUlFixYgOeeew7A+fjM+vp6zJkz\nB0uXLkV/f39eJiuEEFcyWR12YWEhnn32WXz00Ud499138cILL2D//v1oaGhAfX092traUFdXh4aG\nhnzNVwghrliyio7FxcUoLi4GcL5+8Pz589HV1YWtW7di586dAIA1a9bgrrvuGnOnHUdkM51rO+bB\ngwcjtg8//JAee9NNN1nZmKBlEmBsRR3bRrAAFxj3798fsV2430MxCTXMzkRDltVYWlpKx2Sipa1Q\nZNpPJrwxQcu2ATAQT3wyZU8ycZTZ2D02Nfa1FbBtG94C9k2VbRvrAvaCq8u+56IedpyGu2OR1W39\nDrujowOtra1YvHgxent7kUgkAACJRMKpQL0QQojRYeWwjx8/jlWrVmHTpk2RkLUgCL6V9UCEEOJy\nY0SHffbsWaxatQqrV6/GihUrAJz/Vt3T0wMA6O7uNhbyEUIIkZ10Oo1z585lfrKR1WGHYYj169ej\noqICjz76aMa+fPlyNDY2AgAaGxszjlwIIYQbBQUFGD9+fOYnG1l/u3v3brzyyitYuHAhkskkAODp\np5/G448/jgcffBAvv/wyysrK8Nprr1lNLG6WEiNuaUTbhrmmUpPsfFvxyCQUsWPZjWSCmOmGL1y4\nMGJjWZosA9A0z5tvvjliY+LR9OnTIzbWnBYAZsyYEbHZiptfffUVHZOVPWWhqC7NbW2fEfbMm0qh\nmq41HJblyUrDmrAVzkyZuOweMxtr1mv6BmmbZWrKPGWwcru2JV9dhPY4x7mS1WHfcccdxg3avn17\nTiYkhBCCo0xHIYTwBDlsIYTwBDlsIYTwBDlsIYTwhLzWw86FchonVRTgURUsUoI1hwV4VENZWVnE\nxiIDTNEXTN1mNtvIEYDXimaRBbZNRk3ns3myaBTTmKz29g033BCxMTHcFNXAjmVRIrbp0SZYTWrb\neukA3xPb6A0WOWKyX3vttREbe+ZZ82SA7wmr8c1KBbjUF2e4lKJgz2JcH2R7j1yuoya8QghxGSKH\nLYQQniCHLYQQniCHLYQQnpBX0dEljfxSVgCcNGlSxMbSsAFkyswOhYk6tinscTGJN7Zpzy6CFtsn\nJki5NCBmYiAT85igZBIIbY9lqdRnzpyhY9o+n6wBselcF4FzOKY66CwVnN23uXPnRmymxtNffPFF\nxMbS/5kobZqnbdcqF4HPRaC0xVYgzMVnG9A3bCGE8AY5bCGE8AQ5bCGE8AQ5bCGE8IS8io6MuA1z\nXWCiDsvGYsKISbQbqeD4BVxECNs6vi7zZGIkE9SmTp0asbEsONOYTLhj62FiHHC+Hd1w2DptMz8B\nvifsWTh58mTEZhJcWbYiy1xlGaZs7ibYHrNn1iRYsrWzOdnaAN7UuLOzM2Jjz4JLpiMTTF1qV7Pn\nJld1quOgTEchhLgMkcMWQghPyOqwOzs7cffdd6OyshILFizAc889BwB48sknUVpaimQyiWQyiaam\nprxMVgghrmSyvoAtLCzEs88+i5qaGhw/fhy33nor6uvrEQQBNm7ciI0bN+ZrnkIIccWT1WEXFxej\nuLgYwPnsvfnz56OrqwvA6F7eu2Qe2YqOccsYMrHG1mayMzGPiVemrC92PhMtXcpHMkHONrPOdC+Y\nKGTbgHjChAl0TGa3LTFqykpkWXjMxtbJGggDvIkw22N2303CW5zsONN9ty3Vy/aDCcAA32f2LLDz\nTWtnz2IuGt6yexw3a9a2ubhLRqYJ6yeko6MDra2tuP322wEAzz//PKqrq7F+/XrrtFIhhBCjx8ph\nHz9+HA888AA2bdqEa6+9Fhs2bEB7eztSqRRKSkrw2GOP5XqeQghxxTOiwz579ixWrVqFn//851ix\nYgUAoKioCEEQIAgCPPLII2hpacn5RIUQ4nIknU5jcHAw85ONrA47DEOsX78eFRUVePTRRzP27u7u\nzH9v2bIFVVVVMacshBBXJgUFBRg3blzmJxtZRcfdu3fjlVdewcKFC5FMJgEATz31FDZv3oxUKoUg\nCFBeXo4XX3xx7GYvhBCCktVh33HHHVQBve+++0Z1sXzVhHa5PoOp21999RU9lkUBMJtLOjKbJ/tT\niaV3m/6kYvvM6iIzTMo+i3Jh17FV5gEeJcIiGGzrWQNAb29vxMbSyFmzX2YzzdOUwj8cl4gjNiZL\nDXeJlGANcz/55JOI7eDBg/T8o0ePRmy2pQJMkTxxIkJcIs3Y82kb5WGy20Yxxal3nhk39ghCCCHy\nghy2EEJ4ghy2EEJ4ghy2EEJ4Ql7rYdumV5uwbcBpEgzY9Zn4xMQSU0NSJjDa1pRm1wb4nrDUX1a/\n2TQmE8mYaOgijLBr2ab1M+HMNCcmkjFh2DTm1VdfHbGxRslMYGR1rwF7kY3VrjYJ7ewesT1mNpPg\nevjw4Yjts88+i9g+/vjjiG1o+O5QmAhsm64eNzXdRSC09Q1xm/XaCpEmYVX1sIUQ4jJEDlsIITxB\nDlsIITxBDlsIITwhp6Lj8JfxLjViXZpt2mJbD5uJJX19fXRMJgAxIZJdmwlfpmNtj2MCHcDFLybc\nMRsT7QD72tUsI5MJVwDw9ddfR2zffPNNxOaSIctqV7M5mRoDM2ybIruIT+xZMol0tmMeOXIkYmtv\nb4/Yenp6IjbTM28rJrLjmA3ITT3sXGRPMmwzk8eiAbC+YQshhCfIYQshhCfIYQshhCfIYQshhCfk\nNdORYSuq5PP6LGPOJMDYikqslCnLiDSdX1hYGLGxubPjTHZbEdbUszNOqUmXDFcmVJkyOhku5V1t\nzgV4ViObJxN7XYR29iyyUr8dHR10TJbVyERHlhFpErAZcTIVTcfa4pLpaIvpXNvywXGzJ03oG7YQ\nQniCHLYQQnhCVod9+vRpLF68GDU1NaioqMATTzwB4HxsZ319PebMmYOlS5ca/2QWQggxdmR12BMn\nTkRzczNSqRT27t2L5uZm7Nq1Cw0NDaivr0dbWxvq6urQ0NCQr/kKIcQVy4ivRC5kuQ0MDGBwcBBT\npkzB1q1bsWbNGgDAmjVr8M9//jO3sxRCCDFylEg6ncaiRYvw6aefYsOGDaisrERvby8SiQQAIJFI\n0CanQFRpjZtuHjeF1DYqgqUYM7Uf4Co+qw3Mai3PmDGDjjlt2jSrOcWtBW4bPcEiIgAeecJqdLN5\nmpqxsrR+lsLv0uTU9li2x6aoGxYRws5nafVsjQCvac0a3r733ntW1wF4Tesvv/wyYmMlCUzp/yxC\nxyUNnRE3DZ3B/A2zsXWa5p6LebowosMuKChAKpXCsWPHcO+996K5ufmi3wdBMCbhKkIIIbJjHSUy\nefJk3H///fjggw+QSCQyxWK6u7tRVFSUswkKIcTlTDqdxuDgYOYnG1kddl9fXyYC5NSpU9i2bRuS\nySSWL1+OxsZGAEBjYyNWrFgxRlMXQogri4KCAowbNy7zk42sr0S6u7uxZs0apNNppNNprF69GnV1\ndUgmk3jwwQfx8ssvo6ysDK+99tqYLkAIIUSUIMzR23H2Xtsk4DDiCowMJi4wQY2JQiYBhomRrK7y\n/PnzI7Z77rmHjjlz5syIjYklbO6sRjXAhSJ2P2ybEpvOZ/WsmYDjkvrL9tM0JwZbE7u+S6Nk2wbG\nTEg01d1m9/j999+P2Pbt2xexmUTxzz//PGI7cOCA1TxNe2xbC9y2MS+Qm3rY7Nuqbc1z037alqKI\nE2Bx9uxZ8+fDagQhhBCXHDlsIYTwBDlsIYTwBDlsIYTwhLzWw45b+zoX+igTSxgmsYIJGyz7kdUg\nbmtrsx6TNZJlAp1pPUxAst1PU6gRE4rYnGwzzkxj2tbyNona7Fi2H0xoMglvbJ3sOi4CNqu5zgRG\nVuPaBJsTy7pl62TZjwDfJ/bc2TbrzRXsWmye7LPtIjYzcdWl+bKLX9M3bCGE8AQ5bCGE8AQ5bCGE\n8AQ5bCGE8IS8io5xS6HGxVYIYKKKSShiggVbJys7umfPHjomyxasqKiI2FjRLdM8bQUYJtyZCtKw\nvbPNamTCrAk2JxcRlQlIF+q8D4XtkUlsZkITK3HKjjMJb3v37o3YmMDIns8L5Y6HwzJf2XPDGu6y\nZr8A8MUXX0Rsx44di9hcsv1sP+8uzW2Znc2JPd/s+QD488CyRNl+mrI8XYRYfcMWQghPkMMWQghP\nkMMWQghPkMMWQghPkMMWQghPyGmUyHCVlin7JoU0jmoct7Evw5TybVu/mSnRrEEqwBuvsrrGN954\nY8Q2depUOiZLR54yZYrVmMXFxXTMOFE3pia8LKLDtq6xKZrF9hlh1zE9nywKgN23I0eORGyffPIJ\nHfODDz6I2Fh6eElJScRmauzL9tMUSTQcU6o/+xyzvTOldzNsyxe4pHzbRpSw6A12fwH7JtHsOqbI\nqOHPbbZyGfqGLYQQnpDVYZ8+fRqLFy9GTU0NKioq8MQTTwAAnnzySZSWliKZTCKZTKKpqSkvkxVC\niCuZrK9EJk6ciObmZlx99dU4d+4c7rjjDuzatQtBEGDjxo3YuHFjvuYphBBXPCO+ErmQ8TMwMIDB\nwcHMe88ctYIUQghhYETRMZ1OY9GiRfj000+xYcMGVFZW4h//+Aeef/55/PWvf0VtbS3+9Kc/0XrN\nw7EVOy5cdzi2goNL+qutIOWSVspsTKwwiTJMnOjo6IjYDh48aHUdk33atGkR2w9+8IOIzTRPJkqx\n9OxDhw7R8xmTJ0+2Oo4JMyZRhwlI1113XcTG0pFN4ijb+//9738RW1dXV8TG7iXA07vZs9jf3x+x\nsdIHAH9uWekDdm2Wcg3wPWGfbdvUcCA3TbfjNPE1PUtsP9lzY/ocMlzKNIzoQQsKCpBKpXDw4EG8\n9dZb2LFjBzZs2ID29nakUimUlJTgscces76gEEKI0WH9lXfy5Mm4//778f7776OoqAhBECAIAjzy\nyCNoaWmh54RheNGPEEKIizl37hzOnDmT+clGVofd19eX+fPr1KlT2LZtG5LJJHp6ejLHbNmyBVVV\nVfT8C079wo8QQoiLGT9+PCZMmJD5yXpstl92d3djzZo1SKfTSKfTWL16Nerq6vCLX/wCqVQKQRCg\nvLwcL7744pguQAghRJQgzNG7iiAIIkIEE5RcxDzbWsu5aPTp1CiTCDAsO8ylsS/LFmQ2E0wgZMLI\nLbfcErHdeuutdMzq6uqIjX1DYCKZ6b6zfWK1ntk9ZtcBeAZkaWlpxMbm3tvbS8d89913I7bW1taI\njd0j0/PcFWF4AAAHFElEQVTJxCdmc2m+zPaZzcn2swVwMZDtsUtAgClLdTi2zY9NsOeL2VzqyrP9\nYOebmlkPt3d3dxvXpExHIYTwBDlsIYTwBDlsIYTwBDlsIYTwhLyWV43zIh7gwpmtgAK4NQUd7XEA\nF1BsRVSAiyC2YZEuWZ5M0Dpw4EDEZsqiYxmEM2bMiNgmTZoUsbGyoQBw+PDhiG3fvn0RGxMYTRlj\ntkIme5ZMWYlffvml1fkuIhnLNrQVlseiwetwXARCW9HQ9Vo2x+VCHHUJQ7bNdmbPIWAuY8vQN2wh\nhPAEOWwhhPAEOWwhhPAEOWwhhPAEOWwhhPCEnKamD1daWeqvqdgJszOVlVW3MlW8Ykp63EafLsfa\nEifV1ZRSa7v3LPLDpJhfc801ERurM21b4xrgKdYsUoLVdDZFSrA9GanIzkiYopuGw6JhTJEbrP50\ntoasQ8nn8xlnzEtdBM62Jr9pnrZRPy5RM8OfpdOnTys1XQghfEcOWwghPEEOWwghPEEOWwghPCGn\nouNwXEQydiwTHW3FHyBeE1+XNF3b2sImYSOOMGM617YeNpu7Scxja2L3yKUBMbuftiKs6Vli82Ri\nHlun6fmyFWxZU2KX/bRNLXdpbpsLIZJh2+DadCwjF4J+3PNtx7Q9bmBgQKLjWBK3bsK3EdtoBF8w\n1T/xlVw05bjUXG59XvNxj+SwR4Ec9rcf9m3XZ+Swv/3kYz1y2EII4Qly2EII4Qk5Ex3vuusu7Ny5\nMxdDCyHEZcsPf/hD7Nixg/4uZw5bCCHE2KJXIkII4Qly2EII4Qk5d9hNTU2YN28eZs+ejWeeeSbX\nlxtz1q1bh0QigaqqqoztyJEjqK+vx5w5c7B06VLaX/DbSmdnJ+6++25UVlZiwYIFeO655wD4u6bT\np09j8eLFqKmpQUVFBZ544gkA/q5nKIODg0gmk1i2bBkAv9dUVlaGhQsXIplM4nvf+x4Av9fT39+P\nBx54APPnz0dFRQXee++9vKwnpw57cHAQv/71r9HU1IR9+/Zh8+bN2L9/fy4vOeasXbsWTU1NF9ka\nGhpQX1+PtrY21NXVoaGh4RLNzp3CwkI8++yz+Oijj/Duu+/ihRdewP79+71d08SJE9Hc3IxUKoW9\ne/eiubkZu3bt8nY9Q9m0aRMqKioyGXI+rykIAuzYsQOtra1oaWkB4Pd6fvOb3+DHP/4x9u/fj717\n92LevHn5WU+YQ955553w3nvvzfz76aefDp9++ulcXjIntLe3hwsWLMj8e+7cuWFPT08YhmHY3d0d\nzp0791JNLTY//elPw23btl0Wazpx4kRYW1sbfvjhh96vp7OzM6yrqwvffPPN8Cc/+UkYhn4/d2Vl\nZWFfX99FNl/X09/fH5aXl0fs+VhPTr9hd3V1YdasWZl/l5aWoqurK5eXzAu9vb1IJBIAgEQigd7e\n3ks8o9HR0dGB1tZWLF682Os1pdNp1NTUIJFIZF73+LweAPjtb3+LP/7xjxfVR/F5TUEQYMmSJait\nrcVLL70EwN/1tLe3Y/r06Vi7di0WLVqEX/7ylzhx4kRe1pNTh32pu0vkA9ZZxweOHz+OVatWYdOm\nTZEuMb6tqaCgAKlUCgcPHsRbb72F5ubmi37v23r+/e9/o6ioCMlkMmuxJJ/WtHv3brS2tuL111/H\nCy+8gLfffvui3/u0nnPnzmHPnj341a9+hT179uCaa66JvP7I1Xpy6rBnzpyJzs7OzL87OztRWlqa\ny0vmhUQigZ6eHgBAd3c3ioqKLvGM3Dh79ixWrVqF1atXY8WKFQD8XxNwvhXZ/fffjw8++MDr9bzz\nzjvYunUrysvL8fDDD+PNN9/E6tWrvV5TSUkJAGD69OlYuXIlWlpavF1PaWkpSktLcdtttwEAHnjg\nAezZswfFxcU5X09OHXZtbS0OHDiAjo4ODAwM4NVXX8Xy5ctzecm8sHz5cjQ2NgIAGhsbM07PB8Iw\nxPr161FRUYFHH300Y/d1TX19fRk1/tSpU9i2bRuSyaS36wGAp556Cp2dnWhvb8ff//533HPPPfjb\n3/7m7ZpOnjyZKTN74sQJvPHGG6iqqvJ2PcXFxZg1axba2toAANu3b0dlZSWWLVuW+/WM+VvxYfzn\nP/8J58yZE373u98Nn3rqqVxfbsx56KGHwpKSkrCwsDAsLS0N//KXv4SHDx8O6+rqwtmzZ4f19fXh\n0aNHL/U0rXn77bfDIAjC6urqsKamJqypqQlff/11b9e0d+/eMJlMhtXV1WFVVVX4hz/8IQzD0Nv1\nDGfHjh3hsmXLwjD0d02fffZZWF1dHVZXV4eVlZUZP+DresIwDFOpVFhbWxsuXLgwXLlyZdjf35+X\n9Sg1XQghPEGZjkII4Qly2EII4Qly2EII4Qly2EII4Qly2EII4Qly2EII4Qly2EII4Qly2EII4Qn/\nB3WNC0bJYF2yAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x10248cc90>"
       ]
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Why is this image blurry?"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "So what was the subject doing here?\n",
      "\n",
      "https://openfmri.org/dataset/ds000105\n",
      "\n",
      "* 001 House\n",
      "* 002 Scrambled\n",
      "* 003 Cat\n",
      "* 004 Shoe\n",
      "* 005 Bottle\n",
      "* 006 Scissors\n",
      "* 007 Chair\n",
      "* 008 Face"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "event_types = ['house', 'scrambled', 'cat', 'shoe', 'bottle', 'scissors', 'chair', 'face']"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 9
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "You should have the task lists for this subject in your `~/data/ds105` directory.  Let's check:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import os\n",
      "HOME = os.path.expanduser('~')\n",
      "DATA_PATH = os.path.join(HOME, 'data', 'ds105')\n",
      "SUB1_TASK_PATH = os.path.join(DATA_PATH, 'sub001', 'model', 'model001', 'onsets', 'task001_run001')\n",
      "os.path.isdir(SUB1_TASK_PATH)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 10,
       "text": [
        "True"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "SUB1_TASK_PATH"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 16,
       "text": [
        "'/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001'"
       ]
      }
     ],
     "prompt_number": 16
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "open?"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 17
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "If you got \"False\" above - you can download the zip file of the task files you need from:\n",
      "\n",
      "http://nipy.bic.berkeley.edu/practical_neuroimaging/ds105_sub001_tasks.zip\n",
      "\n",
      "Download the file, then `cd ~/data` and `unzip ~/Downloads/ds105_sub001_tasks.zip` (if you downloaded the file to your `Downloads` folder.\n",
      "\n",
      "Recheck:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "os.path.isdir(SUB1_TASK_PATH)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 11,
       "text": [
        "True"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Look in this path:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "os.listdir(SUB1_TASK_PATH)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 12,
       "text": [
        "['cond001.txt',\n",
        " 'cond002.txt',\n",
        " 'cond003.txt',\n",
        " 'cond004.txt',\n",
        " 'cond005.txt',\n",
        " 'cond006.txt',\n",
        " 'cond007.txt',\n",
        " 'cond008.txt']"
       ]
      }
     ],
     "prompt_number": 12
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "What are in these files?"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "fobj = open(COND1_FNAME, 'rt')\n",
      "lines = fobj.readlines()\n",
      "for line in open(COND1_FNAME, 'rt'):\n",
      "    pass\n",
      "file.__iter__"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 42,
       "text": [
        "['156.000\\t0.5\\t1\\n',\n",
        " '158.000\\t0.5\\t1\\n',\n",
        " '160.000\\t0.5\\t1\\n',\n",
        " '162.000\\t0.5\\t1\\n',\n",
        " '164.000\\t0.5\\t1\\n',\n",
        " '166.000\\t0.5\\t1\\n',\n",
        " '168.000\\t0.5\\t1\\n',\n",
        " '170.000\\t0.5\\t1\\n',\n",
        " '172.000\\t0.5\\t1\\n',\n",
        " '174.000\\t0.5\\t1\\n',\n",
        " '176.000\\t0.5\\t1\\n',\n",
        " '178.000\\t0.5\\t1\\n']"
       ]
      }
     ],
     "prompt_number": 42
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "COND1_FNAME = os.path.join(SUB1_TASK_PATH, 'cond001.txt')\n",
      "print(open(COND1_FNAME, 'rt').read())"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "156.000\t0.5\t1\n",
        "158.000\t0.5\t1\n",
        "160.000\t0.5\t1\n",
        "162.000\t0.5\t1\n",
        "164.000\t0.5\t1\n",
        "166.000\t0.5\t1\n",
        "168.000\t0.5\t1\n",
        "170.000\t0.5\t1\n",
        "172.000\t0.5\t1\n",
        "174.000\t0.5\t1\n",
        "176.000\t0.5\t1\n",
        "178.000\t0.5\t1\n",
        "\n"
       ]
      }
     ],
     "prompt_number": 13
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This is a common format for specifying events - any guesses what the columns are?\n"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "I want to read this file to give me something I can use to specify my design.  Something like this:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def evdefs_from_fname(fname):\n",
      "    \"\"\" Get onset, duration, amplitude for event from file\n",
      "    \n",
      "    Parameters\n",
      "    ----------\n",
      "    fname : str\n",
      "        Filename to read data from\n",
      "        \n",
      "    Returns\n",
      "    -------\n",
      "    evdefs : list of lists\n",
      "        list, one element per event, where the elements are lists containing\n",
      "        (onset, duration, amplitude)\n",
      "    \"\"\"\n",
      "    evdefs = []\n",
      "    for line in open(fname, 'rt'):\n",
      "        values = line.split()\n",
      "        float_values = [float(v) for v in values]\n",
      "        evdefs.append(float_values)\n",
      "    return evdefs"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 43
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Of course we need to test\n",
      "evdefs = evdefs_from_fname(COND1_FNAME)\n",
      "assert len(evdefs) == 12\n",
      "assert evdefs[0] == [156.000, 0.5, 1]\n",
      "assert evdefs[-1] == [178.000, 0.5, 1]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 44
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can actually do something neater than that:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "dtype_def = [('onset', 'f'), ('duration', 'f'), ('amplitude', 'f')]\n",
      "evs = np.loadtxt(COND1_FNAME, dtype_def)\n",
      "evs"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 52,
       "text": [
        "array([(156.0, 0.5, 1.0), (158.0, 0.5, 1.0), (160.0, 0.5, 1.0),\n",
        "       (162.0, 0.5, 1.0), (164.0, 0.5, 1.0), (166.0, 0.5, 1.0),\n",
        "       (168.0, 0.5, 1.0), (170.0, 0.5, 1.0), (172.0, 0.5, 1.0),\n",
        "       (174.0, 0.5, 1.0), (176.0, 0.5, 1.0), (178.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')])"
       ]
      }
     ],
     "prompt_number": 52
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "evdefs[2][2]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 53,
       "text": [
        "1.0"
       ]
      }
     ],
     "prompt_number": 53
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This is nice because:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "evs[0]['onset']"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 54,
       "text": [
        "156.0"
       ]
      }
     ],
     "prompt_number": 54
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "evs['onset']"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 55,
       "text": [
        "array([ 156.,  158.,  160.,  162.,  164.,  166.,  168.,  170.,  172.,\n",
        "        174.,  176.,  178.], dtype=float32)"
       ]
      }
     ],
     "prompt_number": 55
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can load our events for all types:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "event_defs = {}\n",
      "for eno in range(8):\n",
      "    fname = os.path.join(os.path.join(SUB1_TASK_PATH, 'cond%03d.txt' % (eno + 1)))\n",
      "    print(fname)\n",
      "    event_name = event_types[eno]\n",
      "    event_defs[event_name] = np.loadtxt(fname, dtype_def)\n",
      "event_defs          "
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond001.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond002.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond003.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond004.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond005.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond006.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond007.txt\n",
        "/Users/mb312/data/ds105/sub001/model/model001/onsets/task001_run001/cond008.txt\n"
       ]
      },
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 56,
       "text": [
        "{'bottle': array([(228.0, 0.5, 1.0), (230.0, 0.5, 1.0), (232.0, 0.5, 1.0),\n",
        "       (234.0, 0.5, 1.0), (236.0, 0.5, 1.0), (238.0, 0.5, 1.0),\n",
        "       (240.0, 0.5, 1.0), (242.0, 0.5, 1.0), (244.0, 0.5, 1.0),\n",
        "       (246.0, 0.5, 1.0), (248.0, 0.5, 1.0), (250.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'cat': array([(84.0, 0.5, 1.0), (86.0, 0.5, 1.0), (88.0, 0.5, 1.0),\n",
        "       (90.0, 0.5, 1.0), (92.0, 0.5, 1.0), (94.0, 0.5, 1.0),\n",
        "       (96.0, 0.5, 1.0), (98.0, 0.5, 1.0), (100.0, 0.5, 1.0),\n",
        "       (102.0, 0.5, 1.0), (104.0, 0.5, 1.0), (106.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'chair': array([(264.0, 0.5, 1.0), (266.0, 0.5, 1.0), (268.0, 0.5, 1.0),\n",
        "       (270.0, 0.5, 1.0), (272.0, 0.5, 1.0), (274.0, 0.5, 1.0),\n",
        "       (276.0, 0.5, 1.0), (278.0, 0.5, 1.0), (280.0, 0.5, 1.0),\n",
        "       (282.0, 0.5, 1.0), (284.0, 0.5, 1.0), (286.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'face': array([(48.0, 0.5, 1.0), (50.0, 0.5, 1.0), (52.0, 0.5, 1.0),\n",
        "       (54.0, 0.5, 1.0), (56.0, 0.5, 1.0), (58.0, 0.5, 1.0),\n",
        "       (60.0, 0.5, 1.0), (62.0, 0.5, 1.0), (64.0, 0.5, 1.0),\n",
        "       (66.0, 0.5, 1.0), (68.0, 0.5, 1.0), (70.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'house': array([(156.0, 0.5, 1.0), (158.0, 0.5, 1.0), (160.0, 0.5, 1.0),\n",
        "       (162.0, 0.5, 1.0), (164.0, 0.5, 1.0), (166.0, 0.5, 1.0),\n",
        "       (168.0, 0.5, 1.0), (170.0, 0.5, 1.0), (172.0, 0.5, 1.0),\n",
        "       (174.0, 0.5, 1.0), (176.0, 0.5, 1.0), (178.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'scissors': array([(12.0, 0.5, 1.0), (14.0, 0.5, 1.0), (16.0, 0.5, 1.0),\n",
        "       (18.0, 0.5, 1.0), (20.0, 0.5, 1.0), (22.0, 0.5, 1.0),\n",
        "       (24.0, 0.5, 1.0), (26.0, 0.5, 1.0), (28.0, 0.5, 1.0),\n",
        "       (30.0, 0.5, 1.0), (32.0, 0.5, 1.0), (34.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'scrambled': array([(192.0, 0.5, 1.0), (194.0, 0.5, 1.0), (196.0, 0.5, 1.0),\n",
        "       (198.0, 0.5, 1.0), (200.0, 0.5, 1.0), (202.0, 0.5, 1.0),\n",
        "       (204.0, 0.5, 1.0), (206.0, 0.5, 1.0), (208.0, 0.5, 1.0),\n",
        "       (210.0, 0.5, 1.0), (212.0, 0.5, 1.0), (214.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')]),\n",
        " 'shoe': array([(120.0, 0.5, 1.0), (122.0, 0.5, 1.0), (124.0, 0.5, 1.0),\n",
        "       (126.0, 0.5, 1.0), (128.0, 0.5, 1.0), (130.0, 0.5, 1.0),\n",
        "       (132.0, 0.5, 1.0), (134.0, 0.5, 1.0), (136.0, 0.5, 1.0),\n",
        "       (138.0, 0.5, 1.0), (140.0, 0.5, 1.0), (142.0, 0.5, 1.0)], \n",
        "      dtype=[('onset', '<f4'), ('duration', '<f4'), ('amplitude', '<f4')])}"
       ]
      }
     ],
     "prompt_number": 56
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Are these onsets in TRs or in seconds?"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "OK - how to make our design?\n",
      "\n",
      "We could do this by hand (right?) - but let's let `nipy` do this for us:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import nipy.modalities.fmri.design as fmrid"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 57
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Try `fmrid.block_design?` at the prompt."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "fmrid.block_design?"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 59
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "So - we need start and stop for our events:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "nipy_dtype = [('start', 'f'), ('end', 'f'), ('type', 'S10')]\n",
      "n_events = sum([len(evs) for e in event_defs])\n",
      "nipy_block_specs = np.zeros((n_events,), dtype=nipy_dtype)\n",
      "event_no = 0\n",
      "for ev_type in event_defs:\n",
      "    ev_def = event_defs[ev_type]\n",
      "    for e in ev_def:\n",
      "        nipy_block_specs[event_no]['start'] = e['onset']\n",
      "        nipy_block_specs[event_no]['end'] = (e['onset'] + e['duration'])\n",
      "        nipy_block_specs[event_no]['type'] = ev_type\n",
      "        event_no += 1\n",
      "nipy_block_defs = np.sort(nipy_block_specs)\n",
      "nipy_block_specs"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 60,
       "text": [
        "array([(156.0, 156.5, 'house'), (158.0, 158.5, 'house'),\n",
        "       (160.0, 160.5, 'house'), (162.0, 162.5, 'house'),\n",
        "       (164.0, 164.5, 'house'), (166.0, 166.5, 'house'),\n",
        "       (168.0, 168.5, 'house'), (170.0, 170.5, 'house'),\n",
        "       (172.0, 172.5, 'house'), (174.0, 174.5, 'house'),\n",
        "       (176.0, 176.5, 'house'), (178.0, 178.5, 'house'),\n",
        "       (84.0, 84.5, 'cat'), (86.0, 86.5, 'cat'), (88.0, 88.5, 'cat'),\n",
        "       (90.0, 90.5, 'cat'), (92.0, 92.5, 'cat'), (94.0, 94.5, 'cat'),\n",
        "       (96.0, 96.5, 'cat'), (98.0, 98.5, 'cat'), (100.0, 100.5, 'cat'),\n",
        "       (102.0, 102.5, 'cat'), (104.0, 104.5, 'cat'), (106.0, 106.5, 'cat'),\n",
        "       (192.0, 192.5, 'scrambled'), (194.0, 194.5, 'scrambled'),\n",
        "       (196.0, 196.5, 'scrambled'), (198.0, 198.5, 'scrambled'),\n",
        "       (200.0, 200.5, 'scrambled'), (202.0, 202.5, 'scrambled'),\n",
        "       (204.0, 204.5, 'scrambled'), (206.0, 206.5, 'scrambled'),\n",
        "       (208.0, 208.5, 'scrambled'), (210.0, 210.5, 'scrambled'),\n",
        "       (212.0, 212.5, 'scrambled'), (214.0, 214.5, 'scrambled'),\n",
        "       (120.0, 120.5, 'shoe'), (122.0, 122.5, 'shoe'),\n",
        "       (124.0, 124.5, 'shoe'), (126.0, 126.5, 'shoe'),\n",
        "       (128.0, 128.5, 'shoe'), (130.0, 130.5, 'shoe'),\n",
        "       (132.0, 132.5, 'shoe'), (134.0, 134.5, 'shoe'),\n",
        "       (136.0, 136.5, 'shoe'), (138.0, 138.5, 'shoe'),\n",
        "       (140.0, 140.5, 'shoe'), (142.0, 142.5, 'shoe'),\n",
        "       (228.0, 228.5, 'bottle'), (230.0, 230.5, 'bottle'),\n",
        "       (232.0, 232.5, 'bottle'), (234.0, 234.5, 'bottle'),\n",
        "       (236.0, 236.5, 'bottle'), (238.0, 238.5, 'bottle'),\n",
        "       (240.0, 240.5, 'bottle'), (242.0, 242.5, 'bottle'),\n",
        "       (244.0, 244.5, 'bottle'), (246.0, 246.5, 'bottle'),\n",
        "       (248.0, 248.5, 'bottle'), (250.0, 250.5, 'bottle'),\n",
        "       (12.0, 12.5, 'scissors'), (14.0, 14.5, 'scissors'),\n",
        "       (16.0, 16.5, 'scissors'), (18.0, 18.5, 'scissors'),\n",
        "       (20.0, 20.5, 'scissors'), (22.0, 22.5, 'scissors'),\n",
        "       (24.0, 24.5, 'scissors'), (26.0, 26.5, 'scissors'),\n",
        "       (28.0, 28.5, 'scissors'), (30.0, 30.5, 'scissors'),\n",
        "       (32.0, 32.5, 'scissors'), (34.0, 34.5, 'scissors'),\n",
        "       (264.0, 264.5, 'chair'), (266.0, 266.5, 'chair'),\n",
        "       (268.0, 268.5, 'chair'), (270.0, 270.5, 'chair'),\n",
        "       (272.0, 272.5, 'chair'), (274.0, 274.5, 'chair'),\n",
        "       (276.0, 276.5, 'chair'), (278.0, 278.5, 'chair'),\n",
        "       (280.0, 280.5, 'chair'), (282.0, 282.5, 'chair'),\n",
        "       (284.0, 284.5, 'chair'), (286.0, 286.5, 'chair'),\n",
        "       (48.0, 48.5, 'face'), (50.0, 50.5, 'face'), (52.0, 52.5, 'face'),\n",
        "       (54.0, 54.5, 'face'), (56.0, 56.5, 'face'), (58.0, 58.5, 'face'),\n",
        "       (60.0, 60.5, 'face'), (62.0, 62.5, 'face'), (64.0, 64.5, 'face'),\n",
        "       (66.0, 66.5, 'face'), (68.0, 68.5, 'face'), (70.0, 70.5, 'face')], \n",
        "      dtype=[('start', '<f4'), ('end', '<f4'), ('type', 'S10')])"
       ]
      }
     ],
     "prompt_number": 60
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Nice - so let's make the design:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "n_trs = img.shape[-1]\n",
      "TR = 2.5\n",
      "t = np.arange(n_trs) * TR\n",
      "X_exper, cons_exper = fmrid.block_design(nipy_block_defs, t, level_contrasts=True)\n",
      "plt.figure(figsize=(10, 20))\n",
      "plt.imshow(X_exper)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 61,
       "text": [
        "<matplotlib.image.AxesImage at 0x109b47790>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAG0AAAR4CAYAAAAMrFhtAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHudJREFUeJzt3X9M1Pfhx/HXMcB20gI6uVOOBqYinD8A54/MzVFmz2Xd\nJLRa1G5KtdU/7B9rs2xx/aPLtqzcZpat3WqyNlbJmqgkTZQ1yjomWrfqTIO1S22LWe9Wfk/lhwqt\nCLy/f3S9bynccRwf+NzrfD0SEuXe9/m89cn7vPfnTnAYYwyESoLdE5DxUzRCikZI0QgpGiFFY2Qs\ndvz4cbNgwQIzb9484/P5RtyemppqAOhjjI/U1NSQf8cOK/dpg4ODWLBgAerq6pCZmYnly5fj4MGD\nyM/PD45xOBwoKSkZdj+/34+cnJywx45kjJXHmsoxo42rr69HqDSWPjyeO3cO8+bNQ3Z2NpKSkrBp\n0yYcPXrUylMILI7W0tKCrKys4O/dbjdaWlqsPIUASLTyYA6HI6Jxfr8/+Ou0tDSkpaWNeZ9IxkQ6\nLtbGAMAXvvCFYX8v4VgaLTMzE01NTcHfNzU1we12jxgXyWP856Wnp1s2LtbGAMA999wz7PeBQCDk\nWEsfHpctW4ZLly4hEAigv78fhw8fRmlpqZWnEFi80hITE/GHP/wB3/rWtzA4OIhHH3102DNHsYal\nT/kjOuEoT/llpCl7yi9TQ9EIKRohRSOkaIQUjZCiEVI0QopGSNEIKRohRSOkaIQUjZCiEVI0QopG\nSNEIKRohRSOkaIQUjZCiEVI0QopGSNEIKRohRSOkaIQUjZCiEVI0QopGSNEIKRohRSOkaIQUjZCi\nEVI0QopGSNEIKRohRSOkaIQs/daB8e6z3+o3lM9+k9LJopVGSNEIKRohRSOkaIQUjZCiEVI0QopG\nSFdExmHHjh1jjnn66acnfR5aaYQUjZCiEVI0QopGSNEIKRohRSOkzfU4rFq1yu4pANBKo6RohBSN\nkKIRUjRCikZI0QgpGiFtrsehr6/P7ikA0EqjpGiEFI2QohFSNEKKRkjRCCkaIW2ux6Gtrc3uKQDQ\nSqOkaIQUjZCiEVI0QopGSNEIKRohRSOkKyLj0NnZafcUAGilUVI0QopGSNEIKRohRSOkaIQUjZA2\n1+Nw8+ZNu6cAQCuNkqIRUjRCikZI0QgpGiFFI6RohLS5HoehoSG7pwBAK42SohFSNEKKRkjRCCka\nIUUjpGiEtLkeB2OM3VMAoJVGSdEIKRohRSOkaIQUjZCiEVI0Qtpcj4NeuZaoKRohRSOkaIQUjZCi\nEVI0QopGSNEI6YrIONTX19s9BQBaaZQUjZCiEVI0QopGSNEIKRohRSOkzfU47N+/f8wxCxYsmPR5\naKURUjRCikZI0QgpGiFFI6RohBSNkDbX49Da2jrmGG2uZVSKRkjRCCkaIUUjpGiEFI2QohHS5noc\nMjIy7J4CAK00SopGSNEIKRohRSOkaIQUjZCiEVI0QnF/RSQ5OXnMMf39/REdy+VyTXQ6lohqpTU1\nNaGkpAQLFy7EokWL8NxzzwEAOjs74fV6kZubi7Vr16K7u9vSyconooqWlJSE3/72t3jnnXdw9uxZ\nPP/883j33Xfh8/ng9XrR2NiINWvWwOfzWT1fQZTRXC4XCgsLAQApKSnIz89HS0sLampqUFFRAQCo\nqKjAkSNHrJupBE34iUggEMD58+excuVKdHR0wOl0AgCcTic6OjomPEEZaUJPRG7cuIH169fj2Wef\nxV133TXsNofDAYfDMer9/H5/8NdpaWlIT0+fyDTiQldXV8TPAaKOduvWLaxfvx5btmxBWVkZgE9W\nV3t7O1wuF9ra2kK+/pSTkxPtaeNWenr6sC/eQCAQcmxUD4/GGDz66KPweDx44okngp8vLS1FVVUV\nAKCqqioYU6zlMFH8ULC///3v+MY3voElS5YEHwIrKyuxYsUKlJeX48MPP0R2djaqq6uRlpY2/IQO\nB0pKSqyZfQSs3KddvXp1zDEzZ86M6Fhjqa+vD/nz2qJ6ePz6178e8nvU19XVRXPISbN+/foxxxw8\neDCiY33+32276DIWIUUjpGiEFI2QohFSNEKKRkjRCMX9K9cbN24cc0ykm+vExNj469JKI6RohBSN\nkKIRUjRCikZI0QgpGqHY2C1OopSUFMuOlZAQG1/jsTELGRdFI6RohBSNkKIRUjRCikZI0QgpGqG4\nvyISyffSZ6OVRkjRCCkaIUUjpGiEFI2QohFSNEJxv7m+fPmy3VOwnFYaIUUjpGiEFI2QohFSNEKK\nRkjRCMX95rq3t9fuKVhOK42QohFSNEKKRkjRCCkaIUUjpGiE4n5zPTAwYPcULKeVRkjRCCkaIUUj\npGiEFI2QohFSNEJxv7mO4mf6xTytNEKKRkjRCCkaIUUjpGiEFI2QohFSNEKKRkjRCCkaIUUjpGiE\nFI2QohFSNEKKRiju327Q1NRk9xQsp5VGSNEIKRohRSOkaIQUjZCiEVI0QnG/uX7hhRfGHJOcnDwF\nM7GOVhohRSOkaIQUjZCiEVI0QopGSNEIxf3m+syZM2OOKS4unoKZWEcrjZCiEVI0QopGSNEIKRoh\nRSOkaIQUjVDcXxG588477Z6C5bTSCCkaIUUjpGiEFI2QohFSNEKKRijuN9ezZ8+2ewqW00ojpGiE\nFI2QohFSNEKKRkjRCCkaoZjcXBcUFEQ07sKFC2OOmTlz5kSnE3O00ggpGiFFI6RohBSNkKIRUjRC\nikYoJjfXO3fujGjc448/PuaYadOmTXQ6MUcrjZCiEVI0QopGSNEIKRohRSOkaIQUjVBMXhHJz8+3\n7FgJCfH3dRl/f6LbgKIRUjRCikZI0QgpGiFFI6RohGJyc93d3W3ZsRwOh2XHihVaaYQUjZCiEVI0\nQopGSNEIKRohRSMUk5vr1tZWu6cQ07TSCCkaIUUjpGiEFI2QohFSNEKKRigmN9c9PT12TyGmTWil\nDQ4OoqioCOvWrQMAdHZ2wuv1Ijc3F2vXrrX0bQPy/yYU7dlnn4XH4wm+D8Pn88Hr9aKxsRFr1qyB\nz+ezZJIyXNTRmpubcezYMTz22GMwxgAAampqUFFRAQCoqKjAkSNHrJmlDBN1tCeffBJ79uwZ9l+J\nOjo64HQ6AQBOpxMdHR0Tn6GMENUTkVdffRUZGRkoKirCyZMnRx3jcDhCvn3N7/cHf52Wlob09PRo\nphFXurq6In4OEFW0N954AzU1NTh27Bg+/vhjXLt2DVu2bIHT6UR7eztcLhfa2tqQkZEx6v1zcnKi\nOW1cS09PH/bFGwgEQo6N6uHxmWeeQVNTE/x+Pw4dOoRvfvOb+NOf/oTS0lJUVVUBAKqqqlBWVhbN\n4WUMlmyuP30Y3L17N/76178iNzcXJ06cwO7du604vHzOhDfXxcXFKC4uBgDMmDEDdXV1E56UhBeT\nV0T6+/vtnkJM07VHQopGSNEIKRohRSOkaIQUjZCiEYrJzfWnr8/J6LTSCCkaIUUjpGiEFI2QohFS\nNEKKRkiba0JaaYQUjZCiEVI0QopGSNEIKRohRSMUk5trvS08PK00QopGSNEIKRohRSOkaIQUjZCi\nEYrJzfUrr7xi9xRimlYaIUUjpGiEFI2QohFSNEKKRkjRCCkaoZi8IlJdXR3RuKVLl07yTGKTVhoh\nRSOkaIQUjZCiEVI0QopGSNEIxeTm+vr163ZPIaZppRFSNEKKRkjRCCkaIUUjpGiEFI1QTG6uZ8+e\nbfcUYppWGiFFI6RohBSNkKIRUjRCikZI0QjF5OZ61qxZdk9h0qWmpkZ9X600QopGSNEIKRohRSOk\naIQUjZCiEVI0QjF5RSQlJcXuKUy6Rx55JOztR44cCXmbVhohRSOkaIQUjZCiEVI0QopGSNEIxeTm\nOjExJqdlqW9/+9tR31crjZCiEVI0QopGSNEIKRohRSOkaIRichfrcDjsnkJM00ojpGiEFI2QohFS\nNEKKRkjRCCkaoZjcXN8OWltbo76vVhohRSOkaIQUjZCiEVI0QopGSNEIKRohXRGxyZUrV6K+r1Ya\nIUUjpGiEFI2QohFSNEKKRkjRCGlzbZOPPvoo6vtqpRFSNEKKRkjRCCkaIUUjpGiEFI2QNtc2GRwc\njPq+WmmEFI2QohFSNEKKRkjRCCkaIUUjpM21TYwxUd9XK42QohFSNEKKRkjRCCkaIUUjpGiEFI2Q\nohFSNEKKRkjRCCkaIUUjpGiEFI2QohHS2w1scuHChajvq5VGSNEIKRohRSOkaIQUjZCiEVI0Qtpc\n2+SPf/xj1PeNeqV1d3djw4YNyM/Ph8fjwT//+U90dnbC6/UiNzcXa9euRXd3d9QTk9CijvaDH/wA\n999/P9599128/fbbyMvLg8/ng9frRWNjI9asWQOfz2flXOV/HCaK/3PT09ODoqIifPDBB8M+n5eX\nh1OnTsHpdKK9vR333nsv3nvvveEndDhQUlIysVnHgba2trC3v/feeyH/O1RUK83v92PWrFnYtm0b\nli5dih07dqC3txcdHR1wOp0AAKfTiY6OjmgOL2OIKtrAwAAaGhqwa9cuNDQ0YPr06SMeCh0OBxwO\nx6j39/v9wY+urq5ophB3ent7cfny5eBHOFE9e3S73XC73Vi+fDkAYMOGDaisrITL5UJ7eztcLhfa\n2tqQkZEx6v1zcnKiOW1cmz59OqZPnx78/dWrV0OOjWqluVwuZGVlobGxEQBQV1eHhQsXYt26daiq\nqgIAVFVVoaysLJrDyxii3qf9/ve/x/e+9z309/dj7ty52L9/PwYHB1FeXo59+/YhOzsb1dXVVs5V\n/ieqZ48TOqGePQIA6uvrxxxj6bNHsZeiEVI0QopGSNEIKRohRSOkaIT0yrVN0tLSwt4e7gVkrTRC\nikZI0QgpGiFFI6RohBSNkKIRUjRCuiICICEhsq/doaEhy845Z86csLfrikicUTRCikZI0QgpGiFF\nI6RohBSNkDbXAO69996Ixp04ccKyc959991R31crjZCiEVI0QopGSNEIKRohRSOkaIS0uQawffv2\niMZZubmeNm1a1PfVSiOkaIQUjZCiEVI0QopGSNEIKRohba4x9lu0J0Okb0Uf9b4WzkOmiKIRUjRC\nikZI0QgpGiFFI6RohBSNkK6IAPjvf/9r9xTGRSuNkKIRUjRCikZI0QgpGiFFI6RohLS5BtDe3m73\nFMZFK42QohFSNEKKRkjRCCkaIUUjpGiEtLkGcO3aNbunMC5aaYQUjZCiEVI0QopGSNEIKRohRSOk\nzTWAwcFBu6cwLlpphBSNkKIRUjRCikZI0QgpGiFFI6TNNaz9oeNTQSuNkKIRUjRCikZI0QgpGiFF\nI6RohBSNkKIRUjRCikZI0QgpGiFFI6RohBSNkKIR0tsNAPT09Ng9hXHRSiOkaIQUjZCiEVI0QopG\nSNEIKRohba4BHDhwwO4pjItWGiFFI6RohBSNkKIRUjRCikZI0Qhpcw3g+PHjEY372te+NskziYxW\nGiFFI6RohBSNkKIRUjRCikZI0QgpGiFFI6RohBSNkKIRUjRCikZI0QgpGiFFI6S3GwCYM2eO3VMY\nF600QopGSNEIKRohRSOkaIQUjZCiEdLmGsCXvvQly46VlZUV0bimpqaoz6GVRkjRCCkaIUUjpGiE\nFI2QohFSNELaXAO44447LDvWjh07Ihr39NNPR30OrTRCikZI0QgpGiFFI6RohBSNkKIRUjRCUV8R\nqaysxMsvv4yEhAQsXrwY+/fvR29vLzZu3Ij//Oc/yM7ORnV1NdLS0qyc76RITLTuwtCqVassO1Yo\nUa20QCCAF198EQ0NDfjXv/6FwcFBHDp0CD6fD16vF42NjVizZg18Pp/V8xVEGe3uu+9GUlIS+vr6\nMDAwgL6+PsyZMwc1NTWoqKgAAFRUVODIkSOWTlY+EVW0GTNm4Ic//CHuuecezJkzB2lpafB6vejo\n6IDT6QQAOJ1OdHR0WDpZ+URUD+b//ve/8bvf/Q6BQACpqal46KGH8PLLLw8b43A44HA4Rr2/3+8P\n/jotLQ3p6enRTCOudHV1obu7O6KxUUV78803sWrVKsycORMA8OCDD+LMmTNwuVxob2+Hy+VCW1sb\nMjIyRr1/Tk5ONKeNa+np6cO+eAOBQMixUT085uXl4ezZs/joo49gjEFdXR08Hg/WrVuHqqoqAEBV\nVRXKysqiObyMIaqVVlBQgK1bt2LZsmVISEjA0qVLsXPnTly/fh3l5eXYt29f8Cm/WM9hjDFTekKH\nAyUlJVN5yjFF+lcQ6t/oz6qrq4voWPfdd1/Y2+vr60POS283QGQxItXX12fZsULRZSxCikZI0Qgp\nGiFFI6RohBSNkKIR0ubaYm1tbZN+Dq00QopGSNEIKRohRSOkaIQUjZCiEdLm2mKdnZ2Tfg6tNEKK\nRkjRCCkaIUUjpGiEFI2QohHS5tpiN2/enPRzaKURUjRCikZI0QgpGiFFI6RohBSNkKIR0hURiw0N\nDU36ObTSCCkaIUUjpGiEFI2QohFSNEKKRkiba4tNxXdi1EojpGiEFI2QohFSNEKKRkjRCCkaIW2u\nLaZXrmVUikZI0QgpGiFFI6RohBSNkKIR0ubaYvX19ZN+Dq00QopGSNEIKRohRSOkaIQUjZCiEVI0\nQroiYrH9+/dHNG7BggVRn0MrjZCiEVI0QopGSNEIKRohRSOkaIS0ubZYa2trROO0ub7NKBohRSOk\naIQUjZCiEVI0QopGSJtri2VkZEz6ObTSCCkaIUUjpGiEFI2QohFSNEKKRsiWzXVycnLY2/v7+6do\nJtZzuVyTfg6tNEKKRkjRCCkaIUUjpGiEFI2QohFSNEK2XBFZv3592NsPHjw4RTOx3l133TXp59BK\nI6RohBSNkKIRUjRCikZI0QgpGiGHmYofgf7ZEzoc6OnpCTumrKxsimZjvUh/6F1CQvj1Ul9fH/Kn\n02ulEVI0QopGSNEIKRohRSOkaIQUjZAtr1ynpKTYcdopMdam2ZJzTPoZxHKKRkjRCCkaIUUjpGiE\nFI2QohGyZXMd6bdBl9FppRFSNEKKRkjRCCkaIUUjpGiEFI2QLZvry5cv23HauKGVRkjRCCkaIUUj\npGiEFI2QohFSNEKKRsiWKyK9vb12nDZuaKURUjRCikZI0QgpGiFFI6RohBSNkC2b64GBATtOGze0\n0ggpGiFFI6RohBSNkKIRUjRCikbIls31FH9X+bgTdqVt374dTqcTixcvDn6us7MTXq8Xubm5WLt2\nLbq7u4O3VVZWYv78+cjLy8Nrr702ebO+zYWNtm3bNtTW1g77nM/ng9frRWNjI9asWQOfzwcAuHjx\nIg4fPoyLFy+itrYWu3btivgHC8j4hI22evVqpKenD/tcTU0NKioqAAAVFRU4cuQIAODo0aPYvHkz\nkpKSkJ2djXnz5uHcuXOTNO3b27ifiHR0dMDpdAIAnE4nOjo6AHzyDV3cbndwnNvtRktLi0XTlM+a\n0BMRh8MBh8MR9vbRHDhwIPjrwsJCFBYWTmQacaGrq2vY84Nwxh3N6XSivb0dLpcLbW1tyMjIAABk\nZmaiqakpOK65uRmZmZmjHuORRx4Z72njXnp6+rB/igKBQMix4354LC0tRVVVFQCgqqoq+BOYSktL\ncejQIfT398Pv9+PSpUtYsWLFeA8vEQi70jZv3oxTp07hypUryMrKws9//nPs3r0b5eXl2LdvH7Kz\ns1FdXQ0A8Hg8KC8vh8fjQWJiIvbu3Rv2oVOiZ8vPTztx4kTYMb/4xS+maDaxSz8/Lc4oGiFFI6Ro\nhBSNkKIRUjRCikZI0QjZ8naDz15YlvHTSiOkaIQUjZCiEVI0QopGSNEIKRohWzbXL7zwQtjbk5OT\np2gmnLTSCCkaIUUjpGiEFI2QohFSNEKKRsiWzfWZM2fC3l5cXDxFM+GklUZI0QgpGiFFI6RohBSN\nkKIRUjRCtmyu77zzTjtOGze00ggpGiFFI6RohBSNkKIRUjRCikZI0QjZckVk9uzZdpw2bmilEVI0\nQopGSNEIKRohRSOkaIQUjZAtm+uZM2dacpyCgoIxx1y4cMGSc8USrTRCikZI0QgpGiFFI6RohBSN\nkKIRsmVzPW3aNEuOs3PnzjHHPP7445acK5ZopRFSNEKKRkjRCCkaIUUjpGiEFI2QLZvrhARrvlby\n8/MtOQ4brTRCikZI0QgpGiFFI6RohBSNkKIRsmVz7XA4LDlOd3e3Jcdho5VGSNEIKRohRSOkaIQU\njZCiEVI0QopGyJYrIlZpbW21ewq20EojpGiEFI2QohFSNEKKRkjRCCkaIerNdU9Pj91TsIVWGiFF\nI6RohBSNkKIRUjRCikZI0QhRb677+/vtnoIttNIIKRohRSOkaIQUjZCiEVI0QopGiHpzbYyxewq2\n0EojpGiEFI2QohFSNEKKRkjRCCkaIUUjpCsihLTSCCkaIUUjpGiEFI2QohFSNEKKRoh6c6338gsN\nRSOkaIQUjZCiEVI0QopGSNEIUW+uX3nlFbunYAutNEKKRkjRCCkaIUUjpGiEFI2QohGi3lxXV1eP\nOWbp0qVTMJOppZVGSNEIKRohRSOkaIQUjZCiEVI0QopGiPqKyPXr1+2egi200ggpGiFFI6RohBSN\nkKIRUjRCikYobLTt27fD6XRi8eLFwc/96Ec/Qn5+PgoKCvDggw8O+6nulZWVmD9/PvLy8vDaa69N\n3qz/Z/bs2WN+xKOw0bZt24ba2tphn1u7di3eeecdXLhwAbm5uaisrAQAXLx4EYcPH8bFixdRW1uL\nXbt2YWhoaPJmfhsLG2316tVIT08f9jmv14uEhE/utnLlSjQ3NwMAjh49is2bNyMpKQnZ2dmYN28e\nzp07N0nTvr1N6N+0l156Cffffz8AoLW1FW63O3ib2+1GS0vLxGYno4r6gvEvf/lLJCcn4+GHHw45\nxuFwjPp5v98f/HVaWtqI1Xw76urqQnd3d0Rjo4p24MABHDt2DH/729+Cn8vMzERTU1Pw983NzcjM\nzBz1/jk5OcN+39XVNWa4SMYAwM2bNzFt2rQJH2sqx3zqs38vgUAg5LhxPzzW1tZiz549OHr0KO64\n447g50tLS3Ho0CH09/fD7/fj0qVLWLFiRUTHjOQrLNKvwps3b07Z+aycd6TjgDFW2ubNm3Hq1Clc\nuXIFWVlZ+NnPfobKykr09/fD6/UCAL761a9i79698Hg8KC8vh8fjQWJiIvbu3Rvy4VEmJmy0gwcP\njvjc9u3bQ45/6qmn8NRTT018VhKemWLFxcUGgD7G+CguLg75d+gw5jb95vbEdO2RkKIRUjRCtker\nra1FXl4e5s+fj1/96lcjbh/tlYbPa2pqQklJCRYuXIhFixbhueeeGzHm448/xsqVK1FYWAiPx4Of\n/OQnIY83ODiIoqIirFu3LuSY7OxsLFmyBEVFRSH3o93d3diwYQPy8/Ph8Xhw9uzZYbe///77KCoq\nCn6kpqaOOvcRpupZ42gGBgbM3Llzjd/vN/39/aagoMBcvHhx2JjXX3/dNDQ0mEWLFoU8Tltbmzl/\n/rwxxpjr16+b3NzcEccxxpje3l5jjDG3bt0yK1euNKdPnx71eL/5zW/Mww8/bNatWxfynNnZ2ebq\n1ath/3xbt241+/btC56zu7s75NjBwUHjcrnMhx9+GPaYxhhj60o7d+4c5s2bh+zsbCQlJWHTpk04\nevTosDGjvdLweS6XC4WFhQCAlJQU5Ofno7W1dcS4L37xiwA++YbVg4ODmDFjxogxzc3NOHbsGB57\n7LExf2pUuNt7enpw+vTp4L42MTERqampIcfX1dVh7ty5yMrKCntOwOaHx5aWlmGTtOKVgUAggPPn\nz2PlypUjbhsaGkJhYSGcTidKSkrg8XhGjHnyySexZ8+e4MtPoTgcDtx3331YtmwZXnzxxRG3+/1+\nzJo1C9u2bcPSpUuxY8cO9PX1hTzeoUOHwl58/yxbo1l9mevGjRvYsGEDnn32WaSkpIy4PSEhAW+9\n9Raam5vx+uuv4+TJk8Nuf/XVV5GRkYGioqIxV9k//vEPnD9/HsePH8fzzz+P06dPD7t9YGAADQ0N\n2LVrFxoaGjB9+nT4fL5Rj9Xf348///nPeOihhyL6c9oa7fOvDDQ1NQ17TW48bt26hfXr1+P73/8+\nysrKwo5NTU3Fd77zHbz55pvDPv/GG2+gpqYGOTk52Lx5M06cOIGtW7eOeoxP38owa9YsPPDAAyNe\n8HW73XC73Vi+fDkAYMOGDWhoaBj1WMePH8dXvvIVzJo1K6I/q61PRG7dumW+/OUvG7/fb27evDnq\nExFjjPH7/WGfiAwNDZktW7aYJ554IuSYy5cvm66uLmOMMX19fWb16tWmrq4u5PiTJ0+a7373u6Pe\n1tvba65du2aMMebGjRtm1apV5i9/+cuIcatXrzbvv/++McaYn/70p+bHP/7xqMfbuHGjOXDgQMi5\nfJ6t0Ywx5tixYyY3N9fMnTvXPPPMMyNu37Rpk5k9e7ZJTk42brfbvPTSSyPGnD592jgcDlNQUGAK\nCwtNYWGhOX78+LAxb7/9tikqKjIFBQVm8eLF5te//nXYeZ08eTLks8cPPvjAFBQUmIKCArNw4cJR\n522MMW+99ZZZtmyZWbJkiXnggQdGffZ448YNM3PmzOAXQSR07ZGQ7ZtrGT9FI6RohBSNkKIRUjRC\nikbo/wAjxt75YX+EQwAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x1098fd2d0>"
       ]
      }
     ],
     "prompt_number": 61
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.set_printoptions(suppress=True)\n",
      "cons_exper"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 62,
       "text": [
        "{'constant_0': array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.]),\n",
        " 'type_bottle_0': array([ 1., -0., -0., -0., -0., -0., -0., -0.]),\n",
        " 'type_cat_0': array([ 0.,  1., -0.,  0., -0., -0., -0., -0.]),\n",
        " 'type_chair_0': array([-0., -0.,  1., -0., -0.,  0., -0.,  0.]),\n",
        " 'type_face_0': array([-0.,  0., -0.,  1., -0., -0., -0., -0.]),\n",
        " 'type_house_0': array([-0., -0.,  0., -0.,  1.,  0., -0., -0.]),\n",
        " 'type_scissors_0': array([-0., -0., -0., -0.,  0.,  1., -0.,  0.]),\n",
        " 'type_scrambled_0': array([ 0., -0., -0.,  0., -0., -0.,  1.,  0.]),\n",
        " 'type_shoe_0': array([-0., -0.,  0.,  0., -0., -0., -0.,  1.])}"
       ]
      }
     ],
     "prompt_number": 62
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next we make the drift terms.  Here we use a spline set:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "drift = fmrid.natural_spline(t, [t.mean()])\n",
      "# Make confounds a little easier to display\n",
      "drift = drift / drift.max(axis=0)\n",
      "plt.figure(figsize=(10, 20))\n",
      "plt.imshow(drift)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 63,
       "text": [
        "<matplotlib.image.AxesImage at 0x109a42150>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAFEAAAR4CAYAAAB+YzrbAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGj5JREFUeJzt3X9sVXf9x/HXZS3+mrKN2FtHMXWFWgqsZXEQl+AuksuS\naUkVbABlDWz6h//Ioslw/1hNpJ3EKHPyj2FLo4mIf1iaBdkkyIYZiAuYmTFTMksspTQOrG6DUdp+\nvn8I90vXW3q4fZWe0z4fSRNoT+/98OS0fZ9zT+9NhRCCMC4zJnsBUwERDYhoQEQDIhoQ0cAecf/+\n/aqqqtL8+fP11FNPjfh4JpNRKpVK7FtewWhgYCBUVFSEzs7O0N/fH2pqasLJkyeHbZPvLr/3ve+N\neJ+kWL7lY90Tjx07pnnz5qm8vFzFxcVat26d9u7d67yLWLJG7O7u1ty5c3N/LysrU3d3t/MuYqnI\neWOjfs94n6amptyfM5mMMpmMcxm3nDXinDlz1NXVlft7V1eXysrKRmx3fcQpwfIT5aorV66Ee+65\nJ3R2dobLly9H/sGSj2LwQyTfWz7WPbGoqEjPPPOMHnroIQ0ODurRRx/VggULnHcRS6mr/+u37g5T\nKUW5y9tuu23MbYaGhhxLuin51s4RiwERDYhoQEQDIhoQ0YCIBkQ0sB6xOM2YMfb/72QM2/mwJxoQ\n0YCIBkQ0IKIBEQ2IaEBEAyIaJPqIJS6Ss9IYI6IBEQ2IaEBEAyIaENGAiAYM2wbJWWmMEdGAiAZE\nNCCiARENiGhARAOGbYPkrDTGiGhARAMiGhDRgIgGRDQgogHDtkFyVhpjRDQgogERDYhoQEQDIhoQ\n0YCIBhyxGCRnpTFGRAMiGhDRgIgGRDQgogERDWI7bEd5EsqoTwQ80c+zyZ5oQEQDIhoQ0YCIBkQ0\nIKIBEQ1iO2xzZnuaIaIBEQ2IaEBEAyIaENGAiAYM2wbJWWmMEdGAiAZENCCiARENiGhARAMiGiT6\niIVrcaYQIhoQ0YCIBkQ0IKIBEQ2IaJDoYTsukrPSGCOiARENiGhARAMiGhDRgIgGiR62o57Znmjs\niQZENCCiARENiGhARAMiGhDRINHDdlwkZ6UxRkQDIhoQ0YCIBkQ0IKIBEQ1iO2w7nz9xorEnGhDR\ngIgGRDQgogERDYhoQEQDIhrE9oiFhwemGSIaENGAiAZENCCiARENiGiQ6GGbhwemECIaENGAiAZE\nNCCiARENiGiQ6GE7LpKz0hgjogERDYhoQEQDIhoQ0YCIBoketqOe2Y6y3Xie0Jw90YCIBkQ0IKIB\nEQ2IaEBEAyIaENFgWhyxTLSC9sSuri6tWLFCCxcu1KJFi/T0009Lki5cuKBsNqvKykqtWrVKfX19\n1sXGVUERi4uL9ZOf/ESvv/66jh49qp///Od644031NLSomw2q46ODq1cuVItLS3u9cZSQRFLS0tV\nW1srSbr99tu1YMECdXd3q729XY2NjZKkxsZGtbW1+VYaY+P+wXL69GmdOHFCy5YtU29vr9LptCQp\nnU6rt7d33AtMgnH9YHnnnXe0Zs0a7dixQx/96EeHfSyVSo36jb+pqSn350wmo0wmM55lTLqCI165\nckVr1qzRxo0bVV9fL+l/e9+5c+dUWlqqnp4elZSU5P3c6yNOBQV9OYcQ9Oijj6q6ulpbtmzJvX/1\n6tVqbW2VJLW2tubiTnWpUMAp3T/96U/63Oc+p3vvvTf3Jdvc3KylS5eqoaFB//znP1VeXq49e/bo\njjvuGH6HqVSks8gLFy4cc5s333wz0nr7+/vH3CZqhnzbFRRxPKJGXLx48ZjbnDp1KtJ9TnREDvsM\niGhARAMiGhDRgIgGRDQgogFntg3YEw2IaEBEAyIaENGAiAZENCCiAcO2AXuiARENiGhARAMiGhDR\ngIgGRDQgokFsj1icr1HFr+omABENiGhARAMiGhDRgIgGRDSI7bAdZUDm4YEphIgGRDQgogERDYho\nQEQDIhrEdth2ntmeaOyJBkQ0IKIBEQ2IaEBEAyIaENEgtsN2XAbpKNgTDYhoQEQDIhoQ0YCIBkQ0\nIKJBbIftW33N9niwJxoQ0YCIBkQ0IKIBEQ2IaEBEAyIaxPaIhQvfpxkiGhDRgIgGRDQgogERDYho\nENthu6ho7KVFeaJKiYcHEoGIBkQ0IKIBEQ2IaEBEAyIaxHbY5sz2NENEAyIaENGAiAZENCCiAREN\nYjts85Qu0wwRDYhoQEQDIhoQ0YCIBkQ0IKJBbI9YnC8Iy7U4CUBEAyIaENGAiAZENCCiARENpsWw\nPdHYEw2IaEBEAyIaENGAiAZENCCiQWyHbeeF75zZTgAiGhDRgIgGRDQgogERDYhoENthO8qF71Gf\n0mWixWMVCUdEAyIaENGAiAZENCCiARENiGgQ2yMWLnyfZohoQEQDIhoQ0YCIBkQ0IKJBbIdtXhB2\nmiGiARENiGhARAMiGhDRgIgGsR22o5zZ5lqcKYSIBkQ0IKIBEQ2IaEBEAyIaJHrYnhJntgcHB7Vk\nyRLV1dVJki5cuKBsNqvKykqtWrVKfX19lkXG3bgi7tixQ9XV1bn/6ZaWFmWzWXV0dGjlypVqaWmx\nLDLuCo545swZ7du3T4899phCCJKk9vZ2NTY2SpIaGxvV1tbmWWXMFRzx8ccf1/bt24d97+rt7VU6\nnZYkpdNp9fb2jn+FCVDQD5bnn39eJSUlWrJkiQ4dOpR3m1QqNeo39KamptyfM5mMMplMIcuIjYIi\nvvLKK2pvb9e+ffv03nvv6b///a82btyodDqtc+fOqbS0VD09PSopKcn7+ddHnAoK+nLetm2burq6\n1NnZqd27d+vzn/+8fvnLX2r16tVqbW2VJLW2tqq+vt662LiyDNvXvmy3bt2qP/zhD6qsrNTBgwe1\ndetWx83HXipc+9F6q+4wlVKUu2xoaBhzmyNHjkS6z/Pnz4+5zaVLlyLdVr61c8RiwLGzARENiGhA\nRAMiGhDRgIgGRDRI9LDNtThTCBENiGhARAMiGhDRgIgGRDRI9LDNme0phIgGRDQgogERDYhoQEQD\nIhoketiOemabYTsBiGhARAMiGhDRgIgGRDQgokFsh23nyysxbCcAEQ2IaEBEAyIaENGAiAZENCCi\nQWyPWJxHGRyxJAARDYhoQEQDIhoQ0YCIBkQ0iO2wHeXhgSjbSAzbiUBEAyIaENGAiAZENCCiAREN\nYjts81um0wwRDYhoQEQDIhoQ0YCIBkQ0SPSwzYXvUwgRDYhoQEQDIhoQ0YCIBkQ0IKJBbI9Yohxl\ncMQyhRDRgIgGRDQgogERDYhoQESD2A7bPC/ONENEAyIaENGAiAZENCCiARENYjtscy3ONENEAyIa\nENGAiAZENCCiARENYjtsR326liiiDuUF3/6E3vo0QUQDIhoQ0YCIBkQ0IKIBEQ2IaBDbI5YoRxk8\nCeUUQkQDIhoQ0YCIBkQ0IKIBEQ0SPWzzJJRTCBENiGhARAMiGhDRgIgGRDRI9LAd9cw21+IkABEN\niGhARAMiGhDRgIgGRDRI9LDNme0phIgGRDQgogERDYhoQEQDIhoketjmmu0phIgGRDQgogERDYho\nQEQDIhoQ0SDRRyxRj0S4FicBiGhARAMiGhDRgIgGRDQgokFsh+0op/6LiqItP7bDdl9fn9auXasF\nCxaourpaf/7zn3XhwgVls1lVVlZq1apV6uvrc641tgqO+K1vfUsPP/yw3njjDb322muqqqpSS0uL\nstmsOjo6tHLlSrW0tDjXGlsFRfzPf/6jw4cPa/PmzZL+92U1a9Ystbe3q7GxUZLU2NiotrY230pj\nrKCInZ2d+vjHP65Nmzbpvvvu09e//nW9++676u3tVTqdliSl02n19vZaFxtXBf1gGRgY0PHjx/XM\nM8/o/vvv15YtW0Z86aZSqVHPsjQ1NeX+nMlklMlkCllGbBQUsaysTGVlZbr//vslSWvXrlVzc7NK\nS0t17tw5lZaWqqenRyUlJXk///qIU0FBX86lpaWaO3euOjo6JEkHDhzQwoULVVdXp9bWVklSa2ur\n6uvrfSuNsYLnxJ/97Gf66le/qv7+flVUVOi5557T4OCgGhoatGvXLpWXl2vPnj3OtcZWwRFramr0\nl7/8ZcT7Dxw4MK4FXeM8sz3ROOwzIKIBEQ2IaEBEAyIaENGAiAaxPbMdZZCOeuG786Wa8mFPNCCi\nARENiGhARAMiGhDRgIgGRDSI7RFLlKOMqNfY8Ku6CUBEAyIaENGAiAZENCCiARENYjtsRxmkow7b\nsb3wHf+PiAZENCCiARENiGhARAMiGiR62ObllaYQIhoQ0YCIBkQ0IKIBEQ2IaJDoYZsz21MIEQ2I\naEBEAyIaENGAiAZENCCiAUcsBuyJBkQ0IKIBEQ2IaEBEAyIaENEgtsO287dMGbYTgIgGRDQgogER\nDYhoQEQDIhrEdtjmwvdphogGRDQgogERDYhoQEQDIhowbF8V5TkWQwj5bz/SKnBDRDQgogERDYho\nQEQDIhoQ0SDRwzaXkUwhRDQgogERDYhoQEQDIhoQ0YCIBok+YnG+ICwPD0wyIhoQ0YCIBkQ0IKIB\nEQ2IaJDoYZvXMp1CiGhARAMiGhDRgIgGRDQgokGih+2oZ7aLisb+Z45nIGdPNCCiARENiGhARAMi\nGhDRgIgGsR22owzSt/q3TEe9/YI/EzlENCCiARENiGhARAMiGhDRgIgGsT1icV6LE/XIplDsiQZE\nNCCiARENiGhARAMiGhDRILbDdpTT9Tw8MIUQ0YCIBkQ0IKIBEQ2IaEBEg9gO285rcVxP6TIa9kQD\nIhoQ0YCIBkQ0IKIBEQ2IaBDbYftWX0YynmeFZ080IKIBEQ2IaEBEAyIaENGAiAZENCg4YnNzsxYu\nXKjFixdrw4YNunz5si5cuKBsNqvKykqtWrVKfX19hS9sxowx32677bZIb1FuK5VKjfk26loL+Qee\nPn1av/jFL3T8+HH97W9/0+DgoHbv3q2WlhZls1l1dHRo5cqVamlpKThikhQU8WMf+5iKi4t18eJF\nDQwM6OLFi7r77rvV3t6uxsZGSVJjY6Pa2tqsi42rgiLedddd+va3v61PfvKTuvvuu3XHHXcom82q\nt7dX6XRakpROp9Xb22tdbFwVdBbnzTff1E9/+lOdPn1as2bN0le+8hX96le/GrbNjb6PNDU15f6c\nyWSUyWQKWcaEGxwc1NDQ0JjbFRTx1Vdf1QMPPKDZs2dLkr785S/ryJEjKi0t1blz51RaWqqenh6V\nlJTk/fzrI8bZtR9M1wwMDOTdrqAv56qqKh09elSXLl1SCEEHDhxQdXW16urq1NraKklqbW1VfX19\nITefOAXtiTU1NXrkkUf0mc98RjNmzNB9992nb3zjG3r77bfV0NCgXbt2qby8XHv27HGvN5ZSYbSX\nu5moO0ylRn2Fnev99re/HXObI0eORLrPjo6OMbc5ePDgmNtc+8p7v0Q/PMC1OFMIEQ2IaEBEAyIa\nENGAiAZENGDYvon7G/VzC/5M5BDRgIgGRDQgogERDYhoQESD2A7bzgGZYTsBiGhARAMiGhDRgIgG\nRDQgokFsh+0ow2+U1yiVGLYTgYgGRDQgogERDYhoQEQDIhoQ0SDRRyzOZ3zniGWSEdGAiAZENCCi\nARENiGhARINpMWwXFxdb7m/Uzy34M5FDRAMiGhDRgIgGRDQgogERDabFsO18qaZ82BMNiGhARAMi\nGhDRgIgGRDQgogHD9k3eVj7siQZENCCiARENiGhARAMiGhDRgIgG0+KIJcqv9HItziQjogERDYho\nQEQDIhoQ0YCIBrEdtqMM0lGfFyfKdlFvKx/2RAMiGhDRgIgGRDQgogERDYhoENthm2txphkiGhDR\ngIgGRDQgogERDYhoENthO5VKjbmN8zIShu1JRkQDIhoQ0YCIBkQ0IKIBEQ2IaBDbIxbnKX2uxUkA\nIhoQ0YCIBkQ0IKIBEQ2IaBDbYdv5grAM2wlARAMiGhDRgIgGRDQgogERDRI9bHNmewohogERDYho\nQEQDIhoQ0YCIBgzbV0V5vdPRsCcaENGAiAZENCCiARENiGhARINED9vOy0gYticZEQ2IaEBEAyIa\nENGAiAZENCCiQWyPWG71b5lyxDLJiGhARAMiGhDRgIgGRDQgokFsh20eHphmiGhARAMiGhDRgIgG\nRDQgokGih+1EnNnevHmz0um0Fi9enHvfhQsXlM1mVVlZqVWrVqmvry/3sebmZs2fP19VVVV68cUX\nC15U0tww4qZNm7R///5h72tpaVE2m1VHR4dWrlyplpYWSdLJkyf1m9/8RidPntT+/fv1zW9+U0ND\nQxO38hi5YcTly5frzjvvHPa+9vZ2NTY2SpIaGxvV1tYmSdq7d6/Wr1+v4uJilZeXa968eTp27NgE\nLTtebvoHS29vr9LptCQpnU6rt7dXknT27FmVlZXltisrK1N3d7dpmfE2rh8sqVTqhi/5MdrHmpqa\ncn/OZDLKZDLjWcaEOX/+vM6fPz/mdjcdMZ1O69y5cyotLVVPT49KSkokSXPmzFFXV1duuzNnzmjO\nnDl5b+P6iHE2e/ZszZ49O/f3U6dO5d3upr+cV69erdbWVklSa2ur6uvrc+/fvXu3+vv71dnZqVOn\nTmnp0qWFrD1xbrgnrl+/Xi+99JLeeustzZ07Vz/4wQ+0detWNTQ0aNeuXSovL9eePXskSdXV1Wpo\naFB1dbWKioq0c+fOSK/uMxWkQgjhlt5hKqUod/n3v/99zG3+8Y9/RLrP67/NjOaFF14Yc5vf/e53\nedfOYZ8BEQ2IaEBEAyIaENGAiAZENCCiwbR4eCDKqf+ZM2dGuq182BMNiGhARAMiGhDRgIgGRDQg\nokGih23nhe8f+MAHIt1WPuyJBkQ0IKIBEQ2IaEBEAyIaENFgWgzbnNlOACIaENGAiAZENCCiAREN\niGgQ22Hb+fyJUYbtD37wg5FuKx/2RAMiGhDRgIgGRDQgogERDYhoQESD2B6x3OqHB7gWZ5IR0YCI\nBkQ0IKIBEQ2IaEBEg9gO21GemMg5bPPwwCQjogERDYhoQEQDIhoQ0YCIBrEdtqNcZ8OwPYUQ0YCI\nBkQ0IKIBEQ2IaEBEg9gO27f6MhKG7UlGRAMiGhDRgIgGRDQgogERDabFsB3l6VoYticZEQ2IaEBE\nAyIaENGAiAZENCCiwbQ4Yony8MCHP/zhSLeVD3uiARENiGhARAMiGhDRgIgGRDSYFsN2lIcHPvSh\nD0W6rXzYEw2IaEBEAyIaENGAiAZENCCiQWyH7Vv9W6ac2Z5kRDQgogERDYhoQEQDIhoQ0SC2w/at\nPrPNsD3JiGhARAMiGhDRgIgGRDQgogERDThiueojH/lIpNvKhz3RgIgGRDQgogERDYhoQEQDIhow\nbF/FsD3JiGhARAMiGhDRgIgGRDQgogHD9lUM25OMiAZENCCiARENiGhARAMiGkyLYTuK22+/veDP\nZU80IKIBEQ2IaEBEAyIaENGAiAZENJgWRyxRbosjlklGRAMiGhDRgIgGRDQgogERHcINbNq0KZSU\nlIRFixbl3ved73wnVFVVhXvvvTd86UtfCn19fbmPbdu2LcybNy98+tOfDi+88ELe2xzjLmNttLXf\n8F/08ssvh+PHjw+L+OKLL4bBwcEQQghPPPFEeOKJJ0IIIbz++uuhpqYm9Pf3h87OzlBRUZHbLspC\nkmC0td/wy3n58uW68847h70vm83mDqOWLVumM2fOSJL27t2r9evXq7i4WOXl5Zo3b56OHTs2AV87\n8TOu74nPPvusHn74YUnS2bNnVVZWlvtYWVmZuru7x7e6hCj4BMQPf/hDzZw5Uxs2bBh1m1Qqlff9\nTU1NuT9nMhllMplClzGhDh06pEOHDo294VjfBzo7O4d9TwwhhOeeey488MAD4dKlS7n3NTc3h+bm\n5tzfH3rooXD06NFI31f++Mc/jrWMSNs4byvfNqPluukv5/3792v79u3au3fvsJc4X716tXbv3q3+\n/n51dnbq1KlTWrp0aeT/ccc2ztuKen/SGF/O69ev10svvaS33npLc+fO1fe//301Nzerv79f2WxW\nkvTZz35WO3fuVHV1tRoaGlRdXa2ioiLt3Llz1C/nqeaGEX/961+PeN/mzZtH3f7JJ5/Uk08+Of5V\nJc2Y3xzMHnzwwSApkW8PPvhg3n9TKoQQJuI/Zzrh2NmAiAZENJj0iPv371dVVZXmz5+vp556asTH\nN2/erHQ6rcWLF496G11dXVqxYoUWLlyoRYsW6emnnx6xzXvvvadly5aptrZW1dXV+u53vzvq7Q0O\nDmrJkiWqq6uL9o+4lT+Z329gYCBUVFSEzs7O0N/fH2pqasLJkyeHbZPvTNL79fT0hBMnToQQQnj7\n7bdDZWXliNsJIYR33303hBDClStXwrJly8Lhw4fz3t6Pf/zjsGHDhlBXVxfp3zGpe+KxY8c0b948\nlZeXq7i4WOvWrdPevXuHbZPvTNL7lZaWqra2VtL/HoRfsGCBzp49O2K7a8+l3d/fr8HBQd11110j\ntjlz5oz27dunxx57TCHi4DKpEbu7uzV37tzc3x1nfk6fPq0TJ05o2bJlIz42NDSk2tpapdNprVix\nQtXV1SO2efzxx7V9+/ZIV01cM6kR3YeF77zzjtauXasdO3bkvSxkxowZ+utf/6ozZ87o5ZdfHnF8\n/Pzzz6ukpERLliyJvBdKkxxxzpw56urqyv29q6tr2DnJm3HlyhWtWbNGX/va11RfX3/DbWfNmqUv\nfOELevXVV4e9/5VXXlF7e7s+9alPaf369Tp48KAeeeSRse880nfOCXLlypVwzz33hM7OznD58uW8\nP1hCyH867npDQ0Nh48aNYcuWLaNu869//Sv8+9//DiGEcPHixbB8+fJw4MCBUbc/dOhQ+OIXvxjp\n3zHpD3js27cvVFZWhoqKirBt27YRH1+3bl34xCc+EWbOnBnKysrCs88+O2Kbw4cPh1QqFWpqakJt\nbW2ora0Nv//974dt89prr4UlS5aEmpqasHjx4vCjH/3ohus6dOhQ5J/OHDsbTPqwPRUQ0YCIBkQ0\nIKIBEQ2IaPB/NxIF0NWaPmMAAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0x105439150>"
       ]
      }
     ],
     "prompt_number": 63
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.plot(drift[:, 3])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 66,
       "text": [
        "[<matplotlib.lines.Line2D at 0x109b0ddd0>]"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEACAYAAAC57G0KAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X1cVVW+x/EPBpVajWZmwsGrAgKGoqmZNVOU42CWTKkV\n2Zijpoxdx7o15VQzpdP4wFjNVNzK7mhWGtKDE5Z2fCqyMsXChymt0MHpiEY+kZYmcNj3jzWRqIBy\nHvY5m+/79dovOLLd+7eKvq3X2muvFWFZloWIiDhOM7sLEBGRwFDAi4g4lAJeRMShFPAiIg6lgBcR\ncSgFvIiIQzUY8KNHj6Zdu3Z069atznMmTpxIQkICqamprF+/3q8FiohI4zQY8KNGjcLtdtf58yVL\nlrB161aKi4t59tlnGT9+vF8LFBGRxmkw4H/2s5/RunXrOn++aNEiRo4cCUDfvn0pLy+nrKzMfxWK\niEij+DwGX1paSmxsbM1nl8vFjh07fL2siIj4yC8PWY9d7SAiIsIflxURER9E+nqBmJgYPB5Pzecd\nO3YQExNz3Hnx8fFs27bN19uJiDQpcXFxbN26tVF/1+cefEZGBi+88AIAa9asoVWrVrRr1+6487Zt\n24ZlWY49HnroIdtrUPvUNrXPGUd1tUW/fhbz51s+dYwb7MHffPPNvPvuu+zZs4fY2FimTJlCZWUl\nAFlZWQwaNIglS5YQHx9Py5Ytee655xpdjIiIwLJlUF4ON90Et9zS+Os0GPC5ubkNXiQnJ6fxFYiI\nSA3LgoceggcfhNNO8+1aepPVT9LS0uwuIaCc3D4ntw3UvnDjdsPBg3DDDb5fK8KyrKBs+BEREUGQ\nbiUiEpYsCy6+GO65B2680fyZL9mpHryISIh44w2oqIBhw/xzPQW8iEgIqK424+5/+hM081MyK+BF\nRELAa69BVBRkZPjvmhqDFxGxmdcL3brBY4/BwIG1f6YxeBGRMJabC61bQ3q6f6+rHryIiI0qKyEp\nCWbPhhPN+FQPXkQkTM2ZA507nzjcfaUevIiITQ4fhoQEWLjQzH8/EfXgRUTC0NNPQ+/edYe7r9SD\nFxGxwcGDpve+fLmZQVMX9eBFRMLMX/8K/fvXH+6+Ug9eRCTIdu82M2fWrTMPWOvjS3Yq4EVEguyu\nu8yaMyez0roCXkQkTHz5JfTsCZ9+Chdc0PD5CngRkTAxejS0bw9Tp57c+b5kp8+bbouIyMn59FN4\n80344ovg3E+zaEREguT++2HSJGjVKjj3Uw9eRCQI3n8fNmyAvLzg3VM9eBGRALMsuPdeePhhOPPM\n4N1XAS8iEmD5+fDtt3DLLcG9r4ZoREQCqKoK7rsPHn0UTjstuPdWD15EJIBmz4boaLj66uDfW/Pg\nRUQC5OBBSEyEN96AXr0adw0tNiYiEoIefRSuuqrx4e4r9eBFRAJg1y5ISYGPP4aOHRt/HS1VICIS\nYsaNg3POgUce8e06WqpARCSEfPKJmRr52Wf21qExeBERP/vd7+CBB6B1a3vrUMCLiPiR2w3/+hf8\n5jd2V6KAFxHxm6oq03v/y1/g9NPtrkYBLyLiN3PmQJs28Mtf2l2JoVk0IiJ+8M035qWmJUvgoov8\nd11NkxQRsdm998LevWZpAn9SwIuI2GjbNujbF/75T7Mdnz8FdKkCt9tNUlISCQkJZGdnH/fzPXv2\nMHDgQHr06EFKSgpz585tVCEiIuHq3nvhrrv8H+6+qrcH7/V6SUxMZMWKFcTExNCnTx9yc3NJTk6u\nOWfy5MkcOXKE6dOns2fPHhITEykrKyMysvY7VOrBi4gTFRTAr38NW7ZA8+b+v37AevCFhYXEx8fT\nsWNHoqKiyMzMJD8/v9Y57du358CBAwAcOHCANm3aHBfuIiJOVFUFd9xhliMIRLj7qt4kLi0tJTY2\ntuazy+Vi7dq1tc4ZO3YsV111FdHR0Rw8eJCXX345MJWKiISYv//dvK06dKjdlZxYvQEfERHR4AWm\nTZtGjx49KCgoYNu2bQwYMICNGzdy9tlnH3fu5MmTa75PS0sjLS3tlAsWEQkF+/fDQw/BsmVwElF5\n0goKCigoKPDLteoN+JiYGDweT81nj8eDy+Wqdc7q1at54IEHAIiLi6NTp058/vnn9O7d+7jrHR3w\nIiLhbPJkGDIEUlP9e91jO79Tpkxp9LXqHYPv3bs3xcXFbN++nYqKCvLy8sjIyKh1TlJSEitWrACg\nrKyMzz//nM6dOze6IBGRUPfJJ5CbCw8/bHcl9au3Bx8ZGUlOTg7p6el4vV7GjBlDcnIys2bNAiAr\nK4v777+fUaNGkZqaSnV1NX/5y18499xzg1K8iEiwWRZMnAgPPgjnnWd3NfXTi04iIqfglVdMz72o\nCIIxYVBvsoqIBMGhQ5CcDC+8AFdcEZx7atNtEZEgmDED+vULXrj7Sj14EZGTsHUrXHIJbNgAx0wm\nDCj14EVEAuiHB6v33hvccPeV1hQQEWnAokVQUgKvv253JadGQzQiIvU4dAi6djXrvPfvH/z7a4hG\nRCRApk0za73bEe6+Ug9eRKQOn38Ol10GGzdCTIw9NagHLyLiZ5YFt98Of/iDfeHuKwW8iMgJLFhg\n9lidMMHuShpPQzQiIsf45hvzxuprr5kXm+ykpQpERPxowgSoqIBnn7W7Et+yU/PgRUSOUlhoeu6f\nfmp3Jb7TGLyIyH9UVUFWltlj1QmrnivgRUT+44knzBrvw4fbXYl/aAxeRAT497+hVy9Yswbi4+2u\n5keaBy8i4oMf5rzfdVdohbuv9JBVRJq8l1+GL7+Ef/zD7kr8S0M0ItKk7dsHKSmwcKFZ7z3UaB68\niEgjjRkDLVrAk0/aXcmJaR68iEgjrFwJy5fDJ5/YXUlg6CGriDRJ330H48bBM8/AOefYXU1gaIhG\nRJqku++GsjKYN8/uSuqnIRoRkVNQWAjz5zt3aOYHGqIRkSblyBEYPRoee8y8tepkCngRaVL+/GeI\ni4Obb7a7ksDTGLyINBnr10N6utmCr317u6s5OVqqQESkAZWVMGoUzJwZPuHuKwW8iDQJ06dDdDTc\neqvdlQSPhmhExPE2bIABA8wQjctldzWnRkM0IiJ1qKiAkSPNJh7hFu6+UsCLiKP9+c/QoUPTGpr5\ngV50EhHH+ugjsxTBhg0QEWF3NcGnHryIONLhw6bX/vjj5uFqU6SHrCLiSHfdBTt2QF5eePfeA/qQ\n1e12k5SUREJCAtnZ2Sc8p6CggJ49e5KSkkJaWlqjChER8ZeCAhPsTz0V3uHuq3p78F6vl8TERFas\nWEFMTAx9+vQhNzeX5OTkmnPKy8u57LLLWLp0KS6Xiz179nDeCRZ4UA9eRILhwAFITYWcHLjmGrur\n8V3AevCFhYXEx8fTsWNHoqKiyMzMJD8/v9Y5L730EkOHDsX1n/lHJwp3EZFgmTjRzHl3Qrj7qt6A\nLy0tJTY2tuazy+WitLS01jnFxcXs27ePK6+8kt69e/Piiy8GplIRkQa89hp88IFZKVIamCYZcRKD\nV5WVlRQVFbFy5UoOHTpEv379uOSSS0hISPBbkSIiDdm5E26/HfLz4ayz7K4mNNQb8DExMXg8nprP\nHo+nZijmB7GxsZx33nk0b96c5s2bc/nll7Nx48YTBvzkyZNrvk9LS9MDWRHxi+pqs8b7b34Dl1xi\ndzW+KSgooKCgwC/Xqvcha1VVFYmJiaxcuZLo6Gguvvji4x6yfvbZZ0yYMIGlS5dy5MgR+vbtS15e\nHl27dq19Iz1kFZEAefxxyM2F996DqCi7q/GvgG3ZFxkZSU5ODunp6Xi9XsaMGUNycjKzZs0CICsr\ni6SkJAYOHEj37t1p1qwZY8eOPS7cRUQCZdMmsxzBmjXOC3df6UUnEQlbhw9Dnz5wzz1mQTEn8iU7\nFfAiErYmTIDdu2HBAue+0BSwIRoRkVC1aBEsXmzWeHdquPtKPXgRCTulpdCrFyxcCJdeanc1gaUN\nP0SkyfB64ZZb4Le/dX64+0oBLyJhZdo0MyTz+9/bXUno0xi8iISNd981K0R+/DGcdprd1YQ+9eBF\nJCzs3m2GZp57rulu4HGq9JBVREJedTVcey107w4zZthdTXDpIauIONrMmVBeDg8/bHcl4UVj8CIS\n0latgr/+Fdat01IEp0o9eBEJWV9/DcOHm3H3o7amkJOkMXgRCUleL1x9tVlrZupUu6uxj8bgRcRx\npkyBykrzVRpHY/AiEnIWL4Y5c8x890ilVKPpH52IhJSSErM708KF0K6d3dWENw3RiEjIOHwYhg2D\n++6Dyy6zu5rwp4esIhISLAtGjYIjR+Cll7QE8A+0HryIhL1nnjFj7mvWKNz9RT14EbHd6tVw3XXm\na3y83dWEFk2TFJGwVVoKN9xgXmZSuPuXAl5EbPP99zBkCPz3f8M119hdjfNoiEZEbGFZcNttcOAA\nvPyyxt3rooesIhJ2cnLMAmKrVyvcA0UBLyJBt3y52Xpv9Wo46yy7q3EuBbyIBNUXX8CvfmWGZTp1\nsrsaZ9NDVhEJmvJyyMgwG3dccYXd1TifHrKKSFBUVcGgQZCUBE88YXc14UPz4EUk5N1xBzRrBo89\nZnclTYfG4EUk4HJyoKDAPFTV8r/Bo3/UIhJQS5aYHZk++AB+8hO7q2laFPAiEjAbNsDIkbBoEXTu\nbHc1TY/G4EUkIEpLzYyZ//1f6NfP7mqaJgW8iPjdgQNmbZnx4+HGG+2upunSNEkR8avKShPunTvD\n009rGQJf+ZKdCngR8RvLMvup7tkD//iHZsz4gxYbE5GQMHkyfPKJmRKpcLdfg2PwbrebpKQkEhIS\nyM7OrvO8devWERkZycKFC/1aoIiEh1mzzF6qixdDy5Z2VyPQQMB7vV4mTJiA2+1m8+bN5ObmsmXL\nlhOeN2nSJAYOHKhhGJEm6PXXYcoUcLvh/PPtrkZ+UG/AFxYWEh8fT8eOHYmKiiIzM5P8/Pzjznvy\nyScZNmwYbdu2DVihIhKa3n8fxo2DN96AuDi7q5Gj1RvwpaWlxMbG1nx2uVyUlpYed05+fj7jx48H\nzAMBEWkaNm6EoUNh/nzo1cvuauRY9Qb8yYT1nXfeyYwZM2qe9GqIRqRp2LbNrA6ZkwMDBthdjZxI\nvc+5Y2Ji8Hg8NZ89Hg8ul6vWOR9//DGZmZkA7Nmzh7feeouoqCgyMjKOu97kyZNrvk9LSyMtLc2H\n0kXELrt2QXo6PPgg3HCD3dU4S0FBAQUFBX65Vr3z4KuqqkhMTGTlypVER0dz8cUXk5ubS3Jy8gnP\nHzVqFIMHD2bIkCHH30jz4EUcYe9es1nH8OFw//12V+N8AZsHHxkZSU5ODunp6Xi9XsaMGUNycjKz\nZs0CICsrq1E3FZHwdPAgXH21eVP1vvvsrkYaojdZReSkHDr0445MWoIgeLRUgYgE1JEjZmXItm3h\n+efhtNPsrqjpUMCLSMBUVsKwYRAVBQsWaAmCYNNaNCISEFVVcMst4PXCK68o3MON/nWJyAlVVcGI\nEWZt99dfh9NPt7siOVUKeBE5jtcLo0aZZX8XLYIzz7S7ImkMBbyI1OL1mjXdd+4068s0b253RdJY\nCngRqeH1wq9//WO4t2hhd0XiC+3JKiKAGXMfORK++krh7hTqwYsIlZXwq19BebkZc9ewjDMo4EWa\nuCNH4KabzPBMfr4eqDqJhmhEmrDDh+H6682bqa+9pnB3GgW8SBN14IBZOKx1a8jL0zx3J1LAizRB\n+/bBz39uFg578UW9oepUCniRJmbnTrOee1qaWRWymVLAsfSvVqQJ2boVfvpTs75MdraW/HU6BbxI\nE7FhA1x+Ofz+9+ZQuDufRt5EmoC334bMTHjqKbP0rzQN6sGLONzLL5twf/llhXtTox68iIM9/jjM\nnAnLl0Nqqt3VSLAp4EUcqLoa7r4bli6FDz6A//ovuysSOyjgRRzm8GG49VbYvduEe+vWdlckdtEY\nvIiDfP01XHWV2T916VKFe1OngBdxiM2b4ZJLYMAAmD8fzjjD7orEbhqiEXGAZcvMcr+PPGKGZ0RA\nAS8S1iwLcnJg6lR49VXzIpPIDxTwImGqogImToT33oPVq6FzZ7srklCjgBcJQ19/bV5aatUKPvwQ\nzjnH7ookFOkhq0iYKSqCPn3MipCvv65wl7qpBy8SRl54wbzA9NRTcMMNdlcjoU4BLxIGKivhrrvA\n7YZ33oGUFLsrknCggBcJcaWlcOON5qWldevMuLvIydAYvEgIW7kSeveGa66BRYsU7nJq1IMXCUFe\nL0yfbsba582D/v3trkjCkQJeJMSUlZm3Uisr4aOPIDra7ookXGmIRiSErFwJF10EffvCihUKd/GN\nevAiIaCyEh56CJ5/3hw//7ndFYkTnFQP3u12k5SUREJCAtnZ2cf9fP78+aSmptK9e3cuu+wyNm3a\n5PdCRZxq2zazhsz69eZQuIu/NBjwXq+XCRMm4Ha72bx5M7m5uWzZsqXWOZ07d2bVqlVs2rSJP/7x\nj4wbNy5gBYs4hWXB3Llmid/MTFi8GM4/3+6qxEkaHKIpLCwkPj6ejh07ApCZmUl+fj7Jyck15/Tr\n16/m+759+7Jjxw7/VyriIHv2wPjxsGWLGXfv3t3uisSJGuzBl5aWEhsbW/PZ5XJRWlpa5/mzZ89m\n0KBB/qlOxIEWLzYbYHfoYGbJKNwlUBrswUdERJz0xd555x3mzJnDBx98cMKfT548ueb7tLQ00tLS\nTvraIuHuwAGzjszy5WbHJf36y4kUFBRQUFDgl2s1GPAxMTF4PJ6azx6PB5fLddx5mzZtYuzYsbjd\nblrXsRHk0QEv0pQsXw633QYDB8KmTVoBUup2bOd3ypQpjb5WgwHfu3dviouL2b59O9HR0eTl5ZGb\nm1vrnC+//JIhQ4Ywb9484uPjG12MiNOUl8Pvfme21Pv73+EXv7C7ImlKGgz4yMhIcnJySE9Px+v1\nMmbMGJKTk5k1axYAWVlZ/OlPf2L//v2MHz8egKioKAoLCwNbuUiIW7QIbr8dMjLgk0/Ua5fgi7As\nywrKjSIiCNKtRGy1a5fZSm/DBtNrv+IKuyuScOZLdmqpAhE/qa6GZ54xs2ISE81Yu8Jd7KSlCkT8\nYP16M689MlIbckjoUA9exAfl5XDHHWZ2zLhxsGqVwl1ChwJepBGqq+G55yA5Gb7/Hj79FEaPhmb6\nL0pCiIZoRE7RmjXmIWqzZvDGG2bHJZFQpP6GyEnyeGDECBg6FH77W1i9WuEuoU0BL9KAgwfhj3+E\nHj3M+jGffWaCXsMxEur0KypSh8pKePpp6NIFtm83M2WmToWzz7a7MpGTozF4kWNYFrz2GjzwgOmx\nL15sttETCTcKeJH/sCyzKNj995tZMk8+qbVjJLwp4EWAd9+FBx+EsjJ4+GHzIFVj7BLuFPDSpK1a\nBVOmwL//bQJ++HDzNqqIE+hXWZocy4K33zY9dY/HjLWPGAFRUXZXJuJfCnhpMqqrzYtJ06fD/v3w\nhz/AzTerxy7OpV9tcbyKCnjpJXjkETjjDLjvPrj+ejjtNLsrEwksBbw41v798H//B088YdaMeewx\nGDAATmGbYZGwpoAXxykuNlMc582DQYPMzkqaxy5NkQJeHKG62ux7+uSTsG6d2eD6n/+EmBi7KxOx\njwJewtrevWbZ3meeMUsITJwIr74KzZvbXZmI/RTwEnYsC957D559Ft5802xqPW8e9O2r8XWRo2nT\nbQkbu3bBiy/CnDnmLdNx48z89TZt7K5MJHB8yU714CWkff+9mbv+wgvw/vtmCYE5c6BfP/XWRRqi\nHryEnOpqs4TA/PmwcCH07AkjR5q562edZXd1IsGlHryEPcuCwkLIy4NXXoHzzjPrwmzYALGxdlcn\nEp4U8GKb6moT6q++atZfP+MMuOkmWLoUuna1uzqR8KeAl6CqqDBL877+OuTnwznnwLBh5nP37hpX\nF/EnBbwE3O7d4Habh6XLl0NiIlx3HaxYAUlJdlcn4lx6yCp+5/Wat0mXLoW33oItW+DKK2HwYLjm\nGrjgArsrFAkfvmSnAl58ZlmwbZvpka9caY6YGLPd3aBB8NOfmvF1ETl1CngJKsuCf/3LTGUsKIB3\n3jG99v79zTFgAERH212liDMo4CWgqqpg40ZYvdq8bPTee+Zh6OWXm6GXtDRISNADUpFAUMCL31iW\n2cZu3TpYu9YcRUXQsSNceqk5fvYz6NRJgS4SDAp4aRTLgu3bzctERUXm+Ogj87M+fcziXZdcYr5v\n1crWUkWaLAW8NGj/fti82ayR/sOxcaNZYjc1FXr1Mpti9OoFLpd65yKhQgEvgHnQ6fHAF1/A55/D\nZ5+ZY8sWOHjQvB3arRukpJivqalmSQARCV0BDXi3282dd96J1+vltttuY9KkScedM3HiRN566y1a\ntGjB3Llz6dmzp1+LlB8dOAD//jeUlPx4bNsGW7ea4ZbzzoMuXcyRlGReKura1aznol65SPgJ2GJj\nXq+XCRMmsGLFCmJiYujTpw8ZGRkkJyfXnLNkyRK2bt1KcXExa9euZfz48axZs6ZRxYSzgoIC0tLS\nGv33LQv27TNrnu/aBaWl5tixw/TKPR4T7FVV0KGDecj5w5GWBvHx0LkztGjhtybV4mv7QpmT2wZq\nX1NWb8AXFhYSHx9Px44dAcjMzCQ/P79WwC9atIiRI0cC0LdvX8rLyykrK6Ndu3aBqzoEHftLduQI\nlJeb0N63z2wtt2eP+bp7tzm+/hrKyszx9dcmnNu3N0dMjDlSUuDqq00PvEMHOPdce3riTv6PyMlt\nA7WvKas34EtLS4k9aq1Wl8vF2rVrGzxnx44dIR/wlmV6w0eO/Hh8/705Dh82Xw8dMt8fOgTffffj\n8e235jh40AyZHDhgHmDm5cE335hgr6oyM0/atDGhfO65ZvikTRto29YMn7RtC+3a/Xiceabd/1RE\nxEnqDfiIk+wqHjs+VNffGzz42L9X9/d1HdXV5jj6e6/3x69VVebrD99XVUFlpTmqqsxqhhUV5nOz\nZuYV+tNPN+F6xhlms+bmzc3nFi3M0bw5tGxpjhYtzMyTtm3N13POMV9feQXuvht+8hNo3dr8HY15\ni4itrHp8+OGHVnp6es3nadOmWTNmzKh1TlZWlpWbm1vzOTEx0frqq6+Ou1ZcXJwF6NChQ4eOUzji\n4uLqi+l61duD7927N8XFxWzfvp3o6Gjy8vLIzc2tdU5GRgY5OTlkZmayZs0aWrVqdcLhma1bt9Z3\nKxER8bN6Az4yMpKcnBzS09Pxer2MGTOG5ORkZs2aBUBWVhaDBg1iyZIlxMfH07JlS5577rmgFC4i\nIvUL2otOIiISXM0CfQO3201SUhIJCQlkZ2cH+nYB5/F4uPLKK7nwwgtJSUnhiSeeAGDfvn0MGDCA\nLl268Itf/ILy8nKbK/WN1+ulZ8+eDP7Pk3Enta+8vJxhw4aRnJxM165dWbt2rWPaN336dC688EK6\ndevG8OHDOXLkSFi3bfTo0bRr145u3brV/Fl97Zk+fToJCQkkJSWxbNkyO0o+JSdq3z333ENycjKp\nqakMGTKEb775puZnp9y+Ro/en4SqqiorLi7OKikpsSoqKqzU1FRr8+bNgbxlwO3atctav369ZVmW\ndfDgQatLly7W5s2brXvuucfKzs62LMuyZsyYYU2aNMnOMn326KOPWsOHD7cGDx5sWZblqPbdeuut\n1uzZsy3LsqzKykqrvLzcEe0rKSmxOnXqZH3//feWZVnWjTfeaM2dOzes27Zq1SqrqKjISklJqfmz\nutrz6aefWqmpqVZFRYVVUlJixcXFWV6v15a6T9aJ2rds2bKauidNmuRT+wIa8KtXr641C2f69OnW\n9OnTA3nLoPvlL39pLV++vNbsoV27dlmJiYk2V9Z4Ho/H6t+/v/X2229b1157rWVZlmPaV15ebnXq\n1Om4P3dC+/bu3Wt16dLF2rdvn1VZWWlde+211rJly8K+bSUlJbUCsK72HDvLLz093frwww+DW2wj\nHNu+oy1cuNC65ZZbLMtqXPsCOkRzopegSktLA3nLoNq+fTvr16+nb9++td7ebdeuHWVlZTZX13j/\n8z//w8yZM2nW7MdfD6e0r6SkhLZt2zJq1Cguuugixo4dy3fffeeI9p177rncfffddOjQgejoaFq1\nasWAAQMc0baj1dWenTt34nK5as5zQt7MmTOHQYMGAY1rX0AD/mRflApH3377LUOHDuXxxx/n7LPP\nrvWziIiIsG37m2++yfnnn0/Pnj3rXOAonNtXVVVFUVERt99+O0VFRbRs2ZIZM2bUOidc27dt2zb+\n9re/sX37dnbu3Mm3337LvHnzap0Trm2rS0PtCee2Tp06ldNPP53hw4fXeU5D7QtowMfExODxeGo+\nezyeWv8HCleVlZUMHTqUESNGcN111wGmJ/HVV18BsGvXLs4//3w7S2y01atXs2jRIjp16sTNN9/M\n22+/zYgRIxzTPpfLhcvlok+fPgAMGzaMoqIiLrjggrBv30cffcSll15KmzZtiIyMZMiQIXz44YeO\naNvR6vpdPDZvduzYQUxMjC01+mru3LksWbKE+fPn1/xZY9oX0IA/+kWpiooK8vLyyMjICOQtA86y\nLMaMGUPXrl258847a/48IyOD559/HoDnn3++JvjDzbRp0/B4PJSUlLBgwQKuuuoqXnzxRce074IL\nLiA2NpYvvvgCgBUrVnDhhRcyePDgsG9fUlISa9as4fDhw1iWxYoVK+jatasj2na0un4XMzIyWLBg\nARUVFZSUlFBcXMzFF19sZ6mN4na7mTlzJvn5+Zx51AJVjWqfn54T1GnJkiVWly5drLi4OGvatGmB\nvl3Avffee1ZERISVmppq9ejRw+rRo4f11ltvWXv37rX69+9vJSQkWAMGDLD2799vd6k+KygoqJlF\n46T2bdiwwerdu7fVvXt36/rrr7fKy8sd077s7Gyra9euVkpKinXrrbdaFRUVYd22zMxMq3379lZU\nVJTlcrmsOXPm1NueqVOnWnFxcVZiYqLldrttrPzkHNu+2bNnW/Hx8VaHDh1q8mX8+PE1559q+/Si\nk4iIQwV8PtHzAAAAMklEQVT8RScREbGHAl5ExKEU8CIiDqWAFxFxKAW8iIhDKeBFRBxKAS8i4lAK\neBERh/p/pxJ1p0VO2e0AAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0x109b7af10>"
       ]
      }
     ],
     "prompt_number": 66
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We combine that with our event design matrix:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "X, cons = fmrid.stack_designs((X_exper, cons_exper),\n",
      "                              (drift, {}))\n",
      "plt.figure(figsize=(20, 10))\n",
      "plt.imshow(X)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 67,
       "text": [
        "<matplotlib.image.AxesImage at 0x1073024d0>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAGIAAAJNCAYAAADKwKDEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGYVJREFUeJztnV9sU9cdx7+mcWJIWg8YsUscFJTEc5xQE0TLVK0ibeZU\nqpQogygCtmIB7aTxslZVC+oTmlRihKaV7s/LlK7eOi3wRKIqi1hUwegYZVUY20hRosoRzj+LLUuW\nEBzH5uwBYdXxdXJzcm3/ru/vI1nFN+eec5JPr78+5957rkkIIcDknHW57gDzCBZBBBZBBBZBBBZB\nBBZBBM1F9PX1weVyobq6GmfOnNG6+vxFaEgsFhOVlZUiGAyKaDQqPB6PGBwc1LKJvEXTI+LGjRuo\nqqpCRUUFzGYzDhw4gO7ubi2byFsKtKxsbGwM5eXlifcOhwOff/55UpmioiJEo1EtmyVNYWEhFhYW\nViynqQiTybRimWg0mpA1MzMDq9W6bPlslbl79y4A4NSpUzh16hS2bdumSXuhUGjZco/RVERZWVlS\nw6FQCA6HI6XczMwMACASiaCoqAgWi0XLbuSUSCSChYUFRCKRVe2naUbs3r0bw8PDGBkZQTQaxfnz\n59HS0pJSzmq1wmq1wmKx5JUEALBYLInfbaWj5utoekQUFBTgF7/4BV5++WXE43EcO3YMNTU1acsX\nFRWtWGc2ywBAQ0MDAGDdutT/Rx8+fKh5e48xCZHdaXCTyZQU6FR4nBGPqaioSCnzWMRqCIVCUPMn\n5pE1EVgEETTNiHxCzVdxLeEjgggsgggsgggsgggc1ml44oknUrbF4/GMtcdHBBFYBBFYBBE4I9LA\nAzqDwiKIwCKIwCKIwGGdBqUBnVKAa3VejY8IIrAIIrAIIrAIIuRVWD/55JMp22ZnZ6Xq4pG1QWER\nRGARRMirjPB4PCnbPvvsM6m6lC65zCR8RBCBRRCBRRCBRRAhr8J669atmtXFYW1QWAQRWAQRWAQR\n8iqsFxcXNatLKaz5VKkBYBFEYBFEyKuMuH//vmZ18YDOoLAIIrAIIrAIIuRVWK92jaTl4LA2KCyC\nCCyCCHmVEbFYTLO61E76adZexmpmVgWLIAKLIAKLIEJehbWWC3aqHdAtDXDZPvARQQQWQQQWQQQW\nQYS8Cmst4dlXg8IiiMAiiJBXGZGLAZ1m7WW1NSYtLIIILIIILIIIeRXWS58TtBZ4dRqDwiKIwCKI\nwCKIkFdhffv27ZRtGzdulKpLKayXe1rjWuEjgggsgggsggh5lRFzc3Mp22QzQmlt8EzCRwQRWAQR\nWAQRWAQR8iqs169fr1ldfKOKQZESEQqF8OKLL6K2thZ1dXX44IMPAABTU1Pwer1wOp1oamrC9PS0\npp3NZ6REmM1m/OxnP8Pt27dx/fp1/PKXv8SXX34Jv98Pr9eLoaEhNDY2wu/3a93fvEUqI+x2O+x2\nOwCgpKQENTU1GBsbQ09PD65cuQIA8Pl8aGho0EyGmqfpFhcXa9IWoMOMGBkZwc2bN7Fnzx6Ew2HY\nbDYAgM1mQzgcXnMHjcKaRMzNzWH//v04d+5cytNMTCZT1s/76hnpr6+Li4vYv38/Xn31VbS2tgJ4\ndBRMTk7CbrdjYmICpaWlivvOzMwk/l1UVASLxSLbDXJEIhEsLCysej+pI0IIgWPHjsHtduONN95I\nbG9paUEgEAAABAKBhKClWK3WxCufJACAxWJJ+v3UInVE/OUvf8HHH3+MZ555BvX19QCAjo4OnDx5\nEu3t7ejs7ERFRQUuXLggU70i1dXVKdvu3LmT9L6wsFCz9rJ97auUiO985ztpTxH29/evqUNGhUfW\nRGARRGARRNDN7KvL5UrZtjSstQxYXnnAoLAIIrAIIugmIwoKVu6qlnNbShmRyUEeHxFEYBFEYBFE\nYBFE0E1Yz8/P57oLGYWPCCKwCCKwCCKwCCJwWKdB6Toq0tc1MdrAIojAIoigm4zQ8sm8alB77atW\nucFHBBFYBBFYBBFYBBF0E9ZarQKjFl731aCwCCKwCCLoJiOyjdoVzLSCjwgisAgisAgisAgicFin\ngZeSMygsgggsgggsggi6CevZ2dmstqe7ZYIYbWARRGARRNBNRty6dStlWyYHXZwRBoVFEIFFEIFF\nEEE3YT0+Pp6yrby8PGPtqQ3mpaG+dC1atfARQQQWQQQWQQTdZITZbM5qe7w6jUFhEURgEURgEUTQ\nTViXlJRktT2efTUoLIIILIIILIIIJMJ6y5YtSe/v3buXUibbT17hlQcMCosgAosgAomM8Hg8Se+V\nnlOU9fsVePbVmLAIIrAIIrAIIpAI682bN69YJtsPoOUBnUFhEURgEURgEUQgEdaRSCTXXUiBT5Ua\nFBZBBBZBBBIZMTc3l+supMDLjRoUFkGENYmIx+Oor69Hc3MzAGBqagperxdOpxNNTU2Ynp7WpJNG\nYE0izp07B7fbnfh+7ff74fV6MTQ0hMbGRvj9fk06aQSkRYyOjqK3txevvfZa4vHyPT098Pl8AACf\nz4eLFy+qqmthYSHpRYF169ZJvaTbk93xzTffxNmzZ5MaD4fDsNlsAACbzYZwOCzdMaMhJeKTTz5B\naWkp6uvrE0fDUkwmU9bPIegZqXHEtWvX0NPTg97eXkQiEfzvf//Dq6++CpvNhsnJSdjtdkxMTKC0\ntFRx/5mZmcS/i4qK5HpOlEgkIvXxKnVEnD59GqFQCMFgEF1dXXjppZfwu9/9Di0tLQgEAgCAQCCA\n1tZWxf2tVmviZbFY8PDhw6QXBWQzwWKxJP1+qtvTotOPP4JOnjyJP/3pT3A6nfj0009x8uRJLao3\nBGue4ti7dy/27t0LANi0aZPixWHMyvDImggsgggkZl/TfQXOJWq/evPlNHkGiyACiyACiyACibCm\niNobVTis8wwWQQQWQQQSGSG7RGcm4ZsZDQqLIAKLIAKLIAKJsB4eHs51F1LI+s2TWW2NSQuLIAKL\nIAKLIAKJsL5z507S+8eXbeYSpWWJMjna5iOCCCyCCCyCCCQyYnFxMdddSIHP0BkUFkEEFkEEFkEE\nEmFdXFyc6y6kwEvJGRQWQQQWQQQSGbFhw4astqfmTlalz36+nMYAsAgisAgisAgikAjrbD8iuba2\ndsUyfIbOoLAIIrAIIrAIIpAI62yvsbp9+/YVy6idaeXZ1zyDRRCBRRCBREZkGzWf63w5jUFhEURg\nEURgEUQwZFireXAIP+zJoLAIIrAIIhgyIx48eLBiGbUDOj5Dl2ewCCKwCCKwCCIYMqxjsdiKZXj2\n1aCwCCKwCCKwCCIYMqzVPNmLw9qgsAgisAgiGDIj1MBn6AwKiyACiyACiyACh3Ua1Ia10k2PUu1p\nUguzZlgEEVgEEVgEEQwZ1vfu3VuxjNpg5tnXPENaxPT0NNra2lBTUwO3243PP/8cU1NT8Hq9cDqd\naGpqwvT0tJZ9zWukRfz4xz/GK6+8gi+//BL/+Mc/4HK54Pf74fV6MTQ0hMbGRvj9fi37mtdIZcTM\nzAyuXr2KQCDwqJKCAlitVvT09ODKlSsAAJ/Ph4aGBpIybt26tWIZXVz7GgwGsWXLFhw5cgS7du3C\n66+/jvv37yMcDicewmGz2RAOhzXppBGQEhGLxTAwMIDjx49jYGAAxcXFKf/nm0ymrD9CTM9IiXA4\nHHA4HHj22WcBAG1tbRgYGIDdbsfk5CQAYGJiAqWlpYr7z8zMJF6RSESy6zSJRCJJv59apETY7XaU\nl5djaGgIANDf34/a2lo0NzcnciMQCKC1tVVxf6vVmnhZLBaZLpDFYrEk/X5qkR7Q/fznP8f3v/99\nRKNRVFZW4je/+Q3i8Tja29vR2dmJiooKXLhwQbb6jDI1NbVimWxfTiMtwuPx4G9/+1vK9v7+/jV1\nyKjwyJoILIIIhpz0U7M2OD8/wqCwCCKwCCKwCCLoOqyXBqUQQtV+JSUlq64bUL9euAx8RBCBRRCB\nRRCBRRBB12FdXl6e9P7u3buq9issLFyxDD/Iw6CwCCKwCCLoOiPq6uqS3qvNCDX3NPAN7waFRRCB\nRRCBRRBB12GtZhZVCdkHeWTy0cp8RBCBRRCBRRCBRRBB12E9Pz+fsbqVQljNjKzsSJuPCCKwCCKw\nCCLoOiPu37+fsbrVDt549jXPYBFEYBFEYBFE0HVYLy4uZqxuXaw8wGgPiyACiyCCrjNCzXMgZOGb\nGQ0KiyACiyACiyCCrsNa7c2LMqg9Q7d0G5+h0zksgggsgggsggi6DutMIns5jXR7mtTCrBkWQQQW\nQQRdZ8TCwkLG6pY9Q8cDOp3DIojAIojAIoig67C+fft2xurmlQcMCosgAosggq4zIhgMJr3funWr\nZnXLnqGTbk+TWpg1wyKIwCKIwCKIoOuwzvblNGoCXPaMHR8RRGARRGARRGARRNB1WMsuJacGtTel\n8OU0eQaLIAKLIIKuM2L9+vUrlnnyySdTts3Ozq64n9Jnf0FB6p+LL6fJM1gEEVgEEVgEEXQd1krh\nuRSPx5Oy7bPPPltxP9nZVw5rnSMtoqOjA7W1tdixYwcOHTqEhYUFTE1Nwev1wul0oqmpCdPT01r2\nNa+REjEyMoJf//rXGBgYwD//+U/E43F0dXXB7/fD6/ViaGgIjY2N8Pv9Wvc3b5ES8dRTT8FsNmN+\nfh6xWAzz8/PYunUrenp64PP5AAA+nw8XL17UtLP5jFRYb9q0CW+99Ra2bduG9evX4+WXX4bX60U4\nHIbNZgMA2Gw2hMNhTTu7FDUzn7LXOsk+Rjmrp0q/+uorvP/++xgZGcH4+Djm5ubw8ccfJ5UxmUya\nXaBrBKSOiC+++ALPP/88Nm/eDADYt28f/vrXv8Jut2NychJ2ux0TExMoLS1V3H9mZibx76KiIlgs\nFplukGR+fl5qFX+pI8LlcuH69et48OABhBDo7++H2+1Gc3MzAoEAACAQCKC1tVVxf6vVmnjlkwQA\n2LBhA775zW8mXmqROiI8Hg8OHz6M3bt3Y926ddi1axd++MMfYnZ2Fu3t7ejs7ERFRQUuXLggU72m\nyC5JqpQHanJD9uNYemT9zjvv4J133knatmnTJvT398tWaWh4ZE0EFkEEFkEEXc++qkH2YR98qtSg\nsAgisAgi5H1GRCIRqf34/giDwiKIwCKIwCKIkPdhHYvFpPaTXXmAw1rnsAgisAgisAgi5H1Yy65O\nwGFtUFgEEVgEEfI+I2RRe3M7Z0SewSKIwCKIwCKIkPdhreWATs02DmudwyKIwCKIwCKIkPdhfffu\nXan91F77yiPrPINFEIFFECHvM0LpgVAbN25ccT/ZM3SyD/bgI4IILIIILIIILIIIeR/Wc3NzKdu0\nDOulgzwOa53DIojAIoiQ9xmhZv1wJfhmRoPCIojAIojAIoiQk7BeGnDxeDxjbRUXF0vtJzv7qma9\ncsX2pPZiNIdFEIFFEIFFECEnYV1dXZ30/s6dOxlrq7CwUGo/nn01KCyCCCyCCDnJCJfLlfQ+kxkh\nOxuq9vkRfDlNnsEiiMAiiMAiiJCTsJadoZRBeq1uHtAZExZBBBZBBBZBhJyEtcxz2rKN7KlSs9ks\n157UXozmsAgisAgicEakQe3sKw/o8gwWQQQWQQQWQYSchLXso42ziVLoqhnkSV++I7UXozksgggs\nggg5yYiHDx/motlVoXYFs6Xb+P4InbOsiKNHj8Jms2HHjh2JbVNTU/B6vXA6nWhqasL09HTiZx0d\nHaiurobL5cKlS5cy1+s8ZFkRR44cQV9fX9I2v98Pr9eLoaEhNDY2wu/3AwAGBwdx/vx5DA4Ooq+v\nD8ePH9fFRxAVlhXxwgsvpCwg0tPTA5/PBwDw+Xy4ePEiAKC7uxsHDx6E2WxGRUUFqqqqcOPGjQx1\nO/9YdbKEw2HYbDYAgM1mQzgcBgCMj4/j29/+dqKcw+HA2NiYRt3MPrKX0+TkDJ3JZFr2uiHZa4qM\nyKqPCJvNhsnJSdjtdkxMTKC0tBQAUFZWhlAolCg3OjqKsrIyxTr+9a9/Jf79eP98YWJiApOTk6ve\nb9VHREtLCwKBAAAgEAigtbU1sb2rqwvRaBTBYBDDw8N47rnnFOuoq6tLvPJNxNNPP436+vrESy3L\nHhEHDx7ElStX8O9//xvl5eX4yU9+gpMnT6K9vR2dnZ2oqKjAhQsXAAButxvt7e1wu90oKCjAr371\nK/5oWgXLivjDH/6guL2/v19x+7vvvot333137b0igOza4Hw5jc5hEURgEUTIyezr7OxsLppdFWoH\ndEszgc/Q6RwWQQQWQQQWQYSchPWtW7eS3steL5pJ1D5Geekgj8Na57AIIrAIIuQkI8bHx5Pel5eX\n56Iby6J20m/pgI4n/XQOiyACiyACiyBCTsJaNtCyiezlNEVFRXLtSe3FaA6LIAKLIAKLIEJOwrqk\npCQXza4K2RtVePZV57AIIrAIIuQkIywWy7I/37JlS8q2e/fuZao7isheTsMDOp3DIojAIojAIohA\n4smMS/F4PCnb0t2TkSlkB3Qc1jqHRRCBRRCBRRAhJ2G90t2mmzdvzlJP0sMja4PCIojAIoiQk4xY\niUgkkusuSGfESjPLaduT2ovRHBZBBBZBBBZBBJJhPTc3l+suSM++cljrHBZBBBZBBJIZsbCwkOsu\nqH7Y09IB3fr166Xa4yOCCCyCCCyCCCyCCCTDmsIq+7IrD/CATuewCCKwCCKwCCKQDGshRK67wCNr\no8IiiMAiiEAyIyigZiUapW0bNmyQao+PCCKwCCKwCCKwCCKQDOt4PJ7rLqiefV16PwQP6HQOiyAC\niyACiyACybAeHh7OdRcUw1pp29LZVx5Z6xwWQQQWQQSSGXHnzp2UbTabLat9UJMHShQXF8u1J7UX\nozksgggsgggsgggkw3pxcTHXXVA9+6pUTqo9TWph1syyIo4ePQqbzYYdO3Yktr399tuoqamBx+PB\nvn37MDMzk/hZR0cHqqur4XK5cOnSpcz1Og9ZVsSRI0fQ19eXtK2pqQm3b9/GrVu34HQ60dHRAQAY\nHBzE+fPnMTg4iL6+Phw/fpzE5fV6YVkRL7zwAjZu3Ji0zev1Jj4X9+zZg9HRUQBAd3c3Dh48CLPZ\njIqKClRVVeHGjRtSnSouLk55ZRuTyZTyeuKJJ1JeZrM56VVSUpL0UsuaMuLDDz/EK6+8AuDRYy8d\nDkfiZw6HA2NjY8vur2Y5oFgstmIZNfWoXXro8uXLWSvzdaS/Nb333nsoLCzEoUOH0pZJt3bf41x5\n/MdZ7i6beDyu+G3l6ywsLKx4p46aMsCjP2BDQ4N0mcuXLydeK9XzdaSOiI8++gi9vb34/e9/n9hW\nVlaGUCiUeD86OoqysjLF/a1WK6xWKywWi/StTlRpaGjAqVOnEv9Vy6pF9PX14ezZs+ju7k76I7a0\ntKCrqwvRaBTBYBDDw8N47rnnVlu9cRHLcODAAfH0008Ls9ksHA6H6OzsFFVVVWLbtm1i586dYufO\nneJHP/pRovx7770nKisrxbe+9S3R19enWOfevXsFAMO89u7du9yfOIFJULgrhOGRNRVYBBFyIqKv\nrw8ulwvV1dU4c+aMYplQKIQXX3wRtbW1qKurwwcffJC2vng8jvr6ejQ3Nyv+fHp6Gm1tbaipqYHb\n7cb169dTynR0dKC2thY7duzAoUOHsLCwoDjFMzU1Ba/XC6fTiaamJkxPT696KkgRVUmiIbFYTFRW\nVopgMCii0ajweDxicHAwpdzExIS4efOmEEKI2dlZ4XQ6FcsJIcRPf/pTcejQIdHc3Kz488OHD4vO\nzk4hhBCLi4tieno66efBYFBs375dRCIRIYQQ7e3t4qOPPhJ//vOfxcDAgKirq0uUffvtt8WZM2eE\nEEL4/X5x4sQJxXKXLl0S8XhcCCHEiRMnxIkTJ5b9u2T9iLhx4waqqqpQUVEBs9mMAwcOoLu7O6Wc\n3W7Hzp07ATx6pGZNTQ3Gx8dTyo2OjqK3txevvfaa4t2oMzMzuHr1Ko4ePQrg0VS21WpNKvPUU0/B\nbDZjfn4esVgM8/PzKCsrU5zi6enpgc/nAwD4fD5cvHhxVVNB6ci6iLGxMZSXlyfeq5kKGRkZwc2b\nN7Fnz56Un7355ps4e/Zs2vMCwWAQW7ZswZEjR7Br1y68/vrrmJ+fTyqzadMmvPXWW9i2bRu2bt2K\nb3zjG/jud7+rWF84HE5cyGCz2RAOh5ftO5A8FZSOrItY6ZEFS5mbm0NbWxvOnTuXMon2ySefoLS0\nFPX19WnvzY7FYhgYGMDx48cxMDCA4uJi+P3+pDJfffUV3n//fYyMjGB8fBxzc3NJswbL/S4r/T5q\npoKAHIhYOhUSCoWSJgu/zuLiIvbv348f/OAHaG1tTfn5tWvX0NPTg+3bt+PgwYP49NNPcfjw4aQy\nDocDDocDzz77LACgra0NAwMDSWW++OILPP/889i8eTMKCgqwb98+XLt2TbFPNpsNk5OTAICJiQmU\nlpam/V2VpoLSkXURu3fvxvDwMEZGRhCNRnH+/Hm0tLSklBNC4NixY3C73XjjjTcU6zp9+jRCoRCC\nwSC6urrw0ksv4be//W1SGbvdjvLycgwNDQF49PSu2trapDIulwvXr1/HgwcPIIRAf38/3G63Ypst\nLS0IBAIAgEAgoPg/CJB+Kigty0Z5hujt7RVOp1NUVlaK06dPK5a5evWqMJlMwuPxJKZT/vjHP6at\n8/Lly2m/Nf39738Xu3fvFs8884z43ve+l/KtSQghzpw5I9xut6irqxOHDx8W0Wg0ZYrnww8/FP/5\nz39EY2OjqK6uFl6vV/z3v/9d9VSQEjzFQQQeWROBRRCBRRCBRRCBRRCBRRCBRRCBRRDh/3GYgA4X\nX7bDAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x109af5a50>"
       ]
      }
     ],
     "prompt_number": 67
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cons"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 68,
       "text": [
        "{'constant_0': array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_bottle_0': array([ 1., -0., -0., -0., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_cat_0': array([ 0.,  1., -0.,  0., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_chair_0': array([-0., -0.,  1., -0., -0.,  0., -0.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_face_0': array([-0.,  0., -0.,  1., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_house_0': array([-0., -0.,  0., -0.,  1.,  0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_scissors_0': array([-0., -0., -0., -0.,  0.,  1., -0.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_scrambled_0': array([ 0., -0., -0.,  0., -0., -0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_shoe_0': array([-0., -0.,  0.,  0., -0., -0., -0.,  1.,  0.,  0.,  0.,  0.,  0.])}"
       ]
      }
     ],
     "prompt_number": 68
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Fit the design, the by-hand way:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Reshape the data for estimation\n",
      "slice_size = np.prod(data.shape[:-1])\n",
      "Y = np.reshape(data, (slice_size, n_trs)).T\n",
      "Y.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 69,
       "text": [
        "(121, 163840)"
       ]
      }
     ],
     "prompt_number": 69
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "B = np.linalg.pinv(X).dot(Y)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 70
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "B.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 71,
       "text": [
        "(13, 163840)"
       ]
      }
     ],
     "prompt_number": 71
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Let's make a contrast map.  But, what are our contrasts?"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cons"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 72,
       "text": [
        "{'constant_0': array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_bottle_0': array([ 1., -0., -0., -0., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_cat_0': array([ 0.,  1., -0.,  0., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_chair_0': array([-0., -0.,  1., -0., -0.,  0., -0.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_face_0': array([-0.,  0., -0.,  1., -0., -0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_house_0': array([-0., -0.,  0., -0.,  1.,  0., -0., -0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_scissors_0': array([-0., -0., -0., -0.,  0.,  1., -0.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_scrambled_0': array([ 0., -0., -0.,  0., -0., -0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.]),\n",
        " 'type_shoe_0': array([-0., -0.,  0.,  0., -0., -0., -0.,  1.,  0.,  0.,  0.,  0.,  0.])}"
       ]
      }
     ],
     "prompt_number": 72
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "contrast = cons['constant_0']\n",
      "contrast"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 78,
       "text": [
        "array([ 1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  0.,  0.,  0.])"
       ]
      }
     ],
     "prompt_number": 78
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "con_data = contrast.dot(B)\n",
      "con_data.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 79,
       "text": [
        "(163840,)"
       ]
      }
     ],
     "prompt_number": 79
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "con_img = np.reshape(con_data, img.shape[:-1])\n",
      "con_img.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 80,
       "text": [
        "(40, 64, 64)"
       ]
      }
     ],
     "prompt_number": 80
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.imshow(con_img[..., 30])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 82,
       "text": [
        "<matplotlib.image.AxesImage at 0x109e6bc50>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAWwAAADwCAYAAAAkTF41AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnW9sleX5x78H7cZUxD/tOa2cauuklP6hPayKL2DialHn\nYBAM02yMAO4Fy5I5fDF9pdsLrVsWg84XxrikmwnTLXMjy+yE2KLgkE16AkUmii0WPG0pf1QQBrbn\n98Jw0p7ne7X30/aU312+n8SkvXie+/9z+fT53td1R9LpdBpCCCH+3zPlQjdACCGEG3LYQgjhCXLY\nQgjhCXLYQgjhCXLYQgjhCXLYQgjhCaN22M3NzSgvL8fMmTPx5JNPjmebhBBCECKj2Yfd39+PWbNm\nYcuWLZgxYwZuvvlmbNy4EbNnz85FG4UQQmCUb9g7d+7ETTfdhJKSEuTl5eG+++7D3/72t/FumxBC\niEFcOpqbDh8+jOLi4szv8Xgcb7/99pBrCgoK0NfXN7bWCSHERUZ+fj6OHDlC/21UDjsSiYx4TV9f\nH1asWAEAaG9vR1VVFcJ8fZkyJfjyz+odGBhwLpPBymR19/f3Z34+358whOm7y/ha11n1jHTt+T6x\n66z2DB6T4a517c9Y7x/cnz179qC6utq833XdWHWzcWa2MH1n959fi4P7wwhTj2vdYdbsaDJcDH6O\nRjPHo2G4MXapZ7j69+7di8rKSgDh1nH2WvzTn/5k1jEqhz1jxgx0dXVlfu/q6kI8Hg9c197eDgDo\n7e1Fb28vCgoKRlOdEEJMWnp7e8036mxG5bDr6urw/vvvo7OzE9dddx1eeuklbNy4MXDd+f97tre3\nIxqNjvn/jkIIMdmIRqOIRqOZ3999913z2lE57EsvvRS//e1vceedd6K/vx9r164ddofI4MZk88UX\nX1D7V7/6VVpvNqdOnQrYwnwWYGV+7WtfC9g+++yzzM/5+fmZzwGuf/qcPXs2YLP+LJ86dWrAxv5s\nYzarTHbt4LEvKChAOp1GXl5e4LpLLrmEltnb2xuwXXnllQEbmw/2OQXgfWfzMdIYx+PxzBpia+nz\nzz93KjPMJ5Fz584FbGw8rfXJnoXzfS8uLs78zObyK1/5Ci3TdY0w2+nTp2mZrJ8MVubg8YzFYpn2\nsfazeti6CfMiyMpk68vyS2xMzvfhqquuwv/+9z8AvD9sLQDcN1iMymEDwN133427777b6drhHLaP\nTMZPO5NtjgoLCy90E8aVoqKiC92EcScWi13oJowr+fn5Oa9DkY5CCOEJcthCCOEJcthCCOEJcthC\nCOEJoxYdXchWqNkOALbLw8LaWZCNpeyzHSHsWqYEW3UzNZhdG2Zzvuu1rmo/4N53ZrPmiEWyMsX9\nxIkTzmWyvjNxiu1GOXnypHOZDLYbxlpLbN7ZHA3eXTTcvQDfwcDKZLterD66BqCx9WHtamDrm607\n13EHeD/ZfLB6wuwKc90NY+2EYTs6XNtp+ZAw46Q3bCGE8AQ5bCGE8AQ5bCGE8AQ5bCGE8IScio7Z\nYoCrqGLBPs6HCaV2hbXTNdOW1aYwGeOYAHQ+5HWk6yxBi9XFBBQmaFntvOqqq5yuZeNx7bXX0jIZ\nTBRi9Vjzzubzsssuc6rHijBk88nEVTYfVjtdxVkm7LL+AHzsXefYWkvMzgRfNu5szVrXsnFi94dJ\nx+AaVh+m7wxWjxXqH8Zf6Q1bCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE8IaeiYzZhhCIm6rBIIdeoQqt+\n16gtq52uAiG738q5y3CNWLP6zkQpdsoFa7vVTiZaMpGL2VikolUmmyMm4Jw5c4aWydrPhDd2nRWR\nyead2Zi4agnYrP1MiGRCvZVT+fLLL3eqn615Jm4CfC275iy3ov1Y/Uy4Y6JfmOedRVu75kYHeN/Z\nM8Pmkq05wF3IBPSGLYQQ3iCHLYQQnjCmTyIlJSW48sorcckllyAvLw87d+4cr3YJIYTIYkwOOxKJ\noLW1Fddcc814tUcIIYTBmEXHsZ6EbokQzO4q3FkiBBMSXO+3IrSY3fWg0DAHpzLYeFhlMsGDiYGu\noh8AHDx4MGCLx+MBW5j/oTOxhgmUYcRmhms6USZIAXzsWSpVJsZZzwy733UtWqlQ2fp2fQ6OHTtG\ny2QCNhsP15TAgHvfw0RksrXM2sRs1rxbQmw2V1xxRcBm9d0SIxlj+oYdiURwxx13oK6uDs8///xY\nihJCCDECY3rD3r59O4qKinDkyBE0NDSgvLwcCxYsGK+2CSGEGMSYHPb5xDgFBQVYtmwZdu7cOcRh\nt7e3Z36ORqMoLCwcS3VCCDHp6OnpQW9vr9O1o3bYn3/+Ofr7+zFt2jScOnUKr732Gh599NEh11RV\nVY22eCGEuCiIxWJDjsHbu3evee2oHXZPTw+WLVsG4EvB4vvf/z4WLVo02uKEEEKMwKgddmlpKZLJ\n5LDXuOwgsfL4uubHdT1I1iqTXcsUd2unBOsjU9yZsm/tamD9ZNcyxTo/P5+WycZ5+vTpAdunn34a\nsB06dIiWyULbXXMgW4o5U/yPHz8esLFwYBaGDdjrweU6q50s5JyFPbN2fvLJJ7RMdj+bd9c85ABf\nn6yfrG5r94XrLijXegD3vPasn1bqBNYmltKAtWnGjBm0zJ6eHqf72RxZsPQDFop0FEIIT5DDFkII\nT5DDFkIIT5DDFkIIT8hpPuxs0WDatGmBa6zwV9cwdGazhA0W9szqYSGt1uGdTBRyFQ2tdl599dVO\nbbKEJga7tqKiImBjoqMl6rC5Y4exMqEnjFDExC9r7BisLtec0pYwzK5l4cwsn7a15ru7uwO2u+66\nK2Dr6Ohwag/Ax4mJXDfccEPAxtIMAO5CPcM6iNY1JzUTca05YmI3W5+sbqs/JSUlARtLvdDX1xew\nWe0Ms5b1hi2EEJ4ghy2EEJ4ghy2EEJ4ghy2EEJ6QU9Ex+yM7i0Szkp64Ro1ZkZIMJiSwepggxQRT\n61om3DGs3LpMnGDRZa6iDMAFVxbByAQUK2qLRVWy3L6snWHyi7vmRbYELbbumBDJ1oc172zu2NiX\nl5cHbIPzRgxmx44dARtbn9FoNGBjgiXA8zIzETdM9CQTTVmfmJBprU8m6jPhjrXdErDZWmLRvQwr\nRzUTd1kub9ZPK+I2TFSk3rCFEMIT5LCFEMIT5LCFEMIT5LCFEMITcio6ZgtlR48eDVxjpZpkggH7\n4M8EBysqkQlFrukaWcQawEUQJpLNnDkzYDt/Yo9LO5lAyKK2WJQlwKPb2HywMq3UsuxwXCYCs/Fk\nIqh1LTupiPWTpWEFeJ+YjQl3VjJ5djgHE9lYOzs7O2mZH3/8ccDmeuisFenIxoQJmaxuJlhabWJr\niUVKWs8RE+lYn1wPk7baydYy66dVpuvzwfpjHexrCdsMvWELIYQnyGELIYQnyGELIYQnjOiw16xZ\ng1gshurq6ozt2LFjaGhoQFlZGRYtWhTqiBshhBCjY0SHvXr1ajQ3Nw+xNTY2oqGhAfv370d9fT0a\nGxtz1kAhhBBfMuIukQULFgRU7U2bNmHr1q0AgFWrVmHhwoXUaWeHYrIQYbbTAHBXg113aViEOdST\nwXakMNWXqeOWaszaxPrOyrQOTmXXslBu1h+rTAbb3cNCcq1dDa550Fk4rxX2zJR9tu5YHnJ22C7A\nd/Kw3QZsPFl7AL57g4U9FxcXB2wFBQW0TNfc6ixk2wp3Z2PP5u29994L2KydJ2x9u6adsHYxsefI\ndeeIVaYVsp4N8wHWoeSuZQKj/Ibd09OTyR0Qi8XoScJCCCHGlzHvw45EImZSkz179mR+jkajoZKc\nCCHExUAqlUIqlXK6dlQOOxaLobu7G4WFhUilUjR7GIAhQqUQQoggRUVFQ4Lo2trazGtH9UlkyZIl\naGpqAgA0NTVh6dKloylGCCFECEZ8w77//vuxdetW9PX1obi4GL/85S/x8MMPY8WKFXjhhRdQUlKC\nl19+md6b/eGe5be1hCL20Z+JEEyosURHFg7NRDImDlgh367CBhMhLFGFaQKu9VhiHoOJeaxMK4yc\nCUiuZVpCJhPPWP5lNnYsP7fVJib4sjKtvN1sfbJ1FyZ9ALMzYZhto7XWvBUKng3rpzWeLKUCywHP\nBFMrx7brIdWuhx8DfD5dnw/rYFw2Tuz5YP20Ph2H2SQxosPeuHEjtW/ZssW5EiGEEGNHkY5CCOEJ\ncthCCOEJcthCCOEJOc2HnS3eMXHAOjjV9dBZFolmCQtMrGFCD7vfyo/rKrwxIdOKqGTtdBVbwoij\nrgeKWvnFmZ2JNUz4CpO3m62RI0eOBGzW4bZsjbA2MfHIOrTVilrLhkWxWVF0TFxl48kihi1x1DVq\nl0U1Wn1n48y29rL+sHkDbOEwG9ZPNh4AH2cmdrNn2xo3S4B3wepjGNFRb9hCCOEJcthCCOEJcthC\nCOEJcthCCOEJEyo6MtHQSgjFIqdYSksmyliClmv6SiZyWZFPLBqMCSNMWGAH6wK876xProIUwNNs\nsnYy8enDDz+kZTLRkdXPxE1LeGMirmvE2sGDB2mZrgIhE5rYXAB83lnUrhXZx2D9ZCIXE6+sKDrX\nlMRhIlyTyWTAdv311wdsLGLX2hDAnjn2vDOstcTqco1stkRHViZ7Dpl4bpVpRXsz9IYthBCeIIct\nhBCeIIcthBCeIIcthBCeIIcthBCekNNdItnKtesBrQDfbcB2SrB8w5bq6hpy7nrAKgDMmjUrYKus\nrAzY2G6DXbt20TJZP9kukxkzZgRsVgg9u/+mm25yKtPa6cDGns0n29VgjScLI2fzyRT3MGH5LOza\nNTwaAI4fPx6wsXVjzQeD7eRh9bC+WzsQmJ3tRmE2a3cN2xXx8ccfB2xsl0k8Hqdlsh1HrH72DFu7\nWVj6AbYrjdVj9Z3tSGH1s7VoHTiu0HQhhJiEyGELIYQnjOiw16xZg1gsNuRA3cceewzxeByJRAKJ\nRALNzc05baQQQggHh7169eqAQ45EIli/fj3a2trQ1taGu+66K2cNFEII8SUjio4LFixAZ2dnwO4S\n7pv94Z0JiZZAyD7Qs4/zTFSxPuIzIYDdz4QvS9CqqakJ2G699daA7S9/+Qu9n8H6zsKJmcDHRCrr\n2t27dwdsN954Y8DGch0PV1c2TKixRDLXsHxWpiVkMvGJrREmOlpCEVu37Jlg9VgimWteZtZOazyt\nsO1sXMcI4M8Hazu7juXIBng7Wbg6W8dWO11FYNfc+wDfOMGEULZmw6QpsBj1N+xnnnkGNTU1WLt2\nLY2bF0IIMb6MymGvW7cOHR0dSCaTKCoqwkMPPTTe7RJCCJHFqPZhD/4T+YEHHsDixYvpdXv27Bly\nzw033DCa6oQQYtLS09NDMxsyRuWwU6kUioqKAACvvPLKkB0kg7HsQgghviQWiw05J3Pwi242Izrs\n+++/H1u3bkVfXx+Ki4vxi1/8Aq2trUgmk4hEIigtLcVzzz1H780WYazDXBmu+XnZdZZgwK5lIgQ7\nKNQSb/75z38GbNu2bQvYWCSYxbRp0wI2Jqx88sknAZslEDKdgbWJ5QcPE63HxLhrrrkmYCspKaH3\ns0g01iYmvFlinquAzcQjdoiuZWdlMmHWGk/rQOpsWJSmtT6ZmMiuZWvOel6ZAM76HkYgZIItGw82\n7pbYzERHViabd9ZHgI+da95uq+/WAcr02pEu2LhxY8C2Zs0a5wqEEEKMD4p0FEIIT5DDFkIIT5DD\nFkIIT8hpetWxwAQPFs3FPthbUV8s+ohFXzKhhqUSBbhAyaK5mKjCUloCvO+sftZP61BjVubRo0ed\nbFaUJ6ufCTBhdgv19vYGbEwQc10LAJ9P18NpLTGPiZ6snUzotspkc+QaUWnNERP+mDjLhF2rnWx9\nM+EvjIDN+snmjaV2DbPJgPXdVYC26mLrjs2HFR2uQ3iFEGISIocthBCeIIcthBCeIIcthBCeIIct\nhBCekNNdItnqK1OSWf5jgCuqTEUPE9bJ1FimhFs5exmpVCpgY/1kB6yGaTtTzFk9lrrN7mc7MpgK\nbqnb7Fq2I+SOO+4I2A4cOEDLZCo+203Dxs4Ky+/r6wvYWLIdFkZuhSi75lVmodTWLibX0HSGdYAw\nq9/aSZRNd3c3tbM0B2xHBrMNzpkxGDaebEcImw+WosHC2k2TjbXzxDUnP6vH2g0SJmWH3rCFEMIT\n5LCFEMIT5LCFEMIT5LCFEMITcio6Zn94Z8KbJWixD/Hsgz/LdcwENsBdFGL5my2hiNXlKiJY+Ztd\ncwszocYqk7WJiU/vvPNOwMbC7wHg5ptvDtiYqMTuLy8vp2Wy+fjggw8CtoMHDwZsbH0BwOzZswM2\nFjb9r3/9K2CzwrMrKysDtrGGuzPBmAlvTNCyckKz+tm1H330Eb2f8d577wVs7Jlh9bCDeQH3UG4m\nQFvPpiXEZsMEbCsPOnsO2XyylBdszQFuB5qfR2/YQgjhCXLYQgjhCcM67K6uLtx+++2orKxEVVUV\nnn76aQBf7lNuaGhAWVkZFi1aRPdlCiGEGF+Gddh5eXl46qmnsHfvXuzYsQPPPvss9u3bh8bGRjQ0\nNGD//v2or69HY2PjRLVXCCEuWoYVHQsLC1FYWAjgSwFh9uzZOHz4MDZt2oStW7cCAFatWoWFCxdS\np50tdB0+fDhwjXUwJYsKYh/82du9VSaDRZcxEeCmm26i9zNhY//+/QEbE5QssYTZWZ+YIGWJJfPm\nzQvY2Hi6ir0Aj1JNJpNO94eJcGURc0y8YtcBXPRkIldBQUHAFkbAZkITwxKZ8vPzAzbWTiYkWn1n\nIh0T6m+44YaAjYmLAFBaWkrt2Vx//fUBmzWerP1sPNmzYYm4TExkNpbD3YqIZHPnamPPK2D7AYaz\nZ+vs7ERbWxvmzZuHnp6ezG6AWCxGw3yFEEKML04O++TJk1i+fDk2bNgQOFUjEomYcfdCCCHGjxH3\nYZ87dw7Lly/HypUrsXTpUgBfvlV3d3ejsLAQqVTKTLrT3t6e+TkajSIej49Ts4UQYnLQ29trxjpk\nM6zDTqfTWLt2LSoqKvDggw9m7EuWLEFTUxN+/vOfo6mpKePIs6mqqgrRbCGEuPiIRqNDXnrfffdd\n89phHfb27dvx4osvYs6cOUgkEgCAJ554Ag8//DBWrFiBF154ASUlJXj55Zfp/dkf3sOkEWRCgms0\nlJWmkgkOrmkdrbYXFRUFbHv37g3Y2OG2VmQeEyGYmMiiGi1B67///W/AxvrObJZGMX/+/ICtrKws\nYGPjaQl0LOUrE9mYUGSlA3VNyclEP+ugZCY2sz65puME+NyxOWY2K60uE+VZ25kQaYlkrrB+WuIo\n6zs71JiNp5UC1zWdMhNCLdGRtZ99EmY+yBLvrTFhDOuw58+fbxa2ZcsW50qEEEKMHUU6CiGEJ8hh\nCyGEJ8hhCyGEJ8hhCyGEJ+Q0H3Y2TMC0djW4quvMxnY6WPUzxZ0pyewgV4Dn/GUKM9tBYB3K6dp3\npvZbIrGrEs7Cjq3c1Sxcnx1ky3ZkzJ07l5Z53XXXOZXJ5thS9pmK7xqOzPI8W21ic+waFm+103Xn\nSZhDeNmOELYbxjowl8GeGba+rdQJbH2znS+s71ZoOruW2dgOsDCh6a6HaYfJe22hN2whhPAEOWwh\nhPAEOWwhhPAEOWwhhPCECRUd2Ud3K9Ofq2DAxAEr9Jfdz0Q611BmgItPTFixhFAGaycL02VY1zGR\njdmYIGWJoyxknYlsLCzfyoddUlISsLExZiKbdQAxWyMsdJjdb4Wmu+YwZgcdW2kO2LxbIefZWMIX\nW3euudktgZDVxcaY5b8/efIkLdM1LJ/5C6vv7Fr2vLNxD7MZgtXDxsPaEGAJnLR+5yuFEEJcUOSw\nhRDCE+SwhRDCE+SwhRDCE3IqOmZ/jGcf1y0BhtmZCMJsVh5fJp4xsYVFeFm5qz/77LOAzTV6Msyh\nnEy8YoKYJZIxMZIJUuwAYSb6AVxsYX1nY8QET4BHZDIhlNVtRaNOnz49YGO5r1k/LZHMVchk65Pl\nqAb4emCCGivTygnNhC7X58Ba8+zkKJZ0P0wOeDafYxENAb6+XSMdLQHbdeNCGHE0DHrDFkIIT5DD\nFkIITxjWYXd1deH2229HZWUlqqqq8PTTTwMAHnvsMcTjcSQSCSQSCTQ3N09IY4UQ4mJm2I8qeXl5\neOqpp1BbW4uTJ0/iG9/4BhoaGhCJRLB+/XqsX79+otophBAXPcM67MLCQhQWFgL4Mo3o7NmzM9FL\n45EqELAjAJnYw8QrJqaFEfNcD1O99tpraZlM0GLCGRPEWGpW61omiDEBh4lpAB9PJlSFibpi4hm7\nnx1UzAQpgK8HNh/MZh2+zAQgK3I1G2veWd9Zm5jQZAnDrO9sfbqmHbVg14Y5iJb1nT0HrumQAT5O\nrJ9hUjQzf+EqBlrtdN0Mwe5n4jkQ7rBj51nu7OxEW1sbbr31VgDAM888g5qaGqxdu9ZUvYUQQowf\nTg775MmTuPfee7FhwwZcccUVWLduHTo6OpBMJlFUVISHHnoo1+0UQoiLnhE3Bp47dw7Lly/HD37w\nAyxduhQAEI1GM//+wAMPYPHixfTe9vb2zM/RaDTzeUUIIcSXpFIpdHd3O107rMNOp9NYu3YtKioq\n8OCDDw6p4Px3yVdeeQXV1dX0/qqqKtc2CyHERUlRUdEQnSeZTJrXDuuwt2/fjhdffBFz5sxBIpEA\nADz++OPYuHEjkskkIpEISktL8dxzz41T04UQQlgM67Dnz59PVdm7777bqfBs9ZaVZeWIZeo8U9FZ\nmK11ICm7lin2x44dC9isA0m/+c1vBmxMif7oo48CNmuHzO7duwM2Juwyxbq4uJiWyUilUgEbGyMr\njJzB2sTUfivXMttpwe5nh+OyA4QBPp9sjbAxtnZ0sDFhOwM+/vhjp7oB95znYQ4gturKhn2utO49\ncOBAwMZSL7BdO9aaZzum2C4m5i+s9BbsWrYjhO3SsHaeuB4WzOqxdvJY7Wco0lEIITxBDlsIITxB\nDlsIITxBDlsIITwhp/mws4UQltPZOoyViTpM2GCikCXAsGuZmFhbWxuwWXvI586dG7AxwZQJE9YB\nq4P3uZ/nzTffDNgOHToUsLW2ttIyXfMys/2gVuisq5jIcgtb+YZZyDqrp6urK2CbOXMmLZMJma45\nti1xlK1FK3e2K6yfTKRjYeTWYdZM/GIh20zgs3JsM4FwcMzFeZjoZ42nay7xMKkTXGECo1UPs7P7\nWdoIazytlAoMvWELIYQnyGELIYQnyGELIYQnyGELIYQnTKjoyHLmWuKT6wGvrgeXsvYAPApux44d\nTu0BuBhYWVkZsLF8K5aw8f777wdsTORiUVfn85VnYwmc2TDhLIyYxkQyFsllRU+yCEYWLciEGksc\nZSIbi+hkIpc172w8WZlMgA6TB72goCBgY8IwE1YBLrQz8Z/lB2f9AfhzyA7m/eCDDwI2K9qPjYnr\nobVhBFfX9W2VyZ5DNvZMcLXERcsHMvSGLYQQniCHLYQQniCHLYQQniCHLYQQnpBT0TH7Az2LrGMf\n8QEezeUqnFkHaDIhgIkgrE1hDjn997//HbD19vYGbGEin5iwwdJfWmkhmdDE5oOJP1aZbEzYGDMB\nx5ojFunI+s7miImTAO8nEyjDCEWs/Sw9qms6TsvORFgmTloCtquYGCYalUWZzpgxI2Bjc2SV6Xpg\nLuunlaKZicisTey6MM87K5PZrHm3Us4y9IYthBCeIIcthBCeMKzDPnPmDObNm4fa2lpUVFTgkUce\nAfDl3uWGhgaUlZVh0aJF9E9OIYQQ48uwDnvq1KloaWlBMpnE7t270dLSgm3btqGxsRENDQ3Yv38/\n6uvr0djYOFHtFUKIi5YRP4mcj5Q6e/Ys+vv7cfXVV2PTpk1YtWoVAGDVqlX461//mttWCiGEGHmX\nyMDAAObOnYsDBw5g3bp1qKysRE9PTyaPdCwWQ09Pj1NlLETZ2oHgquKHUY3Zta75da0yWftZuDsL\nGT9/En02t956a8DGDj5ln6LYrgKAjyfb1XD8+PGAzeq766HKbKeFFfbsukPGNWzZahPbhcTWgqXg\ns7FjOZDZgcpsjAEern/w4MGA7ZZbbgnYOjo6aJk33nijU5lsV4M1xgsWLAjYWF55ViY7+Blw30UV\nZkcFw9qZlY0Vms7uZ9eyvluH7YY55HrEVT9lyhQkk0l88sknuPPOO9HS0jLk3yORiNk5IYQQ44fz\na8r06dNxzz334J133kEsFkN3dzcKCwuRSqXoCSkAsGfPnszP0WgUpaWlY2+xEEJMIrq7u52/Ugz7\nDbuvry/zZ/fp06exefNmJBIJLFmyBE1NTQCApqYmLF26lN5fXV2d+Y/9ySSEEBc7hYWFqKmpyfw3\nHMO+YadSKaxatQoDAwMYGBjAypUrUV9fj0QigRUrVuCFF15ASUkJXn755XHtgBBCiCDDOuzq6mrs\n2rUrYL/mmmuwZcuWEQvPFntYKLUVrukaGsrut76pMzurh4kDVplMPGNlMtHPEjyPHDkSsDEhk4lp\nlijD2snqz8/Pp/czWP0sbJoJs2HEZmZjY8zCmwEuALFxChOOzHIos/VdVFTkZAN42DbrExPoLIGQ\nXcvWAlvf7KBiC3Y/Wx+WwMaEZTbvbN2MVRRn42G1k4nNlpiYTZgc2xaKdBRCCE+QwxZCCE+QwxZC\nCE+QwxZCCE/IaT7sbIEgjEjmGsHomuPaKpOJNaxMK48vEwxcRQyWVxjgwp1r1BcTRQAeXeeaG9jK\nXc1EIXZ/GLGZ3c/6xA50ZmItAHz66acBGxs7FqkY5pBoNm8sDzk7mBfgB+ayQ3iZyGXl7WbjzMaO\nYYlhLGc5EwPZ/WFyq7O2uz5vFq4HLVvr01WwZWveEifDiN16wxZCCE+QwxZCCE+QwxZCCE+QwxZC\nCE/IqeiY/TGefZy3RAgm/LmKGGGiJ5kQ4BoZB3Cxh9XP7rfEEnYQLos6Y0IkOwwV4BFz7FpWd1tb\nGy2TCSssApDNUZjDbVkUHEuPykQ7gI8zGw82xtZaYoIrq5+VyQ7GBXjq0fLy8oCNCdhhRHE2xux+\na32yfrIVncXdAAAKo0lEQVRr2dhZ0b2u0X6s7dazyfwN22TAyrT8Elu3rB5XcTIsesMWQghPkMMW\nQghPkMMWQghPkMMWQghPkMMWQghPyOkukWws5XUs17KdH1aop2sIKtt9YR3eyXZKMBvbVWDtQLAO\nac2mrKwsYLN2ibCw7dtuuy1gY7tmPvjgA1omU9fZwcBMxbfml9XPdjCwusMclByPx53aaeWEZrsN\n2M4TtqPD2i1w3XXXBWzs8OYwhwWzsXNd81aag7EcjsvSBAB83tkzM9ZDeNlulDAHELOxZ2PM2mm1\nXfmwhRBiEjKswz5z5gzmzZuH2tpaVFRU4JFHHgEAPPbYY4jH40gkEkgkEmhubp6QxgohxMXMsJ9E\npk6dipaWFlx22WX44osvMH/+fGzbtg2RSATr16/H+vXrJ6qdQghx0TPiJ5HzUU1nz55Ff39/JgVl\nmO/RQgghxs6IouPAwADmzp2LAwcOYN26daisrMSf//xnPPPMM/j973+Puro6/OY3v6F5f7OdephQ\nVSZssI/7YXJsMyGBffBnoqElaLH2s7FgAgqrB+D/M2RiHhMS33vvPVomC3tm+bhZKLQFE/6YjYnA\nLPc0wNcIu5/NpRXuzuaYiV8sT7QV8s3Gnol0bN1Y+aiZ8Mb6aQliDLaW2Hiy66xnk9XvKm5auObT\nZmNkbQhgZbL7mS1MSgJWD0vRYM2bJe4yRnzDnjJlCpLJJA4dOoQ33ngDra2tWLduHTo6OpBMJlFU\nVISHHnrIuUIhhBCjw/l/1dOnT8c999yD//znP1i4cGHG/sADD2Dx4sX0nvb29szP0WiUvnkKIcTF\nTE9PD3p7e52uHdZh9/X14dJLL8VVV12F06dPY/PmzXj00UfR3d2NwsJCAMArr7yC6upqen9VVdWQ\n38P8iSSEEBcDsVgMsVgs8/vgF91shnXYqVQKq1atwsDAAAYGBrBy5UrU19fjhz/8IZLJJCKRCEpL\nS/Hcc8+NX+uFEEJQIukcbfeIRCL43ve+N8TmGgEIuItKTDCwxBLXA3tZJJoVjcSuZfmCWd1WvmGr\n/dmw3NXWdLK/blg7mQDCBBSAH2jK6meClDWebI6ZqMTm3TosmOWfZmM/+C3nPKyPlp2VyWyWKO4q\nqLE1a42nq+joKvYCPD+56/q2ojxdn2PWTmtDAFsPbC1b64bB+s58FZtjSxTPvra5udk+rNilkUII\nIS48cthCCOEJcthCCOEJcthCCOEJOU2vmi1aMJHMiiRjQgITLMKkemTCBhOPWD2WsMGEBCa2hIkk\nY3W5HjZsCUVMlDp69GjAxtrOxg2wIzWzYelireg0Nh+uc8z6A3Ahld3veoCw1U6G68HRVpvYfIZJ\nr+p6yLTr4bRWmaydrpGwAN984HoQt1Umez5c0/JazxFbYywFL7vfWvNh0Bu2EEJ4ghy2EEJ4ghy2\nEEJ4ghy2EEJ4ghy2EEJ4Qk53iWQrumGSPzGVlam5TN22dl+4HoTLlGRL2Wd1sTa5hsADfKcGU8dZ\n21norNUmVxXeGk8rtN7lfiuUminprqkCLGWflcnC8sPkQWe7g1x3MbFQeQDIz893rj8ba32ycXJ9\njqz1ye533X1hlWntFsuG7Yax0lu4hvqHWUtsLbO15Jo3G3B/jgC9YQshhDfIYQshhCfIYQshhCfI\nYQshhCfkVHQMVBbi8FAmGDDBwjV3tAX74M/KNPPTOobkslDmMDmMWZvC5EVmY8faaYWhM5iwwupn\nokwYgZDZ2LxZ4eJMqGLjGUYcZePkmnrh2muvpWW6Cr5hRDI2764CtCV4jkVUt8KzXQ9vZmVaQqZr\nmgPXA4ABPh+uz4ElLobJxz0hb9iu55X5Ajsx23f6+voudBPGlcm25np6ei50E8adyTZHE+EX5LBH\nwWR02FbiJF+ZbHM02Z4hYPLN0US89OgbthBCeIIcthBCeELODuFduHAhtm7dmouihRBi0nLbbbeh\ntbWV/lvOHLYQQojxRZ9EhBDCE+SwhRDCE3LusJubm1FeXo6ZM2fiySefzHV1486aNWsQi8VQXV2d\nsR07dgwNDQ0oKyvDokWLcOLEiQvYwnB0dXXh9ttvR2VlJaqqqvD0008D8LdPZ86cwbx581BbW4uK\nigo88sgjAPztz2D6+/uRSCSwePFiAH73qaSkBHPmzEEikcAtt9wCwO/+nDhxAvfeey9mz56NiooK\nvP322xPSn5w67P7+fvzkJz9Bc3Mz3n33XWzcuBH79u3LZZXjzurVq9Hc3DzE1tjYiIaGBuzfvx/1\n9fVobGy8QK0LT15eHp566ins3bsXO3bswLPPPot9+/Z526epU6eipaUFyWQSu3fvRktLC7Zt2+Zt\nfwazYcMGVFRUZCL5fO5TJBJBa2sr2trasHPnTgB+9+enP/0pvv3tb2Pfvn3YvXs3ysvLJ6Y/6Rzy\n1ltvpe+8887M70888UT6iSeeyGWVOaGjoyNdVVWV+X3WrFnp7u7udDqdTqdSqfSsWbMuVNPGzHe/\n+9305s2bJ0WfTp06la6rq0u3t7d735+urq50fX19+vXXX09/5zvfSafTfq+7kpKSdF9f3xCbr/05\nceJEurS0NGCfiP7k9A378OHDKC4uzvwej8dx+PDhXFY5IfT09CAWiwEAYrGYt2HDnZ2daGtrw7x5\n87zu08DAAGpraxGLxTKfe3zuDwD87Gc/w69//esh+TR87lMkEsEdd9yBuro6PP/88wD87U9HRwcK\nCgqwevVqzJ07Fz/60Y9w6tSpCelPTh22lZRlMhGJRLzs58mTJ7F8+XJs2LAB06ZNG/JvvvVpypQp\nSCaTOHToEN544w20tLQM+Xff+vP3v/8d0WgUiUTCTDrmW5+2b9+OtrY2vPrqq3j22Wfx5ptvDvl3\nn/rzxRdfYNeuXfjxj3+MXbt24fLLLw98/shVf3LqsGfMmIGurq7M711dXYjH47msckKIxWLo7u4G\nAKRSKUSj0QvconCcO3cOy5cvx8qVK7F06VIA/vcJAKZPn4577rkH77zzjtf9eeutt7Bp0yaUlpbi\n/vvvx+uvv46VK1d63aeioiIAQEFBAZYtW4adO3d62594PI54PI6bb74ZAHDvvfdi165dKCwszHl/\ncuqw6+rq8P7776OzsxNnz57FSy+9hCVLluSyyglhyZIlaGpqAgA0NTVlnJ4PpNNprF27FhUVFXjw\nwQczdl/71NfXl1HjT58+jc2bNyORSHjbHwB4/PHH0dXVhY6ODvzxj3/Et771LfzhD3/wtk+ff/45\nPvvsMwDAqVOn8Nprr6G6utrb/hQWFqK4uBj79+8HAGzZsgWVlZVYvHhx7vsz7l/Fs/jHP/6RLisr\nS3/9619PP/7447mubty577770kVFRem8vLx0PB5P/+53v0sfPXo0XV9fn545c2a6oaEhffz48Qvd\nTGfefPPNdCQSSdfU1KRra2vTtbW16VdffdXbPu3evTudSCTSNTU16erq6vSvfvWrdDqd9rY/2bS2\ntqYXL16cTqf97dOHH36YrqmpSdfU1KQrKyszfsDX/qTT6XQymUzX1dWl58yZk162bFn6xIkTE9If\nhaYLIYQnKNJRCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE8QQ5bCCE84f8A\nHJ4h/8ChMjkAAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0x10a28f510>"
       ]
      }
     ],
     "prompt_number": 82
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Can you find any face activity you believe?"
     ]
    }
   ],
   "metadata": {}
  }
 ]
}