{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "\n",
    "figsize(12,8)\n",
    "import tables\n",
    "from collections import Counter\n",
    "matplotlib.rc(\"image\",cmap=\"gray\")\n",
    "matplotlib.rc(\"image\",interpolation=\"nearest\")\n",
    "from scipy.ndimage import filters,morphology\n",
    "from scipy.ndimage.filters import gaussian_filter as G\n",
    "def crop_black(image):\n",
    "    gray = image\n",
    "    if gray.ndim==3: gray = sum(image,axis=2)\n",
    "    yr = find(sum(gray,axis=1)>0)\n",
    "    y0 = yr[0]\n",
    "    y1 = yr[-1]\n",
    "    xr = find(sum(gray,axis=0)>0)\n",
    "    x0 = xr[0]\n",
    "    x1 = xr[-1]\n",
    "    if image.ndim==3: return image[y0:y1,x0:x1,:]\n",
    "    else: return image[y0:y1,x0:x1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "IOError",
     "evalue": "``1k.h5`` does not exist",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mIOError\u001b[0m                                   Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-2-c2f00917c23a>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mhdf\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtables\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mopenFile\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"1k.h5\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0mimages\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mhdf\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mroot\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0micons\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mimage\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mmean\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcrop_black\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimages\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m77\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mimshow\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/tables/_past.pyc\u001b[0m in \u001b[0;36moldfunc\u001b[0;34m(*args, **kwargs)\u001b[0m\n\u001b[1;32m     33\u001b[0m     \u001b[0;32mdef\u001b[0m \u001b[0moldfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     34\u001b[0m         \u001b[0mwarn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mwarnmsg\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mDeprecationWarning\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mstacklevel\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 35\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0mobj\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     36\u001b[0m     oldfunc.__doc__ = (\n\u001b[1;32m     37\u001b[0m         obj.__doc__ or '') + \"\\n\\n.. warning::\\n\\n    \" + warnmsg + \"\\n\"\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/tables/file.pyc\u001b[0m in \u001b[0;36mopen_file\u001b[0;34m(filename, mode, title, root_uep, filters, **kwargs)\u001b[0m\n\u001b[1;32m    316\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    317\u001b[0m     \u001b[0;31m# Finally, create the File instance, and return it\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 318\u001b[0;31m     \u001b[0;32mreturn\u001b[0m \u001b[0mFile\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmode\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtitle\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mroot_uep\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    319\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    320\u001b[0m \u001b[0mopenFile\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mprevious_api\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mopen_file\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/tables/file.pyc\u001b[0m in \u001b[0;36m__init__\u001b[0;34m(self, filename, mode, title, root_uep, filters, **kwargs)\u001b[0m\n\u001b[1;32m    789\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    790\u001b[0m         \u001b[0;31m# Now, it is time to initialize the File extension\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 791\u001b[0;31m         \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_g_new\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmode\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mparams\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    792\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    793\u001b[0m         \u001b[0;31m# Check filters and set PyTables format version for new files.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/tables/hdf5extension.so\u001b[0m in \u001b[0;36mtables.hdf5extension.File._g_new (tables/hdf5extension.c:3812)\u001b[0;34m()\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/tables/utils.pyc\u001b[0m in \u001b[0;36mcheck_file_access\u001b[0;34m(filename, mode)\u001b[0m\n\u001b[1;32m    155\u001b[0m         \u001b[0;31m# The file should be readable.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    156\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mos\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0maccess\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mos\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mF_OK\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 157\u001b[0;31m             \u001b[0;32mraise\u001b[0m \u001b[0mIOError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"``%s`` does not exist\"\u001b[0m \u001b[0;34m%\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    158\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mos\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpath\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0misfile\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    159\u001b[0m             \u001b[0;32mraise\u001b[0m \u001b[0mIOError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"``%s`` is not a regular file\"\u001b[0m \u001b[0;34m%\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mfilename\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mIOError\u001b[0m: ``1k.h5`` does not exist"
     ]
    }
   ],
   "source": [
    "hdf = tables.openFile(\"1k.h5\")\n",
    "images = hdf.root.icons\n",
    "image = mean(crop_black(images[77]),2)\n",
    "imshow(image)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Harris Corners"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def harris(image,sigma=1.0,area=5.0):\n",
    "    dx = filters.gaussian_filter(image,sigma,(0,1))\n",
    "    dy = filters.gaussian_filter(image,sigma,(1,0))\n",
    "    sxx = filters.gaussian_filter(dx**2,area)\n",
    "    sxy = filters.gaussian_filter(dx*dy,area)\n",
    "    syy = filters.gaussian_filter(dy**2,area)\n",
    "    T = sxx+syy\n",
    "    D = sxx*syy-sxy**2\n",
    "    R = (T**2-4*D)**.5\n",
    "    hi = (T+R)/2\n",
    "    lo = (T-R)/2\n",
    "    return lo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "subplot(121); imshow(image)\n",
    "subplot(122); imshow(harris(image,2.0,6.0)**2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Corners from Median"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "smooth = filters.gaussian_filter(image,3.0)\n",
    "subplot(121); imshow(smooth)\n",
    "subplot(122); imshow(abs(smooth-filters.median_filter(smooth,10))**2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Corners from Level Curves"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'image' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-3-33d85c2061c0>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mI\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m3.0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0ms\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m2.0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mIx\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mIy\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mIxx\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'image' is not defined"
     ]
    }
   ],
   "source": [
    "I = filters.gaussian_filter(image,3.0)\n",
    "s = 2.0\n",
    "Ix = filters.gaussian_filter(image,s,(0,1))\n",
    "Iy = filters.gaussian_filter(image,s,(1,0))\n",
    "Ixx = filters.gaussian_filter(image,s,(0,2))\n",
    "Ixy = filters.gaussian_filter(image,s,(1,1))\n",
    "Iyy = filters.gaussian_filter(image,s,(2,0))\n",
    "k = Ix**2*Iyy + Iy**2*Ixx - 2*Ix*Iy*Ixy\n",
    "subplot(121); imshow(smooth)\n",
    "subplot(122); imshow(abs(k))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Orientation Histograms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'image' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-4-e82b995efd37>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mI\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m1.0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0ms\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m1.0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mIx\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mIy\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgaussian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'image' is not defined"
     ]
    }
   ],
   "source": [
    "I = filters.gaussian_filter(image,1.0)\n",
    "s = 1.0\n",
    "Ix = filters.gaussian_filter(image,s,(0,1))\n",
    "Iy = filters.gaussian_filter(image,s,(1,0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7f3a45565b10>]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAscAAAHfCAYAAAC8m3VFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl41eWd9/HPCQFkXyVAEmVJyMISkLAIaoOKKFOxFhew\nVQetZVymtdqO87TPPIPTKUqdqVVxOmir4obUrVjFVEGDKPsaIAHCpiEIGiAFAQlJfs8fdwJi9uSc\nc/+W9+u6zkUCx3M+V/u7wof7fH/3HXIcxxEAAAAAxdgOAAAAALgF5RgAAACoRDkGAAAAKlGOAQAA\ngEqUYwAAAKAS5RgAAACoVG85vu222xQXF6fBgwfX+pyf/OQnSk5OVkZGhtavXx/WgAAAAEC01FuO\np02bpuzs7Fr/fOHChdqxY4cKCgr01FNP6c477wxrQAAAACBa6i3HF198sbp06VLrn7/11lu69dZb\nJUmjRo1SSUmJDhw4EL6EAAAAQJQ0e+a4qKhIiYmJp79PSEjQ3r17m/uyAAAAQNTFhuNFvn0CdSgU\nqvacpKQk7dy5MxxvBwAAANSqf//+2rFjR5P+22aX4/j4eBUWFp7+fu/evYqPj6/2vJ07d1Yr0cCM\nGTM0Y8YM2zHgMlwXqAnXBWrCdYGa1LRQ21DNHquYNGmSnn/+eUnSihUr1LlzZ8XFxTX3ZQEAAICo\nq3fleOrUqVqyZImKi4uVmJioBx98UKdOnZIkTZ8+XRMnTtTChQuVlJSkdu3a6dlnn414aAAAACAS\n6i3H8+bNq/dFZs+eHZYwCJ6srCzbEeBCXBeoCdcFasJ1gXALOVEaBA6FQswcAwAAIOKa0zs5PhoA\nAACoRDkGAAAAKlGOAQAAgEqUYwAAAKAS5RgAAACoRDkGAAAAKlGOAQAAgEqUYwAAAKAS5RgAAACo\nRDkGAAAAKlGOAQAAgEqUYwAAAKAS5RgAAACoRDkGAAAAKlGOAQAAgEqUYwAAAKAS5RgAAACoRDkG\nAAAAKlGOAQAAgEqUYwAAAKAS5RgAAACoRDkGAAAAKlGOAQAAgEqUYwAAAKBSrO0AQFCUl0u7dkmf\nfy7t33/248ABaehQ6R//UUpOtp0UAIDgohwDUZCTI917r3TokHTeeVLPnmceY8ZI3btLH30kXXSR\nNGCANG2adP31UocOtpMDABAsIcdxnKi8USikKL0V4Bq7dkm/+IW0dq30yCPSdddJoVDtzz91Slq4\nUHr2WWnJEumaa6Qf/ciUZgAA0DDN6Z3MHAMRcOSI9K//Ko0cKQ0fLuXnm5XguoqxJLVsaQrxX/4i\nbdsmDRki/fCH0v/7fxL/tgQAIPJYOQbCyHGk556TfvUr6YorpJkzpd69m/eaX3whTZpkxi3++Eep\nVauwRAUAwLea0zspx0AY/frX0vz5piBnZobvdY8fl266STp6VHrjDalTp/C9NgAAfsNYBeACf/iD\nNHeutHhxeIuxJLVtK73+upSebuaPCwvD+/oAAMCgHANh8Oc/S7/5jfTee1JcXGTeo0UL6fHHzXZv\nY8ZIGzdG5n0AAAgyxiqAZnrvPenmm6X33zc30EXDq69Kd98tvfyydPnl0XlPAAC8gpljwJKVK6Wr\nrzZzwNHebu2jj6TJk6XVq6U+faL73gAAuBnlGLAgP18aN07605+kf/gHOxkeeUT661+lDz80YxcA\nAIAb8oCo++wzacIEU05tFWNJuu8+KSZG+q//spcBAAA/YeUYaKSyMmnECHM4x/33205jinpmpvS3\nv0nDhtlOAwCAfawcA1H0P/8jde1qVm3d4LzzpEcflX7wA+nECdtpAADwNlaOgUbYv18aPNjcDJeW\nZjvNGY4jTZ0q9ehhtnsDACDIuCEPiJJbbpF69ZJmzbKdpLrDh6WMDOnpp808NAAAQUU5BqJg6VJz\nhHN+vtS+ve00NVu82BT43FypWzfbaQAAsIOZYyDCysrMoRv//d/uLcaSdNll0o03Sj/+sRm1AAAA\njUM5BhrgySfNPO/119tOUr+ZM6Xt26X5820nAQDAexirAOrh1pvw6rJokVnpzsvjcBAAQPAwVgFE\n0L/8i3Tbbd4pxpIZr+jaVXr1VdtJAADwFlaOgTp44Sa82ixcaIp9bq45RQ8AgKBg5RiIAK/chFeb\nq66SzjlH+stfbCcBAMA7KMdALf74R+ncc71xE15NQiHp//5f6T//k50rAABoKMoxUIOKCul3v5Me\nfNCUTK+aNMmsgL/zju0kAAB4A+UYqEF2thmlGDvWdpLmiYkxq8e//jWrxwAANATlGKjB449LP/2p\nt1eNq0yeLB05YrZ3AwAAdWO3CuBbtm6VvvMd6dNPzQ1tfvDii9JTT5m9mgEA8Dt2qwDC6IknzPHL\nfinGkjRlirRvn7Rkie0kAAC4GyvHwDeUlEh9+0pbtki9e9tOE17PPCO9/DLjFQAA/2PlGAiTZ581\n+wP7rRhL0s03Szt2SMuX204CAIB7sXIMVCovl5KTpXnzpFGjbKeJjP/9X+mtt8zpeQAA+BUrx0AY\nvPOOOfTDr8VYkqZNkzZsMGMjAACgOsoxUOnxx6Wf/MR2ishq3dqMV7zwgu0kAAC4E2MVgKTNm6Ur\nrpD27JFatbKdJrI2bzZz1Z9+ag4JAQDAbxirAJrpiSekf/on/xdjSRo0SOreXcrJsZ0EAAD3YeUY\ngXfokNS/vzn8Iy7Odpro+N3vpE2bzO4cAAD4TXN6J+UYgffb35ob1ObOtZ0kevbvl9LSpKIiqW1b\n22kAAAgvxiqAJiork5580v834n1bz57S6NHSX/5iOwkAAO5COUagffCBGaUYPtx2kuhj1woAAKqj\nHCPQXntNuvFG2yns+N73pBUrzIgFAAAwKMcIrLIy6c03pcmTbSexo21bU5Bfftl2EgAA3INyjMBa\nskTq08c8gorRCgAAzkY5RmC99pp0/fW2U9iVlSUVF5uDQQAAAOUYAVVeHuyRiioxMdIPf8jqMQAA\nVSjHCKSPP5Z69zaHfwTdzTdLL75o/sEAAEDQUY4RSK+9Jl13ne0U7pCebvY9/vBD20kAALCPcozA\nqaiQXn+dcvxN3JgHAIBBOUbgLF8ude8uDRhgO4l7TJ0qvfWWdOyY7SQAANhFOUbgvPoqq8bfFhcn\njRljblIEACDIKMcIlKqRiqBv4VaTqhvzAAAIMsoxAmXVKqljRyktzXYS9/nud6VPPpGOHrWdBAAA\neyjHCBR2qahd+/bSqFHSBx/YTgIAgD2UYwSG41CO6zNxovTuu7ZTAABgD+UYgbF2rdS6tTRokO0k\n7nXVVaYcO47tJAAA2FFvOc7OzlZqaqqSk5M1a9asan9eXFysK6+8UkOHDtWgQYP03HPPRSIn0GxV\nq8ahkO0k7pWaav73yc+3nQQAADtCjlP7GlF5eblSUlK0aNEixcfHa8SIEZo3b57SvnE304wZM3Ty\n5Ek99NBDKi4uVkpKig4cOKDY2Niz3ygUUh1vBUSU40jJyaYgDx1qO4273XmnOVb75z+3nQQAgKZp\nTu+sc+V41apVSkpKUp8+fdSyZUtNmTJFCxYsOOs5vXr10pEjRyRJR44cUbdu3aoVY8C2DRvMrxkZ\ndnN4QdVoBQAAQVRnOS4qKlJiYuLp7xMSElRUVHTWc+644w5t2bJFvXv3VkZGhh577LHIJAWagZGK\nhrv0UrPlHVu6AQCCqM5yHGpAk5g5c6aGDh2qffv2acOGDbr77rt1lL9V4TJvvilNnmw7hTe0by+N\nHs2WbgCAYKpz/iE+Pl6FhYWnvy8sLFRCQsJZz1m2bJl+9atfSZL69++vvn37atu2bcrMzKz2ejNm\nzDj9dVZWlrKyspoRHWiYvXulL7+Uhg+3ncQ7qkYrrrnGdhIAAOqXk5OjnJycsLxWnTfklZWVKSUl\nRYsXL1bv3r01cuTIajfk3XffferUqZP+/d//XQcOHNDw4cOVm5urrl27nv1G3JAHS+bOld55R/rz\nn20n8Y78fGnCBOnTTxlFAQB4T3N6Z50rx7GxsZo9e7YmTJig8vJy3X777UpLS9OcOXMkSdOnT9cv\nf/lLTZs2TRkZGaqoqNBvf/vbasUYsGnxYumyy2yn8JbUVKlFCykvTxo40HYaAACip86V47C+ESvH\nsMBxpPh4aelSsz0ZGu6uu6R+/djSDQDgPRHbyg3wuq1bpVatTMlD47ClGwAgiCjH8LXFi6XLL2du\ntinY0g0AEESUY/jaokXMGzdVu3ZmS7fFi20nAQAgeijH8K2yMmnJErMCiqaZOJHRCgBAsFCO4Vtr\n10qJiVJcnO0k3lU1d8y9tACAoKAcw7fYwq35UlLObOkGAEAQUI7hW5Tj5guFzOrxwoW2kwAAEB2U\nY/jSiRNmp4XvfMd2Eu9jSzcAQJBQjuFLn3wiDRkidehgO4n3XXqptHo1W7oBAIKBcgxfYqQifNq1\nky68kC3dAADBQDmGLy1aZA7/QHhMmCC9/77tFAAARB7lGL5z+LC0bZs5wALhcfHF0scf204BAEDk\nUY7hOx9+KI0ZI7VqZTuJfwwbJu3aJZWU2E4CAEBkUY7hO4sXM1IRbi1bSiNHSsuX204CAEBkUY7h\nO9yMFxkXXSQtXWo7BQAAkUU5hq/s3SsVF0sZGbaT+M9FFzF3DADwP8oxfGXxYrMvbwxXdtiNHi2t\nWyedPGk7CQAAkUOFgK8sWsRIRaR06CClpkpr1thOAgBA5FCO4RuOw7xxpDFaAQDwO8oxfGPrVrN9\nW//+tpP4F+UYAOB3lGP4RtVIRShkO4l/XXSR9MknUkWF7SQAAEQG5Ri+sXSp9J3v2E7hbz17St26\nSXl5tpMAABAZlGP4xooV0oUX2k7hf4xWAAD8jHIMXygqko4fl5KSbCfxv4sv5jAQAIB/UY7hCytW\nmH14mTeOPFaOAQB+RjmGLyxfzkhFtCQnS19/LX32me0kAACEH+UYvsC8cfSEQqweAwD8i3IMzyst\nlTZskEaMsJ0kOCjHAAC/ohzD8zZsMAd/dOhgO0lwUI4BAH5FOYbnVd2Mh+gZNkzas0c6fNh2EgAA\nwotyDM/jZrzoi42VRo6Uli2znQQAgPCiHMPzuBnPDkYrAAB+RDmGp+3fL/3972Z7MUTXRRdxGAgA\nwH8ox/C0qnnjGK7kqBs9Wlq/3ux5DACAX1Ap4GnLl3Mzni3t20vp6dKaNbaTAAAQPpRjeBo349nF\naAUAwG8ox/CsU6ekdevMrgmw4+KLuSkPAOAvlGN4Vm6u1KeP1KmT7STBNXas2c6tosJ2EgAAwoNy\nDM9iCzf74uKkc8+VtmyxnQQAgPCgHMOzuBnPHUaNklavtp0CAIDwoBzDs1g5dofMTHasAAD4B+UY\nnvTFF9LBg1Jqqu0kyMxk5RgA4B+UY3jSihVmlwoO/7Bv6FAzc3zypO0kAAA0H9UCnsRIhXu0ayf1\n7y9t3mw7CQAAzUc5hidxM567jBjB3DEAwB8ox/CcsjJTxEaNsp0EVbgpDwDgF5RjeM7mzVJiotSl\ni+0kqEI5BgD4BeUYnrNiBSMVbjNkiLRtm3TihO0kAAA0D+UYnrN8OTfjuc0555ht9TZutJ0EAIDm\noRzDc7gZz50YrQAA+AHlGJ5y8KB04ICUnm47Cb6NHSsAAH5AOYanrFplVihbtLCdBN/GyjEAwA8o\nx/CUdeuk4cNtp0BNBg6Udu+WvvrKdhIAAJqOcgxPWbdOuuAC2ylQk1atpEGDpPXrbScBAKDpKMfw\nFMqxuzFaAQDwOsoxPOPQIXNDXlKS7SSoDeUYAOB1lGN4xvr10tChUgxXrWuxYwUAwOuoGfCM9esZ\nqXC71FSpqEgqKbGdBACApqEcwzOYN3a/2Fizur9une0kAAA0DeUYnkE59gbmjgEAXkY5hiccPSoV\nFpqP7eFulGMAgJdRjuEJGzdKgwebj+3hbpRjAICXUY7hCYxUeMeAAWbLveJi20kAAGg8yjE8gXLs\nHTEx5v+rtWttJwEAoPEox/AEyrG3MFoBAPAqyjFc78QJaccOaeBA20nQUJRjAIBXUY7heps2mV0q\nWre2nQQNRTkGAHgV5Riux0iF9/TrJx07Ju3fbzsJAACNQzmG61GOvScUYvUYAOBNlGO4HuXYmyjH\nAAAvohzD1UpLpbw8acgQ20nQWJRjAIAXUY7hanl5Ut++Utu2tpOgsYYPpxwDALyHcgxXY6TCu847\nz6z8c1MeAMBLKMdwNcqxd4VC0tCh0saNtpMAANBwlGO4GuXY2zIypA0bbKcAAKDhKMdwrfJyKTfX\nrD7Cm1g5BgB4DeUYrrVtm9Srl9Spk+0kaKqhQ1k5BgB4C+UYrsVIhfelpUm7d0snTthOAgBAw1CO\n4VqUY+9r1UpKSZE2b7adBACAhqEcw7Uox/7ATXkAAC+hHMOVKiqk9eulYcNsJ0FzcVMeAMBL6i3H\n2dnZSk1NVXJysmbNmlXjc3JycjRs2DANGjRIWVlZ4c6IANq1S+rcWere3XYSNBc35QEAvCS2rj8s\nLy/XPffco0WLFik+Pl4jRozQpEmTlJaWdvo5JSUluvvuu/W3v/1NCQkJKi4ujnho+B8jFf6RkWG2\n5KuokGL4rAoA4HJ1/lW1atUqJSUlqU+fPmrZsqWmTJmiBQsWnPWcl19+WZMnT1ZCQoIkqTtLfQgD\nyrF/dO1qPgXYvdt2EgAA6ldnOS4qKlJiYuLp7xMSElRUVHTWcwoKCnTo0CGNGzdOmZmZeuGFFyKT\nFIFCOfYXbsoDAHhFnWMVoVCo3hc4deqU1q1bp8WLF+v48eO68MILNXr0aCUnJ1d77owZM05/nZWV\nxXwyauQ4phxzM55/VN2UN3my7SQAAD/KyclRTk5OWF6rznIcHx+vwsLC098XFhaeHp+okpiYqO7d\nu6tNmzZq06aNLrnkEm3cuLHecgzUZt8+M5vaq5ftJAiXoUOluXNtpwAA+NW3F10ffPDBJr9WnWMV\nmZmZKigo0J49e1RaWqr58+dr0qRJZz3nmmuu0ccff6zy8nIdP35cK1euVHp6epMDARs3SkOGSA34\n4AIewVgFAMAr6lw5jo2N1ezZszVhwgSVl5fr9ttvV1pamubMmSNJmj59ulJTU3XllVdqyJAhiomJ\n0R133EE5RrPk5poyBf/o108qKZEOHTI36AEA4FYhx3GcqLxRKKQovRU8bupU6aqrpFtusZ0E4TR2\nrPSf/ymNG2c7CQDA75rTO9l1FK6Tm2vGKuAvnJQHAPACyjFc5euvzel43zhnBj7BSXkAAC+gHMNV\ntmyRkpKk1q1tJ0G4ZWSwcgwAcD/KMVyFm/H8a9AgaetWqbTUdhIAAGpHOYarbNxIOfartm2lvn2l\n/HzbSQAAqB3lGK7CzXj+xmgFAMDtKMdwDcdh5djvuCkPAOB2lGO4RlGR1KKFFBdnOwkihZVjAIDb\nUY7hGlU343FstH9VrRxzHhAAwK0ox3ANRir8r2dPqWVLae9e20kAAKgZ5Riuwc14wcBoBQDAzSjH\ncA1WjoOBm/IAAG5GOYYrnDgh7d4tpabaToJIY+UYAOBmlGO4Ql6elJzMsdFBwMoxAMDNKMdwBUYq\ngmPAAGnfPunoUdtJAACojnIMV+BmvOCIjZXS06VNm2wnAQCgOsoxXIGV42BhtAIA4FaUY1hXdWw0\nK8fBkZFBOQYAuBPlGNYVFZmDIXr2tJ0E0ZKRwVgFAMCdKMewjpGK4Bk8WNq8WaqosJ0EAICzUY5h\nHTfjBU/nzlLXrmZvawAA3IRyDOtYOQ6mwYPNP4wAAHATyjGsoxwH05AhlGMAgPtQjmHViRPSnj0c\nGx1ElGMAgBtRjmHVli3mxLRWrWwnQbRRjgEAbkQ5hlXcjBdcAwaYbfyOHbOdBACAMyjHsIp54+CK\njTXjNFu22E4CAMAZlGNYRTkONkYrAABuQzmGNY7DWEXQUY4BAG5DOYY1e/eaG/Hi4mwngS2UYwCA\n21COYQ2rxqgqx45jOwkAAAblGNYwb4wePcynB0VFtpMAAGBQjmENK8eQGK0AALgL5RjW5Oaycgxp\n8GDKMQDAPSjHsOLrr6Xduzk2GqwcAwDchXIMK/LypORkjo0G5RgA4C6UY1jBvDGqpKVJO3dKJ0/a\nTgIAAOUYllCOUeWcc6R+/aStW20nAQCAcgxLKMf4JkYrAABuQTlG1DmO2eOYcowqlGMAgFtQjhF1\nBw6Ygtyrl+0kcAvKMQDALSjHiLqqkYpQyHYSuAXlGADgFpRjRB3zxvi2hASz9/UXX9hOAgAIOsox\noo5yjG8Lhcw1sWmT7SQAgKCjHCPqKMeoCaMVAAA3oBwjqk6dkrZtk9LTbSeB2wweTDkGANhHOUZU\nbdsmnX++1Lat7SRwG1aOAQBuQDlGVDFSgdoMGiTl50tlZbaTAACCjHKMqKIcozbt20u9e0sFBbaT\nAACCjHKMqKIcoy7sWAEAsI1yjKiiHKMuzB0DAGyjHCNqDh6Ujh41N+QBNaEcAwBsoxwjajZtMtt1\ncWw0akM5BgDYRjlG1DBSgfr06ycVF0t//7vtJACAoKIcI2oox6hPTIzZ0o3VYwCALZRjRA3lGA3B\naAUAwCbKMaKivFzassWsCgJ1yciQNm60nQIAEFSUY0TFzp1SXJzUsaPtJHA7yjEAwCbKMaKCkQo0\n1ODB0ubN5tMGAACijXKMqKAco6E6dZJ69JB27LCdBAAQRJRjRAXlGI3BaAUAwBbKMaKCcozGoBwD\nAGyhHCPijhyRDhyQ+ve3nQRekZHBdm4AADsox4i4zZulgQOlFi1sJ4FXsHIMALCFcoyIy801ZQdo\nqL59pZIS6dAh20kAAEFDOUbEMW+MxoqJMVu6MVoBAIg2yjEijnKMpmC0AgBgA+UYEeU40qZNZhUQ\naAzKMQDABsoxIurTT82R0V272k4Cr6EcAwBsoBwjonJzWTVG0wweLOXnS2VltpMAAIKEcoyI2rBB\nGjrUdgp4Ubt2UkKCtG2b7SQAgCChHCOiKMdojiFDGK0AAEQX5RgRRTlGc3BSHgAg2ijHiJiSEunL\nLzk2Gk3HTXkAgGijHCNiqm7G49hoNBXlGAAQbZRjRAwjFWiu886TTpwwn0AAABANlGNEDOUYzRUK\ncVMeACC6KMeIGMoxwoHRCgBANFGOERGlpdLWrdKgQbaTwOsoxwCAaKIcIyK2bpX69JHatrWdBF5H\nOQYARFO95Tg7O1upqalKTk7WrFmzan3e6tWrFRsbqzfeeCOsAeFNjFQgXAYOlLZvN59GAAAQaXWW\n4/Lyct1zzz3Kzs5WXl6e5s2bp/z8/Bqf98ADD+jKK6+U4zgRCwvvoBwjXNq0kfr2NZ9GAAAQaXWW\n41WrVikpKUl9+vRRy5YtNWXKFC1YsKDa85544gldd911OvfccyMWFN5COUY4MVoBAIiWOstxUVGR\nEhMTT3+fkJCgoqKias9ZsGCB7rzzTklSKBSKQEx4ieNQjhFebOcGAIiW2Lr+sCFF995779XDDz+s\nUCgkx3HqHKuYMWPG6a+zsrKUlZXV4KDwjs8+Mx+F9+hhOwn8IiNDevRR2ykAAG6Vk5OjnJycsLxW\nneU4Pj5ehYWFp78vLCxUQkLCWc9Zu3atpkyZIkkqLi7Wu+++q5YtW2rSpEnVXu+b5Rj+xaoxwq1q\nrMJxzMEgAAB807cXXR988MEmv1ad5TgzM1MFBQXas2ePevfurfnz52vevHlnPWfXrl2nv542bZqu\nvvrqGosxgoNyjHDr3VuqqJD275d69bKdBgDgZ3XOHMfGxmr27NmaMGGC0tPTdeONNyotLU1z5szR\nnDlzopURHkM5RriFQtyUBwCIjpATpb3XqmaS4X99+0rvvSclJ9tOAj+57z4pLk564AHbSQAAbtec\n3skJeQirkhKpuFjq3992EvgNK8cAgGigHCOsNm40227FcGUhzCjHAIBooMIgrJg3RqSkpUm7d0sn\nTthOAgDwM8oxwopyjEhp3VpKTZVyc20nAQD4GeUYYUU5RiRdcIG0bp3tFAAAP6McI2xKS6Vt26RB\ng2wngV8NHy6tXWs7BQDAzyjHCJv8fLONW5s2tpPAry64gHIMAIgsyjHChpEKRNqQIebTiZMnbScB\nAPgV5RhhQzlGpLVpYw6X2bTJdhIAgF9RjhE2lGNEA6MVAIBIohwjLBzHlOOMDNtJ4HfDh7NjBQAg\ncijHCIvPPpPatpV69LCdBH7HyjEAIJIoxwgLRioQLUOHSnl5ZutAAADCjXKMsKAcI1ratpX69ZO2\nbLGdBADgR5RjhAXlGNHEaAUAIFIoxwgLyjGiiZvyAACRQjlGsx06JB08KPXvbzsJgoKVYwBApFCO\n0Wxr1piyEsPVhCgZNkzavFk6dcp2EgCA31Bn0GyrV0sjRthOgSBp31467zwpP992EgCA31CO0WyU\nY9jAaAUAIBIox2g2yjFs4KY8AEAkUI7RLPv2SSdPSn362E6CoGHlGAAQCZRjNEvVqnEoZDsJgmbY\nMCk3Vyors50EAOAnlGM0CyMVsKVTJ6l3b2nbNttJAAB+QjlGs1COYROjFQCAcKMco8kcx+xxTDmG\nLcOHU44BAOFFOUaT7doltW0r9expOwmC6oIL2LECABBelGM0GSMVsO2CC6QNG6TycttJAAB+QTlG\nk1GOYVuXLlKPHtL27baTAAD8gnKMJlu9WsrMtJ0CQcdoBQAgnCjHaJLycmn9esox7OOmPABAOFGO\n0ST5+eZGvC5dbCdB0LFyDAAIJ8oxmoR5Y7jF8OHmU4yKCttJAAB+QDlGk1CO4RbduplPMHbssJ0E\nAOAHlGM0CeUYbsJoBQAgXCjHaLSTJ6UtW6Rhw2wnAQxuygMAhAvlGI2WmyslJUnt2tlOAhiZmebT\nDAAAmotyjEZjpAJuM2qUWTk+dcp2EgCA11GO0WiUY7hN587S+eebTzUAAGgOyjEabc0ayjHc58IL\npeXLbacAAHgd5RiNcuyYtGuXNHiw7STA2caMkZYts50CAOB1lGM0yrp10qBBUqtWtpMAZ2PlGAAQ\nDpRjNArzxnCrAQOkI0ekzz+3nQQA4GWUYzQK5RhuFRMjjR7N6jEAoHkox2gUyjHcjLljAEBzUY7R\nYIcOSV84ybeVAAAXv0lEQVR8IaWk2E4C1Iy5YwBAc1GO0WBr1kgXXCC1aGE7CVCzkSOlDRvMEecA\nADQF5RgNxkgF3K59e3Nj3vr1tpMAALyKcowGoxzDC5g7BgA0B+UYDeI40iefmOIBuBlzxwCA5qAc\no0Hy86UOHaSEBNtJgLpVrRw7ju0kAAAvohyjQZYulS6+2HYKoH59+0plZVJhoe0kAAAvohyjQSjH\n8IpQyIxWMHcMAGgKyjEa5KOPpEsusZ0CaJgxY5g7BgA0DeUY9fr0U7NvbHKy7SRAw7ByDABoKsox\n6lU1UhEK2U4CNExmppSXJx0/bjsJAMBrKMeoFyMV8Jo2baSBA82pjgAANAblGPXiZjx4EXPHAICm\noByjTl98IX3+uTRkiO0kQOMwdwwAaArKMer08cemZLRoYTsJ0DhVK8ccBgIAaAzKMeq0dCnzxvCm\nxESpVStp507bSQAAXkI5Rp2YN4aXMXcMAGgsyjFqdfSotHWrNGKE7SRA0zB3DABoLMoxarVsmTR8\nuNS6te0kQNOwcgwAaCzKMWrFSAW8btgwqaDAfAoCAEBDUI5RK8oxvK5VK1OQV660nQQA4BWUY9To\n5Elp7Vozswl42cUXS0uW2E4BAPAKyjFqtHq1lJoqdexoOwnQPJdfLi1aZDsFAMArKMeoESMV8Iux\nY6XNm6W//912EgCAF1COUaOPPqIcwx/OOceMB334oe0kAAAvoByjmvJys/0V5Rh+wWgFAKChKMeo\nJjdX6tVLOvdc20mA8Bg/nnIMAGgYyjGqYaQCfpORIRUXS4WFtpMAANyOcoxquBkPfhMTI112GavH\nAID6UY5xFscx5fiSS2wnAcKL0QoAQENQjnGWggKpdWvp/PNtJwHCq+qmvIoK20kAAG5GOcZZlixh\npAL+1KeP1KGD2fMYAIDaUI5xluxsacIE2ymAyGC0AgBQH8oxTistlRYvlq680nYSIDIuv1x6/33b\nKQAAbkY5xmmffCINGCD16GE7CRAZ48aZ6/zkSdtJAABu1aBynJ2drdTUVCUnJ2vWrFnV/vyll15S\nRkaGhgwZorFjxyo3NzfsQRF5CxdKEyfaTgFETteuUmqqtGKF7SQAALeqtxyXl5frnnvuUXZ2tvLy\n8jRv3jzl5+ef9Zx+/frpo48+Um5urv7t3/5NP/7xjyMWGJFDOUYQMFoBAKhLveV41apVSkpKUp8+\nfdSyZUtNmTJFCxYsOOs5F154oTp16iRJGjVqlPbu3RuZtIiYPXukL7+UMjNtJwEiq2pLNwAAalJv\nOS4qKlJiYuLp7xMSElRUVFTr8//0pz9pIsuPnvPuu+ZGvBim0OFzY8ZIW7ZIhw/bTgIAcKN6q1Ao\nFGrwi3344Yd65plnapxLhrsxUoGgOOccU5BzcmwnAQC4UWx9T4iPj1dhYeHp7wsLC5WQkFDtebm5\nubrjjjuUnZ2tLl261PhaM2bMOP11VlaWsrKyGp8YYff11+bwj7lzbScBoqNqtOLaa20nAQCEQ05O\njnLCtOoRchzHqesJZWVlSklJ0eLFi9W7d2+NHDlS8+bNU1pa2unnfPbZZ7r00kv14osvavTo0TW/\nUSiket4Klvztb9Kvfy19/LHtJEB0bNgg3XCDtH277SQAgEhoTu+sd+U4NjZWs2fP1oQJE1ReXq7b\nb79daWlpmjNnjiRp+vTp+o//+A8dPnxYd955pySpZcuWWrVqVZMCIfoYqUDQDBliZo4//VQ6/3zb\naQAAblLvynHY3oiVY9dKTpZefVUaOtR2EiB6pkyRrrhCuu0220kAAOHWnN7J3gQBV1AgHTsmZWTY\nTgJE1/jx7HcMAKiOchxwVSMVjdiUBPCFyy+XFi+WKipsJwEAuAnlOOCYN0ZQnX++1K2bxO0RAIBv\nohwH2LFj0rJlZgUNCKLrrpNee812CgCAm1COA+yDD6QRI6SOHW0nAeyoKsfcKwwAqEI5DjBGKhB0\nQ4ZIrVpJa9faTgIAcAvKcUA5DuUYCIUYrQAAnI1yHFB5eaYYfOOgQyCQrr/e7PPNaAUAQKIcBxZb\nuAFG1eE3GzbYzQEAcAfKcUAxUgEYjFYAAL6J46MDqKRESkyU9u+X2rWznQawb80a6aabpG3b+DQF\nAPyA46PRKK+/bo7OpRgDxvDhUmmptGmT7SQAANsoxwH04ovSD39oOwXgHoxWAACqMFYRMJ99Jg0b\nJu3bJ7VubTsN4B4rV0r/+I9ndnIBAHgXYxVosJdeMitkFGPgbCNHmiPV8/JsJwEA2EQ5DhDHkV54\nQbr5ZttJAPcJhaTJkxmtAICgoxwHyPr10okT0tixtpMA7lR1IAgAILgoxwFSdSMe85RAzUaPlg4f\nlvLzbScBANhCOQ6IsjJp3jx2qQDqEhNjRitef912EgCALZTjgFi0yBz8kZJiOwngboxWAECwUY4D\n4sUXuREPaIgxY6QvvpC2b7edBABgA+U4AL76Snr7bWnKFNtJAPdr0UL6/vcZrQCAoKIcB8Cbb0oX\nXSSde67tJIA33HCD2ROcc4sAIHgoxwHA3sZA41xyibmJdelS20kAANHG8dE+t2+fNHCg+bVNG9tp\nAO944gnp44+l+fNtJwEANBbHR6NW8+ZJ115LMQYa65ZbpPfeM/+wBAAEB+XY5xipAJqmUydp6lTp\nqadsJwEARBNjFT62aZM0caL06afmcAMAjbNlizR+vLRnj9Sqle00AICGYqwCNXrxRekHP6AYA001\ncKA5OOfNN20nAQBEC7XJp77+Wnr+eTM3CaDp7r5bevJJ2ykAANFCOfap55+Xhg2T0tNtJwG87Zpr\npF27pNxc20kAANFAOfah8nLpkUekBx6wnQTwvpYtpenTWT0GgKCgHPvQm29K3bqZgwwANN8dd0h/\n/rNUUmI7CQAg0ijHPuM40qxZZtU4FLKdBvCHnj2lK6+UnnvOdhIAQKRRjn3mww+lo0fNnCSA8Lnn\nHul//keqqLCdBAAQSZRjn5k1S/rFL9i+DQi3MWOktm2l99+3nQQAEElUKB9Zv17avFn64Q9tJwH8\nJxQyq8fcmAcA/sYJeT4ydao0fLj085/bTgL40/Hj0nnnSatXS3372k4DAKhNc3on5dgndu2SRo40\nv3bsaDsN4F8//7l08qT0xBO2kwAAakM5hu6+W+rUSZo503YSwN++/FJKS5NWrJCSkmynAQDUhHIc\ncF98IaWmSvn5Ulyc7TSA/82cKW3YYPY+BgC4T3N6Jzfk+cATT0g33EAxBqLl3nulZcuklSttJwEA\nhBsrxx539KjUr5+0fDkf8QLR9Mwz5lCQJUs4cAcA3IaV4wD7wx+kceMoxkC03XqrdOiQ9Ne/2k4C\nAAgnVo49bO9eaehQ8/HugAG20wDBs3ChdP/90qZNUmys7TQAgCqsHAfUz34m3XUXxRiw5aqrpF69\nzIgFAMAfWDn2qOxss33b5s1Smza20wDBtWaNNGmStH271L697TQAAImV48A5ccIcYzt7NsUYsC0z\nU8rKkn73O9tJAADhwMqxB82YYVaMX3vNdhIAkrR7tynJW7ZIPXvaTgMA4BCQACkokMaMkdavlxIS\nbKcBUOW++8ynOn/4g+0kAADKcUA4jjRhgnTlleYvYgDucfCgOany/ffNLjIAAHuYOQ6IV1+V9u+X\n/vmfbScB8G3dukmPPipNnSodO2Y7DQCgqVg59ogjR6T0dGn+fGnsWNtpANTmlluk1q2lp5+2nQQA\ngouxigC4917pq6+kP/7RdhIAdTl6VLrgAuk3v5FuuMF2GgAIJsqxzy1bJl17rbkTvnt322kA1Gft\nWnNAyKpVUp8+ttMAQPAwc+xje/dK119vTuCiGAPeMHy49MAD0k03SadO2U4DAGgMyrGLnTghfe97\n0k9/Kv3DP9hOA6AxfvYzqWNH6cEHbScBADQGYxUu5TjSD34gxcRIL7wghUK2EwForAMHpGHDpJde\nksaNs50GAIKDsQofmjXLHPjx9NMUY8Cr4uKkZ581O1gUF9tOAwBoCFaOXejtt6Xp06WVKzkFD/CD\nX/xCysuTFiyQYmNtpwEA/2O3Ch/Jy5OysqS33pJGj7adBkA4lJZK11wjdelixqRatLCdCAD8jbEK\nnzh0yPwF+sgjFGPAT1q1kt54Q/ryS2naNKm83HYiAEBtWDl2iePHpauvljIypN/9znYaAJFw/LjZ\neaZfP3M/QQzLEwAQEawce9zBg9Lll0u9e0u//a3tNAAipW1b6a9/lbZvl+66y+xKAwBwF8qxZZ99\nJl18sTR2rDR3LjfrAH7Xvr20cKGUmyv98z9TkAHAbSjHFm3ebErxj35k5oz5iBUIhg4dpHffNcdL\n33cfBRkA3ISZY0s+/liaPNnMF//gB7bTALDh8GEzUjVypPT730utW9tOBAD+wMyxxyxYIF17rdnS\niWIMBFeXLtLixeYkvbFjpV27bCcCAFCOo6iiwqwO/dM/mZnDK66wnQiAbZ07S6+/bk7RGz3afA0A\nsIexiijZuVO6/Xbp5EmzYpyUZDsRALdZvVq68Uaz3dt//RdjFgDQVIxVuFhFhfTYY9KoUdKkSWbW\nmGIMoCYjRkjr1kn79pkxi507bScCgOChHEfQ9u3SJZdIr70mLV9u7krn2FgAdenc2fzMuPVW6cIL\npdmzzfHTAIDooBxHQFmZ9N//bVZ+brxRWrJESk62nQqAV4RCZg/kxYuld96R0tOlV19lyzcAiAZm\njsPo5Enp+eelWbOkPn2kp54yx8QCQHMsWiT9y79ILVuaPdEvucR2IgBwt+b0TspxGBw/Lj39tLmB\nZuBA6Ve/MqfeAUC4VFRIr7wi/fKX0pAh0kMPmZ83AIDquCHPkr//3fwF1a+fGZ14800pO5tiDCD8\nYmKkm26Stm6VsrKkcePMASJ//rP51AoAEB6sHDdSWZmZA3zlFXOYx8SJ0v/5P6zgAIiur782/yB/\n6ilpyxazT/Idd0gpKbaTAYB9jFVEWEWFtGyZNG+euYu8Tx9p6lTphhuk3r1tpwMQdAUF0h//KM2d\na8rxLbdIV13FzycAwUU5joDDh6WlS6UPPpDeeEPq2NEU4ilTpP79bacDgOpKS6W33jKjFu+/b/4h\nf9VV5nHhhVJsrO2EABAdlOMwqCrDOTnmUVBg/jL5znfM4R2DB9tOCAANV1YmrVghvfuueezebWaU\nL7nEHDaSkSG1aWM7JQBEBuW4ERxH+uwzKTf3zGPjRqmoyJThrCzzyMyUWrWynRYAwuPzz6X33jMj\nYqtXmxv7UlJMUc7MNI+UFKldO9tJAaD5IlqOs7Ozde+996q8vFw/+tGP9MADD1R7zk9+8hO9++67\natu2rZ577jkNGzYsrCEbq6LCHL+6e/fZj507pU2bpLZtzVZIGRnm1yFDzF8KlOHoy8nJUVZWlu0Y\ncBmui8j7+muzMLBmjSnLa9dKO3ZI3bpJAwaceSQnm1GyhAQzXmYT1wVqwnWBmjSnd9Y5gVZeXq57\n7rlHixYtUnx8vEaMGKFJkyYpLS3t9HMWLlyoHTt2qKCgQCtXrtSdd96pFStWNClMXUpLzehD1aO4\nWNq/3zwOHDjz9f790t69UpcuUt++Zx4XX2yOYx00SDr33LDHQxPxQw014bqIvHPOkUaNMo8q5eXm\n5+f27WceixaZhYWiInNyX+/eUnz8mUdcnNS9uynV3/y1Y0fz/HDiukBNuC4QbnWW41WrVikpKUl9\n+vSRJE2ZMkULFiw4qxy/9dZbuvXWWyVJo0aNUklJiQ4cOKC4uLhqr/f662a14sQJ86j6+tgx6auv\nzjyOHj3za1UZLi01hbfq0b271LOn+cGckmLm6Kq+T0w0q8MAgIZr0UI6/3zzGD/+7D9zHOnIEVOS\nv/nYs8esOhcXSwcPmkdxsfnZ3rHjmUeHDmd/3bZtzY9zzpFat675UVws7dplTgr89iM21jzCXcgB\nBE+d5bioqEiJiYmnv09ISNDKlSvrfc7evXtrLMcvvWRuAGnTxvwArPq6Wzfzw7h9++qPqjLcrh0/\n9ADAllBI6tTJPNLT639+aalZ4DhypPrj6FFzsmjV4/Bh8+uxY2bR5OTJmh/790vvvCOdOnX2o7TU\nrHqXl5ucVUW5RYuzHzEx1b8OhczX33xU/V4oVP3r+h5V/1vV9H3V19/837S23/v21/X9WW3Pa8yf\nhfO/aaqmvNe2beYfaAie++8394mFW53lONTAq/TbMx01/Xf9+/fXm2/SblHdgw8+aDsCXIjrAjU5\ndKju68JxzpRmBEdBAT8vgujtt2v/s/7N2He3znIcHx+vwsLC098XFhYqISGhzufs3btX8fHx1V5r\nx44dTQ4JAAAARENMXX+YmZmpgoIC7dmzR6WlpZo/f74mTZp01nMmTZqk559/XpK0YsUKde7cucaR\nCgAAAMDt6lw5jo2N1ezZszVhwgSVl5fr9ttvV1pamubMmSNJmj59uiZOnKiFCxcqKSlJ7dq107PP\nPhuV4AAAAEC4Re0QEAAAAMDt6hyrCIfs7GylpqYqOTlZs2bNivTbwaUKCws1btw4DRw4UIMGDdLj\njz8uSTp06JDGjx+vAQMG6IorrlBJSYnlpLChvLxcw4YN09VXXy2J6wJSSUmJrrvuOqWlpSk9PV0r\nV67kuoAeeughDRw4UIMHD9ZNN92kkydPcl0E0G233aa4uDgNHjz49O/VdR089NBDSk5OVmpqqt57\n7716Xz+i5bjqEJHs7Gzl5eVp3rx5ys/Pj+RbwqVatmypRx99VFu2bNGKFSv05JNPKj8/Xw8//LDG\njx+v7du367LLLtPDDz9sOyoseOyxx5Senn56pxuuC/z0pz/VxIkTlZ+fr9zcXKWmpnJdBNyePXv0\n9NNPa926ddq0aZPKy8v1yiuvcF0E0LRp05SdnX3W79V2HeTl5Wn+/PnKy8tTdna27rrrLlVUVNT9\nBk4ELVu2zJkwYcLp7x966CHnoYceiuRbwiOuueYa5/3333dSUlKc/fv3O47jOJ9//rmTkpJiORmi\nrbCw0LnsssucDz74wPnud7/rOI7DdRFwJSUlTt++fav9PtdFsB08eNAZMGCAc+jQIefUqVPOd7/7\nXee9997jugio3bt3O4MGDTr9fW3XwcyZM52HH3749PMmTJjgLF++vM7XjujKcU0HhBQVFUXyLeEB\ne/bs0fr16zVq1KizTlOMi4vTgQMHLKdDtP3sZz/TI488opiYMz+OuC6Cbffu3Tr33HM1bdo0XXDB\nBbrjjjt07NgxrouA69q1q+6//36dd9556t27tzp37qzx48dzXUBS7X9v7Nu376xtiBvSRSNajht6\niAiC46uvvtLkyZP12GOPqUOHDmf9WSgU4poJmLfffls9evTQsGHDqh0mVIXrInjKysq0bt063XXX\nXVq3bp3atWtX7aNyrovg2blzp37/+99rz5492rdvn7766iu9+OKLZz2H6wJS/ddBfddIRMtxQw4R\nQXCcOnVKkydP1s0336zvfe97ksy/7vbv3y9J+vzzz9WjRw+bERFly5Yt01tvvaW+fftq6tSp+uCD\nD3TzzTdzXQRcQkKCEhISNGLECEnSddddp3Xr1qlnz55cFwG2Zs0ajRkzRt26dVNsbKy+//3va/ny\n5VwXkFR7n2joYXXfFNFy3JBDRBAMjuPo9ttvV3p6uu69997Tvz9p0iTNnTtXkjR37tzTpRnBMHPm\nTBUWFmr37t165ZVXdOmll+qFF17gugi4nj17KjExUdu3b5ckLVq0SAMHDtTVV1/NdRFgqampWrFi\nhU6cOCHHcbRo0SKlp6dzXUBS7X1i0qRJeuWVV1RaWqrdu3eroKBAI0eOrPvFwj0g/W0LFy50BgwY\n4PTv39+ZOXNmpN8OLrV06VInFAo5GRkZztChQ52hQ4c67777rnPw4EHnsssuc5KTk53x48c7hw8f\nth0VluTk5DhXX3214zgO1wWcDRs2OJmZmc6QIUOca6+91ikpKeG6gDNr1iwnPT3dGTRokHPLLbc4\npaWlXBcBNGXKFKdXr15Oy5YtnYSEBOeZZ56p8zr4zW9+4/Tv399JSUlxsrOz6319DgEBAAAAKkX8\nEBAAAADAKyjHAAAAQCXKMQAAAFCJcgwAAABUohwDAAAAlSjHAAAAQCXKMQAAAFDp/wPzmWnb3n3u\n9wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3a455cab50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def between(x,a,b,r=0.03):\n",
    "    return 1.0/(1+exp((a-x)/r))/(1+exp((x-b)/r))\n",
    "plot(between(linspace(0.0,1.0,100),0.2,0.4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'Ix' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-6-e7a9bf748be4>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0mbins\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mi\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ma\u001b[0m \u001b[0;32min\u001b[0m \u001b[0menumerate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mlinspace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mpi\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m8\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m     \u001b[0mp\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mabs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcos\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mIx\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0msin\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mIy\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      4\u001b[0m     \u001b[0mp\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mbetween\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0.6\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m1.0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0.05\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m     \u001b[0mf\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmaximum_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m10\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'Ix' is not defined"
     ]
    }
   ],
   "source": [
    "bins = []\n",
    "for i,a in enumerate(linspace(0,pi,8)[:-1]):\n",
    "    p = abs(cos(a)*Ix+sin(a)*Iy)\n",
    "    p = between(p,0.6,1.0,0.05)\n",
    "    f = filters.maximum_filter(p,10)\n",
    "    f /= amax(f)\n",
    "    bins.append(f)\n",
    "subplot(131); imshow(bins[0])\n",
    "subplot(132); imshow(bins[2])\n",
    "subplot(133); imshow(minimum(bins[0],bins[3]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "axes don't match array",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-7-2b887e7c41b7>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mohist\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtranspose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbins\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0mplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mohist\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m5\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m80\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/numpy/core/fromnumeric.pyc\u001b[0m in \u001b[0;36mtranspose\u001b[0;34m(a, axes)\u001b[0m\n\u001b[1;32m    531\u001b[0m         \u001b[0mtranspose\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtranspose\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    532\u001b[0m     \u001b[0;32mexcept\u001b[0m \u001b[0mAttributeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 533\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0m_wrapit\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'transpose'\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0maxes\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    534\u001b[0m     \u001b[0;32mreturn\u001b[0m \u001b[0mtranspose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0maxes\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    535\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/usr/lib/python2.7/dist-packages/numpy/core/fromnumeric.pyc\u001b[0m in \u001b[0;36m_wrapit\u001b[0;34m(obj, method, *args, **kwds)\u001b[0m\n\u001b[1;32m     41\u001b[0m     \u001b[0;32mexcept\u001b[0m \u001b[0mAttributeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     42\u001b[0m         \u001b[0mwrap\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mNone\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 43\u001b[0;31m     \u001b[0mresult\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mgetattr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0masarray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mobj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmethod\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     44\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mwrap\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     45\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mresult\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmu\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mndarray\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mValueError\u001b[0m: axes don't match array"
     ]
    }
   ],
   "source": [
    "ohist = transpose(bins,[1,2,0])\n",
    "plot(ohist[5,80])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Morphological Corner Detection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f3a454b02d0>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArUAAAFMCAYAAADV1UQNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFTBJREFUeJzt3V9sVvX9B/DPYWUXiGNiKhDajAVltPxrFUayyHzYnFtm\nMBK7BJjDgHjjzealV5Qbp9NlgWzZhdlMvAGTXQw1QCZZHmAa05jCtsRETISsoC5zCSogKS3nd8Ev\nVWS0T/+e8+15vZImrT05eSM87757ep6nWZ7neQAAQMJmFB0AAADGy6gFACB5Ri0AAMkzagEASJ5R\nCwBA8prGe4JarRZHjhyZiCwAU+6ee+6Jer1edIwplWVZ0REAxuxGL9w17iu1R44ciTzPJ+xt586d\nE3o+ueSSS67h3nxTDjA9uP0AAIDkGbUAACSvdKO2VqsVHeF/kmt05BoduUanrLkAKE6W5/m4fk1u\nlmUxzlMAFKaKHeaJYkDKbtTZpbtSCwAAo2XUAgCQPKMWAIDkGbUAACTPqAUAIHlGLQAAyTNqAQBI\nnlELAEDyjFoAAJI34qg9dOhQLF26NO6444545plnpiITAOOgt4EqGvbX5A4ODsa3vvWtOHz4cCxc\nuDDWrFkTe/fujba2ts9PUMFfMQlMH9OtwxrtbYBUjenX5Pb09MTtt98eixYtipkzZ8amTZti//79\nkxIQgPHT20BVDTtqz549G62trUMft7S0xNmzZyc9FABjo7eBqmoa7pON/oiqu7t76P1arRa1Wm08\nmQAmTb1ej3q9XnSMSePWAqCqhh21CxcujL6+vqGP+/r6oqWl5brjvjhqAcrsy99479q1q7gwk6DR\n3gaYboa9/WD16tXx7rvvxunTp6O/vz9eeumleOCBB6YqGwCjpLeBqhr2Sm1TU1P89re/jR/+8Icx\nODgYjz766DXPoAWgXPQ2UFXDvqRXQyeYZi+HA1RLFTvMfbdAysb0kl4AAJACoxYAgOQZtQAAJM+o\nBQAgeUYtAADJM2oBAEieUQsAQPKMWgAAkmfUAgCQPKMWAIDkGbUAACTPqAUAIHlGLQAAyTNqAQBI\nnlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5RCwBA8oxaAACSZ9QC\nAJA8oxYAgOQZtQAAJM+oBQAgeUYtAADJG3HUbt++PebNmxcrVqyYijwAjIPOBqpqxFG7bdu2OHTo\n0FRkAWCcdDZQVSOO2nXr1sUtt9wyFVkAGCedDVSVe2oBAEieUQsAQPKaJuIk3d3dQ+/XarWo1WoT\ncVqACVev16NerxcdA4AJluV5no900OnTp2PDhg3xz3/+8/oTZFk0cAqAUpqOHTZcZ0dc/TMDpOpG\nnT3i7QebN2+O73znO3Hy5MlobW2NF154YcLDATAxdDZQVQ1dqR32BNPwKgdQHVXsMFdqgZSN+Uot\nAACUnVELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5RCwBA8oxaAACS\nZ9QCAJA8oxYAgOQZtQAAJK+p6ABVk2VZ0REAAKYdV2oBAEieUQsAQPKMWgAAkmfUAgCQPKMWAIDk\nGbUAACTPqAUAIHlGLQAAyTNqAQBInlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkb8RR29fXF+vX\nr49ly5bF8uXLY8+ePVORC4Ax0NlAVWV5nufDHfDhhx/Ghx9+GB0dHXH+/Pm466674s9//nO0tbVd\nPUGWxQin4AuyLCs6AvAl06nDRursCD0EpO1GnT3ildr58+dHR0dHRETMnj072tra4v3335/YdABM\nCJ0NVNWo7qk9ffp0HD9+PNauXTtZeQCYIDobqJKGR+358+ejq6srdu/eHbNnz57MTACMk84Gqqap\nkYMuX74cDz30UDz88MPx4IMPXvf57u7uofdrtVrUarWJygfAKI3U2QDT0YhPFMvzPB555JG49dZb\n4ze/+c31J/BEsVHxBA0on+nUYSN1doQeAtJ2o84ecdT+7W9/i+9+97uxcuXKoSL85S9/GT/60Y+u\nnsCoHRVfTKB8plOHjdTZEXoISNuYR+1IjNrR8cUEyqdqHaaHgJSN+SW9AACg7IxaAACSZ9QCAJA8\noxYAgOQZtQAAJM+oBQAgeUYtAADJM2oBAEieUQsAQPKMWgAAkmfUAgCQPKMWAIDkGbUAACTPqAUA\nIHlGLQAAyTNqAQBInlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5R\nCwBA8oxaAACSZ9QCAJA8oxYAgOQZtQAAJG/EUXvp0qVYu3ZtdHR0RHt7ezz55JNTkQuAMdDZQFVl\neZ7nIx108eLFmDVrVgwMDMTdd98dzz33XNx9991XT5Bl0cAp+H9ZlhUdAfiS6dZhw3V2hB4C0naj\nzm7o9oNZs2ZFRER/f38MDg7G3LlzJy4ZABNKZwNV1NCovXLlSnR0dMS8efNi/fr10d7ePtm5ABgj\nnQ1UUUOjdsaMGXHixIk4c+ZMHD16NOr1+iTHAmCsdDZQRU2jOXjOnDlx//33x1tvvRW1Wm3ov3d3\ndw+9X6vVrvkcAMW4UWcDTEcjjtqPPvoompqa4utf/3p89tln8dprr8XOnTuvOeaLoxaA4jTS2QDT\n0Yij9oMPPohHHnkkrly5EleuXImf/exn8f3vf38qsgEwSjobqKqGXtJr2BN4Sa9R8VI6UD5V6zA9\nBKRsXC/pBQAAZWbUAgCQPKMWAIDkGbUAACTPqAUAIHlGLQAAyTNqAQBInlELAEDyjFoAAJJn1AIA\nkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5RCwBA8oxaAACSZ9QCAJA8oxYAgOQZtQAAJM+o\nBQAgeUYtAADJM2oBAEieUQsAQPKMWgAAkmfUAgCQPKMWAIDkGbUAACTPqAUAIHkNjdrBwcHo7OyM\nDRs2THYeAMZJZwNV1NCo3b17d7S3t0eWZZOdB4Bx0tlAFY04as+cORMHDhyIHTt2RJ7nU5EJgDHS\n2UBVjThqn3jiiXj22Wdjxgy33wKUnc4GqmrY1nv11Vfjtttui87OTt/xA5SczgaqrGm4T77xxhvx\n8ssvx4EDB+LSpUvxySefxNatW+PFF1+85rju7u6h92u1WtRqtcnICsAwGu1sgOkoyxv8dv7IkSPx\n3HPPxSuvvHLtCbLMFYFR8MQNKJ/p2GE36uwIPQSk7UadPaqbrhQhQDp0NlAlDV+pveEJXKkdFV9k\noHyq1mF6CEjZhFypBQCAMjJqAQBInlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0A\nAMkzagEASJ5RCwBA8oxaAACSZ9QCAJA8oxYAgOQ1FR0AAGAq5HledATGKcuyG37OlVoAAJJn1AIA\nkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5RCwBA8oxaAACSZ9QCAJA8oxYAgOQZtQAAJM+o\nBQAgeUYtAADJa2rkoEWLFsXXvva1+MpXvhIzZ86Mnp6eyc4FwBjpbKCKGhq1WZZFvV6PuXPnTnYe\nAMZJZwNV1PDtB3meT2YOACaQzgaqpqFRm2VZ3HvvvbF69ep4/vnnJzsTAOOgs4Eqauj2g9dffz0W\nLFgQ//nPf+IHP/hBLF26NNatWzfZ2QAYA50NVFFDo3bBggUREdHc3BwbN26Mnp6eawqyu7t76P1a\nrRa1Wm1CQwLQuJE6GyAV9Xo96vV6Q8dm+Qg3Xl28eDEGBwfj5ptvjgsXLsR9990XO3fujPvuu+/q\nCbLMvVujkGVZ0RGAL5lOHTZSZ0foIaprOj3Wq2q43Tnildp///vfsXHjxoiIGBgYiJ/+9KfXlCMA\n5aGzgaoa8UrtiCdwpXZUXCGB8qlah+khqqpqj/XpaLjd6TeKAQCQPKMWAIDkGbUAACTPqAUAIHlG\nLQAAyTNqAQBInlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5RCwBA\n8oxaAACSZ9QCAJA8oxYAgOQZtQAAJM+oBQAgeUYtAADJM2oBAEieUQsAQPKMWgAAkmfUAgCQPKMW\nAIDkGbUAACTPqAUAIHlGLQAAyRtx1J47dy66urqira0t2tvb480335yKXACMkd4GqqhppAN+/vOf\nx49//OP405/+FAMDA3HhwoWpyAXAGOltoIqyPM/zG33y448/js7OznjvvfdufIIsi2FOwZdkWVZ0\nBOBLplOHNdrbUEXT6bFeVcPtzmFvPzh16lQ0NzfHtm3b4s4774zHHnssLl68OCkhARg/vQ1U1bCj\ndmBgIHp7e+Pxxx+P3t7euOmmm+Lpp5+eqmwAjJLeBqpq2HtqW1paoqWlJdasWRMREV1dXf+zHLu7\nu4fer9VqUavVJjQkAI1ptLcBUlCv16Nerzd07LCjdv78+dHa2honT56MJUuWxOHDh2PZsmXXHffF\nUQtAcRrtbYAUfPli6a5du2547LBPFIuI+Pvf/x47duyI/v7+WLx4cbzwwgsxZ86cz0/giWKj4gka\nUD7TrcMa6W2ooun2WK+i4XbniKN2PCfner6YQPlUrcP0EFVVtcf6dDTmVz8AAIAUGLUAACTPqAUA\nIHlGLQAAyTNqAQBInlELAEDyjFoAAJJn1AIAkDyjFgCA5Bm1AAAkz6gFACB5Ri0AAMkzagEASJ5R\nCwBA8oxaAACSZ9QCAJA8oxYAgOQZtQAAJM+oBQAgeUYtAADJM2oBAEieUQsAQPKMWgAAkmfUAgCQ\nPKMWAIDkGbUAACTPqAUAIHlGLQAAyTNqAQBInlELAEDyRhy177zzTnR2dg69zZkzJ/bs2TMV2QAY\nJZ0NVFWW53ne6MFXrlyJhQsXRk9PT7S2tl49QZbFKE5ReVmWFR0B+JLp2mH/q7Mj9BDVNV0f61Uy\n3O4c1e0Hhw8fjsWLF19TjgCUk84GqmRUo3bfvn2xZcuWycoCwATS2UCVNHz7QX9/fyxcuDDefvvt\naG5u/vwEWRY7d+4c+rhWq0WtVpvwoNOFH/tB+UzHH0neqLMj9BDVNR0f69NdvV6Per0+9PGuXbtu\n+PfY8Kjdv39//P73v49Dhw5dewL31I6KLyZQPtOxw27U2RF6iOqajo/1qpmQe2r37t0bmzdvnrBQ\nAEwenQ1UTUNXai9cuBDf+MY34tSpU3HzzTdfewJXakfFFRIon+nWYcN1doQeorqm22O9iobbnaN6\nSa/Rnpzr+WIC5VO1DtNDVFXVHuvT0YS9pBcAAJSRUQsAQPKMWgAAkmfUAgCQPKMWAIDkGbUAACTP\nqAUAIHlGLQAAyTNqAQBInlELAEDySjdq6/V60RH+p7LmAgCmVlk3QdVzGbUNKmsuAGBqlXUTVD1X\n6UYtAACMllELAEDysjzP8/GcoFarxZEjRyYqD8CUuueee0r7I7vJoreBVA3X2eMetQAAUDS3HwAA\nkDyjFgCA5Bm1AAAkrzSj9tChQ7F06dK444474plnnik6TkREbN++PebNmxcrVqwoOso1+vr6Yv36\n9bFs2bJYvnx57Nmzp+hIERFx6dKlWLt2bXR0dER7e3s8+eSTRUe6xuDgYHR2dsaGDRuKjjJk0aJF\nsXLlyujs7Ixvf/vbRccZcu7cuejq6oq2trZob2+PN998s+hI8c4770RnZ+fQ25w5c0rzb7+KytjZ\nEXp7tMrc22Xs7Ai9PRpT3tt5CQwMDOSLFy/OT506lff39+erVq3K33777aJj5UePHs17e3vz5cuX\nFx3lGh988EF+/PjxPM/z/NNPP82XLFlSiv9feZ7nFy5cyPM8zy9fvpyvXbs2P3bsWMGJPvfrX/86\n37JlS75hw4aiowxZtGhR/t///rfoGNfZunVr/oc//CHP86t/l+fOnSs40bUGBwfz+fPn5//617+K\njlJJZe3sPNfbY1HW3i5jZ+e53h6rqejtUlyp7enpidtvvz0WLVoUM2fOjE2bNsX+/fuLjhXr1q2L\nW265pegY15k/f350dHRERMTs2bOjra0t3n///YJTXTVr1qyIiOjv74/BwcGYO3duwYmuOnPmTBw4\ncCB27NgRecle8KNseT7++OM4duxYbN++PSIimpqaYs6cOQWnutbhw4dj8eLF0draWnSUSiprZ0fo\n7bEoY2+XubMj9PZYTEVvl2LUnj179po/ZEtLS5w9e7bAROk4ffp0HD9+PNauXVt0lIiIuHLlSnR0\ndMS8efNi/fr10d7eXnSkiIh44okn4tlnn40ZM0rxT35IlmVx7733xurVq+P5558vOk5ERJw6dSqa\nm5tj27Ztceedd8Zjjz0WFy9eLDrWNfbt2xdbtmwpOkZl6ezx0dsjK2tnR+jtsZqK3i7Fv5Ysy4qO\nkKTz589HV1dX7N69O2bPnl10nIiImDFjRpw4cSLOnDkTR48eLcWL2r/66qtx2223RWdnZ+m+u379\n9dfj+PHjcfDgwfjd734Xx44dKzpSDAwMRG9vbzz++OPR29sbN910Uzz99NNFxxrS398fr7zySvzk\nJz8pOkpl6eyx09sjK3NnR+jtsZiq3i7FqF24cGH09fUNfdzX1xctLS0FJiq/y5cvx0MPPRQPP/xw\nPPjgg0XHuc6cOXPi/vvvj7feeqvoKPHGG2/Eyy+/HN/85jdj8+bN8de//jW2bt1adKyIiFiwYEFE\nRDQ3N8fGjRujp6en4ERXr7q1tLTEmjVrIiKiq6srent7C071uYMHD8Zdd90Vzc3NRUepLJ09Nnq7\nMWXu7Ai9PRZT1dulGLWrV6+Od999N06fPh39/f3x0ksvxQMPPFB0rNLK8zweffTRaG9vj1/84hdF\nxxny0Ucfxblz5yIi4rPPPovXXnstOjs7C04V8dRTT0VfX1+cOnUq9u3bF9/73vfixRdfLDpWXLx4\nMT799NOIiLhw4UL85S9/KcUztufPnx+tra1x8uTJiLh6H9SyZcsKTvW5vXv3xubNm4uOUWk6e/T0\nduPK2tkRenuspqy3J+0paKN04MCBfMmSJfnixYvzp556qug4eZ7n+aZNm/IFCxbkX/3qV/OWlpb8\nj3/8Y9GR8jzP82PHjuVZluWrVq3KOzo68o6OjvzgwYNFx8r/8Y9/5J2dnfmqVavyFStW5L/61a+K\njnSder1emmfSvvfee/mqVavyVatW5cuWLSvNv/s8z/MTJ07kq1evzleuXJlv3LixNM+iPX/+fH7r\nrbfmn3zySdFRKq+MnZ3nenu0yt7bZersPNfbYzGVvZ3leQlvWAEAgFEoxe0HAAAwHkYtAADJM2oB\nAEieUQsAQPKMWgAAkmfUAgCQPKMWAIDk/R+TzRKydTzkjAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f3a4555b1d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "outside = zeros((8,8),'B')\n",
    "inside = zeros((8,8),'B')\n",
    "outside[:2,:] = 1; outside[:,:2] = 1\n",
    "inside[-3:,-3:] = 1\n",
    "subplot(121); imshow(outside)\n",
    "subplot(122); imshow(inside)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'image' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m\n\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-9-d5bc24bcfc07>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0msmooth\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfilters\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmedian_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mimage\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m5\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0mimshow\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msmooth\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'image' is not defined"
     ]
    }
   ],
   "source": [
    "smooth = filters.median_filter(image,5)\n",
    "imshow(smooth)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "oimage = morphology.grey_opening(smooth,footprint=outside)\n",
    "iimage = morphology.grey_closing(smooth,footprint=inside)\n",
    "subplot(121); imshow(oimage)\n",
    "subplot(122); imshow(iimage)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "imshow(maximum(0,oimage-image)**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {},
 "nbformat": 4,
 "nbformat_minor": 0
}