gistfile1.txt

{
 "metadata": {
  "name": "RANSAC similarity transform"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "source": [
      "RANSAC similarity transformation: Oct 21, 2013"
     ]
    },
    {
     "cell_type": "markdown",
     "source": [
      "In this notebook, we combine RANSAC with least squares similarity\n",
      "transformation to learn a similarity transform that maps a set of 2D\n",
      "points to another set of 2D points.\n",
      "\n",
      "RANSAC is based on [the implementation in SciPy's\n",
      "cookbook](http://wiki.scipy.org/Cookbook/RANSAC) while the similarity\n",
      "transformation is based on scikit-image's representation.\n",
      "\n",
      "Consider the following scenario: We have a bunch of points in some\n",
      "space. We learn an embedding of these points, which corrupts our\n",
      "points in three ways:\n",
      "\n",
      "- Many points are completely incorrect.\n",
      "- Correct points are distorted by a rigid-body (similarity) transform,\n",
      "- Correct points are also displaced by some noise due to the\n",
      "  Voroni-like nature of the triplet embedding. Here, we model that as\n",
      "  Gaussian noise.\n",
      "\n",
      "Our goal is to find the similarity matrix that maps the correct points\n",
      "to the correctly-mapped points while ignoring incorrect points."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "actual_data = np.random.randn(50, 2)\n",
      "# Our embedding learns a noisy version of a mapped version of the data...\n",
      "mtx = np.array([[ 0.61900171, -0.42348186, 3],\n",
      "                [ 0.42348186,  0.61900171, 2]])\n",
      "embedded_data = mtx.dot(np.hstack([actual_data,np.ones((50,1))]).T).T\n",
      "embedded_data = embedded_data + 0.1*np.random.randn(50,2)\n",
      "\n",
      "# ...but many datapoints are incorrect!\n",
      "N_CORRECT=25\n",
      "embedded_data[N_CORRECT:,:] = 2*np.random.randn((50-N_CORRECT),2)\n",
      "\n",
      "fig,(ax1,ax2,ax3) = subplots(1,3, figsize=(12,4),sharex=True,sharey=True)\n",
      "ax1.scatter(actual_data[:,0], actual_data[:,1])\n",
      "ax1.set_title(\"Actual Data\")\n",
      "ax2.scatter(embedded_data[:N_CORRECT,0], embedded_data[:N_CORRECT,1],color='g')\n",
      "ax2.scatter(embedded_data[N_CORRECT:,0], embedded_data[N_CORRECT:,1],color='r')\n",
      "ax2.set_title(\"Learned Embedding\")\n",
      "\n",
      "# So, learn a model:\n",
      "model = learn_ransac_similarity_xform(\n",
      "    actual_data,\n",
      "    embedded_data,\n",
      "    n=3, # Min required to fit the model\n",
      "    k=100, # Max number of iterations allowed in the algorithm\n",
      "    t=0.5, # Threshold value for determining when a data point fits a model\n",
      "    d=10, # Number of close data values to assert that the model fits the data\n",
      ") \n",
      "ax3.scatter(model.transform(actual_data)[:,0], model.transform(actual_data)[:,1])\n",
      "ax3.set_title(\"Transformed Data\") \n",
      "print \"Actual matrix was\"\n",
      "print mtx\n",
      "print \"RANSAC found this matrix:\"\n",
      "print model._matrix\n",
      "print \"Difference:\", np.mean(np.abs(model._matrix-mtx)),\"which should be about 0.1\""
     ],
     "language": "python",
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Actual matrix was\n",
        "[[ 0.61900171 -0.42348186  3.        ]\n",
        " [ 0.42348186  0.61900171  2.        ]]\n",
        "RANSAC found this matrix:\n",
        "[[ 0.60684115 -0.42554073  2.99253736]\n",
        " [ 0.42554073  0.60684115  2.00898689]]\n",
        "Difference: 0.00748139615944 which should be about 0.1\n"
       ]
      },
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAr8AAAEKCAYAAAAb9j21AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3XdcleX/x/HXOezDRgQRFHeK5EhLzYW5TTMX7m1p2dDS\n6qc5szRHlpWlfs09Wq7EVSpaLjRNzb0QBdzIOsDhnHP//kCPkqAo49zI5/l48Kj7Pvd430c+cHGf\n674ujaIoCkIIIYQQQhQBWmsHEEIIIYQQoqBI41cIIYQQQhQZ0vgVQgghhBBFhjR+hRBCCCFEkSGN\nXyGEEEIIUWRI41cIIYQQQhQZ0vgtQsaPH0/v3r2tHUOIIqVfv36MGTMm388TGRmJVqvFbDbnyfHK\nlCnD1q1bs3wtPDycUqVKWZaDg4PZuXNnnpxXiIJw9epVGjVqhJubGyNHjrR2nEzyupbFg6TxW4BC\nQkLw8vLCYDDkaPuFCxfSsGHDPDu/RqPJ9rXw8HC0Wi2urq64urpSqlQpunbtyoEDB3J8fGlci/z2\nsAaZWmk0mmxrb+HChdjY2FjqztXVFTc3N65cuVLAKR/0sNz/9e+//9KoUaN8TiSeFi4uLpbvd61W\ni06nsyyvWLGiQDLMnTsXHx8fEhISmDZtWoGcM6+UKVMGnU6Hm5sbnp6e1K9fnzlz5pDTaRukcS2N\n3wITGRlJREQEPj4+rFu3zioZHlUY/v7+JCYmkpiYyN69e6lcuTINGzZk27ZtBZRQiId7nAbZk8qP\nXwgPq7369etb6i4xMZGEhARKlCiR5xmEUIukpCTL93tgYCDr16+3LHfv3t2yndFozLcMFy9epEqV\nKk+0b37mygmNRsP69etJSEggKiqKjz76iM8//5yBAwc+1nGK8hxn0vgtIIsXL6ZZs2b07t2bRYsW\nZXrt0qVLdOzYER8fH7y9vXn77bc5efIkQ4YMYc+ePbi6uuLl5QVk3D2eP3++Zd//3h1+9913KV26\nNO7u7tSuXZu//vrrifL6+/szYcIEBg0axIcffvjI42/atInJkyfz448/4urqSs2aNQFYsGABQUFB\nuLm5Ub58eebOnftEeYR4GEVRmDJlChUqVMDb25uuXbsSFxdneb1Lly74+fnh4eFB48aNOX78uOW1\nfv368cYbb9CmTRtcXFzYvn07ZcqUYcaMGVSvXh0PDw+6detGWlqaZZ/169dTo0YNy12Xo0ePWl47\ndOgQzz33HG5ubnTr1o3U1NRHZs9OmTJlmD59OtWqVcPV1ZWBAwdy9epVWrdujbu7O82bN+f27duZ\n9pk/fz7+/v6ULFmSGTNm5Pg9WrJkCYGBgXh7e/PZZ59lOmZKSgr9+vXDy8uLqlWrsn///gdy3v0j\nefz48YSGhtK3b1/c3NwIDg7m77//tmx78OBBatasiZubG6GhoXTt2rVAuoUI9QsPDycgIICpU6fi\n5+fHwIEDuX37Nm3btsXHxwcvLy/atWtHdHS0ZZ+QkBDGjh1LgwYNcHNzo2XLlty8eROA1NRUevXq\nhbe3N56enrzwwgtcu3aNfv36sXjxYqZOnYqrqyvbtm3DYDAwbNgw/P398ff3Z/jw4ZZPaf+ba8CA\nAUyYMIEuXbrQu3dv3NzcqFatGmfOnGHy5Mn4+voSGBjI77//bskZHx/PwIEDKVmyJAEBAYwZM8by\nh7bZbGbEiBEUL16c8uXLExYWluP3zNXVlXbt2vHjjz+yaNEijh07BkBYWBg1a9bE3d2d0qVLM2HC\nBMs+dz+l8fDwwNXVlX379nHu3DleeuklvL29KV68OL169SI+Pv4J/yXVTxq/BWTx4sV07dqV0NBQ\nNm/ezLVr1wAwmUy0bduWsmXLcvHiRaKjo+nevTuVK1dmzpw51KtXj8TERG7dugU8+s7XCy+8wOHD\nh4mLi6NHjx506dIlx90sstKhQwcOHjxISkrKQ4/fqlUrRo0aRbdu3UhMTOTQoUMA+Pr6EhYWRkJC\nAgsWLGD48OGW14TIK7NmzWLdunXs3LmT2NhYPD09GTp0qOX1l19+mbNnz3L9+nWee+45evbsmWn/\nFStWMGbMGJKSkmjQoAEajYaff/6ZzZs3c+HCBY4cOcLChQuBjMbtwIEDmTdvHrdu3WLw4MG88sor\npKenYzAYePXVV+nbty9xcXF06dKFX3/99YnvVms0GlatWsXWrVs5deoU69evp3Xr1kyZMoVr165h\nNpuZNWtWpn3Cw8M5e/YsW7Zs4fPPP7d0E3nYe3T8+HHefPNNli1bRkxMDDdv3uTy5cuWY06YMIEL\nFy5w/vx5Nm/ezKJFizJd03+v77fffqN79+7Ex8fzyiuv8NZbbwFgMBjo0KEDAwYMIC4uju7du7Nm\nzZp8v5svCo+rV68SFxdHVFQUc+bMwWw2M3DgQKKiooiKisLJycny/XTXihUrWLhwIdeuXcNgMDB9\n+nQAFi1aREJCApcvX+bWrVvMmTMHJycnFi5cSM+ePfnwww9JTEzkpZdeYtKkSURERHD48GEOHz5M\nREQEkyZNyjLX3LlzURSF9evX06dPH+Li4qhZsybNmzcHICYmhjFjxjB48GDL/v369cPe3p5z585x\n6NAhtmzZwv/+9z8gowtGWFgY//zzDwcOHOCXX3557Jp4/vnnCQgIsNyQcnFxYenSpcTHxxMWFsZ3\n333H2rVrAfjzzz+BjAZ5YmIiderUAWD06NHExsZy4sQJLl26xPjx4x8rQ6GiiHz3559/Ko6OjkpC\nQoKiKIpSvXp1ZebMmYqiKMru3buV4sWLKyaT6YH9FixYoDRo0CDTupCQEGX+/PkP3eZ+np6eypEj\nRxRFUZRx48YpvXr1ynK77du3KwEBAQ+sP3HihKLRaJSYmJhcHf+uV199Vfnqq68euo0Q2SlTpoyy\ndevWB9ZXqVIl0/qYmBjFzs4uy7qKi4tTNBqNpR779u2r9O3b94HzLFu2zLL8wQcfKEOGDFEURVGG\nDBmijBkzJtP2zzzzjLJjxw5lx44dSsmSJTO99uKLLz6w/V0LFixQbG1tFQ8PD8tXhQoVMuVYvny5\nZblTp07Km2++aVn++uuvlVdffVVRFEW5cOGCotFolFOnTmXKPXDgQEVRFKVy5cpZvkdGo1GZMGGC\n0r17d8trycnJir29vWX7cuXKKZs3b7a8Pnfu3Ew/L+7/dxk3bpzSvHlzy2vHjh1TnJycFEVRlB07\ndij+/v6Z3oMGDRpk+/6Ip9/93zvbt29X7O3tlbS0tGy3P3TokOLp6WlZDgkJUT799FPL8uzZs5VW\nrVopiqIoP/zwg/Liiy9afkfdr1+/fsrHH39sWS5fvryyceNGy/LmzZuVMmXKZJtr3LhxSosWLSzL\n69atU1xcXBSz2awoiqIkJCQoGo1GiY+PV65cuaI4ODgoKSkplu2XL1+uNGnSRFEURWnSpIkyZ84c\ny2tbtmxRNBpNlj+//vue3a9u3brKZ599luU+7777rjJ8+HBFUe79rMju+IqiKKtXr1Zq1qyZ7euF\nndz5LQCLFi2iRYsWuLq6Ahkfwd7t+nDp0iUCAwPRavPmn2L69OkEBQXh4eGBp6cn8fHx3Lhx44mP\nFx0djUajwcPD44mOv3HjRurWrUuxYsXw9PRkw4YNlo+khMgrkZGRdOjQAU9PTzw9PQkKCsLW1par\nV69iMpn46KOPqFChAu7u7pQtWxbA8n2r0WgyjVxw1/39bp2cnEhKSgIy+grOmDHDci5PT08uX75M\nbGwsMTEx+Pv7ZzpOYGDgQ7s21K1bl7i4OMvXmTNnMr3u6+ubKcf9y46OjpZcd91/LaVLlyYmJsaS\nO7v3KDY2loCAAMt+Op2OYsWKWZZjYmIeOO7D3J9Rp9ORmpqK2WzO8v0pVapUke57KDIrXrw49vb2\nlmW9Xs/gwYMpU6YM7u7uNG7cmPj4+EzfM9nVau/evWnZsiXdunXD39+fDz/8MNv+ujExMQQGBlqW\n76+drHIB+Pj4ZDqvt7e35Y6tk5MTkNG/+eLFi6Snp+Pn52epvyFDhnD9+nUAYmNjH6u+shMdHW3p\nIrlv3z6aNGmCj48PHh4ezJkz56G/e69evUq3bt0ICAjA3d2d3r17P9W/q6Xxm89SUlL46aef2LZt\nG35+fvj5+TFjxgwOHz7MkSNHKFWqFFFRUZhMpgf2zepjD2dnZ5KTky3L9z8V/ueffzJt2jR+/vln\nbt++TVxcHO7u7rn6xbJ69Wpq1aqFk5PTI4//37xpaWl06tSJDz74gGvXrhEXF0ebNm3kF53Ic6VL\nl2bTpk2ZGpF6vR4/Pz+WL1/OunXr2Lp1K/Hx8Vy4cAF4/Ic97n5/ly5dmtGjR2c6V1JSEl27dsXP\nzy9Tf0TIaHTm5cf6j8odFRWV6f/vNjaze49KliyJn58fly5dsuyn1+sz/eLz8/N74LhPIqv3Jyoq\nSro9CIv/fi/MmDGD06dPExERQXx8PDt27EBRlBzVr62tLWPHjuXYsWPs3r2b9evXs3jx4iy3LVmy\nJJGRkZblqKgoSpYsmW2ux/meLVWqFA4ODty8edNSe/Hx8ZZnBfKivvbv3090dDQNGjQAoEePHrz6\n6qtcvnyZ27dvM2TIEEsf46yyjxo1ChsbG/7991/i4+NZsmTJUz0ahDR+89maNWuwtbXlxIkTlr5E\nJ06coGHDhixevJg6derg5+fHRx99hF6vJzU1ld27dwMZd08uX75Menq65Xg1atRg1apVpKSkcPbs\nWebPn2/5Rk5MTMTW1hZvb28MBgMTJ04kISHhsTMrikJ0dDQTJkxg/vz5lodfHnX8EiVKEBkZafmh\nZDAYMBgMeHt7o9Vq2bhxI1u2bHni91IIyPi+Sk1NtXwZjUaGDBnCqFGjLL80rl+/bhlVJSkpCQcH\nB7y8vEhOTmbUqFGZjpfTRvDd7V577TW+//57IiIiUBSF5ORkwsLCSEpK4sUXX8TW1pZZs2aRnp7O\nqlWrHng4LL9NmjSJlJQUjh07xsKFC+natSvAQ9+jzp07s379enbt2oXBYGDs2LGZfvGFhoYyefJk\nbt++zeXLl/n666+fKFu9evWwsbHhm2++wWg0snbt2gJ/f0ThkpSUhJOTE+7u7ty6dSvTg1t3ZVfD\n27dv5+jRo5hMJlxdXbGzs8PGxibLfbp3786kSZO4ceMGN27cYOLEiQ8duvNx/nj28/OjRYsWvPfe\neyQmJmI2mzl37pxlbOzQ0FBmzZpFdHQ0cXFxTJky5ZHHvHv+hIQE1q9fT/fu3enduzdVq1YFMt43\nT09P7O3tiYiIYPny5Za2QvHixdFqtZw7d85yvKSkJJydnXFzcyM6OrrQDf/2uKTxm88WL17MgAED\nCAgIwMfHBx8fH3x9fXnrrbdYvnw5kPFwyNmzZyldujSlSpXip59+AqBp06ZUrVqVEiVKWD5eGT58\nOPb29vj6+tK/f3969eplOVerVq1o1aoVlSpVokyZMjg5OWX6+ORhD8tpNBpiYmIsYy2+8MILHDt2\njB07dtCsWbMcHb9Lly4AFCtWjNq1a+Pq6sqsWbMIDQ3Fy8uLFStW0L59+zx8d0VR1KZNG3Q6neVr\n4sSJvPvuu7zyyiu0aNECNzc36tWrR0REBAB9+vQhMDAQf39/goODqVev3gMPaz3qLs7929SqVYt5\n8+bx1ltv4eXlRcWKFS13k+zs7Fi1ahULFy6kWLFi/PTTT3Tq1Omhx707osv9X/ePjpDVPtll12g0\nNG7cmAoVKtCsWTNGjhxpqd+HvUdBQUF8++239OjRg5IlS+Ll5ZXpY9hx48YRGBhI2bJladWqFX36\n9Hnoz5Ls7pLZ29uzatUq5s+fj6enJ8uWLaNt27YPfJwsiq7/fu8MGzaMlJQUvL29efHFF2nduvVD\n78Le//139epVunTpgru7O0FBQYSEhFgatP/9Pv3444+pXbs21apVo1q1atSuXZuPP/4421wP+z7P\nannx4sUYDAaCgoLw8vKiS5culk9uX3vtNVq2bEn16tWpXbs2nTp1euTPpHbt2uHm5kbp0qWZPHky\n77//PgsWLLC8Pnv2bMaOHYubmxuffPKJ5Y9gyOiKNHr0aOrXr4+XlxcRERGMGzeOgwcP4u7uTrt2\n7XKUoTDTKHnwGbTJZKJ27doEBATw22+/5UUuIYQQRUCdOnV488036du3r7WjCCGKiDy58/vVV18R\nFBT0VP+VIIQQIvd27tzJlStXMBqNLFq0iH///ZdWrVpZO5YQogjJdeP38uXLbNiwgUGDBsmDTEII\nIR7q1KlTlglCZs6cyS+//JJpdAghhMhvue720KVLF0aNGkVCQgLTp0+Xbg9CCCGEEEK1cnXnd/36\n9fj4+FCzZk256yuEEEIIIVQvV3d+R40axZIlS7C1tSU1NZWEhAQ6deqUaRy9ChUqZBpOQ4iirnz5\n8pw9e9baMbIk9SpEZmquV5CaFeK/clSzeTVVXHh4uNK2bdsH1ufhKZ7YuHHjrB3B6hmsfX7JcI8a\naiI7asimhn8jyaCODNY+v6KooyYextr51PBvJBkkw/1yUhN5Os6vjPYghBBCCCHUzDavDtS4cWMa\nN26cV4cTQgghhBAizxWJGd5CQkKsHcHqGax9fskgckoN/0aSQR0ZrH1+8Whq+DeSDJLhceXJDG8P\nPYFGIyNBCHEfNdeEmrMJYQ1qrwm15xOioOWkJorEnV8hhBBCCCFAGr9CCCGEEKIIkcavEEIIIUQR\nkZKSQr9+b1C8eBkqVnyOzZs3WztSgZM+v0IUMDXXhJqzCWENaq8JtecT6tOjx0BWr75Fauo04DQ6\nXT927/6d6tWrWztanpA+v0IIIYQQwmLdurWkpn4HVADaYDD0YtOmTdaOVaCk8SuEEEIIUUQ4ObkA\nlyzLdnaXcHFxsV4gK5DGrxBCCCFEETF9+ifodB2AT7C3703x4v/Sq1cva8cqUNLnV4gCpuaaUHM2\nIaxB7TWh9nxCncLDw9m4cQvFinnw+uuv4+HhYe1IeSYnNSGNXyEKmJprQs3ZhLAGtdeE2vMJUdDk\ngTchhBBCCCHuI41fIYQQQghRZEjjVwghhBBCFBnS+BVCCCGEEEWGNH6FEEIIIUSRIY1fIYQQQoin\nzOnTp6lVqzFubr48//xLnDt3ztqRVEOGOhOigKm5JtScTQhrUHtNqD2fsA69Xk/ZslW5fv09FKUT\nWu1KfH1nc/78vzg6Olo7Xr6Soc6EEELkqQtxF/hy75d8G/Et15OvWzuOECILx44dIzXVA0V5GyiJ\n2fweSUkOnDp1ytrRVMHW2gGEEEIUDodiD9FoYSMMJgNajZbxO8bzz+B/8Hfzt3Y0IcR93NzcSE+/\nCiQDzkAi6enXcXd3z3L7uLg4Ll26RGBgYLbbPE3kzq8QQogceW/LeyQZkjCYDKQaU4lLiWPSzknW\njiWE+I9KlSrx6qutcXZuAozD2TmErl07UaZMmQe2XbZsBSVLlqNBg+6ULFmOsLANBZ63oEmfXyEK\nmJprQs3ZhPUFzw7m2PVjmdZ1qNyBVV1XWSlR/lN7Tag9n7Aes9nMypUrOXHiJMHBVQkNDUWj0WTa\nJjo6mooVq5OSEg4EA3twdm5HTMx53NzcrBE716TPrxBCiDzzauVX0dnpLMs6Ox0dqnSwYiIhRHa0\nWi09evTgk08m0rVr1wcavgBnz57F3r4yGQ1fgHpotd5ERUXl+vwxMTE0bNgaFxdvnnmmFgcOHMj1\nMfNKrhq/qamp1KlThxo1ahAUFMT//d//5VUuIYQQKjM+ZDx9qvVBZ6fD1d6V0Q1H0+vZXtaOJYR4\nQmXLlsVgOAncHQbtMEbjNQICAnJ1XEVReOmlduzZU4vk5GOcPj2Cl156mWvXruU6c17IdbcHvV6P\nTqfDaDTSoEEDpk+fToMGDe6dQD6SESITNdeEmrMJ9fv32r98+MeH3NTfJLRqKMPqDkOrKdwfMKq9\nJtSeT6jfnDn/Y9iwD3FwqITBcJoFC76na9cuuTrmlStXKFs2mNTU60DGHWc3t9YsXfom7dq1y4PU\n2ctJTeR6tAedLuMjMIPBgMlkwsvLK7eHFEIIUchciLtAvfn1SDYko6Bw9NpRbuhv8FnTz6wdTQhB\nRjvt1KlTuLi4UKZMGUs3iMGDB/Hyy62IjIykdOnSJCYmcvHiRQIDA5/4XK6urphMKcBVoASQjtl8\nUTUjSeT6T3Kz2UyNGjXw9fWlSZMmBAUF5UUuIYQQKmNWzPwd8zc7L+4kMS0x02u/HP+FNGMaChl3\nXPTper7d/y1bzm3hxPUT1ogrRJGmKApXrlwhISGBqKgoKlasTv36XahatS49egzEbDZbtg0ICCAw\nMJBGjVpTt24HKld+np49B2Xa5nE4Ozvz0Uf/h07XCI1mNM7OTalbt1KmngHWlOs7v1qtln/++Yf4\n+HhatmxJeHg4ISEhmbYZP3685f9DQkIeeF2Ip1l4eDjh4eHWjpFjUq8iK0azkZeXv8yuqF3YaG1w\nsHFg14BdVCxWEcj4qFFD5gdqEtISCP05FIPJwPC6w/m06afWiP5YClu9gtSseND169dp1qw9p0+f\nxmRKpXjxUly92gOTaQyQzLp1zVi8eDH9+vWz7NOr1xAuX+5m2Wbt2qYsXbqUPn365OicSUlJ2NjY\n4OTkBMDEiR9Tt+5z7N9/gMDAAfTq1QutNu+7QT1JzebpUGeffPIJTk5OjBgx4t4JpD+SEJmouSbU\nnE1Y15wDc3hvy3vo0/UAaDVaXij5AnsG7QHgUvwlgr8LJsmQhFl58G6Rzk7HrgG7qFGiRoHmzi21\n14Ta8wnrePnlUH7/PYD09BlAHFAO+Bsof2eLybi7f83LL7/Mt99Ow8PDg+LFy3DjxtZM2wwbdouZ\nM6c99FwpKSl07NibP/7YACj07TuQuXNn5UtDNyfyfaizGzducPv2bSDj4n///Xdq1qyZm0MKIYRQ\noRM3TlgavpDRBeJs3FnLcin3Uhx47QChQaE0LNUQO61dpv1ttbacjztfYHmFKMr27z9AevobZDxs\n5gUEoNH8cufVVGAt8fGdWblyNyVLBvPOOyMpV648Wu0ayzY63UaCgys/8lwjRnxMeLgGozEOo/EK\nK1Yc4ptvvsuX68oruWr8xsbG8tJLL1GjRg3q1KlDu3btaNq0aV5lE0IIoRK1/GrhbOdsWbbR2FDd\nt3qmbSoWq8iKzivY3m87rg6umV5LN6UT7BOMECL/lS4diEaz7c6SCQeHYri6folOFwwE3Plaj9nc\nmZSUxcydexYnJwd8fGbj5lYLna4STZuWzNQtIjvh4XtITX0HcADc0etfZ9u2Pfl1aXkiV31+n332\nWQ4ePJhXWYQQQqhUr2q92B65neVHl2OrtcXP1Y/FHRZnua2N1oYNPTbQellrDCYDRrORr9t8TaVi\nlQo4tRBF06JF39CwYQtMpl8xm2OpXr0kYWEnWblyJe+99wV6fU8gAZgAQFpafXbt8uLSpXNERkbi\n4uJC1apVs5wY47/KlAng5MldmM0NAQV7+134+3tz8OBBSpcujbe3d75e65OQ6Y2FKGBqrgk1ZxPq\ncCXpCkmGJMp4lMFW+/D7J2nGNKLio/B18cXN4cmnSr2adJXfTv+GBg2vPPMKxZ2LP/GxHpfaa0Lt\n+YT13Lx5k7179+Li4kKDBg2wsbEhPT2dunWbcvSomfR0E7CbjK4Ridja+hAXdx0XF5fHOs+FCxd4\n4YXGpKUFA8k4O0eTkHALW9tSGAxRfPfdLPr1650PV5i1nNSENH6FKGBqrgk1ZxNPh/Wn1zNz70xs\nNbZ8UP8DmpZ7eFe5c7fO8fy850kzpQEZD84dfP0gpdxLFURc1deE2vMJ9dHr9Xz66RS++moeqant\nMJkao9PNo2PHiixZMu+B7RVFITY2lvDwcNav34aXlxsjR76baRzgW7duER4eTnp6Ov37v0FKykag\nDnACJ6eGnD79T65njcspafwKoUJqrgk1ZxOF37pT6+j+a3fLg3NOtk6E9QijSdkm2e7T8ceOrD21\n1jKChI3Ghu7B3VnScUmBZFZ7Tag9n1CvW7duMW7cZ5w/f5mQkDq899472NjYZNomLi6OFi068M8/\nBzEaHYEJaLWRuLkt4d9/9+Pv759p++PHj1O3bgcSE09Z1rm7N2T16ok0aZJ9neelApnhTQghhMiJ\nabumZRoxIsWYwsy9Mx/a+I1JjMk0dJpJMRGTGJOvOYUoCry8vPj66+kP3WbIkPc4cqQyRuMVYCFQ\nF7MZkpISWLhwEaNHj8q0fUBAAEbjdeAg8BxwFoPhBOXLl3/g2NZUuCddF0IIUXhk8ezMox6oaVup\nLTo7nWVZZ6ej7TNt8zqZEIXWvn37mDt3Ln/88UeefwoQEXEQg+F1IA24NzWxyeRBWprhge3d3NxY\nunQ+Ol0L3N1fwMmpDl9++TmlS5fO01y5Jd0ehChgaq4JNWcThd/qE6sJ/SUUo9kIgIONA5t6bSKk\nTEi2+5jMJt7c8CYLDi1Ao9EwpPYQZraciVZTMPdu1F4Tas8n8teMGV8xdux0oAUazS66d2/BvHmz\n8uz4TZu2Jzy8PmZzIrAVmApcRKcbxp49W6lWrVqW+12/fp2zZ88SGBhIyZIl8yxPTkifXyFUSM01\noeZsqpeWBvb2kIOhgYqqtze8zbyD8ywPrzloHdg5YCcv+L/wyH3vfl/mZOilvKT2mlB7PpF/4uPj\n8fEJwGA4BpQGEtHpgti9ez3Vq1d/1O6ZHD9+nC5d+nPu3AnKl6/Czz8vICgoiPPnz1Ov3kukpASS\nmnoerdbIs88G8cUX42nYsGG+XFduSZ9fIYTIb7Gx8PLLcPgwODjA7NmQg4Hh1UJRFEyK6ZHDlj0u\nk9nEj8d+5ELcBZ7ze47WFVuz8thKS8MXwGA2sPrE6hw1fgu60SuEWh07doywsDD0ej22tp4YDHe7\nFLhiZ1eJK1euPFbjV6/X07hxa27eHI2ihHLixE80btyayMjjlCtXjjNnjrBnzx4cHR158cUXsbOz\ne/RBVU4av0IIkRuvvgpHjoDZDCkpMHQoBAdD7drWTvZIiw4vYmjYUFKMKTxf8nnWdV+Hj7NPro+r\nKAodfuwxlnLIAAAgAElEQVTAtgvbSDGm4GTrxLt13s3Udxcypjz+70xwQojsbd++nbZtQzEYemJj\nE0N6ehwwDxgI/IHReIQaNWo81jFPnDiBweCBorwOgKK8jsEwmxMnTlC7dm3c3Nxo2bJlnl+LNckD\nb0II8aQUBf7+G0yme+tMJti923qZcmh/9H7eWP8GyenJmBUzf8f+TaefOuXJsSOiI9h2YZvl2Mnp\nyUzfM50JIRPQ2WY0gG21tng6ejLouUF5ck4hioK33x6FXj8Ho/FL0tJ+QqttgYfHRDQaO7y9B7F+\n/c/4+vo+1jE9PT0xGGKB+DtrEjAYYvDy8srz/Gohd36FEOJh4uNh9WpITYU2beD+p5Y1GvDwgJs3\n762ztYUSJQo+52P6M+pPTMq9RrvRbGTf5X15cuy41LgHulHYamxpWrYpG3puYPXJ1Xg4evBG7Tfy\n5E6zEEVFXFwccG+acKOxDj17BjJz5udP3B2hXLly9OnTnWXLGpKa2gIHh0306NGNcuXK5VFq9ZHG\nrxBCZOfGDahRA27fzujWMHIk7NwJNWve22bhQggNBa02ozFcqxZ07Gi1yDlVwqUEdlo7DKZ7wxV5\nOHrkybFrl8zc5cNGY4Oviy8lXUtSyr0Ujcs0zpPzCFHUtGvXisWLPyIlZR4Qi073Le3azc11P9xZ\ns6ayf39jjh6dg9HoSETEAW7fvo2HR978TFAb6fYghBDZ+fxzuHYNkpMz+vMmJcHbb2fepm1bOHQI\nZs2CZctg69aMu78q1yWoC9V8q+Fi74LOTofOTscP7X/I9XEVReFK0hVmtJhBec/yONk6UcuvFuH9\nwrHR2jz6AEKIbH311ed07uyPs/OzFCv2KjNnjslRf1yz2fzQERCmTZvJyZNeGI03MRi2cPRoIkFB\n9Vm58kcuXLjA559/ztSpU4mKisrLy7EaGepMiAKm5ppQczar6NkTli/PvK5CBThzxjp58pjRbGTt\nybXcTLlJw9INqVK8Sq6OZ1bMdPulG2FnwrDV2mKntWNn/50kpCWw5dwWvJy86FejHy72Lnl0BflP\n7TWh9nzCuhITE+nadQBbtqzD3l7HlCmTeOedoQ9s17ZtN8LC2gDPA42A/wN8cXT8CEVJwmzuCZhx\nclpFRMQOnnnmmYK9kMcg4/wKoUJqrgk1Z7OK5cvhtddAf2dKXicnGDQo4y5vEXQh7gKf/vkpt1Ju\n0S24G6FVQzO9vvzocl7/7XWS05MB0KChlFspruuvk2ZKw8HGAX9Xfw4NOVRoGsBqrwm15xPW1a3b\nANasMZCWNheIRqdrwa+/fkerVq0s25jNZkqXrkx0dABQAfABJt15tSXQAngfAI1mCl26nOLHHxcU\n6HU8jpzUhHR7EEKI7HTvDh99lNHotbOD9u1h2jRrpypwcSlxvLToJcrNKsf8Q/NZfXI1/df25+uI\nrzNtd+rGKUvDF0BB4VLCJVKMKZgVMynGFKITo1lyeElBX4IQRdLWrdtISxsP6ICK6PWD2Lo1PNM2\nO3fuJCbmKnAWWAbc3z1JD1S0LClKeW7cuJ3PqfOfNH6FECI7Gg2MGZPR5zctDVasyJjIoojp8GMH\ndl7cmWmdPl3Ppzs/zbTuWd9ncbZztixrNdoHJqdIN6UTnxaPECL/ubm5A3OA3YAZR8fDlChRPNM2\no0Z9hqIMBS4C24EvgLnAWuzszmNvPxo4BZxAp/uELl3aZNo/KSmJb775hokTP2HfvrwZMSa/qf+p\nDCGEsLYiPLuYyWzir6i/Mg2Ldv9r9+tUpRN/nP+Dhf8sxM7GDg8HD6r6VCU8Mtwys5u9jT3NyzUv\nkOxCFGXr1q0jOvoy8DewEltbe0qVcmPw4MwPtkZGXgS+BDTAC8AAbG1H8fzzz9O//0QuXrzM7NlN\n0Wg0DBs2lMGD743NnZSURM2a9YmOrkBaWkWmTGnPwoVfExrapeAu9AlIn18hCpiaa0LN2YR1KIqC\ny2cu6I36TOvtbewZXnc4U5pNeWCf6IRo4tPiqeBVAYPJQP81/dl8bjOuDq7MbjOb9pXbF1T8XFN7\nTag9n7AORVFwdfUmOXkDUAdIxsGhOqtXf029evWIjY0lMDAQnU7HSy+9ws6dz2MyDSfjLvFXDB7c\nhu+///6R5/nuu+94//0tpKSsvrNmF76+fbhy5Vz+XdwjyANvQqiQmmtCzdmE9Xx/4Hve3/I+Kekp\naDQa7LX2fNTwI8Y0GoNW83T3nlN7Tag9n8hfer2epUuXsnXrNmJiEqhQIZCxYz+gePHieHoWx2jU\nk3FHF1xcetK1qzPLlv2EnZ0PEM9vv/1E2bJlqVevGbGxN1CUekAtnJwWMGvWBAYNGvDQ80+ePJmx\nY29iNE6/s+YqOl0Qyck3H7pffpLGrxAqpOaaUHM2YV3hkeHsiNxBCZcS9K3RF0dbR2tHKhBqrwm1\n5xP5R6/XU6tWI86dSyM9XQ+MQ6M5i7v7//j33/3Ur9+CqKh3UJTBwEkcHRugKJCW9hdQGfgdN7de\nXLsWxZIlS3jnnZ9ISdlMRmP5KG5uzYiPv/rQDH///TeNGrVBr/8RqISDw/u0bWvHL78szu/Lz1ZO\nakL6/AohhHjA7dTbrDu1jnRTOm0qtiGkTAghZUKsHUsIcceSJUuIivIjPf0IEAYEoyig119l+fLl\nbNr0K82bt+f69VGAgSFDXmPBguOkpVW+c4TmGI32xMbGkpKSgtlcnrt3iSEQvT7hkRlq1arF8uVz\nefPNwSQm3qZly1YsWPBtvlxvXpLGrxBCiEyuJl2l5pyaJKQloKDw/pb32TtoL5W9Kz96ZyFEgYiL\ni8NgeAY4CNz7JMZsdsRkMlG5cmWiok5y/fp1PDw8OH/+PHPmhAAxQEngAJCMr68vzZo1w8ZmItAO\nCMbBYTRNm7bLUY727dvTvn3h6ccPuRzq7NKlSzRp0oSqVasSHBzMrCI68LsQQjxNPtn5Cdf110lO\nT0africhLYF3N75b4DkUReFywmUib0fKR/tC/EfTpk2xs1tKxiQUPYHfgdk4OCynU6dOQEYXAB8f\nH+zt7alcuTJjx47E0bE67u6N0OlasXTpDzg5OVGlShXWrl1BuXKj8PCoR7t2WlaunG/Fq8tfuerz\ne+XKFa5cuUKNGjVISkqiVq1arFmzhipV7k2RKf2RhMhMzTWh5myi4LRf2Z51p9ZlWhdcPJijbx4t\nsAwGk4EOP3Zg24VtaNAQ7BPMH33+wM3BrcAygPprQu35RP76+edfePPNEdy+fQOdrjg1a1Zl5sxP\nqFmzZrb7nD9/nosXL1K5cmX8/PwKMG3ByPc+vyVKlKBEiRIAuLi4UKVKFWJiYjI1fkXhpygKa9as\n4fTp01StWpWXX375gYHrhRBPj5crvswf5/9An54xvJmTrROtK7bOt/MZzUYURcHOxs6ybuquqWy/\nsJ1UYyoAR64eYfjm4cx/5em9GyXE4+rSpTNdunR+rH3KlStHuXLlnuh8JpOJffv2ERsby8mTJzGb\nFV55pR1Xr17l228XYWtrw4gRb1CvXr0nOn5BybM+v5GRkRw6dIg6derk1SGFSgwYMJSff95NWloz\nHBw+YODAP/nqq8+tHUsIkU9ee+41zt46y1d7v8KMmQ5VOjDppUl5fh6zYmZo2FD+d+h/KIpCp6BO\nLOmwBHsbe/Ze3kuKMcWybZopjf3R+/M8gxAiZwwGA82atefvv8+SknIdRWmNVhvAp5++hKIoGAxT\ngVQ2bXqFrVt/o27dutaOnK08afwmJSXRuXNnvvrqK1xcXB54ffz48Zb/DwkJISQkJC9OKwrAqVOn\n+PHHNaSknAZcMBpHM3duBUaOfJuAgABrxysUwsPDCQ8Pt3aMHJN6FRqNhqnNp/J5s89RUPJtLN9Z\n+2ax+MhijGYjAL+d+o2Pt33M1OZTqVq8Kn+c/8MyM5yd1q5AHrgrbPUKUrOiYMyZM4cDByAlpTtw\nHfgOsxnS0l4EPgYyZn7T683MmPE9P/9cMI3fJ6nZXI/zm56eTtu2bWndujXDhg178ATSH6lQ2717\nN61bDyMhIcKyztW1Crt2/cSzzz5rxWSFl5prQs3ZxNOnzbI2bDy7MdO6miVqcnDwQZIMSTRc0JCz\nt86i1WjxcvRiz6A9lHApUaAZ1V4Tas8n1Gnv3r306TOUK1cu88IL9VixYh7Fixd/6D5Dhw5n9mx/\nIJqM0SJG3nnlMNABOH9n+Xs6dNjNqlXWGes3JzWRqz/nFUVh4MCBBAUFZdnwFeqkKArffPMdL730\nKl279ufs2bPZbvvss89iaxsDLARuo9HMRqdLo2LFigUVVwjxKIqS8VXIlHYvja323geQWo2WUm6l\nAHCxdyFiUARbem1hfff1HB96vMAbvkI8jaKjo2ne/BXOnPmIxMRD7NxZjlatOj1yvzp1auLs/BPQ\nGJgF7AMisbcfjq1tIrAKWI6T03jeeqtffl5CruXqzu9ff/1Fo0aNqFatmuUBqMmTJ9OqVat7J5C/\nSlXn448nMHPmWvT6UWi1p3B1/YZjxw7g7++f5fZHjhwhNHQAkZGnqFgxmF9+WcgzzzxTwKmfHmqu\nCTVnE1kwGGDQIFi5EmxtYeRIGD8eCskDqdeSr/HcnOeIT4sHBRxsHYh4LYJynk/2ME5+UHtNqD2f\nUJ8ff/yR1177icTEX++sMWNr68rNm7G4uWU/moqiKLz++jssXrwYRdGiKBp0Oke6detCo0Z1+Oab\nxdjZ2TJq1FuZ2oEFTaY3Fllyc/MlMXE3UB4AB4eBfP55Nd59t+DH8SyK1FwTas4msjBiBMyeDSl3\nHgzT6TKW+/a1bq5sKIrC1F1Tmbl3JgoKw+oM483n32TT2U2YFBMtyrfAW+dt7ZiZqL0m1J5PqM/G\njRsJDf2YpKQIwAaIxta2Anp9AvHx8fTrN5SIiP2ULh3IwoVfExwcnGn/GzduoNfrCQgIQKvNn+cB\nckOmNxbZUMjc40UrPzyFKIzCwu41fAH0eli/XrWN33kH5zFx50TLEGqT/pyEp5MnQ2oPsXIyIYqO\n5s2b8+yzX3D4cCv0+nrodMsZNWoCtra2NG/+KseO1SI9fTM3bmynYcMWnDlzBG/ve3+U3v//hZX6\nmuwi3w0Z8jo6XVdgHRrNNBwc1llmgxFCFCK+vpm7ONjaQsmS1svzCCuOrrA0fAH06Xre2vAW3+7/\n1oqphChabG1tCQ8P48svQxkzRsOvv37D6NEfcOPGDU6cOEZ6+kygIoryOmZzNfbs2ZNp/yNHjtC6\ndRdeeKE506bNxGw2W+dCckHu/BZBU6ZMpEQJH3799Xt8fLyYMmUnpUqVsnYsIcTjmjULGjQAozGj\nEezmBqNGWTvVAy4nXCbZkIy7ozsaNCjc+6TJpJj44PcPCHQPpG2ltlZMKUTRYW9vz2uvvZZpnZOT\nE2azAYgDigEmzObYTEPYnj9/nvr1m5GcPAZFqcCxY+O4eTOOKVMmFmj+3JI+v0IUMDXXhJqziWxc\nupTR1cHeHjp2BE9PayeyMCtm+qzuw68nfsVGY4OXkxc3U25muvt712vPvcbcdnOtkPLh1F4Tas8n\nCpfhwz9i3ryNJCf3wMlpJ+XK3cDV1YmjR8/g5ORAjRqV2LGjPOnpdz+tOYu7eyNu346xau77SZ9f\nIYTIjb174YcfwM4Ohg6FoCBrJ3pQqVLwxhvWTpGlRf8sYvXJ1ZYpitOMadQNqMuZW2e4mnzVsp2d\n1g4fZx9rxRRC3PHFF5OpW7cmu3fvx2wux/ff/4nR6A78RHKyC+HhbVCU+ye4MmBjY2OtuE9MGr9C\nCJGVbdugbduMB8o0Gli0CPbsAZncJccOxBzIdJfXqBg5c+sMv4T+QsulLUk3pWOrtcXD0YNhdWWs\neCGsTaPR0LVrV7p27UqjRm0xGp8DegL1ATAav0CrHYRWWxqzuQI63We8//7bVs38JKTxK4QQWRkz\n5t5ICoqSMZLC1KmwZIl1cxUiVX2qorPTWRrAWo2Wil4VaVC6AX+//jdhp8NwsnOie3B3PJ3U011D\nCAFpaQYgFrh031oDwcHPUq3aZW7c+JcuXUbQv786R5d5GGn8qpxer+fddz9i27a/KFmyBN99N/WB\nMfeEEPng/iHEIKMBnJxsnSyF1Ou1XmftqbXsitqFjdYGnZ2ORR0WAVDZuzKVvStbOaEQIjtDhvQk\nImIQMAdIBFyAr2jTZiiTJ0+2brhckgfeVK5t21C2boXU1JFoNAdwdZ3AyZOH8PPzy9H+JpOJzz6b\nxtq1v+Pr68306eOpUqVKPqcWD6PmmlBztkdKTIQ1ayA1FVq2hNKlc3e8777LmERCf+dje50Ofv4Z\n2rTJfdYixKyYOXL1CMmGZGqUqIGzvbO1Iz0WtdeE2vOJwktRFOzsnDCZdgPrgHQcHQ/wzTehDBw4\n0NrxsiUzvBVyBoMBnc4VkykecATAxSWU7757hUaNGrFu3TpsbW3p3LlztoNOv/XW+yxYEIFe/zEa\nzTFcXady7NgBAgICstxe5D8114Sasz3UrVtQsybcvJlxh9bGBnbuhBo1nvyYigLffJPxZWsLY8dC\n1655l7kI+yvqL345/gtuDm4MqT2Ekq7qHZtY7TWh9nyicBsyZBhLlvyNXj8SG5vDuLvP4eTJQxQv\nXtza0bIljd9CzmQy4ejojNF4EfAFFFxcmjN+fGsmTJhKenpbNJoUXFx2cejQbvz9/R84hpOTO6mp\nJ4GMO8WOjv2ZPr02Q4cOLdBrEfeouSbUnO2h/u//YMYMSE+/t65u3YwH1ISqrDm5hp6/9kRv1GOr\ntcXdwZ0jbxxRbQNY7TWh9nyicDOZTHz++ReEhW3H39+HqVPHU6ZMGWvHeqic1ITM8KZiNjY2vPfe\nSJydWwLfYW8/EF/fa4SF7SApaQypqfNJSVlOXFw3Jk78PMtjaLU2gMGyrNEUzmFJhHio6OjMDV+A\nq1ez3lZY1YgtI9AbM7qSGM1G4tPimXNgjpVTCSGyYmNjw6hRI9m1awM//bRQ9Q3fnJLGr8pNmTKR\nb799j169DjFiRCkOHNjJzZu3UZR7440ajVW5cuVmlvu/887b6HQdgBXY2IxBp9tJx44dCyi9EAWk\ndWtwvq8vqZNTRr9foTopxswPEhrNRhINiVZKI4QoiqTxq3IajYa+ffuwZMlcPv10Ah4eHrzySnN0\nuk+Aq8AFdLrptG/fPMv9P/tsPNOnD6FFi1X06nWDgwd3kZaWRr16zXF2LkblyrU5ePBggV6TEHmu\nWzcYORIcHDL657ZuDV98Ye1UIgvdg7ujs9NZlp1snegc1NmKiYQQRY30+S2EjEYjQ4e+z+LFC7Gx\nseX994czfvxoNBrNI/c1mUw888xzREZ2xmQaDPyOu/sIzp49mu1DcyJvqbkm1JwtRxQl40srf9er\nldFs5KM/PmLlvyvR2emY1nwa7Su3t3asbKm9JtSeT4iCJg+8iQdERUVRpUpd9PpoIKOx7O7ejJUr\nR9CqVSvrhisi1FwTas6WL65fh99/z5i+uHVrcHGxdiKhMmqvCbXnE6Kg5aQmZJKLIsbNzQ2jMRG4\nDvgAaRiNUXh4eDywbUJCAl9//Q3R0ddo2bIJ7dur9+6MEI/tzBmoU+feg3JeXnDwIBQrZt1ceWHl\nSliwANzcMmaqq1bN2omEEEI15M5vETRq1HhmzVpJSkpHnJx2EBLiz2+//Zip20RycjLVq9fj0qVq\nGAw10OnmMG7cG3zwwXtWTP50UHNNqDlbnmvVKuOur9mcsWxnB2+8AV99Zd1cuTVvHgwbdm9yDmdn\n2L8fcju5zfnzcPlyxnFUPMZnXlN7Tag9nxAFTbo9iGyFhYVx8OBBypYtS/fu3R8Y/mzp0qUMGbKM\n5OQNZHSPOI+jY030+ts56lsssqfmmlBztjwXHAzHjmVe1759xixxhVn58hkN1bs0Ghg+PGMc5Cf1\nySfw2WcZDxQajfDrr0VmNA2114Ta8wlR0GScX5Gtl19+mTFjxhAYGMgPP/zArl27Mr2u1+tRFF/u\n9gsGX9LTU+WHrHh6NG2aMSTaXTodNM961JRC5e6d7LvuPgT4pA4fhilTMqaNjo+H5GTo3DmjESyE\nEIWQNH6LsA8+GEOrVn0ZNmwPLVr0ZOzYSZbXmjdvjkYTBiwDjuHoOIA2bTqglafoxdPi888z7l7a\n2GQMj9a/f0a3h8Ju2LCMhvxdOh306/fkxzt7NuP9uZ/RmDGVtBDCahRFkRtST0i6PRRRFy5cICjo\nhTtTHxcDruLgUIVz545apknet28fgweP4Nq1azRv3oTZs2fgfP9EAuKJqLkm1Jwt3xgMGUOj/beB\nV1gpCsyfn/HAm4sLTJyY8WDfkzp+HGrXhpT7Jqfw8IAbNzL+cHjKqb0m1J5P5I8vv/ya0aPHYTDo\nadOmA8uX/y9Hv58VRWHZsmXs3LmPihUDeeutoTjd9wmY0Wjkt99+49atWzRs2JBKlSrl52XkC+nz\nK7K1d+9eWrV6m/j4/ZZ1bm7PEh6+mJo1a1ox2dNPzTWh5mzCir79Ft5/H+ztM/oQh4VBgwbWTlUg\n1F4Tas8n8l5YWBihoe+g128E/HB0HERoqCeLFn2f5fZpaWnY29uj0Wh4552RzJ+/Fb2+D46OOwgK\nusnevVuxs7MjPT2dxo3bcPRoEopSCUXZwKpVS2lZyPr3F0if3wEDBuDr68uzzz6b20OJAlSlShUU\n5RKwBjADn2IyRePp6WnlZEII1Rk6NGOkh127ICYmTxu+cSlx/H7ud/Zd3odZMT96ByGKuM2bt6HX\nvw5UAlxJTR3Pli3bHtjuypUrPPdcI3Q6V5yc3Pi//xvF7NnfoNf/AQwjNfVXTp9OZseOHQCsXLmS\nI0eMJCXtIjl5EXr9Svr3f6tAr62g5Lrx279/fzZt2pQXWUQBcnd3Z/PmNfj4DAdc0Wi+B2pSvXpd\nDhw4YO14Qgi18faGZ5/NGDotj/x77V/KzypP558703RxU1ovbY3RLA/SCfEwfn7FcXA4et+aoxTP\nYvjBjh37cPToi5jNcaSl1WDKlB8wmbSA+50ttGi1PujvDIt45coV0tJqcK9pWItbt67k45VYT64b\nvw0bNpS7hflo69atTJkyheXLl2MymfL02HXr1uWLLybh7FwdRTlLcvJWEhK+pEeP1/P0PEIIkZXe\nq3tzO/U2CWkJJKcn89elv1j4z0JrxxJCdc6fP0/Tpu0pW7Y6+/f/i5/fIZydX8bRcRDOzkP57rup\nD+yzf/+fGI0fA18CxYHLQB1gMHAMjeZbbGyOUL9+fQAuXozGaFwEnACM2NiMp06dRgV1iQXqKXnC\n4+k0efJ0Jk2ajcHQCQeHdSxZ8ithYT/n6YgLkZGRpKQ0AhzurGnG5cv3nnjfv38/x48fp1KlStSr\nVy/PziuEEBdvX0ThXt88fbqes7fOWjGREOqTkJBAvXpNuXHjDczmccTEzKVqVS/GjeuCXq+nZcv/\no3z58g/s5+VVgmvXDgCngbZkNPl+BXphY9OI55+vxQ8//EGxYsW4cOECP/ywDBgF1AX0mM06/ve/\ne58EX7x4kU8+mcb163F07tyGXr16FNpx/wuk8Tt+/HjL/4eEhBASElIQpy3UUlJSGDduHOnpp4AA\njMbdbN8eytKlS+nTp0+enadWrVo4Ob1DcvJwwAcbm++pVq0WAJMmTWXy5K/RaBqjKON4993+fPbZ\nuDw7d1ERHh5OeHi4tWPkmNTrfdLTITIyY+rjp2Ha4wKiKArTdk/jm4hvsNHaMLrhaAY9N+iB7aqX\nqM6fF//EpGR8quVs50wtv1oFHTeTwlavIDX7tNuzZw+pqYGYzR8AYDDM5sQJX1q0aEHJkiWz3W/h\nwtl07tyF9PRA0tOjgJ6AO/b2ZenYsQQrVsy3bHvx4kXM5jLACOB9IB1HxyDS70z/HhsbS82aLxIf\n3x+z+Xn++GMSMTFX+PDD9/PrsnPsiWpWyQMXLlxQgoODs3wtj05R5Fy9elVxcPBUwKzAhwqUVqCV\nYm/vpSxcuDhPzzVmzCeKnZ2z4uTkq5QvX025ePGiEhsbqzg4eCgQc2eE/OuKo6O3cv78+Tw9d1Gk\n5ppQc7YCd+qUovj5KYqzs6LY2yvK2LHWTlRofPj7h4p2glZhPArjUZwmOSk/H/v5ge2iE6KVSl9X\nUnSTdIr9J/bK2xveVsxmsxUSZ0/tNaH2fCL3wsPDFReXagqY7vw+TlDs7V2VmzdvPnLf06dPK3Pn\nzlWee66B4uRUQnF2LqMEB9dRbty4kWm7jRs3KqBT4MCdc/yhgJOSkJCgKIqifPHFF4q9/QDl3qw5\nxxRPT/98ud7cyklNSLcHlSpevDiBgWU4e/ZNzObfgKOAJwbDCQYPrkOXLp3Q3T+QfS5MnPgxI0a8\nQ2JiIn5+fmi1Wg4fPoyDgz9paX53tvLGwaEcsbGxlC1bNk/OK4SqtW8PV67cmx1txgwICYEmTawa\nS+2i4qOYtntappEbUowpLD68mM5BnTNtW9K1JMffPM7lhMu4Orji5eRV0HGFUL369evzzDMeHDvW\nhdTUpuh0y+jcuQdeXo+ul4oVK1KxYkUGDRpEZGQkBoOBChUqYPOfMbrj4uJwdKxBampzwBVIwdY2\nY9xfyPivotw3IyZOmEwZr5lMJubN+x8HDhwlOLgiQ4e+iZ2dXR5dff7IdefR7t278+KLL3L69GlK\nlSrFggUL8iJXkafRaNi6dR3ly+8GSgN3Hyqsgo2NCzfzeHYlNzc3/P39Lf2JK1SogFZ7C1gFKMBG\nzOZIgoKC8vS8QqjWmTOZpwU2GjOm+hUPteHMhizXu9q7ZrneRmtDoEegNHyFyIatrS1//rmJMWOe\np3fvw8yY0ZcFC2Y/1jE0Gg1ly5blmWeeeaDhC1C1alW02gvAHmA78C3u7sXw8PAAoGPHjjg4/AR8\nD2xFp+vFgAH9AOjRYyDvv7+U+fMrMmrUBlq37oT5v9Osq4xMcqFyUVFRVKlSC71+A/A8sBJv7w+I\njRbHMwYAACAASURBVD2P7X9mpDKZTHz66VTWr9+Gn19xpk+fQMWKFZ/43BEREbRrF8rNm1dwd/dm\n7dqVNCgiA9vnJzXXhJqzFbjSpeHSpXvLzs6wYgW0a2e9TIXAgkMLeHPDm6QaUzOt//eNf6nqU9VK\nqZ6c2mtC7flE4WAymfjww9F8/fV32Nv7YWubwKZNq6lz3+yQhw4dYsSICdy8eZtOnVoxatRIoqOj\neeaZWqSmXgR0QDrOzpX5889frDZhlszw9pRYs2YtPXv2w2TKGJ9306bVWX5Tvf76OyxbdgS9/kO0\n2iO4uX3FiRMHKVGixBOfW1EU9Ho9Op2u0D7VqTZqrgk1Zytw+/ZBixYZM5qlp/9/e/cdHkW1/3H8\nvckmJJvQQ+DSjDSJKKGjWOhdVJCLFFEUAQvwowiIhXalXhAUFBRBBIEoNq4IiApBvJeiIqj0KkWC\ndEI2bTfz+2MhEEMJyWZnknxez5PHzGRmzmeJ3+zJ5Mw58Mgj8MEHnm25pnOJ57j97duJjYsllVRs\n2OhUrRPRHaPNjpYlVq8Jq+cT64uPj6dx4wfYvv0IAOHhgXz//deULVv2hufu3LmTOnVaEx+/H/D8\nbCxUqA4rVrxJgwYNcjL2Nanzm4e4XC7OnDlD8eLFM0x1duTIEZYuXUr//gNJTf0VqApAcHA33nij\nEb169TIhsVyLlWvCytlMceoUbNkCJUp4FnhQxzdTenzeg4W/LcRleMYEOuwOfn32VyoWyzgdk9VZ\nvSasnk+sb+jQV5g+fR+JiR8CfgQG9qVzZzczZvybkJCQ606v6nK5qFatHvv3N8Pl6o6//38oWXIe\ne/Zs9dpzSTfLJ8sbi2/Y7XZKlCiR4X/CnTt3Uq1aHV544UdSU9sAzQDPn2ptthSvzgksku8ULw5N\nm0L16ur43oQlO5akdXwBXKkulu5aamIiEbmWrVt3kZj4EOAP2EhOrkd09CcUK1aS0NBiREd/fM1z\n7XY733+/ghYt/qBMmU40bLiJ//3vW9M6vpml2R5yucGDRxIXNxTDGHRxzzCgF/7+tXE4NvDQQzc3\nKF5EJLvsfunfWvz8/Aj0DzQpjYhcT5061Vi3bgkJCR3x3BMdTnLyS8BAXK5feeqp5kRF3UlkZORV\nzy9ZsiRfffWRLyNnm24L5nLHj5/EMO64Yk8NihffR/fup9m8+b+EhYWZlk1E8o6ziWeZvnE649eN\nZ0vsluse+9K9L+EI8Nz5sdvshAaG8mi1R30RU0Ru0iuvDKNOHScOx604HLcCJ4GBeMbwRuHv34yf\nfvrp+hfJZTTmN5cbPXo8kyatwOn8GEjG4XiYyZOf4dlne5sdTa7ByjVh5WxintMJp4maFcVJ50lS\n3CkUsBfgs06f0bJSy2ues+i3RXy+43NKhpZk+L3DKVOojA8Te4/Va8Lq+cQ8sbGxzJ8/n4SERDp0\naM+dd955zWMNw2D37t2kpKRQv/59OJ2rgZpAAiEhtVm6dDpNmzb1Wfbs0ANv+YDb7aZfvxeYO3cO\nfn7+/N//9WfcuFE3NTPDsmXLGDRoFPHxF+jU6WEmTfqX5Seozs2sXBNWzibmmfDDBEbGjCTZnZy2\nr3Kxyuzut9vEVL5h9Zqwej4xx9GjR4mKuovz51vhdhcjKGguK1Z8yv3333/Dc5cs+YQePZ7Dz68p\nsJU2beoTHT0318z4pM6v3NDGjRtp3PhBEhLmAeVwOPrTs2cd3nxzktnR8iwr14SVs4l5hqwawuT1\nk9PtCw8J5/gLx01K5DtWrwmr5xNzDBw4lOnT3bjdUy7uWUzt2rP56afVmTp/586d/PTTT5QuXZrG\njRvnmo4vaLYHyYTPPltKQsKzQGvgDpzOt/noo0/NjiUiFtK2Stu0MbwAwfZgHqjygImJROR6Tp8+\nj9sdccWeWzl79lymz69atSqPPfYYTZo0yVUd38xS5zefK1QolICAY1fs+ROHI8S0PCJiPY0iGjGz\n7UzCQ8IJDQyl4+0deavNW2bHEpFr6NSpHQ7HFGATsBeHYxidOml1yks07CGfO378OHfeWY8zZ9rg\ncpXD4ZjBvHlv8s9/djQ7Wp5l5ZqwcjYRM1i9JqyeT8wze/YcRoyYQHJyEo891pkpU8Zht+f9GW41\n5jef2rRpEz169CM29ih3392A+fNnUrx48WseHxsby8yZ73Du3AU6dGiXYUB8amoqo0ePY+7cxQQF\nBfHaa8N49NFOOf0y8iwr14SVs4mYweo1YfV8Ir6mzm8+dPToUSIjaxEXNw1oQEDAZGrV2s6GDd9l\nONYwDObNm8+77y7G4Qhi9OgXuPfeezMcN2bMeCZO/AKncyZwhuDgx/nPfz6gWbNmOf+C8iAr14SV\ns4mYweo1YfV8Ir6mzm8+FB0dTe/eS4iLu/TQmhu7vSBnzvxFaGhoumNnznyXF16YgtM5CTiNwzGU\ntWtXUKdOnXTHVa5ch7173wQaAOeBJ4iKOsLMmW9y9913++BV5S1WrgkrZxMxg9Vrwur5RHxNsz3k\nQ4UKFcIwDgOpF/fEAgYFChRIO8YwDJYvX87Aga/idM4CHgKexOkcxHvvLchwzZAQB3AcT8e3AeDH\n1q3NadasAx99dO01v0VERESsRp3fPKZ58+ZERoYQHPwAMAaHozEjRoxMW7TiwoUL1KnTkLZtu5GU\nlASMAJIAsNlc+Ptn/F9i4sSXcTieAboCFYFPgXE4nZ8wYMBLvnlhIiIiIl6gYQ95UFJSEnPmzOHw\n4aPcc8/dPPDA5fk4n3tuELNn/4nLtRBwA+2AUOBeQkLGsX796qsugbhx40YGDRrC+vV1MYxLk2b/\nSWhoFHFxJ3zwqvIOK9eElbOJmMHqNWH1fCK+pjG/kkGdOs34+echQMuLe5YAg4Dz+PklUKpUBB9/\n/D733HNPhnN//vln7ruvDQkJi4FKBAUN5qGHChIdPdd3LyAPsHJNWDmbiBmsXhNWzyfiaxrzKxnc\nfntFAgK+AoyLH18A57DZRpCaeo4//5xM69YdOH4847KltWvXJjp6NuXK9aNw4bto374wc+fOuGo7\nmzZton79ZlSoUJP+/YdeHGIhIiIiYi7d+c1nTp48yd13NyM21kZSkpOgoPMkJblJTv4r7ZjChVuy\naNH/0aZNmyy1sXfvXmrUuJv4+ClAJMHBI+nYsRzz57/jpVeRu1m5JqycTcQMVq8Jq+cT8TXd+ZUM\nwsLC+P33jaxcOYO1a+dx4MDvGIYTOHrxiARcrr2UKFEiy2189dVXuFyPAI8DdUlImMeSJdHpjjEM\ng7fffod69ZrTtOnDrF+/PsvtiYiIiGRW3l/nTtIxDINx4/7N5Mmv43an0L17D1599RUmTLgHt7sd\ndvsPtG17X4a5fm9GUFAQfn5nrthzmoCAoHTHTJ48jdGj5xIfPwE4RrNmD7J+/XdUr149y+2KiIiI\n3Iju/OYz77//AZMnf4zTuZmkpH0sWrSdxMQkli//gPHjK/DhhyOJjn4fm82W5TY6depEkSI/ExDw\nPPAWDseDjBgxPN0x06fPIT5+DtAWeBqn81nmz1+UrdcmudjKlVC2LISEQKtWcObMjc8RERHJgmx3\nfleuXEnVqlWpXLkyEydO9EYmyUH/+c+3OJ2DgAggHKfzVb788lsaNmzIwIEDefjhhzl37hwtWrQn\nICCYokVLs2jR4ptqo2jRomzdup6BA4vw+OO/MX/+eF54YUC6Y/z9/YHktG2bLfmqcwxLPrBjBzzy\nCBw9Ck4nrFnj2RYREckB2Rr24Ha76du3L99++y1lypShbt26PPjgg0RGRnorn9yE1NRUpkx5g+jo\nLylSpBCTJr1K7dq10x3zj3+EYbdvw+XybNts2wgPL57umM6de7J2bRgu11+cPbuLXr3aUalSRerV\nq5fpLCVKlGDixLHX/PpLL/VnwIAncDpHYrMdIyRkLj17/pD5Fyt5x5o1cOXDCcnJ8P33kJoKfvqF\nSEREvCtbnd9NmzZRqVIlIiIiAOjcuTNLly5V59ckI0e+xuuvf4nT+RrwBw0btuKnn9ZRtWrVtGNe\nfXUon33WgLi4wxiGg4CAr5g69dt011m79luSk/cDBYE6JCd3Zc2aNTfV+b2RXr16UrhwIT744FMK\nFQrhlVfWUKVKFa9dX3KRokXB3z/9vuBgdXxFRCRHZOvd5ejRo5QrVy5tu2zZshw9evQ6Z+RdBw4c\noE+f/nTq9CRLly41JcOsWXNxOj/As4BFbxITn+Djj5ekO6Z06dJs3/4zU6c2ZvLk2vz++48ZVnQr\nXDgM2HZxyyAwcDthYWFez9up0z/56qtoFi+eQ7Vq1bx+fcklOnSAihXB4fB0gh0OeOMNs1OJiEge\nla07v5l9KGrUqFFpnzdq1IhGjRplp1nLOXToEDVrNiAuriepqdX46qsBTJt2kl69erJx40befXc+\ndrs/zz3Xk6ioqBzL4e9vBy4vJmGzJeLvXyzDccWLF6dPnz7XvM67706lc+eOpKZ2wm7fRYUKF+jW\nrVuWc6WmpnL8+HGKFClCcHBwlq+TW8XExBATE2N2jEzzeb0WKAAbNsCCBXDiBDRsCFdZYTBfWLYM\n3nkHgoLgxRfhb8OWJOfltnqFvP8eK3I9WanZbC1ysWHDBkaNGsXKlSsBGD9+PH5+fgwbNuxyA/lg\nAu4xY/7FmDF/4XZPv7hnPWXL9mThwlm0avUICQlDgRRCQqby/fdfU6tWrRzJMW3adF5+eQZO58v4\n+f1BaOhb/PbbJsqXL3/T1/r1119ZvXo1xYoVo1OnTgQFBd34pKvYs2cPTZu248SJU6SmJvD665N5\n/vlnsnStvMLKNWHlbHnekiXQo4fnoT/w3AH/4QeoWdPUWPmd1WvC6vlEfC0zNZGtzq/L5eK2227j\nu+++o3Tp0tSrV4/FixenG/ObHwrzlVdGMG6cC8MYd3HP75Qs+TCRkdWIiWkP9Li4fxodO25hyZJ5\nOZblww8XER39JcWKFWLEiCFUqlQpx9rKjNtuq82ePU9gGP2B/Tgc9xET8wV169Y1NZeZrFwTVs6W\n59WoAVu3pt/35JMwd645eQSwfk1YPZ+Ir2WmJrI17MFutzNjxgxatmyJ2+2mZ8+e+fJht06dOjJ1\nalOczkigPA7HUJ5+ujvffbceKHrFkUVISEi6xlW847HHuvLYY11ztI3Mcrlc7NmzBcPYdHFPBaA1\nmzdvztedX5Grcrszt09ERLIl2yu8tW7dmtatW3sjS65VvXp1vv76c4YM+RdxcRfo0uURhg9/gcqV\nF/Lrry/gdIYCyTgcI+jde4bZcX3GbrdTrFhpTp36HmgMOPHz20j58prDVSSDgQOhX7/Lwx6Cg+GZ\n/D1ESEQkJ2Rr2EOmGsjnf5KZPXsOU6e+h7+/Py+91JcuXTqbHcmnvvnmG9q374q/f33c7p08+GBD\nFi58L1sryOV2Vq4JK2fLFz74AGbN8jwEOGIENGlidqJ8z+o1YfV8Ir6W42N+vRVC8rbDhw+zefNm\nwsPDueuuu/J1xxesXRNWzuZVsbGwfTuULw8mj4sXa7N6TVg9n4ivqfMrYkFWrgkrZ/OaZcvg0Uch\nIMCzmtyLL3rusopchdVrwur5RHxNnV8RC7JyTVg5m1ckJ3tWlLs0rhY8U4pt2AB/W+xFBKxfE1bP\nJ+JrmakJrR+aj+gHpOR7J0/C3+vAboe9e83JIyIiPqfObz6wdetWKlSojt0eSIUK1dn697lERfKL\n8HDPw2RXcrlAy2uLiOQb6vzmcRcuXKBJk7YcODCE1NQLHDgwhCZN2nLhwgWzo4n4nt3uGfNbuDCE\nhnqWEZ42DapUMTuZiIj4iDq/edzOnTtxucKA7kABoDsuVxi7du3K1nVPnTrFnj17SElJ8UZMEd+5\n5x44dgx+/NHz3169zE4kIiI+pM5vHhcWFkZy8lHgzMU9Z0hOPkpYWFiWrzly5FhKl76VmjVbUL58\n1Wx3pEV8LjgYqlaFIkXMTiIiIj6mzm8eFxERQZ8+PQkJuZsCBfoSEnIXffr05JZbbsnS9WJiYpg8\neQ7JyXuIjz/A8eODeeihbl5OLSIiIpIzNNVZPrFq1Sp27NhBZGQkLVq0uOZxmzZtYtasedhsNp57\n7ilq166d7uvTpk3jxRf3kZQ0/eKeRPz8CuJyJef7xSsyy8o1YeVskgP27oUxYzyzYHTqBE88Aarj\ndKxeE1bPJ+JrmakJu4+ySCYkJycTEBCQI53IFi1aXLfTC/DDDz/QokV7EhKGAW6io1uxevUy6tev\nn3ZMpUqVsNvnkpQUD4QAyylduqI6viK5zZEjUKcOxMVBaiqsXQt//QVDh5qdTEQkR2nYgwWcOnWK\ne+5pQXBwKA5HYd5++x1TcowZM5WEhAnAC8AwnM5RjB37Rrpj2rZtyyOPNMDhiKRw4YYUKvQcn3zy\ngSl5RW6aYcCJE57FLvK7RYs8i32kpnq2nU6YNClr10pI8CwXfeKE9/KJiOQQdX4toEuXp/nxx9tI\nTXWSmPgzQ4aMZe3atT7PkZCQBBS+uPU9MIOVK1fRvXtv4uPjAc+fE+bNm8n69cv47LOR7N+/Ld2d\nYRHL2r0bIiKgXDkoVAhmz86ZdpYuhXr1PHdVP/wwZ9rwBrc744IflzrCN+OXX6BsWbjrLs+/7YQJ\n3sknIpJDNObXAkJCiuN07gDCAfDze5HRo0N55ZVXvN7WkSNH2LdvHxUrVqRs2bLpvrZ4cTRPP/0S\nTuerwGBgNnAHQUEjaNXKzuefL/R6nvzIyjVh5WzZVrky7Nt3ucPncMB//ws1anivjRUroGPHy8sn\nOxzw3nvQpYv32vCWvXuhZk24NOe3wwH9+8P48Td3nTJl4M8/L287HLBmjecXgDzA6jVh9Xwivqbl\njXOJsLBSwM8Xt1IJCtrMP/7xD6+38/7786lSpQYPPfQyVarUYM6ceem+3qVLZ956axT/+Mdr2GwP\nAo8At5GYOJtlyz7TD1jJvZKSYP/+9Hc6bTb4+edrn5MVb711ueMLns9nzPBuG95SqRKsWwctWnju\nUo8aBWPH3tw1EhMhNjb9PpsNfv/dazFFRLxND7xZwPvvT6ddu0ex2VqTmvorwcFnWb16Pbfddhv3\n3nuvV9o4ceIEzz8/gISEtSQkLAKO06vXQAIC/Hj88cfTjuvR43EMw03fvp9e8R5+lKCgUD3UJrlX\nYKBnRbfz59PvL1/eu+38felkgIAA77bhTTVqwNdfZ/38AgWgaFE4dSr9/sqVs5dLRCQH6c6vBTRp\n0oStW9czePCtpKYe4tSpPixadActWz7CggULmDdvHl988UW2VlM7dOgQAQHlgQ+BjcBnGMYievce\nzMaNG9Md+89//pNSpQ5RoEA3YDwORxvGjx+TnZcoYi6bDRYv9vxJvlAhCAmBhx6CZs28287QoZ42\nLgkOhhwYvmQZNht8/jkULOhZMjo4GHr3hvvuMzuZiMg1acyvhXTp0pPo6GrAoIt7XsbP702Cgx/C\nZtvH7bcHs27dSgIDA2/62qdPn6ZMmQokJtqA/wGRANhsoxgyJJmJE8elO/78+fPMnDmL48dP0qJF\nE1q1asWePXvYvn07FSpU4M4778zWa83PrFwTVs7mFQcPeoY6lCoFDRrkzJy2mzbBtGngckHfvnD/\n/d5vw2pOnoRt2zz/rrfdZnYar7J6TVg9n4ivaZ7fXCYxMRkodMWeaFJTo4mPbwuk8vvvrViwYAE9\ne/a86WsXLFiQoKBQEhMDgGNc6vz6+x+jYMGMf/otVKgQw4Zdnu9zzpx59Os3lICAerhcmxk6tB8j\nRw6/6RwipoqI8HzkpHr1PNOI5SdhYdCwodkpREQyRcMeLKR3764EB48EvgS+wdNJrXPxq34kJtbm\n2LFjWbr24cOHSUnxByYD3YBxQG8cjv/Qq9fTbN++nc6dn6JNm0dZtCg63bnnzp2jb98BJCT8wPnz\ny3A6f2HixGns3r07i69URERExBy682shrVu35sMPpzNmzGTcbjcuVzX27RtLSsrrwB8EBS3i3nuz\ntqBEsWLFSEk5A9wFfAR8it3+MZ9++glxcXHUr9+I+PjBGEYZ1q59lTNnzvL8888AEBsbi90eBlS5\neLWSBAZGcvjwYapUqXL1BkVEREQsSHd+LaZDhw5s2bKW3377gXXrVlC79g78/BwEBkYxceKLNGrU\nKEvXLVKkCKNHj8ThaEBw8EJCQr6le/duNGvWjHnz5uN0PoFhDAMew+mcz8SJ09POLV++PP7+8cBX\nF/dsIiVlG5GRkdl9uSIiIiI+leU7v0uWLGHUqFHs3LmTH3/8kVq1ankzlwBhYWGsX/8NSUlJBAQE\n4OeXvd9VXnxxMPfffzdbtmyhYsUOtGjRAgCXy41hXDkdUyButzttKzg4mOXLP6Vt244kJaVisyWz\nePEHlC5dOlt5RERERHwty7M97Ny5Ez8/P/r06cOUKVOu2fnVk6jWt23btovDHl4DyuBwvMSrrz7O\niy++kO44l8vF8ePHKVGiRJZmnBAPK9eElbOJmMHqNWH1fCK+lqMrvFWtWlXjPfOIatWqsXr1VzRr\n9jX16r3JpEnPMWzY4AzH2e12ypQpo46vyN8lJMDhw57pzXLKqlVQty7cfjtMmZJ+tToREck0PfAm\nANSrV49vvvnM7Bgiuc+CBZ6FHfz8PAtcfP01eHsY2Pr10L795aWTR470dLSHDfNuOyIi+cB1O7/N\nmzcn9u/rtgPjxo2jXbt2mW5k1KhRaZ83atQoyw9t5XYrVqzg6acHcObMX9x3X2MWL36PYsWKmR1L\nclhMTAwxMTFmx8g01etN2LMH+vSBxETPttMJrVpBbKynM+wtCxZwxXrjEB8P77yjzm8OyG31CqpZ\nyd+yUrPZXuGtcePGGvObCTt27KBOnftxOqOBKAIDX6VBg8OsWbPM7GjiY1auCStns6RPP4WnnoLz\n5y/vK1AADh2C8HDvtTNwILzxRvqhDlWrwo4d3mtDrsrqNWH1fCK+lqNjfq+kwruxNWvWYBjtgaZA\nGMnJU1m3bhWpqanZuu7mzZt58MGuNGnyMAsWLPRKVhHJpIiIjON8/fygaFHvtvP88xAaenk55uBg\nGD3au22IiOQTWR7z+/nnn9O/f39OnjxJ27ZtqVmzJitWrPBmtjylaNGi+PntBgzABuwhOLhQtqYv\n27ZtG/ff35L4+BFASTZufIW4uAs891wfL6UWkeuqXRueeQZmzQK73dMRnj8fAgJufO7NqFQJfvwR\nXn8dLlyAJ56Ai1MViojIzcn2sIcbNqA/yQCQlJTEXXc1ZffuUJKSqlOgwEJmzBjHk08+kaXrxcbG\n0rXrE6xZUxd47eLeH7j11r7s37/Fa7nF+6xcE1bOZmlbt3qGOlSvDrfcYnYa8SKr14TV84n4WmZq\nQrM9+EiBAgVYv/5bFi5cyF9//UXDhkto0KBBlq518OBBatW6h/Pni+O5i3yJfgiKmCIqyvMhIiKW\np86vDwUFBdGzZ89sX+fll1/j7NmnMIxHgcZAKaAkDsfLDB48INvXFxEREcmrNOwhl4mLi6NUqUic\nzinAo8BPwFAKFtzP9On/4oknupucUG7EyjVh5WwiZrB6TVg9n4iv+Wy2B/GdyZOnkpRUFpgIHARK\n4ed3juHD+6rjK2Jlv/4Ky5bBH3+YnUREJF/TsIdcZt++I7jdjwN/AXUAF6GhoQwbNsjkZCJyTUOG\nwNtve2aBSEnxLFrRoYPZqURE8iUNe8gFnE4ngwa9xJo1/8Pf383Bg5CQ8A1QkAIFnqRz54LMmzfT\n7JiSSVauCStny7V+/hnuvz/9Cm0OB5w96/0p0cTrrF4TVs8n4mua7SGP6NDhMdauDSAx8Q1stv8S\nGDgef//S2Gx+NGjQlBkzZpkdUUSu5eBBzxzAV0pNhdOnoWRJUyKJiORnuvNrcfHx8RQpEobLdRYo\nAEDBgu2YNasLDz74IKGhoeYGlJtm5ZqwcrZca/duqFEDEhIu7ytRAmJjPavBiaVZvSasnk/E1/TA\nWx5gt9svfhMv/cnUAM4TEhKijq9YS0ICfPQRzJ2rh7quVKWKZwW4oCDPcIfixWHFitzT8VXHSkTy\nGN35zQWef34Q8+b9F6ezF4GB6ylf/me2bv0fDofD7GiSBVauiSxnu3AB6taFI0cud5a++w7q1/du\nwNwsIQFOnIDSpTMOg7Ci996DgQM9uZs0gSVLoHBhs1P5nJXrFayfT8TXMlMT6vzmAoZh8M47s1m9\nej233lqG4cNfoEiRImbHkiyyck1kOdu//w0jRkBi4uV9d97pmd5Lcp/vv4fWrS8/pBcYCK1awdKl\n5uYygZXrFayfT8TX9MBbHmGz2Xjmmd4880xvs6OIXN2RI+k7vgDHj5uTRbJv9er0Y5STk2HNGvPy\niIh4US4ZdCYiltasmWc86yUFCkDjxublkewJD/eMUb5SsWLmZBER8TJ1fnOZuLg42rfvhsNRlJIl\nK/DRRx9nOMbtdtO372AcjiKEhhZnxIh/6c9ikrPatfMMewgMBH9/z7y2s2ebnUqyqkcPqFABQkIu\nP6in76eI5BEa85vLdOjQneXLDZKSpgJ7cTg6sHr1F9S/4sGi0aPHM2nSCpzOaCCJkJD2TJ3aj169\nepqWWy6zck1kO5thgNudOx7okutLTIRPPoFz5zwPvEVGmp3IFFauV7B+PhFf0wNveVBoaBjx8b8D\npQDw9x/G6NGFePnll9OOqV27CZs3DweaX9yzkDZtvuSrr6J9nlcysnJNWDmbiBmsXhNWzyfia5rn\nNw8qVKgYsOvilkGBArsp9rexeOHhxbHZtqdt+/tvp1Sp4r4LKSIiImJRuvObyyxdupQuXXqTktKN\nwMC9lC17hM2b1xESEpJ2zLZt27j77iakpDyAzZaEw7GWX375H+XKlTMxuVxi5ZqwcjYRM1i9Jqye\nT8TXNOwhj9q8eTPffvstRYoUoVu3buk6vpccPnyYL774An9/fzp27Eh4eLgJSeVqrFwTVs4mm4bz\nIAAAC3VJREFUYgar14TV84n4mjq/IhZk5ZqwcjYRM1i9JqyeT8TXNOZXREREROQK6vyKiIiISL6h\nzq+IiIiI5BtZ7vwOGTKEyMhIoqKi6NChA+fOnfNmLhER8bY//4R334U5c+D0abPTiIiYIssPvH3z\nzTc0bdoUPz8/XnzxRQAmTJiQsQENxhdJx8o1YeVskk07dsDdd0NyMthsEBoKv/wCpUubnczSrF4T\nVs8n4ms5+sBb8+bN8fPznF6/fn2OHDmS1UuJDxmGwdmzZ3G5XGZHERFfGjgQzp+HhARwOj13fkeP\nNjuViIjPeWXM79y5c2nTpo03LiU56ODBg1SpUpPw8HKEhBRh9uy5ZkcSEV+JjYUr74a4XHD0qHl5\nRERMYr/eF5s3b05sbGyG/ePGjaNdu3YAjB07lsDAQLp27XrN64waNSrt80aNGtGoUaOspZVseeCB\nzuzf34XU1KHAHgYMaEStWlHUrl3b7Gh5WkxMDDExMWbHyLRcXa9xcfDHH1C2LBQpYnYaa2nbFvbs\n8dz1BXA44IEHzM1kQbmtXiGX16xINmWlZrO1yMW8efOYPXs23333HUFBQVdvQOORLMHtdhMQEIhh\nJHHpd57g4N5MmVKTZ5991txw+YyVa8LK2W7o66/hkUfAzw9SUuC996BbN7NTWUdKCvTpAx9+6Pk3\n6t8fJk70jP+Va7J6TVg9n4iv5egKbytXrmTw4MGsXbuWsLCwbIUQ3yhWrAxnzkQD9wFJhIbexaJF\nY9Lu4otvWLkmrJztui5cgFKlID7+8r7gYM+dzjJlzMtlRZe+v+r0ZorVa8Lq+UR8LUcfeOvXrx8X\nLlygefPm1KxZk+eeey6rlxIfWbjwPRyODhQs+E9CQ2vRpElV2rZta3Yskew7fNhzN/NKgYGwe7c5\neazMZlPHV0TytWwNe8hUA/qt1FL279/Pxo0bKVmyJI0bN8amN0Gfs3JNWDnbdZ0/77nzm5BweV9w\nsGd6r1tuMS+X5HpWrwmr5xPxtRwd9uDNECL5iZVrwsrZbuijj+CppyAgwDOX7aRJ0Lev2akkl7N6\nTVg9n4ivqfMrYkFWrgkrZ8uUY8c843wjIqB8ebPTSB5g9Zqwej4RX1PnV8SCrFwTVs4mYgar14TV\n84n4Wo4+8CYiIhb05ZfQubNnWrP9+81OIyJiObrzK+JjVq4JK2eTTHj/fc84Z6fTM/tFaCj8+qse\n+ssGq9eE1fOJ+Jru/IqI5CejRl1ewS011TPv8VwtYy4iciV1fkVE8oqUlPTbbjckJZmTRUTEotT5\nFRHJK55+GhyOy9sOh2f8r4iIpLGbHUBERLxk1CgICoKFC6FgQZg4EWrUMDuViIil6IE3ER+zck1Y\nOZuIGaxeE1bPJ+JreuBNREREROQK6vyKiIiISL6hzq+ISE7YuBE6doSHHoKvvzY7jYiIXKQxvyI+\nZuWasHK2XGXTJmjc+PKcuw4HREdDu3bm5pKbZvWasHo+EV/TmF8RETNMnXq54wuez197zbw8IiKS\nRp1fERFvc7ky7ktN9X0OERHJQJ1fERFve/75jItNDBxoXh4REUmjMb8iPmblmrBytlxn1SrPUIeU\nFOjXD7p2NTuRZIHVa8Lq+UR8LTM1oc6viI9ZuSasnE3EDFavCavnE/E1PfAmIiIiInIFdX5FRERE\nJN9Q51dERERE8o0sd35fffVVoqKiqFGjBk2bNuXw4cPezCUiIiIi4nVZ7vwOHTqUrVu3smXLFh5+\n+GFGjx7tzVxeFRMTY3YE0zOY3b4ySGZZ4XukDNbIYHb7cmNW+B4pgzLcrCx3fgsWLJj2+YULFwgL\nC/NKoJxghW+G2RnMbl8ZJLOs8D1SBmtkMLt9uTErfI+UQRlulj07J7/88sssWLAAh8PBhg0bvJVJ\nRERERCRHXPfOb/PmzbnzzjszfHz55ZcAjB07lkOHDtGjRw8GavUiEREREbE4ryxycejQIdq0acPv\nv/+e4WuVKlVi37592W1CJM+oWLEie/fuNTvGValeRdKzcr2Calbk7zJTs1ke9rBnzx4qV64MwNKl\nS6lZs+ZVj7PyDw0RSU/1KpK7qGZFbl6W7/x27NiRXbt24e/vT8WKFZk5cybh4eHeziciIiIi4jVe\nGfYgIiIiIpIb+GSFt+nTpxMZGckdd9zBsGHDfNHkVU2ZMgU/Pz9Onz7t87aHDBlCZGQkUVFRdOjQ\ngXPnzvms7ZUrV1K1alUqV67MxIkTfdbuJYcPH6Zx48ZUq1aNO+64gzfffNPnGQDcbjc1a9akXbt2\nprR/9uxZOnbsSGRkJLfffrulZ0hRzZpXs6rXy1SzmaN61XusFWo2V9WrkcNWr15tNGvWzEhOTjYM\nwzD++uuvnG7yqg4dOmS0bNnSiIiIME6dOuXz9letWmW43W7DMAxj2LBhxrBhw3zSrsvlMipWrGgc\nOHDASE5ONqKioozt27f7pO1Ljh07Zvzyyy+GYRhGXFycUaVKFZ9nMAzDmDJlitG1a1ejXbt2Pm/b\nMAzj8ccfN+bMmWMYhmGkpKQYZ8+eNSXHjahmPcyoWdVreqrZG1O9eug91vyazU31muN3fmfOnMnw\n4cMJCAgAoESJEjnd5FUNGjSISZMmmdI2eKaN8/Pz/HPXr1+fI0eO+KTdTZs2UalSJSIiIggICKBz\n584sXbrUJ21fUqpUKWrUqAFAaGgokZGR/Pnnnz7NcOTIEZYvX87TTz+NYcJIn3PnzrFu3Tqeeuop\nAOx2O4ULF/Z5jsxQzXqYUbOq18tUs5mjevXQe6zeY2+mXnO887tnzx6+//577rrrLho1asRPP/2U\n001msHTpUsqWLUv16tV93vbVzJ07lzZt2vikraNHj1KuXLm07bJly3L06FGftH01Bw8e5JdffqF+\n/fo+bXfgwIH8+9//Tvvh6GsHDhygRIkSPPnkk9SqVYtevXrhdDpNyXIjqtmMfFWzqtfLVLOZo3rN\nSO+xeo+9Ub1ma4W3S5o3b05sbGyG/WPHjsXlcnHmzBk2bNjAjz/+SKdOndi/f783ms10hvHjx7Nq\n1aq0fTn1W8m1MowbNy5tDMzYsWMJDAyka9euOZLh72w2m0/ayYwLFy7QsWNH3njjDUJDQ33W7rJl\nywgPD6dmzZqmLb3ocrnYvHkzM2bMoG7dugwYMIAJEyYwZswYU/KoZq+fwayaVb16qGbTU71eP4Pe\nYz30HnsT9ZrTYzBatWplxMTEpG1XrFjROHnyZE43m+a3334zwsPDjYiICCMiIsKw2+3GLbfcYhw/\nftxnGS55//33jQYNGhgJCQk+a3P9+vVGy5Yt07bHjRtnTJgwwWftX5KcnGy0aNHCmDp1qs/bHj58\nuFG2bFkjIiLCKFWqlOFwOIzu3bv7NMOxY8eMiIiItO1169YZbdu29WmGzFLNXubrmlW9eqhmM0/1\nepneY/Uee8mN6jXHO7+zZs0yRowYYRiGYezatcsoV65cTjd5XWYNxl+xYoVx++23GydOnPBpuykp\nKUaFChWMAwcOGElJSaYMxk9NTTW6d+9uDBgwwKftXk1MTIzxwAMPmNL2fffdZ+zatcswDMMYOXKk\nMXToUFNy3Ihq1sOMmlW9ZqSavT7Vq4feY61Rs7mlXnO885ucnGw89thjxh133GHUqlXLWLNmTU43\neV233nqrKYVZqVIlo3z58kaNGjWMGjVqGM8++6zP2l6+fLlRpUoVo2LFisa4ceN81u4l69atM2w2\nmxEVFZX2+lesWOHzHIbhKUyznkTdsmWLUadOHaN69epG+/btLfnkuGGoZi8xq2ZVr+mpZq9P9eqh\n91hr1GxuqVctciEiIiIi+YY5j+WJiIiIiJhAnV8RERERyTfU+RURERGRfEOdXxERERHJN9T5FRER\nEZF8Q51fEREREck31PkVERERkXxDnV8RERERyTf+HwrrfS1F8/HuAAAAAElFTkSuQmCC\n"
      }
     ],
     "prompt_number": 414
    },
    {
     "cell_type": "heading",
     "level": 1,
     "source": [
      "Implementation of similarity transform fitting"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from sklearn.metrics import euclidean_distances\n",
      "\n",
      "# Taken from /usr/lib/python-2.7/site-packages/skimage/transform/_geometric.py:499\n",
      "class SimilarityTransform:\n",
      "    def transform(self, points):\n",
      "        assert points.ndim==2\n",
      "        assert self._matrix is not None, \"Please run fit() before estimating points\"\n",
      "        return self._matrix.dot(np.hstack([points, np.ones((len(points),1))]).T).T\n",
      "\n",
      "    def get_error(self, points, test_points):\n",
      "        \"\"\"For each of test_points, returns the minimum Euclidean\n",
      "        distance between it and the rest of the data.\"\"\"\n",
      "        points = self.transform(points)\n",
      "        return np.sqrt(np.sum((test_points-points)*(test_points-points), axis=1))\n",
      "        \n",
      "\n",
      "    def fit(self, src, dst):\n",
      "        \"\"\"Learns a 2x3 2D affine transformation matrix.\n",
      "\n",
      "        Number of source and destination coordinates must match.\n",
      "\n",
      "        The transformation is defined as::\n",
      "\n",
      "            X = a0*x - b0*y + a1\n",
      "            Y = b0*x + a0*y + b1\n",
      "\n",
      "        These equations can be transformed to the following form::\n",
      "\n",
      "            0 = a0*x - b0*y + a1 - X\n",
      "            0 = b0*x + a0*y + b1 - Y\n",
      "\n",
      "        which exist for each set of corresponding points, so we have a set of\n",
      "        N * 2 equations. The coefficients appear linearly so we can write\n",
      "        A x = 0, where::\n",
      "\n",
      "            A   = [[x 1 -y 0 -X]\n",
      "                   [y 0  x 1 -Y]\n",
      "                    ...\n",
      "                    ...\n",
      "                  ]\n",
      "            x.T = [a0 a1 b0 b1 c3]\n",
      "\n",
      "        In case of total least-squares the solution of this homogeneous system\n",
      "        of equations is the right singular vector of A which corresponds to the\n",
      "        smallest singular value normed by the coefficient c3.\n",
      "\n",
      "        Parameters\n",
      "        ----------\n",
      "        src : (N, 2) array\n",
      "            Source coordinates.\n",
      "        dst : (N, 2) array\n",
      "            Destination coordinates.\n",
      "\n",
      "        \"\"\"\n",
      "        # print \"Fitting with\", src,\" to \",dst\n",
      "        xs = src[:, 0]\n",
      "        ys = src[:, 1]\n",
      "        xd = dst[:, 0]\n",
      "        yd = dst[:, 1]\n",
      "        rows = src.shape[0]\n",
      "\n",
      "        # params: a0, a1, b0, b1\n",
      "        A = np.zeros((rows * 2, 5))\n",
      "        A[:rows, 0] = xs\n",
      "        A[:rows, 2] = - ys\n",
      "        A[:rows, 1] = 1\n",
      "        A[rows:, 2] = xs\n",
      "        A[rows:, 0] = ys\n",
      "        A[rows:, 3] = 1\n",
      "        A[:rows, 4] = xd\n",
      "        A[rows:, 4] = yd\n",
      "\n",
      "        _, _, V = np.linalg.svd(A)\n",
      "\n",
      "        # solution is right singular vector that corresponds to smallest\n",
      "        # singular value\n",
      "        a0, a1, b0, b1 = - V[-1, :-1] / V[-1, -1]\n",
      "\n",
      "        self._matrix = np.array([[a0, -b0, a1],\n",
      "                                 [b0,  a0, b1]])\n",
      "        return self\n"
     ],
     "language": "python",
     "outputs": [],
     "prompt_number": 159
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "st = SimilarityTransform()\n",
      "st.fit(np.array([[0,0], [0,1], [1,0]]),\n",
      "       np.array([[1,1], [1,3], [3,1]]))\n",
      "st.get_error(np.array([[0,0], [0,1], [1,0]]),\n",
      "             np.array([[1,1], [1,3], [3,1]]))\n",
      "# st.transform(np.array([[0,0], [0,1]]))"
     ],
     "language": "python",
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 160,
       "text": [
        "array([  1.29473141e-15,   1.79018084e-15,   4.96506831e-16])"
       ]
      }
     ],
     "prompt_number": 160
    },
    {
     "cell_type": "heading",
     "level": 1,
     "source": [
      "Implementation of RANSAC"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def ransac(data,fit_model,get_error,n,k,t,d,debug=False,return_all=False):\n",
      "    \"\"\"Fit model parameters to data using the RANSAC algorithm\n",
      "    \n",
      "This implementation written from pseudocode found at\n",
      "http://en.wikipedia.org/w/index.php?title=RANSAC&oldid=116358182\n",
      "\n",
      "{{{\n",
      "Given:\n",
      "    data - a set of observed data points\n",
      "    fit_model - a function that returns the model fit to the points\n",
      "    get_error - a function that returns, for each point, the error of that point\n",
      "    n - the minimum number of data values required to fit the model\n",
      "    k - the maximum number of iterations allowed in the algorithm\n",
      "    t - a threshold value for determining when a data point fits a model\n",
      "    d - the number of close data values required to assert that a model fits well to data\n",
      "Return:\n",
      "    bestfit - model parameters which best fit the data (or nil if no good model is found)\n",
      "}}}\n",
      "\"\"\"\n",
      "    iterations = 0\n",
      "    bestfit = None\n",
      "    besterr = numpy.inf\n",
      "    best_inlier_idxs = None\n",
      "    while iterations < k:\n",
      "        maybe_idxs, test_idxs = random_partition(n,data.shape[0])\n",
      "        maybeinliers = data[maybe_idxs,:]\n",
      "        test_points = data[test_idxs]\n",
      "        maybemodel = fit_model(maybeinliers)\n",
      "        test_err = get_error(maybemodel, test_points)\n",
      "        also_idxs = test_idxs[test_err < t] # select indices of rows with accepted points\n",
      "        alsoinliers = data[also_idxs,:]\n",
      "        if debug:\n",
      "            print 'test_err.min()',test_err.min()\n",
      "            print 'test_err.max()',test_err.max()\n",
      "            print 'numpy.mean(test_err)',numpy.mean(test_err)\n",
      "            print 'test_err',test_err\n",
      "            print 'iteration %d:len(alsoinliers) = %d'%(\n",
      "                iterations,len(alsoinliers))\n",
      "        if len(alsoinliers) > d:\n",
      "            betterdata = numpy.concatenate( (maybeinliers, alsoinliers) )\n",
      "            bettermodel = fit_model(betterdata)\n",
      "            better_errs = get_error(bettermodel, betterdata)\n",
      "            thiserr = numpy.mean( better_errs )\n",
      "            if thiserr < besterr:\n",
      "                bestfit = bettermodel\n",
      "                besterr = thiserr\n",
      "                best_inlier_idxs = numpy.concatenate( (maybe_idxs, also_idxs) )\n",
      "        iterations+=1\n",
      "    if bestfit is None:\n",
      "        raise ValueError(\"did not meet fit acceptance criteria\")\n",
      "    if return_all:\n",
      "        return bestfit, {'inliers':best_inlier_idxs}\n",
      "    else:\n",
      "        return bestfit\n",
      "\n",
      "def random_partition(n,n_data):\n",
      "    \"\"\"return n random rows of data (and also the other len(data)-n rows)\"\"\"\n",
      "    all_idxs = numpy.arange( n_data )\n",
      "    numpy.random.shuffle(all_idxs)\n",
      "    idxs1 = all_idxs[:n]\n",
      "    idxs2 = all_idxs[n:]\n",
      "    return idxs1, idxs2"
     ],
     "language": "python",
     "outputs": [],
     "prompt_number": 325
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def learn_ransac_similarity_xform(src, dst, **kws):\n",
      "    \"\"\"Learns a similarity transformation between two sets of\n",
      "    corresponding points, if possible, using RANSAC.\n",
      "\n",
      "    Arguments:\n",
      "        src - The original points\n",
      "        dst - The mapped points\n",
      "        n - the minimum number of data values required to fit the model\n",
      "        k - the maximum number of iterations allowed in the algorithm\n",
      "        t - a threshold value for determining when a data point fits a model\n",
      "        d - the number of close data values required to assert that a model fits well to data\n",
      "\n",
      "\n",
      "    If fitting is impossible, throws an exception. There are three\n",
      "    ways to fix this: Increase k, increase t, or decrease d. \n",
      "    \"\"\" \n",
      "    m = SimilarityTransform()\n",
      "    point_correspondences = np.rollaxis(np.array([src, dst]), 1,0)\n",
      "    assert len(src) == len(dst), \"You must have corresponding points\"\n",
      "    return ransac(point_correspondences,\n",
      "                  fit_model=lambda pairs: SimilarityTransform().fit(pairs[:,0,:], pairs[:,1,:]),\n",
      "                  get_error=lambda model,pairs: model.get_error(pairs[:,0,:], pairs[:,1,:]),\n",
      "                  **kws\n",
      "                  )\n"
     ],
     "language": "python",
     "outputs": [],
     "prompt_number": 411
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "np.array([[cos(0.6),-sin(0.6)],[sin(0.6),cos(0.6)]]).dot(0.75*np.eye(2))"
     ],
     "language": "python",
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 345,
       "text": [
        "array([[ 0.61900171, -0.42348186],\n",
        "       [ 0.42348186,  0.61900171]])"
       ]
      }
     ],
     "prompt_number": 345
    }
   ]
  }
 ]
}