Produce 300W results.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "import tensorflow as tf\n",
    "import menpo.io as mio\n",
    "import utils\n",
    "import pandas\n",
    "\n",
    "from menpofit.visualize import plot_cumulative_error_distribution\n",
    "from menpofit.error import compute_cumulative_error\n",
    "from scipy.integrate import simps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# Download pretrained model.\n",
    "!wget https://www.doc.ic.ac.uk/~gt108/theano_mdm.pb --no-check-certificate\n",
    "\n",
    "# Download 300W dataset from http://ibug.doc.ic.ac.uk/resources/300-W/.\n",
    "images = mio.import_images('PATH_TO_300w')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "images = mio.import_images('/Users/gtrigeo/db/300w')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "MDM_MODEL_PATH = 'theano_mdm.pb'\n",
    "bounding_boxes_path = Path('bbs/300w/')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# the image to fit (rgb image of HWC) where H: height, W: weight and C\n",
    "# the number of channels (=3).\n",
    "image = tf.placeholder(tf.float32, shape=(None, None, 3), name='images')\n",
    "\n",
    "# we only use the upper-left (x0, y0) and lower-down (x1, y1) points\n",
    "# of the bounding box as a vector (x0, y0, x1, y1).\n",
    "initial_bb = tf.placeholder(tf.float32, shape=(4), name='inits')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "with open(MDM_MODEL_PATH, 'rb') as f:\n",
    "    graph_def = tf.GraphDef.FromString(f.read())\n",
    "    pred,  = tf.import_graph_def(graph_def,\n",
    "                                 input_map={\"image\": image, \"bounding_box\": initial_bb},\n",
    "                                 return_elements=['prediction:0'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Start TensorFlow session.\n",
    "sess = tf.InteractiveSession()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Fitting images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def outer_eye_corner_68(prediction, ground_truth):\n",
    "    \"\"\"Computes the outer-eye corner normalised rmse error.\n",
    "    \n",
    "    Args:\n",
    "      prediction: A menpo `PointCloud' with the predicted landmarks.\n",
    "      ground_truth: A menpo `PointCloud' with the ground truth landmarks.\n",
    "    Returns:\n",
    "      A scalar with the normalized error.\n",
    "    \"\"\"\n",
    "    \n",
    "    num_landmarks = 68\n",
    "    assert prediction.points.shape[0] == num_landmarks\n",
    "    assert ground_truth.points.shape[0] == num_landmarks\n",
    "\n",
    "    interocular_dist = np.linalg.norm(ground_truth.points[36] - ground_truth.points[45])    \n",
    "    rmse_error = np.sum(np.linalg.norm(prediction.points - ground_truth.points, axis=-1))\n",
    "\n",
    "    return rmse_error / (num_landmarks * interocular_dist)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def fit_image(im):\n",
    "    \"\"\"Estimates the landmarks usign the provided image and bounding box.\n",
    "    \n",
    "    Args:\n",
    "        A menpo `Image' with the bounding_box attached.\n",
    "    Returns:\n",
    "        A menpo `PointCloud' with the predicted landmarks.\n",
    "    \"\"\"\n",
    "        \n",
    "    # menpo stores images CHW instead of HWC that TensorFlow uses.\n",
    "    pixels = im.pixels.transpose(1, 2, 0)\n",
    "\n",
    "    # the model expects rgb images.\n",
    "    if im.n_channels == 1:\n",
    "        pixels = np.dstack([pixels]*3)\n",
    "\n",
    "    bounding_box = im.landmarks['bounding_box'].lms\n",
    "\n",
    "    prediction, = sess.run(pred, feed_dict={\n",
    "        image: pixels,\n",
    "        # grab the upper-left and lower-down points of the bounding box.\n",
    "        initial_bb: bounding_box.points[[0, 2]].ravel()}\n",
    "    )\n",
    "    \n",
    "    return menpo.shape.PointCloud(prediction)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Loading and fitting images: [==========] 100% (600/600) - done.                 \n"
     ]
    }
   ],
   "source": [
    "errors = []\n",
    "\n",
    "images_iter = menpo.visualize.print_progress(images, prefix='Loading and fitting images')\n",
    "\n",
    "for im in images_iter:\n",
    "    bb_name = im.path.stem + '.pts'\n",
    "    bounding_box = mio.import_landmark_file(bounding_boxes_path / bb_name)[None]\n",
    "    im.landmarks['bounding_box'] = bounding_box\n",
    "\n",
    "    # Temporary fix. \n",
    "    # There is a bug that causes the model to output\n",
    "    # NaNs sometimes on a new image. If we fit the\n",
    "    # same image a few times it stops occuring.\n",
    "    _ = [fit_image(im) for i in range(3)]\n",
    "\n",
    "    mirrored_fit = utils.mirror_landmarks_68(\n",
    "        fit_image(utils.mirror_image_bb(im)), im.shape)\n",
    "\n",
    "    prediction = menpo.shape.mean_pointcloud((fit_image(im), mirrored_fit))\n",
    "    error = outer_eye_corner_68(prediction, im.landmarks['PTS'].lms)\n",
    "    errors.append(error)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualising results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdkAAAFaCAYAAABSYyd0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYFOW5/vHvwyDiAhoXPOeIYWQzStzQoMSFieQXUaIm\niopLREXx6NFIYhJQDgyLSkyiMcYVxfVEMZooRo0rDCpqJCyihnEBB1BEIUQzIuI48/z+qGpsmu6Z\nnqFrqqvn/lxXX3RVv1V9d830PLy1vGXujoiIiBReu7gDiIiIlCoVWRERkYioyIqIiERERVZERCQi\nKrIiIiIRUZEVERGJiIqsNMnMdjGz+8zsbTObY2aPmlnPiN9zppn1baLNxWbWMW36UTPrXID3rjGz\nV81svpnNM7NrN3edjbxXNzP7zMzmmtk/zOxlMxuW9voxZvaLRpbf18yOauT1A1L5zazSzH7azHyR\nbGORtqJ93AEkER4C7nD3UwDMbG9gF+CdWFPBSOAe4HMAd/9+gdbbAFS4+79yNTCzMnevzzWd73Kh\nd9z9gPD1cuAhM8Pd73L3vwB/aWSV+wEHAn/N8V5zgblN5WpEVNtYpE1QT1YaZWbfAb5w91tT89z9\nNXefbWYDzOwvaW1/b2ZnhM/fNbMrw97gK2a2v5k9EfaGR4Rtci6fkeHGcB2vmVllOO8i4L+AmWb2\nbNp77mBmk83sgrTlN/TgzOxn4boWpNaV7WOT5bsR9q5/a2avAD82szvM7CYzexm4ysy+ZmYPhb3g\nF83sm2nvf7eZvQDc3dj2dvca4KfAxeGyw8zs9+HzE8NtMN/MqsxsC2AicFLY4z4x870ytzGwX5jt\nTTM7p7GfQ2PbOHz+0zDPQjNL5e0W9sinmNnr4c98y8Y+s0gpU5GVpnyTxntCjQ0ZVuPu+wMvAHcA\nxwP9CQpDPsunXObu/YB9gQoz+6a7/x54n6DHOTBjXfcDJ6UtfxJwv5n9P6BXuK79gQPN7NAc7zkj\nbXfxxWnzt3D3fu7+23B6V3c/2N1/BkwA5rn7vsAYgh5gyp7AEe5+Wh6fdx6wR9p06nONBb4XbtNj\n3b0OGAfc7+593f2BHO+Vvo33BiqAbwPjzOw/srQJZgTbeAVZtnG4K38Y8C2Cn+m5ZrZv2KYn8Ht3\n/ybwCXBCHp9ZpCRpd7FEKdU7eg3Yxt0/Az4zs8+beVxvqJmdS/D7+h/AXsDrBD1OS2tnAO6+wMx2\nDgtIF2CNu79vZiOB/2dm88K22wC9CP4TkCnX7uL7M6YfSHt+KMF/JHD3mWGvetvwtUfc/Ys8P6/l\nmP8CcJeZ/RH4cyPLN/Ze08PX/mlmM4B+BIWwuXkOBR5y988BzOzPwGEEP/N33f21sN1coLyJ9YuU\nLBVZacobwJAcr33JxntDOma8vj78tyHteWq6fR7Lp45RXgIc4O7/NrM7srXL4gHgRIKinCqMBkxO\n3/XdiFyFbm0j0431yjOXa0xfYFHmTHe/wMy+BXwfmGu5Twxr7L3SM1o4/SVQljY/n+3bmPSfdX0B\n1ieSWNpdLI1y9xlAh9TxOwhOfDKzQ4ClwF5mtoWZbQ8MzLWeDKkCls/ynYFPgVoz2wVIP5P23+Hr\n2fwRGEqwqzLV23wSONvMtgk/x3+Z2c55Zs7H88Dp4borgNXu/mkey20o6OF/Kn4NXLdJI7Pu7j7H\n3SuBj4DdgFpyb4NsjjOzDma2IzAAmEPwc9gzx88hcxunsj4P/MDMOobb84fhvI0+j0hbp56s5OOH\nwO/MbDSwDqgBRrr7e+Guy9eBdwmOJaY01qtzgHyWd/eFZraAoGe3nI137d4KPGFm74fHDDe8p7v/\nw8w6Ae+5+4fhvKfN7BvAS2YGQYE6HViVJd9MM0udBbzQ3c/M8pkypycAt5vZqwS9yU1O4sqhu5nN\nBbYiKGrXuvs9Wdr92sx6hc+fDbfNcmB0uAt8cpZMmRYCVcCOwER3XwnQyM8h6zZ29/lmdidBkXZg\niru/ambd8sgg0mZYlLe6M7OpBLu2PnT3fXK0uY6gd7IWONPdF0QWSEREpBVFvbv4DuDIXC9acBF9\nD3fvBZwH3BxxHhERkVYTaZF19xeAnBf0A8cRXjfo7n8DtguPu4mIiCRe3Cc+7UpwnC3l/XCeiIhI\n4iXmxCcz08kUIiIt4O464zsmcfdk3ye4DCGlazgvK3cv+seAAQNiz6CcypnUjMpZ+IfEqzWKbOao\nPOkeIbzMwcwOBj728HILERGRpIt0d7GZ3UswTuqOZrYMqAQ6AO7uU9z9cTM72szeIbiE56wo87SG\n8vLyuCPkRTkLKwk5k5ARlFNKS6RF1t1PzaPNhVFmaG1nnnlm3BHyopyFlYScScgIyimlJdLBKArJ\nzDwpWUVEioUF9ybWiU8xifvEJxERkZKlIisiIhIRFVkREZGIqMgWWFVVVdwR8qKchZWEnEnICMop\npUVFtsCS8sVTzsJKQs4kZATlTLKtttpqpZl5W3xstdVWK7Ntk8QMqygiIsXt888/36WtXgWS6+Y2\n6smKiIhEREVWREQkItpdLCIiBWXWvLEvUruYW7pcMVNPtsAqKirijpAX5SysJORMQkZQTiktKrIF\nlpQvnnIWVhJyJiEjKGcpaelt+Db39n3t2rXjjDPO2DBdX1/PzjvvzLHHHgvAXXfdRZcuXTjggAPo\n3bs3Rx11FC+99NKG9meeeSbbbLMNa9eu3TBv5MiRtGvXjjVr1jRrG6jIiohISdlmm214/fXXWb9+\nPQBPP/00u+2220Zthg4dyty5c3nrrbcYNWoUxx9/PG+++SYQ7Lbu1asX06dPB4KiP3PmTLp27drs\nLCqyIiIlwMyyPtqqo48+msceewyA++67j1NOOSVn24qKCkaMGMGUKVM2zBs6dCj3338/EFwTfcgh\nh9C+ffNPY1KRFRGRWBX6PwRmxtChQ7nvvvtYv349Cxcu5KCDDmp0mb59+1JdXb1hulevXqxatYqP\nP/64ySLdGBVZEZES0tzjl6Xqm9/8JjU1Ndx3330MHjy4yW2R+bqZcfzxxzNt2jReeeUVDj300BZt\nTxXZAkvKUGvKWVhJyJmEjKCcbVFU/yE49thj+fnPf55XL3T+/PnsueeeG8076aSTGDt2LN/73vda\n3NNWkS2wpHzxlLOwkpAzCRlBOWXzpQr22WefTWVlJX369MnZBmDWrFnceuutjBgxYqM2X//617ny\nyis5//zzW5xFg1GIiEhJSfU6d911Vy688MKsbf74xz8ye/Zs1q5dS/fu3fnzn/9M7969N1oe4Nxz\nz91kvc0ReZE1s0HAtQS95qnuflXG69sDtwM9gHXA2e7+j6hziYiUiuZeu9laWrqLdXNPgvr3v/+9\nybwBAwYwYMAAAIYNG8awYcNyLn/77bdnnb9kyZJmZ4m0yJpZO+B6YCCwAphjZtPdvTqt2WXAfHc/\n3sz2AG4AvhtlLhGRJPvwww957rnneO6555g1axavvfZa3JEkh6h7sv2At919KYCZTQOOA9KL7F7A\nZAB3f9PMys1sZ3dfFXE2EZFEeO+99zYU1Oeee26jS00Attxyyw0DLxSDlp7EVIpnQ0ddZHcFlqdN\nv0dQeNO9ChwPzDazfsDXga5AIotsUoZaU87CSkLOJGQE5XR3ampqmDVr1oaimrmbcuutt+bb3/42\nhx9+OAMGDKBfv35stdVWkeSRzVMMJz79Evidmc0DXgPmA/XZGqb/UpeXl1NeXk5FRUXWX/aqqqqs\nZ/+1Rvts4syTK9P48eOLJo+2Z+u0Ty1TLHmyta+oqCiqPI21Ty2zOet3d+655x4eeOABli5dSk1N\nDbW1tRu16dy5M4ceeuiGolpbW8vs2bOpr69nxowZzJgxY0PbtjzKUzGyKLvnZnYwMN7dB4XTowHP\nPPkpY5l3gb3d/dOM+V6KuxJEpLgV+vZrDQ0NvPHGGxt6qc899xwffvjhRm122GEHDjvsMAYMGMDh\nhx/OfvvtR1lZWYtzunurVN62/HfazLJu56h7snOAnmbWDfgAGApsdFWwmW0HfObudWZ2LjArs8CK\niCRVfX09CxYs2HBM9fnnn9/kbOAuXbpsOPv18MMPp0+fPrRr17xhDHIVt9bs2Xbs2PFDM9ul1d6w\niHTs2PHDbPMjLbLuXm9mFwJP8dUlPIvM7LzgZZ8C7AncZWYNwBvA8CgziYi0RFM9tFQxq6urY+7c\nuRuOqc6ePXuTS0q6du26oaAOGDCA3r17l8Ru3nXr1v1H3BmKTaS7iwupLe+GEJH4pIpfvkV26623\n5rPPPtvote7du28oqIcffji77757qxXVXLsxpXVoWMUCS8pQa8pZWEnImYSMkPycn332GXvssQcj\nRozg//7v/1i2bBmLFy/mjjvu4Mwzz6R79+4l0WuV/KjIFljS/0AUG+UsnCRkhOTnXLlyJdXV1dxy\nyy2cdtppm9wsXNoWFVkRkQLaZZc2ed6P5KAiKyKSw+LFi+OOIAmnIisikmHp0qWcc8457LHHHnFH\nkYQrhhGfRESKwvvvv88VV1zBbbfdRl1dXbOvVRXJpN+gAmvr464WmnIWThIyQjw5V65cyciRI+nR\nowc33XQTX375JaeddhqLFi3a0MbMNnpMmDBho2mRbHSdrIi0WatWreJXv/oVN9xwA+vWrQPgxBNP\nZPz48ey1115A4YdVbG26TjZe2l0sIm3OmjVruPrqq7nuuuv49NNgFNfjjjuOCRMmsO+++27UttiK\npiSLiqyItBmffPIJ1157Lddcc82GoQ6PPvpoJk6cyAEHHBBzOilFKrIiUvI+/fRTfv/73/PrX/+a\nf/3rXwB897vfZeLEifTv3z/mdFLKVGRFpGR99tln3HjjjVx11VWsXr0agMMPP5xJkyZx+OGHx5xO\n2gKdXVxgSR8SrtgoZ+EkISMUJufnn3/OddddR48ePfj5z3/O6tWr6d+/P8888wxVVVUFKbBJ2Z4S\nLxXZAkvKF085CysJOZOQETYv5xdffMFNN91Ez549ufjii1m5ciUHHnggjz/+OLNnz2bgwIEFu9wm\nKdtT4qXdxSKSeHV1ddx9991MmjSJpUuXArDPPvswceJEjj32WF3HKrFRkRWRxKqvr+fee+9lwoQJ\nG8YZ3nPPPZkwYQInnHCCRmyS2KnIikjiNDQ08MADDzB+/Hiqq6sB6NWrF+PHj+fkk0+mrKws5oQi\nARVZEYlFY7twJ0yYsMk8d8fdeeihh6isrOT1118HYPfdd2fcuHGcfvrptG+vP2lSXPQbWWAaH7aw\nlLNwkpCxMY8++ijjxo1j/vz5AOy222787//+L2eddRZbbLFFq+dJ+vaU1hH52MVmNgi4luBM5qnu\nflXG652B/wO+DpQBV7v7nVnWo7GLRUpIqifb1Pc6s8f7n//5n4wZM4ZzzjmHLbfcMrJ8pUJjF8cr\n0iJrZu2At4CBwApgDjDU3avT2lwKdHb3S81sJ+BNYBd3/zJjXSqyIiWkuUW2S5cujB49mv/+7/9m\nq622ijxfqVCRjVfUu4v7AW+7+1IAM5sGHAdUp7VxoFP4vBPwz8wCKyKyZMkSttlmm7hjiDRL1Oe3\n7wosT5t+L5yX7npgLzNbAbwKXBxxJhFJIBVYSaJiOPHpSGC+ux9hZj2Ap81sH3f/NLNh+okG5eXl\nlJeXU1FRkfUEhKqqqqwjsqi92qt98bTPJVv78ePHF13+Ymyfel2KROq0+CgewMHAE2nTo4FRGW0e\nBQ5Jm34WODDLujwJZs6cGXeEvChnYSUhZ7FlJDhUtMn8zJy52sWt2LZnLuG2i/RvvR65H1HvLp4D\n9DSzbmbWARgKPJLRZinwXQAz2wXoDSyJOFdkkvI/SOUsrCTkTEJGUE4pLZHuLnb3ejO7EHiKry7h\nWWRm5wUv+xTgcuBOM1sYLvYLd18TZS4REZHWEPkxWXd/AtgjY94tac8/IDguKyJtxPTp0zc8zzby\nU7YRn0SSSKNni0iruvnmmzn++OPjjiHSKlRkRaRVuDtjxozh/PPPp6GhgYkTJ9LQ0LDJiSKVlZW5\nTqQUSZxiuISnpCRlPFPlLKwk5IwzY11dHeeccw533303ZWVlTJkyhbPPPjtr2yRsS0hOTolX5GMX\nF4qGVRRJptraWoYMGcJTTz3F1ltvzYMPPshRRx0Vd6w2Q8Mqxks9WRGJzMqVKzn66KOZP38+Xbp0\n4bHHHuPAAw+MO5ZIq1GRFZFIvPnmmwwaNIiamhp69uzJE088QY8ePeKOJdKqdOKTiBTcSy+9xLe/\n/W1qamro168fL774ogqstEkqsiJSUA8//DBHHHEEa9as4fvf/z4zZsxg5513jjuWSCxUZAssKUOt\nKWdhJSFna2S86aabOOGEE/j8888599xzeeihh5p995wkbEtITk6Jl4psgSXli6echZWEnFFmTF0D\ne8EFF2y4BvaWW26hffvmn/aRhG0Jyckp8dKJTyKyWZpzDaxIW6MiKyItpmtgRRqnIisiLZJ+DezO\nO+/MY489xre+9a24Y4kUFRVZEWk2XQMrkh+d+FRgSRnPVDkLKwk5C5Ux6mtgk7AtITk5JV4au1hE\n8vbwww9zyimn8PnnnzN48GDuv//+Zl+iI61LYxfHSz1ZEclL5jWwDz/8sAqsSBNUZEWkUZnXwE6Y\nMKHF18CKtDX6lohITpnXwN5yyy0MHz487lgiiRF5T9bMBplZtZm9ZWajsrz+MzObb2bzzOw1M/vS\nzLaPOpeINK62tpbvf//73H333Wy99dY88sgjKrAizRRpkTWzdsD1wJFAH+AUM/tGeht3/4277+/u\nfYFLgSp3/zjKXFFKylBryllYScjZnIwrV65kwIABPPXUU+y8885UVVVx9NFHRxcuTRK2JSQnp8Qr\n6p5sP+Btd1/q7nXANOC4RtqfAtwXcaZIJeWLp5yFlYSc+WZ888036d+/P/Pnz6dnz5689NJLrTrI\nRBK2JSQnp8Qr6iK7K7A8bfq9cN4mzGwrYBDwp4gziUgOug+sSGEV04lPxwAvNLarOP3i7/LycsrL\ny6moqMh6UXhVVVXW/2lG3b6mpiZr9rjy5GpfU1PD+PHjiyaPtmf07VPzcrW/+eabefDBB6mvr6dX\nr14MHDiQN954o2jyF1v7prZnXPlTr0uRcPfIHsDBwBNp06OBUTna/hkY2si6PAkqKyvjjpAX5Sys\nJORsLOONN97o7dq1c8DPPfdcr6ura71gGZKwLd2TkzP82xnp33o9cj+i3l08B+hpZt3MrAMwFHgk\ns5GZbQcMAKZHnEdE0rjrGliRKEX6TXL3ejO7EHiK4PjvVHdfZGbnBS/7lLDpD4An3X1dlHlaQ1LG\nM1XOwiqmnGa5R9CbMGFC1vnFdA1sMW3LxiQlp8RLYxeLlJjGimw2W2+9NQ888ECrXaIjrUtjF8dL\nRVakxKSKbFPfl1S7V155RfeBLWEqsvFSkRUpMc0tsvpelTYV2XjpBgEiIiIRUZEVERGJiIpsgSXl\nInDlLKwk5ExCRlBOKS0qsgWWlC+echZWEnImISMop5QWFVkREZGI5DUYhZntCnRLb+/uz0UVSkRE\npBQ0WWTN7CrgZOAfQH042wEVWZEilm1QilwjPolINPLpyf4A2MPd10cdRkREpJTkc0x2CbBF1EFK\nRVLGM1XOwiqmnO7OqlWr6NSpExCcoOPuzJw5M9cdropKMW3LxiQlp8SryRGfzOxPwL7As8CG3qy7\n/zjaaJvk0IhPInm65JJLuOaaazjyyCN54okn4o4jMdKIT/HKp8gOyzbf3e+KJFHuHCqyInlYvnw5\nvXr1Yv369cydO5e+ffvGHUlipCIbryaPybr7XeG9YHuHs95097poY4lIS02YMIH169dz8sknq8CK\nxCyfnmwFcBdQAxiwGzCstS/hUU9WpGnV1dX06dMHM2PRokX06tUr7kgSM/Vk45XP2cVXA99z9zcB\nzKw3cB9wQJTBRKT5xo4dS0NDAyNGjFCBFSkC+ZxdvEWqwAK4+1vobOOckjLUmnIWVjHk/Pvf/86D\nDz5Ix44dGTdu3CavF0PGfCinlJJ8iuzfzew2M6sIH7cCf486WFIl5YunnIVVDDkvu+wyAC666CJ2\n3XXXTV4vhoz5UE4pJfnsLj4f+B8gdcnO88CNkSUSkWabMWMGTz/9NJ07d2bUqFFxxxGRUD5nF68H\nrgkfzWZmg4BrCXrNU939qixtKoDfEuyGXuXu32nJe4m0Re7OpZdeCsAvfvELdtxxx5gTiUhKziJr\nZn9095PM7DWCsYo34u77NLVyM2sHXA8MBFYAc8xsurtXp7XZDriB4OSq981spxZ8DpE26+GHH+aV\nV16hS5cuXHzxxXHHEZE0jfVkU9/W72/G+vsBb7v7UgAzmwYcB1SntTkV+JO7vw/g7qs34/1E2pT6\n+nrGjBkDBGcWb7vttjEnEpF0OU98cvcPwqcXuPvS9AdwQZ7r3xVYnjb9XjgvXW9gBzObaWZzzOxH\n+YYvRkkZz1Q5CyuunPfccw+LFi2ivLycESNGNNpW27KwkpJT4pXPiU//D8g8k+KoLPM2J0Nf4Ahg\nG+AlM3vJ3d/JbJj+S11eXk55eTkVFRVZf9mrqqqynv3XGu2ziTNPrkzjx48vmjzans1v/+WXX3Lr\nrbcCMHHiRDp06NDk+lPvEUWeQrWvqKgoqjyNtU8tUyx5Um115nPxyDnik5mdT9Bj7QGkF7xOwGx3\nP73JlZsdDIx390Hh9GjA009+MrNRQEd3nxBO3wb81d3/lLEujfgkkuZ3v/sdI0eOpE+fPrz66quU\nlZXFHUmKkEZ8ildjRXY74GvAZGB02ku17r4mr5WblQFvEpz49AHwCnCKuy9Ka/MN4PfAIGBL4G/A\nye7+j4x1qciKhGpra+nRowerVq1i+vTpHHvssXFHkiKlIhuvnLuL3f0TM/sU2D914lJzuXu9mV0I\nPMVXl/AsMrPzgpd9irtXm9mTwEKgHpiSWWBFZGO//e1vWbVqFf379+eYY46JO46I5JDPDQKmAxe5\n+7LWiZQzh3qyIsDq1avp3r07tbW1VFVVMWDAgLgjSRFTTzZe+Qyr+DXgDTN71sweST2iDpZUSTnh\nQDkLqzVzTp48mdraWo488shmFVhty8JKSk6JVz5FdizBtbITCe7Ik3pIFkn54ilnYbVWzuXLl3PD\nDTcAcOWVVzZrWW3LwkpKTolXPsMqzjKzXYBvhbNecfePoo0lItnohuwiydJkT9bMTiI4K/hE4CTg\nb2Y2JOpgIrKx6upq7rjjDsrKypg0aVLccUQkD/kMRjEG+Faq92pmOwPPAA9GGUxENqYbsoskTz7H\nZNtl7B7+Z57LiUiBNHVDdhEpTvn0ZJ8Ir2O9L5w+GXg8ukjJlpTxTJWzsKLO2dQN2fOhbVlYSckp\n8WryOlkAMzseODScfN7dH4o0VfYMuk5W2qQZM2YwcOBAOnfuzJIlS3S/WGkWXScbr3x6sgAvEozG\n1ADMiS6OiKTTDdlFki2fs4vPITi7+IfAEOBlMzs76mAiohuyiyRdPj3ZnxOMX/xPADPbkaBne3uU\nwUTaOt2QXST58jlL+J9Abdp0bThPRCLUnBuyi0hxyqfIvkMwAMV4M6sEXgbeMrOfmtlPo42XPEkZ\nak05C6vQOdevX09lZSWw6Q3ZW6qtbsuoJCWnxCufIrsYeBhIndo7HXiX4ObtnSLKlVhJ+eIpZ2EV\nOufNN9/MsmXL6NOnD6eeempB1tlWt2VUkpJT4pXP2MUTAMxs23D606hDibRltbW1XHHFFUBwE4Cy\nsrKYE4lIS+VzdvE3zWw+8AbBLe/mmlmf6KOJtE26IbtI6chnd/EU4Kfu3s3duwGXALdGG0ukbVq9\nejW/+c1vgOC+sWYaQ0AkyfIpstu4+8zUhLtXAdtElkikDWvpDdlFpDjlc53sEjMbC9wTTp8OLIku\nUrIlZTxT5SysQuTcnBuy56MtbcvWkJScEq8mxy42s68BEwjGLnbgeWCCu/8rrzcwGwRcS9Brnuru\nV2W8PoDgjOVU4f6zu1+eZT0au1hK2jnnnMPUqVM5+eSTmTZtWtxxpERo7OJ4NVpkzawMuMrdf9ai\nlZu1A94CBgIrCMY9Huru1WltBgCXuPuxTaxLRVZKVnV1NX369MHMWLRoke4XKwWjIhuvRo/Juns9\nX919pyX6AW+7+1J3rwOmAcdlaadfAGnTUjdkHz58uAqsSAnJ55jsfDN7BHgAWJua6e5/zmPZXYHl\nadPvERTeTP3NbAHwPvBzd/9HHusWKQm6IbtI6cqnyHYkGKv4iLR5DuRTZPMxF/i6u39mZkcRjC7V\nO1vD9BMNysvLKS8vp6KiIusJCFVVVVlHZFF7tS+29vfcE5xT2LdvX95+++2sN2Uv5vxqX1ztU69L\nkXD3Rh/ATk21aWTZg4En0qZHA6OaWOZdYIcs8z0JZs6cGXeEvChnYbU057PPPuuAd+7c2VevXl3Y\nUBlKfVu2tqTkDP92tuhvuB6b/8h5TNbMjjGzVcBCM3vPzL7dgho+B+hpZt3MrAMwFHgk4312SXve\nj+BkrDUteK+ikJT/QSpnYbUkp3vr3pC9lLdlHJKSU+LV2O7iK4DD3L3azA4CfgU06+p4d683swuB\np/jqEp5FZnZe8LJPAYaY2flAHbAOOLklH0QkaXRDdpHS11iR/dLDS23c/W9m1qI77rj7E8AeGfNu\nSXt+A3BDS9YtklS6IbtI29BYke2Scb/Yjabd/ZroYomUNt2QXaRtaKzI3srG94vNnBaRFojihuwi\nUpxyFlkP7yMrzZOU8UyVs7CakzOKG7LnoxS3ZZySklPi1eTYxcVCwypKKaitraVHjx6sWrWK6dOn\nc+yxjY4mKrLZNKxivPK51Z2IFIhuyC7StqgnK9JKVq9eTffu3amtraWqqkr3i5VWoZ5svJrsyZrZ\nxWbW2QJTzWyemX2vNcKJlBLdkF2k7cnnfrKvuvu+ZnYkcB4wFrjH3fu2RsC0HOrJSmItX76cXr16\nsX79eubOnUvfvq369ZE2TD3ZeOVzTDb1wzmaoLi+gW5Nl1NShlpTzsJqKueECRNYv349J598cmwF\ntlS2ZbGkYeyNAAAgAElEQVRISk6JVz5Fdq6ZPUVQZJ8MR35qiDZWciXli6echdVYzurqau644w7K\nysqYNGlS64XKUArbspgkJafEK59b3Q0H9gOWeHA7uh2Bs6KNJVI6UjdkHzFihG7ILtLG5NOTdWAv\n4Mfh9DYE95gVkSbohuwibVs+RfZGoD9wSjhdiwb0F8nLZZddBsBFF12U9WbsIlLa8tldfJC79zWz\n+QDu/q/w3rAi0ogZM2bw9NNP07lzZ0aNGhV3HBGJQT492TozKyPYbYyZ7YxOfMopKeOZKmdhZeZs\n7Ruy5yOp27JYJSWnxCuf62RPI7iRel/gLmAI8L/u/kD08TbKoetkJTEeeughjj/+eLp06cLixYt1\nv1iJja6TjVeTu4vd/Q9mNhcYSHB97A/cfVHkyUQSSjdkF5GUfHqyO2SZXevuddFEyplDPVlJhDvv\nvJOzzjqL8vJy3nzzTd0vVmKlnmy88imyNcBuwL8IerLbAyuBD4Fz3X1uxBlTOVRkpeitX7+e3r17\ns2zZMu6++25+9KMfxR1J2jgV2Xjlc+LT08DR7r6Tu+8IHAU8ClxAcHlPo8xskJlVm9lbZpbzFEsz\n+5aZ1ZnZ8fmGFyk2cd2QXUSKUz5F9mB3fzI14e5PAf3d/WVgy8YWNLN2wPXAkUAf4BQz+0aOdr8E\nnsx8LWmSMtSachZWVVUVtbW1XHHFFQBceeWVlJWVxZxqY0nalkmQlJwSr3yK7AdmNsrMuoWPXwAf\nhpf1NHUpTz/gbXdfGh7DnQYcl6XdRcCDwEfNCV+MkvLFU86WMbOsj+985zt07tyZVatWAXDcccdh\nVlx76IptW+ainFJK8imypwJdgYfDx9fDeWXASU0suyuwPG36vXDeBmb2XwRnLN+E7u4jIiIlJJ9L\neFYT9DSzeacAGa4F0o/V5iy06Rd/l5eXU15eTkVFRdaLwquqqrL+TzPq9jU1NVmzx5UnV/uamhrG\njx9fNHmSsj1TMk/CO/PMMykvL98wPWHChA3rKZbtmZpXLHmS3r5Yt2fqdSkS7t7oA9gZ+DXwODAj\n9WhquXDZg4En0qZHA6My2iwJH+8SjIu8Ejg2y7o8CSorK+OOkBflbBmCkc82mZ+ZM1e7OBXbtsxF\nOQsr/D1s8u+1HtE88tld/AegGtgdmADUAHPyrOFzgJ7hsdwOwFDgkYwi3z187E5wXPYCd38ky7pE\nREQSJZ8iu6O7TwXq3H2Wu58NHJHPyt29HrgQeAp4A5jm7ovM7DwzG5FtkXyDF6ukjGeqnIWVhJxJ\nyAjKKaUln8EoXnb3g83sSeA6YAXwoLv3aI2AaTm8qawiUUudMZzH9yavdiJR02AU8crnVneXm9l2\nwCXA74HOwE8iTSUiIlICmuzJFgv1ZKUYqCcrSaOebLya7Mma2e4El/CUp7d392OjiyUiIpJ8+ewu\nfhiYCvwF3axdBKDoRnMSkeKUz9nFn7v7de4+Mzy7eJa7z4o8WUIl5SJw5Wx7krItlVNKST5F9ndm\nVmlm/c2sb+oRebKESsoXTzlbpr6+nn322QeAa665ZsMF55WVlbkGZCkaxbYtc1FOKSX57C7eG/gR\nwbWxqd3FTp7XyoqUkmnTprFw4UJ22203zj///LjjiEiRy6fIngh0d/cvog4jUszq6uoYN24cAJWV\nlXTs2DHmRCJS7PLZXfw6sH3UQUSK3e23387ixYvp3bs3w4YNizuOiCRAPj3Z7YFqM5sDrE/N1CU8\n0pasW7eOiRMnAjBp0iTat8/nqyMibV0+fykqI09RQpIynqlyNs/111/PihUr2H///RkyZMgmrxdL\nzsYkISMop5QWjfgk0oRPPvmE7t27s2bNGh5//HGOOuqouCOJ5E0jPsUrZ0/WzGrJflccI7g/YefI\nUokUkauvvpo1a9Zw2GGHMWjQoLjjiEiCqCcr0oiPPvqI7t27s3btWl544QUOOeSQuCOJNIt6svHK\n5+xikTbryiuvZO3atQwePFgFVkSaTT1ZkRyWLVtGr169+OKLL5g/fz777bdf3JFEmk092XipJ1tg\nSRlqTTmbNmHCBL744guGDh3aZIFNwvZMQkZQTiktKrIFlpQvnnI2rrq6mjvvvJOysrIN18c2Jgnb\nMwkZQTmltEReZM1skJlVm9lbZjYqy+vHmtmrZjbfzF4xMx34ktiNHTuWhoYGhg8fTq9eveKOIyIJ\nFemwNWbWDrgeGAisAOaY2XR3r05r9oy7PxK23xv4I7BnlLlEGjN37lwefPBBttxyS8aOHRt3HBFJ\nsKh7sv2At919qbvXAdOA49IbuPtnaZPbohvDS8zGjBkDwIUXXkjXrl1jTiMiSRZ1kd0VWJ42/V44\nbyNm9gMzWwT8BTg74kwiOc2aNYsnn3ySTp06MXr06LjjiEjCFcUo5+7+MPCwmR0KXA78v2zt0scK\nLS8vp7y8nIqKiqxjiFZVVWU9MSHq9ttvn/2GRXHlydV+++23Z/z48UWTpxi2p7tz2WWXAXDCCSdw\n/fXX573+JGzPmpoaqqqqiiZPrvYVFRVFlSdX+2LdnqnXpThEep2smR0MjHf3QeH0aIIhGa9qZJnF\nwLfcfU3GfF0nK5F69NFHOeaYY9hpp51YsmQJnTp1ijuSyGbTdbLxinp38Rygp5l1M7MOwFDgkfQG\nZtYj7XlfoENmgRWJWkNDw4ZjsZdddpkKrIgURKS7i9293swuBJ4iKOhT3X2RmZ0XvOxTgBPM7Azg\nC2AdcFKUmUSymTZtGgsXLmS33Xbj/PPPjzuOiJQIDasobV5dXR177rknixcv5rbbbmP48OFxRxIp\nGO0ujpdGfJI27/bbb2fx4sX07t2bYcOGxR1HREqIimyBJeWsPuUMrFu3bsOwiZMmTaJ9+5YdQUnC\n9kxCRlBOKS0qsgWWlC+ecgauv/56VqxYwf7778+QIUNavJ4kbM8kZATllNKiIitt1ieffMIvf/lL\nAK644gratdPXQUQKS39VpM26+uqrWbNmDYcddhiDBg2KO46IlCAVWWmTPvroI6655hoAJk+ejJlO\nvhSRwlORlTZp8uTJrF27lsGDB3PIIbq7oohEQ0W2wLKNM1qM2nLOZcuWceONNwJw+eWXF2SdSdie\nScgIyimlRYNRSJszfPhwbr/9doYOHcp9990XdxyRSGkwinipyEqbUl1dTZ8+fTAzFi1aRK9eveKO\nJBIpFdl4aXextCnjxo2joaGB4cOHq8CKSOTUk5U2Y968eRxwwAFsueWWvPPOO3Tt2jXuSCKRU082\nXurJSpuRuiH7hRdeqAIrIq1CRbbAkjLUWlvLOWvWLJ588kk6derE6NGjC7LOdEnYnknICMoppUVF\ntsCS8sVrSzndfUMv9mc/+xk77bTTZq8zUxK2ZxIygnJKaVGRlZL32GOP8eKLL7LTTjvxk5/8JO44\nItKGqMhKSWtoaGDMmDFAcEy2U6dOMScSkbZERVZK2rRp01i4cCFdu3bl/PPPjzuOiLQxKrJSsurq\n6hg3bhwAlZWVdOzYMeZEItLWRF5kzWyQmVWb2VtmNirL66ea2avh4wUz2zvqTFFKynimbSHn7bff\nzuLFi+nduzdnnnlmwTJlk4TtmYSMoJxSWiIdjMLM2gFvAQOBFcAcYKi7V6e1ORhY5O6fmNkgYLy7\nH5xlXRqMQvK2bt06evbsyYoVK7j//vs56aST4o4kEgsNRhGvqHuy/YC33X2pu9cB04Dj0hu4+8vu\n/kk4+TKwa8SZpA244YYbWLFiBfvvvz9DhgyJO46ItFFRF9ldgeVp0+/ReBE9B/hrpImk5H3yySdM\nnjwZgCuuuIJ27XTqgYjEo33cAVLM7DvAWcChudqkHwMpLy+nvLycioqKrMdGqqqqsl4srval337G\njBmsWbOGww47jEGDBsWeR+3VvjXbp16X4hD1MdmDCY6xDgqnRwPu7ldltNsH+BMwyN0X51iXjslK\nkz766CO6d+/O2rVreeGFFzjkkEPijiQSKx2TjVfU+9HmAD3NrJuZdQCGAo+kNzCzrxMU2B/lKrBJ\nkpT/QZZqzsmTJ7N27VoGDx7cqgU2CdszCRlBOaW0RFpk3b0euBB4CngDmObui8zsPDMbETYbC+wA\n3Ghm883slSgzRS0pX7xSzLls2TJuvPFGAC6//PKIEmWXhO2ZhIygnFJaIj8m6+5PAHtkzLsl7fm5\nwLlR55DSN2HCBL744guGDh3KfvvtF3ccERGN+CSlobq6mjvvvJOysjImTpwYdxwREUBFVkrEuHHj\naGhoYPjw4fTq1SvuOCIigIqslIB58+bxwAMPsOWWWzJ27Ni444iIbKAiW2BJGc+0lHKmbsh+4YUX\n0rVr14gTZZeE7ZmEjKCcUloivU62kHSdrGQza9YsKioq6NSpE0uWLGGnnXaKO5JIUdF1svFST1YS\ny9039GIvueQSFVgRKTrqyUpiPfrooxxzzDHstNNOLFmyhE6dOsUdSaToqCcbL/VkJZEaGhoYM2YM\nEByTVYEVkWKknqwk0n333cepp55K165defvtt+nYsWPckUSKknqy8VJPtsCSMtRaknPW1dVtuFSn\nsrKyKApsErZnEjKCckppUZEtsKR88ZKc8/bbb2fx4sX07t2bM888s9UzZZOE7ZmEjKCcUlpUZCVR\n1q1bt2HYxEmTJtG+fdHcEllEZBMqspIoN9xwAytWrGD//fdnyJAhcccREWmUiqwkxieffMLkyZMB\nuOKKK2jXTr++IlLc9FdKEuPqq69mzZo1HHbYYQwaNCjuOCIiTVKRLbCkjGeatJwfffQR11xzDQCT\nJ0/GrLiuSEjC9kxCRlBOKS26TlYS4Sc/+QnXXnstgwcP5tFHH407jkhi6DrZeKnIStFbtmwZvXr1\n4osvvmD+/Pnst99+cUcSSQwV2XipyEpRaO7uX/0uiORHRTZekR+TNbNBZlZtZm+Z2agsr+9hZi+a\n2edm9tOo84iIiLSWSK/kN7N2wPXAQGAFMMfMprt7dVqzfwIXAT+IMoskQ1M91GI74UlEpDFR92T7\nAW+7+1J3rwOmAcelN3D31e4+F/gy4iytIilDrSlnYSUhZxIygnJKaYm6yO4KLE+bfi+cV7KS8sVT\nzsJKQs4kZATllNKSqIFf069LKy8vp7y8nIqKiqzXq1VVVWX9EkTdvqamJmv2uPLkal9TU8P48eOL\nJk9KZiZtz8K1T80rljxJb1+s2zP1uhQJd4/sARwMPJE2PRoYlaNtJfDTRtblSVBZWRl3hLwUW07A\ns/2MM3Pmahe3Ytue2SQho7tyFlr4fYn0b70euR9R7y6eA/Q0s25m1gEYCjzSSHud1SIiIiUj0t3F\n7l5vZhcCTxEc/53q7ovM7LzgZZ9iZrsAfwc6AQ1mdjGwl7t/GmU2KQ4NDQ3ceuutcccQEYlE5Mdk\n3f0JYI+MebekPf8Q2C3qHK0lKeOZFkPOt956ixEjRjBr1qycbYohZz6SkDMJGUE5pbRoxCdpdXV1\ndVx99dWMHz+e9evX06VLFz766KNmrUO/CyL50YhP8dJdeKRVzZs3j379+nHppZeyfv16zjjjDP7x\nj3/EHUtEJBIqstIq1q1bx6hRo+jXrx8LFiygvLycJ598krvuuosdd9yxJWeui4gUvURdJyvJVFVV\nxbnnnss777yDmTFy5EgmTZrEtttuG3c0EZFIqchKZD7++GN+8YtfbDh7uE+fPkydOpWDDjoo5mQi\nIq1Du4sLLCkjrUSd8+GHH2avvfbi1ltvZYsttmDChAnMmzev2QVW27NwkpARlFNKi4psgSXlixdV\nzpUrV3LiiSfywx/+kA8++ID+/fuzYMECxo0bR4cOHZq9vra+PQspCRlBOaW0qMhKQbg7d9xxB3vt\ntRcPPvgg22yzDddddx3PP/88e+21V9zxRERioWOystmWLFnCeeedxzPPPAPAoEGDuPnmm+nWrVvM\nyURE4qWerLRYfX0911xzDXvvvTfPPPMMO+64I/fccw+PP/64CqyICOrJSgstXLiQc845hzlz5gBw\nyimncO2119KlS5eYk4mIFA/1ZAssKeOZtjTn+vXrGTt2LAcccABz5syha9euPProo9x7772RFNhS\n356tKQkZQTmltGjsYsnb7NmzOeecc6iurgbgggsuYPLkyXTu3DnmZCKSi8Yujpd2F0uTamtrufTS\nS7nxxhtxd/bYYw9uu+02Dj300LijiYgUNe0ulkY99thj7LXXXtxwww2UlZUxZswYFixYoAIrIpIH\n9WQlq1WrVjFy5EjuvfdeAA488ECmTp3KPvvsE3MyEZHkUE9WNuLu/OEPf2DPPffk3nvvZauttuI3\nv/kNL730kgqsiEgzJarImtkmj2LK0tgjKTkHDx7M6aefzj//+U+OOOIIXnvtNS655BLat49np0dS\nhq5LQs4kZATllNKSqCIr0fvrX//K9ttvz9SpU3nmmWfo0aNHrHmS8ocsCTmTkBGUU0pL5EXWzAaZ\nWbWZvWVmo3K0uc7M3jazBWa2X651FeNNuzNvJj5s2LBE5xwyZAiLFi3i7LPPjrUHnlJTUxN3hLwk\nIWcSMoJySmmJdB+gmbUDrgcGAiuAOWY23d2r09ocBfRw915mdhBwM3BwlLmilJQvXq6cDzzwQOsG\naULSt2cxSUJGUE4pLVEfaOsHvO3uSwHMbBpwHFCd1uY44G4Ad/+bmW1nZru4+4cRZytqDQ0N1NfX\nb/Jo7nwREYlP1EV2V2B52vR7BIW3sTbvh/PyKrIHHHDA5uRrNeXl5c0qmCIiknyJuk422zHCefPm\nxZDkK9kyZZu3dOnS1oiTU745i+E4bKZizJRNEnImISMop5SOqIvs+8DX06a7hvMy2+zWRBuNvSki\nIokT9dnFc4CeZtbNzDoAQ4FHMto8ApwBYGYHAx+39eOxIiJSGiLtybp7vZldCDxFUNCnuvsiMzsv\neNmnuPvjZna0mb0DrAXOijKTiIhIa0nMre5ERESSpihGfNqcASvyWTamjPunzZ9qZh+a2cKo8m1G\nzv3CeV3NbIaZvWFmr5nZj4s055Zm9jczmx/mrCzGnGmvtTOzeWaWeZgk7pzpv581ZvZquE1fKbKc\n6d/17czsATNbFP6eHlRsOc2sd7gd54X/fhLVd2kzt+VPzOx1M1toZn8ID+dJFDJHAmrtB0Ghfwfo\nBmwBLAC+kdHmKOCx8PlBwMv5Lht3xnD6UGA/YGERb8v/APYLn28LvBnFtizQ9tw6/LcMeBnoV4w5\nw3k/Af4PeKQYf+7h9BLga1H+bhYo553AWeHz9kDnYsyZsZ4VwG7FlBH4r/Bn3iGcvh84I+qff1t9\nFENPdsOAFe5eB6QGrEi30YAVwHZmtkuey8adEXd/AfhXBLkKltPdV7r7gnD+p8AiguuViypnOP1Z\n2GZLgj+2UR3z2KycZtYVOBq4LaJ8BckJGK2zV6vFOc2sM3CYu98Rvvalu/+72HJmtPkusNjdl1N4\nm5uxDNjGzNoDWxP8Z0AiUAxFNtuAFZl/3HO1yWfZuDK+n6VN1AqS08zKCXrefyt4wuwZmpUz3AU7\nH1gJPO3uc4oxJ/Bb4OdE95+AXBmam9OBp81sjpmdG1nKzcu5O7DazO4Id8VOMbOtijBnupOB+wqe\nLvv7553R3VcAVwPLwnkfu/szEeVs84qhyLaErpmNiJltCzwIXBz2aIuOuze4+/4E11QfZGZ7xZ0p\nk5kNBj4M9w4Yxf07e4i79yXodf+PmR0ad6As2gN9gRvCrJ8Bo+ONlJuZbQEcCxTXYOCAmW1P0Mvt\nRrDreFszOzXeVKWrGIrs5gxYkc+ycWdsTZuVM9x19CBwj7tPL9acKeHuwpnAoAgypjK0NOchwLFm\ntoSgN/MdM7u7CHPi7h+E/64CHmLToU+LIed7wHJ3/3s4/0GColtsOVOOAuaG2zQKm5Pxu8ASd1/j\n7vXAn4FvR5RT4j4oTHBsIHUAvwPBAfw9M9oczVcH8A/mqwP4TS4bd8a018uB14p1W4bTdwPXFPnP\nfCdgu/D5VsBzwNHFljOjzQCiPfFpc7bn1sC24fNtgNnA94otZzg9C+gdPq8ErirGnOG8+4BhRfoz\n7we8BnQk2MNyJ/A/UWVt64/YA4Q/9EEEZ7O+DYwO550HjEhrc334S/Uq0LexZYsw470EJxasJzgO\nclYR5dw/nHcIUB9+WecD84BBRZSzbzhv7zDbAmAhMKZYfzfTXo+0yG7m9tw97Wf+WpTfoc3dnsC+\nBKPILSDofW1XpDm3BlYBnYp4W1YSnNy4ELgL2CLKrG35ocEoREREIlIMx2RFRERKkoqsiIhIRFRk\nRUREIqIiKyIiEhEVWRERkYioyIqIiERERVY2YWYNZvbrtOlLzGxcK2e4w8yOD5/fambfKMA6a3PM\nrw/Hw33NzO43s45NrOeFPN7r4lzrMbN9zeyo/FJvtNzM8NZmC8zseTPr1UT7KU1tNzM7LlcbM6s0\ns/fSbts2LxyoX0TypCIr2awHjjezHVqysJmVFTKMu5/r7tWFWFWO+Wvdva+77w3UAf/dRJ58xvYd\nSTAoQTb7EYzG0xKnuPt+BKNz/aaxhu4+Io/t9gOgTyOvXxNum/3Dfze6803mz7o5P3szK+bxnEUK\nQkVWsvkSmAL8NPMFM+tmZs+Gvamnw9u5pXqeN5nZS8BVYS/oTjN7zszeNbMfmtlV4U2iH0/9MTaz\nsRbchH2hmd2cLUzYg+sb3n3njrDtq2Z2cfh6dzP7a3gXmVlm1jucX25mL4ZtJ+X52Z8HeobL/zTs\n3S5MvVc4vzb8d0CYLXUj8XvC+RcRDLw+08yezfgsWwATgZPCnuGJZvY1M3sozPmime3dSL5UYXoO\n6BGuc2C4rlfN7LbwPTZst1RmM7s8/Lm9aGY7m1l/gkHsfxUuv3sj75f+GYaZ2fTwsz0TbofnzGw6\n8EaubRf+7lSb2V1m9hrBWLoipS3uIaf0KL4H8G+CG7e/C3QCLgHGha89ApwePj8LeCh8fgdpQwcS\nDNv2HMF/5PYB1hKOiUswJN6x4fPt05a5Gxictr7jw+czCQaD7ws8lda+c/jvM0CP8Hk/4Nnw+XTg\ntPD5BcC/c3ze2vDf9sDDBEPT9SUYiq4jwZi+rwP7prZP+O8AgvsE/ydBMXoR+Hb4Ws4boQPDgOvS\npq8DxobPvwPMz7HcTL4aDvHnBOPjbkkwVGfq898F/DhL+wbCMZ6Bq4DLMrdzlverJBiYfx7BsIvP\npuVfxldjSA8AaoGvh9NZtx3BOLtfAt+K+3dcDz1a66GerGTlwW3u7gIuznipP1/dI/MegjGPUzJv\n6/VXd28gGBO3nbs/Fc5/jeCGCQADzexlM1tIUGAa23W5BNjdzH5nZkcCtWa2DcEdRB6w4B6ztwCp\nG1MfQnAz61TWXLYys3nAK0ANMBU4lOA/EJ+7+1qC/xgclmXZV9z9A3d3gjF1U5+rObe3OzSVz91n\nAjtYcMvBbP4QZu0P/AzYg+COKovD1+8CDs+y3Hp3fzx8PjctZ1PSdxcPTJv/tLt/kjb9irsvS/s8\nubbdUo/u/r8iRad93AGkqP2OoBdzR9q8xga7XpsxvR7A3d3M6tLmNwDtzWxL4AaC3tYKM6sk6P1k\n5e4fm9m+wJEEx01PBH4C/MuDe4xuskha3sYK3meZyzfjcOH6tOf1ZPlOmdkPCHqFDpyTI+cmzOwJ\noAvwd3cfEc4+1d3np7XZkfyKefr2z5qzmTJ/1pnT+S4nUtLUk5VsDMDd/wX8ERie9tqLwCnh89MJ\njmHmvc4MHQkKzD/DntuQRlcQFJQyd38I+F+C4lwLvGtmQ9La7RM+nZ2W9bRmZnse+IGZdQx7yz8k\n2P2dq32mfwOdAdz9Yf/qxKF5BLtW08/SfZ5gW2JmFcBqd//U3QeFy4xIa5v53m8C3cysezj9I6Aq\nz89Iliz5LteYbNsu9Xuik52kTVGRlWzSe1ZXAzumzfsxcJaZLSAoXBdnWaapdQYzgt2NtxGcLPNX\ngt212dqnnu8KVIW7he8BRofzTweGhyf1vE5wMg8EZ/j+j5m9SnDctDnZ5hPcZ3MO8BIwxd0X5mqf\nZf6twBOZJz6FZgJ7pU58AsYDB4Q5ryQ45plvzvUEx8YfDJevJ9hlntk+V+ZpwM/NbG6OE59G2saX\n8Hw9S5vMTNm23atN5BApSbrVnYiISETUkxUREYmIiqyIiEhEVGRFREQioiIrrS4867TKAv9pZn/M\nY5lLWyNbI+//rrVwmMks6+oWjniEmR1gZtc20faUXK8XmpntZGZ/ba33Eyl1KrISh7OBP3ngA3c/\nKY9lLmvum5hZIX+/m3WGYB7v7QDuPtfdRzbSbnfg1Oa89+Zw99XAinDIRRHZTCqyEofTCIY8zOzV\nDTOzP1kwDvGbZvbLcP5kwlGZ0sYHPs2CMY/nWTBmsoXza83sN+FlPgc3Mq7vt8xsdnjZz8tmtk34\n/r9PhTSzv5hZavQkS5v/kAXjJL9mZuekzd/ovdM/cNhjXRC+9j9p8weY2V/SnqculZkbXmM6GTg0\nnHdxuL2eM7O/h4+D05bdZBzlRj5rOzP7VbgNF5jZuWlxpxNetysimynucR31aFsPYAtgRdp0N2Bh\n+HwY8A7BuMlbEgxxuGv42r/TlvkGwRjKZeH0DXw1nnIDcEL4POu4vmGGxXw1ru+2QBmbjin8F+Dw\n8Pm7wA7h8+3DfzsSDBH5tcz3zvK5XwUOCZ//Ku0zDyAc8zn8TP3D51sT/Cd4w+tp79khfN4TmJO2\nnk3GUW7ks57LV+MXdyC4prVbOP1fqXx66KHH5j3Uk5XWthPwcSOvP+vBaEfrgX8QFGHYeKSggQSD\n0M8Je4ZHEOxWhWAwhj+Hz3ON67sHQaGfB8E4ze5e34zPMDIcjONlgjvJpO7r+mXae29gZtsRDKY/\nO5yVaxzl2cBvLbiLz9c8GPc5UwfgNgvGen4A2DPttWzjKOf6rN8Dzgi339+AHdI+x0c0PniHiORJ\nY4w6k+0AAAH1SURBVBdLa1tHI+MTk8dYwAQF9y53H5Nt/WGRSW+bTbb5X7LxIZRNcprZAIKifpC7\nrzezmWntPs9472Zx96vM7FFgMDDbzL6XpdlPgJXuvo8Ftwtcl/Zarm2X7bMacJG7P53ltY4Z6xWR\nFlJPVlqVu38MlJlZh2Yu+oV9dUPwZ4EhZrYzgAX3Y90tfC29oOQa1/dN4D/M7IBw+W3DddcA+4Vn\nPe9GcNu8TNsR3JBgvZl9g42PvWYt6B4MH/mxmX07nJX1eKeZdXf3N9z9VwS7b7/BpmMLbwd8ED4/\ng2DXb2NyfdYngQvMrH04v5eZbRUu05vg9nQispnUk5U4PEVwO7QZTbRL7xVOAV4zs7nu/iMzGws8\nFZ7F+wXByUTL05cJC2FqXN8ygsJ1i7vXmdnJwPVhYfkM+K67zzazGoKxlBcR3BIuM8sTwH+b2RsE\nBeylHHkznQ3cbmYN4efPZqSZfYegF5oaz9mB+nC37p0Ex5//bGZnhFly3dUmdfZy1s9KMGZ0OTAv\nPGnsI+AH4bLfAR5r5LOISJ40drG0OjPbHxjp7rkGwpcYmVkVcJxvfL9YEWkB7S6WVufBXVpmpi67\nkeJhZjsR3KhdBVakANSTFRERiYh6siIiIhFRkRUREYmIiqyIiEhEVGRFREQioiIrIiISkf8PRfl5\n8uJCYN8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x11bbf1150>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.yticks(np.linspace(0, 1, 11))\n",
    "\n",
    "plot_cumulative_error_distribution(\n",
    "    [errors],\n",
    "    error_range=(0., 0.09, .01), \n",
    "    legend_entries=['MDM'], \n",
    "    line_colour=[ 'k'], \n",
    "    marker_edge_colour=['k'],\n",
    "    marker_style=['s'],\n",
    "    figure_size=(6, 5),\n",
    "    marker_size=10,\n",
    "    x_label='Normalised Point-to-Point Error\\n(interocular distance)',\n",
    ");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def AUC(errors, max_error, step_error=0.0001):\n",
    "    x_axis = list(np.arange(0., max_error + step_error, step_error))\n",
    "    ced = np.array(compute_cumulative_error(errors, x_axis))\n",
    "    return simps(ced, x=x_axis) / max_error, 1. - ced[-1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The AUC @ 0.08 is 41.01 and 8.67% of the errors are failures\n"
     ]
    }
   ],
   "source": [
    "auc, failures = AUC(errors, max_error=0.08)\n",
    "print('The AUC @ 0.08 is {:.2f} and {:.2f}% of the errors are failures'.format(auc * 100, failures * 100))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}