control_refine.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:07.804451Z",
     "start_time": "2020-08-06T18:56:07.802729Z"
    }
   },
   "outputs": [],
   "source": [
    "# default_exp control_refine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This contains code for refining the following control points:\n",
    "\n",
    "* Checker\n",
    "* Ellipse\n",
    "\n",
    "The term \"refine\" is used here because an initial guess for the location is required."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Import"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:08.802501Z",
     "start_time": "2020-08-06T18:56:07.805888Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "import math\n",
    "\n",
    "import numpy as np\n",
    "import skimage.draw\n",
    "import torch\n",
    "from torch.distributions.multivariate_normal import MultivariateNormal\n",
    "\n",
    "from camera_calib.utils import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:08.806618Z",
     "start_time": "2020-08-06T18:56:08.803807Z"
    }
   },
   "outputs": [],
   "source": [
    "from pathlib import Path\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.core.debugger import set_trace\n",
    "\n",
    "from camera_calib.image import File16bitImg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Base class"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I tried to make a generic control flow loop for refining control points which you can inherit and modify/implement"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:08.816143Z",
     "start_time": "2020-08-06T18:56:08.807954Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "class CPRefiner:\n",
    "    def __init__(self, cutoff_it, cutoff_norm):\n",
    "        self.cutoff_it   = cutoff_it\n",
    "        self.cutoff_norm = cutoff_norm\n",
    "\n",
    "    def proc_arr(self, arr):            return (arr,)\n",
    "    def it_preproc(self, p, b):         pass\n",
    "    def get_bb(self, p, b):             raise NotImplementedError('Please implement get_bb')\n",
    "    def get_W(self, p, b, bb):          return None    \n",
    "    def refine_point(self, arrs, p, W): raise NotImplementedError('Please implement refine_point')\n",
    "\n",
    "    def __call__(self, arr, ps, bs):\n",
    "        arrs = self.proc_arr(arr)\n",
    "        bb_arr = array_bb(arr)\n",
    "        ps_refined = []\n",
    "        for idx, (p, b) in enumerate(zip(ps, bs)):\n",
    "            self.it_preproc(p, b)\n",
    "            b_init = b\n",
    "            for it in torch.arange(self.cutoff_it):\n",
    "                p_prev = p\n",
    "                bb = self.get_bb(p, b)\n",
    "                if not is_bb_in_bb(bb, bb_arr): p = arr.new_full(2, math.nan); break\n",
    "                W = self.get_W(p, b, bb)\n",
    "                p = self.refine_point(tuple(bb_array(arr, bb) for arr in arrs), p-bb[0], W)+bb[0] \n",
    "                if torch.any(torch.isnan(p)): break\n",
    "                if not is_p_in_b(p, b_init):    p = arr.new_full(2, math.nan); break\n",
    "                if torch.norm(p-p_prev) < self.cutoff_norm: break\n",
    "                b = b-p_prev+p\n",
    "            ps_refined.append(p)\n",
    "        return stackify(tuple(ps_refined))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Checker"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:08.824881Z",
     "start_time": "2020-08-06T18:56:08.817116Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "class CheckerRefiner(CPRefiner):\n",
    "    def __init__(self, hw_min, hw_max, cutoff_it, cutoff_norm):\n",
    "        super().__init__(cutoff_it, cutoff_norm)\n",
    "        assert_allclose(type(hw_min), int)\n",
    "        assert_allclose(type(hw_max), int)\n",
    "        self.hw_min, self.hw_max = hw_min, hw_max\n",
    "\n",
    "    def it_preproc(self, p, b):\n",
    "        d_min = min([pld(p, l) for l in b_ls(b)])\n",
    "        hw = torch.floor(d_min/math.sqrt(2))\n",
    "        hw = hw.clamp(self.hw_min, self.hw_max)\n",
    "        self.hw = hw\n",
    "        self.bb = stackify(((-hw, -hw),\n",
    "                            ( hw,  hw)))\n",
    "\n",
    "    def get_bb(self, p, b): return self.bb+p.round()\n",
    "\n",
    "    def get_W(self, p, b, bb):\n",
    "        zero = p.new_tensor(0)\n",
    "\n",
    "        sigma = self.hw/2\n",
    "        cov = stackify(((sigma**2,     zero),\n",
    "                        (    zero, sigma**2)))\n",
    "        d = MultivariateNormal(p, covariance_matrix=cov)\n",
    "        return torch.exp(d.log_prob(grid2ps(*bb_grid(bb)))).reshape(tuple(bb_sz(bb).long()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Opencv method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:08.833785Z",
     "start_time": "2020-08-06T18:56:08.826226Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "@numpyify\n",
    "def checker_opencv(arr_dx, arr_dy, W=None):    \n",
    "    # Condition array points\n",
    "    ps_cond, T = condition(array_ps(arr_dx))\n",
    "        \n",
    "    # Form linear system \n",
    "    A = grid2ps(arr_dx, arr_dy)\n",
    "    b = (A*ps_cond).sum(dim=1)\n",
    "    \n",
    "    # Get weighted least squares estimate\n",
    "    p = wlstsq(A, b, W)\n",
    "    \n",
    "    # Convert back to unconditioned points\n",
    "    return pmm(p, torch.inverse(T), aug=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.146356Z",
     "start_time": "2020-08-06T18:56:08.834962Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAACECAYAAACJbXCEAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAF75JREFUeJztnc+vXMWVx7/n+QdgP9vPxsZYBkOE0EgILyKhIcOAFBGQmGzYjcJIIxYgNoyUSFmEzPwD2ZDdCAkpiCyijEZKJFhEihgr0WjEKOItohl+iGB+GrANtvFvwDy7ZvHaPed+/e6pW923b/ft9/1Ilrv63q57bp3b9bq+deqUpZQghBCifyxM2wAhhBCjoQ5cCCF6ijpwIYToKerAhRCip6gDF0KInqIOXAghespYHbiZPWJmb5vZYTN7pi2jxHSRX+cX+Xa+sFHjwM1sA4C/AHgYwMcAXgPwWErpzfbME10jv84v8u38sXGMz/41gMMppfcAwMz+DcCjAGofBjNLZrbmMf5DUnfeuHC9vrxhw4bKsRtuuKFS3rZt2/D19ddfXzm2cWO1KX1dXG9k05UrV8Jz666x1mfr2vDIkSM4efJkXQMX+3Xr1q1paWmpkc1sU/QDouQZ4HN927Aft27dWin7drt48WLl2DfffFMpe3svX75cOcblkh9H3l7+HNdb1y6nT5/GhQsXokYr8u2WLVtq/Vrix0mRu2bJ8+OfgZJ7WVioihjRNXP1RsePHz9+IqW0h98fpwPfD+CIK38M4N7oA2aGTZs2rXmMH1JuGF/mG+Wy/zJwg3K9mzdvHr7esWNH5djBgwcr5QcffHD4+o477qgc2717d6Xs67rxxhsb2/vll19WjkVfFLb3q6++qpT9HxX/uYceeggBxX5dWlrC008/PSxHf4T4D50/lz/Hf6C877hduF7fNnfffXfl2L33Vm/Ht/ny8nLl2Kefflop+3b84osvKsdOnz5dKftnmu+Nn8Pt27ev+TkAOHPmTKXM93qV5557bs33HUW+XVpawpNPPjkse5vZfn6mo2dgnM7fn5v7seOvw/byZ7/++uvh60uXLoXn+vbfsmVL7TGg6svcffofC/wMPPvssx+u9ZmJT2Ka2VNmtmxmy1q2Pz94v164cGHa5oiW8H7l0YiYPcbpwD8BcKsr3zJ4r0JK6fmU0j0ppXsmJYuIVin2K0sSYmbJ+tb7lX9ditljHAnlNQB3mtm3sPoQ/ADAP5RU4Dv0OmnlKn4Yk5NF6q6x1rl+VJCTZrwNPFTi63g549y5c5VjKysrlbK/dy/prHUdP8zywz6uBwDOnz+/pn2ZoWexX82sdi6Bh4J8794fubkC74/I50AszbBffbuxH7mNvf18L2x/dD88TOe6InvrrtHgB1Kxb+vqH2cuI/LdODpx5I+chOLL7BvG1xU9Hzn7ou9Crj+8ysgdeEppxcz+CcDvAWwA8EJK6Y1R6xOzgfw6v8i388c4v8CRUvodgN+1ZIuYEeTX+UW+nS/G6sBHIYok8ETD/NzwOboGf9ZLFHxuFBbG0gYPh3zoGk8GXXfddbX1RpEkbBOHMn7++eeVstemm7b7NCkZlpfIXbmhaxTyyc+A9w/XEz0TUeQRH89FQbDfr9K2X80sDG/08Pcqko9K5JaoH8hFvvh2jMJMgar93L6Rr1hCiUJ5+VniENWmsolHS+mFEKKnqAMXQoieog5cCCF6SqcauJnVriJjIg0tp5VG4U6RxpbTVb2W7UP0gGuXa/vVfayp8WpLf5z1z0hbPHnyZGjvsWPHhq937tw5fF2yXL8JKaXKtaP0BKwDcj1RueRcf498TdYtvY38fHL7RyFmUajaOHpnNI9TGEZYBH9fozDO6DuZ06qjkM8S+Drez7n5Cv+MNE1dsNa50fMdzb/xZzn0uLbORmcJIYSYOdSBCyFET1EHLoQQPaVzDdwvFfd6F2ulrDt57ZG1o5LltZGmnItB9Trm2bNnK8c4M6DXuTnZU6SXM5F2yu0QZfrzmn3bGnhUZy5rndf9xtFw+Tre76xbR8m3cilhfb25tQDRXAzjr8PpFLher+36+ZO2k8Vt3rwZt91227Ds24ZtinTt3LlRrDSf6+viuYxoXoR9zustou9EZH/uO+ht4D6O14N4+N7q0C9wIYToKerAhRCip3QqoWzevBkHDhxY81hu+bCHhzSRzFASmsZDVx6W+9CeKHk/UB0ucbpVXvLu4WFVFELGbeQ3BQCqw0bfJtxebVAXtpdbwhztcBOFnPI9RKFfvCHCkSNHas89ceJE5RgPtX07liyNzslzvszPAJf9d8Vfs20JZXFxEffdd9+a1+V2idqf7zWSSfi7zf2ClxvZr++//36l/NFHHw1fv/fee5VjHIIbEWUYzEkonlw6Bf/ZpjKnfoELIURPUQcuhBA9RR24EEL0lE418G3btuGBBx4Ylr2eFYVKAVU9kY9FIWQ5rTTSt3gJ/PHjx4evOUVstFky18M2eY0zpxlHx1jD92Wvv7W912FKqXJP/nVuhyFPbmNrTy7lgNcQP/jgg8oxn2IAqPqHtUcO8fRtzs8O2+/9yvZFu92zDdGy+0mEhF6FwwgjvT1aWh+FjnK9JRo4bzjNZT+fwfMevCG171Mi3zC5dNVRO0Q7QzVNtaBf4EII0VPUgQshRE9RBy6EED2lUw18cXER999//7Ds46pz+mcUW8zL2qOYZ9Z/fb283Ja10nfffXf4OorjBWKtl/G635YtWyrHIg2c9eXouLcnWko+CimlSjtGW89F5La080RzF0BV0/TxwMC1sd67du0avuZ1ChyD7Z+t3G733v4ozQFQ9XNO1/Z1TSKm/yqXLl3Chx9+OCz77yTfT0n65+i7nZujKkkR69Nb7Nmzp3KMvzvRVnlROoVcOtkoXXWUJiB69j36BS6EED1FHbgQQvSUTiWUTZs24ZZbbhmWfYhQlKGMj+fkCi8R8FCEw8J8mXfBYGnGD3FYXomGUrlMi77eKHSOz80t/fd4G9oOIzSz2iyJLEHwkD8KOcyFlEXH/JD4k08+qRw7fPhwpbx3797ha5awbrrppko52pWG8feWW04e7f6ey6g5Kc6fP49XX311WPZZNLmdcpkZPeznKISS28J/XzmdBV9z//79w9f8HHI/EPkq6os4pDl6Rhl+Bvxn+V6Wl5fXrEO/wIUQoqeoAxdCiJ6iDlwIIXpK5zvyeK00WkK+uLhYKfsdZVhHZU3N63OsSbG2dOrUqdpzI12Pz2Wb/P1w6s1o2Xcu1MjblNs9pm6n+KYhSk0xs4ovox1L2EY/z5Bbdu/bPApFA6o6P2ul3udA9TkcJ01qybxNBD+jXJ7kTvSey5cvV/zj25TbMAoXze1+4z/L9fD9+e8Zh8Py3I7X7G+++ebQhiiEL/Iz28fPd/S9K/keSgMXQog5I9uBm9kLZvaZmb3u3ttlZq+Y2TuD/3dO1kzRNvLr/CLfrh+a/AJ/EcAj9N4zAA6llO4EcGhQFv3iRciv88qLkG/XBVkNPKX0n2Z2O739KIDvDl7/EsAfAfwke7GNG7G0tDQsez082roKqMbj+rSuwLW7vHvNircz41jvnTv//4cIa2pRvC1rX5FWx0QpSnPLyUdNH8pxxm361cwqbRXpxNGSa24z1kPrrgFcqz1G239x2cfy5nTKyAYu+2c4t/Q/mgcpSVE6eK8V36aUKm0T7ZQexbWzH/lc9ocn0vhzS+lLtmX0bVyyDWO0RSBfh8/NrQ9pwqga+N6U0tHB62MA9kYni94gv84v8u0cMvYkZlr9k1P7J8rMnjKzZTNb5plrMbuU+LXt5FhiskS+9X71kV9iNhk1jPC4me1LKR01s30APqs7MaX0PIDnAeDgwYPJDxv8UIrDBhk/FOHlzTwc8kMRHvax3OKP8xCMd9Lxx3OhXtHy+EgmKVmeHckRXJe3JxgSjuTX/fv3p7r7zWVt9CGfUYoBIL53Lvs2zS1j9/WWSFY5GScaekftkpNm6q6R2ZW+kW+9Xw8cOJD8taNQuyjUMYdv/0imGdg3fJ2TEyNZqqSNo7BfhuuJMoZm/NWIUX+Bvwzg8cHrxwG8NLYlYhaQX+cX+XYOaRJG+GsA/w3gr8zsYzN7AsDPADxsZu8AeGhQFj1Cfp1f5Nv1Q5MolMdqDn2vZVtEh8iv84t8u37odCn9xo0bK7ufeN2MNTTWmbx2ndOzotSpnFayJNTIn8u6WMkS2ojcuSW7fXsbSzTJcYnalG308wo5rTqaK4jmFXKhgSW6aqSVjjO3ET1bJTa0SUqpNn0B2xiF8JWkL+bw4Ugnjnbr4c9GfUQpkX+iYyWafWNbij8hhBBiJlAHLoQQPUUduBBC9JRONfArV65UdCqvR7MmxVqYX1TAsZVRTK1fug9cu1x++/btw9e80CiKFS3Rs/jeOBVtiaYWxTdHS3MbxoGPhJnVatlsPy+Xj86NYmpzcya+LtZGo3mG3JZeo5J7XqLrtJWmtpQrV65U0j745yu3DNy3cW6bsagfiFIQRClsmZLd43NzDCXL+yOitQxNnzv9AhdCiJ6iDlwIIXpKpxLKhg0bKtkB/TCXh1k8vPDDrNxyZ38uD9mj7GG50Cgv6+TCkKLQtKZLo4Fr7Y9Co5ou8Z10GJr3D7cht1uUyoDtjHbOido4SrUAxL6KJJXcUnp/rzlppsQ/JSGH41InRXGbsuQZyS1MtNNVSToCpiTkc9RshDkbSjJoKoxQCCHWEerAhRCip6gDF0KIntKpBg5UdR+vfeXSvpYsC49C7Vh3irT1SGvP6drjLKmN6m16LLpO22GEKaXaMKzcjiUlmnJkdy69QnSsICVrkX0laUdHrWeSXL58GWfOnBmWvc7NqZbZxpLl856SuSU+N0pvMc4cRPS85GyIrjnNdLJCCCGmjDpwIYToKerAhRCip3S+lL4uPeWePXsq57L+7DW13NZbUWw3fzZaos+xrb6ci8GOltKXLNGPdijPpZ6N9OW28Xb7eQWO7fZLs4H4HiJ9MUoxwMc5hXC0Hdg4scUlRNt0TUvzZhYWFq6Zi7pKydZzfG6UIiH3vfLPRBR7DlTXULSZyqAkpXPJXMwoKW71C1wIIXqKOnAhhOgpnUooPCTzcsrZs2cr5/LQzQ+zeNjCy2/9cn0vkayFD4fiYS0PvaMdsyOpoyTD3ThDvZLds9skpVQZOvpl7lFKBKA65OS0AX7HeoaHqlFmNx5qM9Fu6yWSSUl4YiQPRZkV+fgkfbywsFDxl8/kee7cucq5LD9u27Zt+LpEZoikMD6XpdOo3pJd6HPhfaOGnZbY0NSv+gUuhBA9RR24EEL0FHXgQgjRUzrVwM2sdqkppwdlvMbGu+owXhMvWe6c2+XFa+s5PSsKDWT902t3uTCqknSyJSGH4+Lr93MFOQ3c2xzNMfC5OZ04Wu7MoY1R6gWmRNduWg9Q9WUu9DWqp0047Nd/B3mnK56H8vMZrI9Hz3/u++qvw37k60SpXMfZeSnS1qN6c9eInoE69AtcCCF6ijpwIYToKerAhRCip3SeTtbjNZ8dO3ZUjrHO7bUljhFnHcrH/XJsMRNtT8VxyF7Dyu2IHcX1RrGuuTS1kaYW6WbR8t+28VpktCUcEMd6s64aLaWP9PPcVm1Rm0blEk0zlyI20msnnQYhom4ugds0ur+SLdVyvop2j8/NYTWtdxytOpcSOWIiW6qZ2a1m9gcze9PM3jCzHw7e32Vmr5jZO4P/dxZfXUwN+XU+kV/XF00klBUAP04p3QXgOwCeNrO7ADwD4FBK6U4AhwZl0R/k1/lEfl1HZCWUlNJRAEcHr8+Z2VsA9gN4FMB3B6f9EsAfAfwkV180VPTwMS875LIG+mW+LLdwNrxo92welkchhzyUijKhjbNzSzSEZDLSQKt+rcsAmUt7sLi4OHwd+RGoyi18LktukXxRsnsPUyKhlOzOEkkOJZny2vQrh/36Z55tiiQV/h6Nkw4iCjnM1TUqkdyVe5aic6Pv/UTCCM3sdgDfBvAnAHsHDwsAHAOwt6QuMTvIr/OJ/Dr/NO7AzWwRwG8A/CilVMk8lVb/tKz559DMnjKzZTNbPnHixFjGivZpw6+5hVWie+TX9UGjDtzMNmH1YfhVSum3g7ePm9m+wfF9AD5b67MppedTSveklO7ZvXt3GzaLlmjLr36Fqpg+8uv6IauB26qo8wsAb6WUfu4OvQzgcQA/G/z/UoO6KvpelEaSNU6vH+V2A/EpYiMtncu55dneplwYlb9Pvmak3bH+GWntJcu+eVl9m35NKdVqw9yG3G5eE889A75toh2QcjZEGnPu2Yp8N4627smlVI2Ote1Xfv6ukttlyrcbz3vwZ/33NZd6NqonCu2NdmFi+F5yIcNNKQlPbEoTy/4WwD8C+F8z+/PgvX/G6oPw72b2BIAPAfx98dXFNJFf5xP5dR3RJArlvwDU/Wn4XrvmiK6QX+cT+XV9oaX0QgjRUzpdSs9bb3mdjJdUR3G9rF/xcm0f+80xqCUpJ702lyOnyY5az6gxyrl626au/ih1KJe5jshXufspuV9/bi7t7qj+yMU3R2mCmbbim3OwBu5t5O9ctMM6z2XwuV4jz82ZsH0etilaHh+lK8ilMojSGo+TBiGai6lDv8CFEKKnqAMXQoie0nk2wrqdanj3+GjoHYUhAcDFixeHr3kHGMYP2Xio5+vhc0tCvXI78kQhTSVEw8RJZyD0RL7KhZ95eEgc+Yrbv0Ruidopsi8nZZQslR51p/lJyilmdo2MdZXcM+v9HKWZAKoSSi6Ez9ebS5FQksogklCYkp1zSvqFaBep2vobnSWEEGLmUAcuhBA9RR24EEL0lKnuyON1TNaKeJd6r8WxfsX6uYfTx7J2HYUnss7q9bic9uXryu1cXaKTeRtyeu2kQwc9dfeQWxodhYmxP6LUv5FWmtM0o3CzklSubTHNHXg8Zlb7zPPcEvvKh+/mnsNobin6ruTCE0vCQ0vmK0ZNExyFWvJnm85Z6Re4EEL0FHXgQgjRU9SBCyFET+lUAzez2h3jWVfiXMRnzpypPZdjVb3uzTGnUXw5b7/GGtuo2znlYsZLtkmLbIhiR0fR15rCWmmkVfOchL93bm9OQ+rLOQ3ct2Pkx7XKnkgTH2fZ9DhL6btiZWUFp06dGpb9egvWwEt04ijFczSXAVTbiXX3tuaAcs9H020hmUnMr+gXuBBC9BR14EII0VM6lVBWVlZw8uTJYdkPw3JZ6/y50U7VQDWEieUVHnb5YQwP2aMwvZJl1CW7jDDR8uBcvU13ih+XhYWFivzkJQtuUyZqpyisM9qxictr7FpTa09uGO7rzaVIKFnOPyolS8BLSSlVJEj/XYp2igLicFeWNes+txa+HaNUC0BzSZHJZY4cNeRw6rvSCyGEmB3UgQshRE9RBy6EED3FukwzamafY3VD1d0ATnR24TyzZg8wWZtuSyntaasy+bUI+XV8Zs0eYPI2renbTjvw4UXNllNK93R+4RpmzR5gNm3KMWs2z5o9wGzalGPWbJ41e4Dp2SQJRQgheoo6cCGE6CnT6sCfn9J165g1e4DZtCnHrNk8a/YAs2lTjlmzedbsAaZk01Q0cCGEEOMjCUUIIXpKpx24mT1iZm+b2WEze6bLazsbXjCzz8zsdffeLjN7xczeGfy/s0N7bjWzP5jZm2b2hpn9cNo2lSK/rmmP/NqODfJrQGcduJltAPCvAP4OwF0AHjOzu7q6vuNFAI/Qe88AOJRSuhPAoUG5K1YA/DildBeA7wB4etAu07SpMfJrLfJrO7wI+bWelFIn/wD8DYDfu/JPAfy0q+uTLbcDeN2V3wawb/B6H4C3p2HX4PovAXh4lmySX+VX+XU2/dqlhLIfwBFX/njw3iywN6V0dPD6GIC90zDCzG4H8G0Af5oVmxogv2aQX1tnJtpwFvyqSUwirf4J7Tw0x8wWAfwGwI9SSmdnwaZ5Qn6dT9a7X7vswD8BcKsr3zJ4bxY4bmb7AGDw/2ddXtzMNmH1YfhVSum3s2BTAfJrDfLrxJBfB3TZgb8G4E4z+5aZbQbwAwAvd3j9iJcBPD54/ThWda1OsNUs/L8A8FZK6eezYFMh8usayK8TRX69SseC//cB/AXAuwD+ZUqTDr8GcBTAN1jV9Z4AcCNWZ47fAfAfAHZ1aM/9WB1u/Q+APw/+fX+aNsmv8qv82g+/aiWmEEL0FE1iCiFET1EHLoQQPUUduBBC9BR14EII0VPUgQshRE9RBy6EED1FHbgQQvQUdeBCCNFT/g+RwE86nrVzdAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "file_img = Path('data/dot_vision_checker/SERIAL_16276942_DATETIME_2019-06-07-00:38:48-109736_CAM_3_FRAMEID_0_COUNTER_2.png')\n",
    "bb = torch.FloatTensor([[540, 965], [566, 991]])\n",
    "arr = bb_array(File16bitImg(file_img).array_gs(torch.float), bb)\n",
    "arr_dx, arr_dy = grad_array(arr)\n",
    "fig, axs = plt.subplots(1,3)\n",
    "for a, ax in zip([arr, arr_dx, arr_dy], axs): ax.imshow(a, cmap='gray')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.151946Z",
     "start_time": "2020-08-06T18:56:09.147392Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([11.9810,  6.2994])"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "p = checker_opencv(arr_dx, arr_dy)\n",
    "p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.264583Z",
     "start_time": "2020-08-06T18:56:09.152888Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fa8174610f0>]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPsAAAD4CAYAAAAq5pAIAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAEj9JREFUeJzt3V+IHed5BvDnWf1bafXHK68shGU3qfCNqalSFlGIXFzHDY5u5BAw0UVQwaBcxJBALmrSi/jSlCahFyWg1CJqSR0CibEulDaqCJjcBK+NYst2W/9BRhLalbSqrV1JK2m1by92VDbOnu8ZnW/OOVN9zw+WPTvfmZnvzDnvzu683/cOIwJmdvcbGnQHzKw/HOxmhXCwmxXCwW5WCAe7WSFW9nNnJINk1+urzEHOtttCvQbVvmLFimT72rVrk+0bNmxItg8PDyfbV65Mf6RU/1S7oo7PwsJC1vbrUK9B9SHnc3z69GlMT08vu4GsYCf5JIB/ALACwD9FxAvi+Vi1alXX+7t161ayfWgo/YeKale/TFR7nQ+qeiNVH1evXp1s37RpU7L9kUceSbY//vjjyfYdO3Yk28fGxpLtqn/33ntvsj33Pbh27VqyHdDvkeqDeo1zc3PJdvULM7X/J554omNb13/Gk1wB4B8BfAnAwwD2kXy42+2ZWW/l/M++C8D7EfFhRNwA8FMAe5vplpk1LSfY7wdwesnPZ6plv4fkAZITJCc8Ws9scHp+gS4iDgI4CABDQ0OOdrMByTmznwXwwJKft1fLzKyFcoL9NQAPkfwsydUAvgrgSDPdMrOmdf1nfETMk3wWwL9jMfV2KCLezumMSnnkpO2A/PymSovVkbuP3PSgalfHSKWF1OtTaaeZmZlk+/z8fLJdfUZU6hLQr/HmzZvJ9uvXryfbVR9nZ2eT7aljnHr/sv5nj4ijAI7mbMPM+sPDZc0K4WA3K4SD3awQDnazQjjYzQrhYDcrRF/nswPpPHLufPTcucq5efQ66+fm2XPz3GqasMrDq/2rPLiaT3/16tVk+5o1a5Ltqv8qzw/o16iOoZrzf+HChWT7yMhIsr3bGPKZ3awQDnazQjjYzQrhYDcrhIPdrBAOdrNCONjNCtHvuvEyh5mj1zXHm6hbn1sXXo0lyJ3PrvLcaq61yqOrUs4qR527/o0bN5LtgB7roD5n09PTyXb1HkxOTibbR0dHO7alPh8+s5sVwsFuVggHu1khHOxmhXCwmxXCwW5WCAe7WSH6nmdP1e1WOWSV31Q5apVjzblVbl1qvndujjc3T69qol++fDnZrm5XrPLkV65cSbbn5vHrUMdA1X1XnyPVrt6j1FgH59nNzMFuVgoHu1khHOxmhXCwmxXCwW5WCAe7WSH6mmdfvXo1Hnzwwa7XV3ny3LrxKgeem38F8uebK+r+4yqPr+6P/vHHHyfbVU11NU5A1UxXNdcVVXceyK/trz6HGzduTLarsQapz1nqM5oV7CRPAZgBcAvAfESM52zPzHqniTP7X0bExQa2Y2Y95P/ZzQqRG+wB4FckXyd5oIkOmVlv5P4ZvzsizpK8D8Axkv8ZEa8ufUL1S+AAUO8Clpn1RtaZPSLOVt/PA3gZwK5lnnMwIsYjYryXlWXNLK3rYCc5QnLD7ccAvgjgZFMdM7Nm5ZxqtwJ4ucpJrgTwrxHxb6kVNmzYgEcffbRju5qLrPLg169fT7arfyPU+iq/qvoH6Fy9ylPX2UeKqqs+NTWVbFd15dUxUuMAVP/U8VN59Dq1/XPvP6DWV8dAtafGYqTen66DPSI+BPCn3a5vZv3l1JtZIRzsZoXw5fHCvHj0KEaXuTbxP2vW4Jk9ewbQI+sXn9kLs1ygp5bb3cPBblYIB7tZIRzsZoXo6wW69evXY/fu3R3bVeEENeBEDehQ66sBF+oGCXXG/qtBIWrQinoNqvABfvnLjk3r1q3DBx98kFxdFWZQ7b0eMq0GHa1bt05uI3dQjRoUk7t+6him3n+f2QvzcYcRap2W293DqbfCPPuVrwy6CzYgPrObFcLBblYIB7tZIRzsZoVwsJsVoq9X41etWoXt27d3bM+9CYTKk/c6Ryxz3NB5clXAQ7WrsQpqrIDKU09OTibb1evLvYmEeo9V/9VYjDrUPnJv1KGkjlFqnIbP7GaFcLCbFcLBblYIB7tZIRzsZoVwsJsVwsFuVoi+5tlJJud85xb4X79+fbJ9dnY22a7mmqv8aZ250mosgcrhXrp0KWv7aqyA2r/avjqG6j2cm5tLtqsctRorkbrBQt1tqGOoxhKoPqj21DFM7dtndrNCONjNCuFgNyuEg92sEA52s0I42M0K4WA3K0Rf8+wrV67EPffc07Fd1V1XOVi1/n333Zdsn5qaSravXbs22a7yswAwMjKSbFfz0UdHR5Ptak597lxqlSdX7SoPr6j5/Lnz4QGdJ6/zPvdS6jWmcvTynSd5iOR5kieXLNtM8hjJ96rv6U+gmQ1cnV/zPwbw5KeWPQfgeEQ8BOB49bOZtZgM9oh4FcCnx2juBXC4enwYwFPNdsvMmtbtP3BbI+Jc9XgSwNZOTyR5gOQEyQk1rtvMeif7anwsXhHoeFUgIg5GxHhEjG/evDl3d2bWpW6DfYrkNgCovp9vrktm1gvdBvsRAPurx/sBvNJMd8ysV2SeneRLAB4DMEbyDIDvAngBwM9IPgPgIwBP191hKkeo8ptqvrqi5mqrPLxav04O9/r168l2lctX66u51sPi1sy5891z29U4AHWM1WeoibrxaqxA7v3dVR9TxzCVZ5fBHhH7OjR9Qa1rZu3h4bJmhXCwmxXCwW5WCAe7WSEc7GaFcLCbFaLv89lTQ2Zzc7RqnrHKYefU6wbq5VdT988GdO383LnUucdQ5YBzj6HKk+dqYvvqNahjkDsWoNs8vs/sZoVwsJsVwsFuVggHu1khHOxmhXCwmxXCwW5WiL7m2RcWFpLzsVWOWeUXVd14dX92df91lYOukwNP1c0HdN33jRs3JttVnb/cuu+5eX61vnqP1f3hm5ivrraRW1c+t7Z9t/PZfWY3K4SD3awQDnazQjjYzQrhYDcrhIPdrBAOdrNC9DXPvmLFiuT9yVV+U+VYVf5S5fFz86dq+0B+zXE1V1odQ1UXXo1VyK2JrvqvjnGd2vw5+wf0e5Q7Zz93rEJq+6k2n9nNCuFgNyuEg92sEA52s0I42M0K4WA3K4SD3awQfc2zA+kcpappnntv817V476tTk1ylWNVOdzcsQK5Ywly8+C5881zc9R15Obyc9fPeY1Z89lJHiJ5nuTJJcueJ3mW5Inqa0/XvTOzvqjza/bHAJ5cZvkPImJn9XW02W6ZWdNksEfEqwDStY7MrPVy/oF6luSb1Z/5o52eRPIAyQmSExcvXszYnZnl6DbYfwhgB4CdAM4B+F6nJ0bEwYgYj4jxsbGxLndnZrm6CvaImIqIWxGxAOBHAHY12y0za1pXwU5y25IfvwzgZKfnmlk7yDw7yZcAPAZgjOQZAN8F8BjJnQACwCkAX6+zs4WFheScdJWD3bJlS7Jd5ZDVvdFz6nUD9eZKq9eo+pBb+17NV1ftuXO1e11TvYm55Lk1A3LvAZ87X74TGewRsW+ZxS92tTczGxgPlzUrhIPdrBAOdrNCONjNCuFgNyuEg92sEH2dzz40NJScc67qwl++fDnZruazqxy0ym+q/qVq4t+m8uTK8PBwsl3lgNV8eFUzQMnNUefO91eamA+f24fcY9D1fnuyVTNrHQe7WSEc7GaFcLCbFcLBblYIB7tZIRzsZoXoa56dZFbd8Lm5uaz9qzz7lStXsrav8vCAztHm3r9czblX66uxArnzyXPnm6v57moueZ254Lm5/tzX2Kv58j6zmxXCwW5WCAe7WSEc7GaFcLCbFcLBblYIB7tZIfp+f/YUlZ/ctGlTsl3lyVV+U82HV/1TNdcB4ObNm/I5Kbk52nXr1iXbVQ555cq8j0zuXG6VY25irrh6n3PHSuTeoz7Vv6z7s5vZ3cHBblYIB7tZIRzsZoVwsJsVwsFuVggHu1kh+p5nT+UYu73vdN31VY42997mMzMzyXZA5/KvXbuWbFfzuVW7qgufO19e5ZDVOIPcueRN1FzPnY+e+znOzeN3Is/sJB8g+WuS75B8m+Q3q+WbSR4j+V71fbSrHphZX9T5M34ewLcj4mEAfw7gGyQfBvAcgOMR8RCA49XPZtZSMtgj4lxEvFE9ngHwLoD7AewFcLh62mEAT/Woj2bWgDu6QEfyMwA+B+C3ALZGxLmqaRLA1g7rHCA5QXLi4sWLOX01swy1g53kegA/B/CtiPi9OyzG4hWDZa8aRMTBiBiPiPGxsbGszppZ92oFO8lVWAz0n0TEL6rFUyS3Ve3bAJzvTRfNrAl1rsYTwIsA3o2I7y9pOgJgf/V4P4BXmu+emTWlTp798wC+BuAtkieqZd8B8AKAn5F8BsBHAJ7uSQ/NrBEy2CPiNwA6jRL4wp3sjGRy0EfuoBg16EUNRlDbV4MphoeHk+1A/mvM7aMa9KKOoSpeoY5x7k0euh1QUnf/QP7AodyBPTnFL1LHz8NlzQrhYDcrhIPdrBAOdrNCONjNCuFgNyuEg92sEH0tXhERyTygyk+qGxyom0SoHLfKr6rCC6owBaCLR6gCGirHq15jnbEAOdR72ERxiRxN3CQiV+5YgW6Poc/sZoVwsJsVwsFuVggHu1khHOxmhXCwmxXCwW5WiL7fJCKVw1T5w9nZ2WR77lxptb7KUV+9ejXZDgBr1qyRz0nJnY+u+qi2r46hWl9RnwH1HqmxEm2g8uwqz9+zm0SY2d3BwW5WCAe7WSEc7GaFcLCbFcLBblYIB7tZIfqeZ09ROWKVw52bm0u2q7niKn+p8vx1XLt2Ldmu8ti5c/Jz68Ln5tFV/9Try62LX4faR27t/NzX2C2f2c0K4WA3K4SD3awQDnazQjjYzQrhYDcrhIPdrBB9zbOTTOZ5VV12lUMdGRlJtn/yySdZ21d5epVDB/RrzJ2Tr2rXqzx7r++P3uv58ipH3URN+NxjVOce8d1uP9Umz+wkHyD5a5LvkHyb5Der5c+TPEvyRPW1p5uOm1l/1DmzzwP4dkS8QXIDgNdJHqvafhARf9+77plZU2SwR8Q5AOeqxzMk3wVwf687ZmbNuqMLdCQ/A+BzAH5bLXqW5JskD5Ec7bDOAZITJCcuXLiQ11sz61rtYCe5HsDPAXwrIi4D+CGAHQB2YvHM/73l1ouIgxExHhHjW7Zsye+xmXWlVrCTXIXFQP9JRPwCACJiKiJuRcQCgB8B2NW7bppZrjpX4wngRQDvRsT3lyzftuRpXwZwsvnumVlT6lyN/zyArwF4i+SJatl3AOwjuRNAADgF4OtqQ/Pz85ienu7Yrmqq5+ag1fZVjlZtX917HdC5ejXfW+Vwb9y4kbW+ykP3ui686p96D5qg+qCOQe5ryNl/at06V+N/A2C5T8BRta6ZtYeHy5oVwsFuVggHu1khHOxmhXCwmxXCwW5WCObOT76jnZEXAHy0ZNEYgIt968Cda3v/gPb30f3Lc6f9+6OIWHZcel+D/Q92Tk5ExPjAOiC0vX9A+/vo/uVpsn/+M96sEA52s0IMOtgPDnj/Stv7B7S/j+5fnsb6N9D/2c2sfwZ9ZjezPnGwmxViIMFO8kmS/0XyfZLPDaIPCslTJN+qymRPtKA/h0ieJ3lyybLNJI+RfK/6vmwdwAH3sRUlxxMl0VtzDHtdtr3v/7OTXAHgvwH8FYAzAF4DsC8i3ulrRwSSpwCMR0QrBlyQ/AsAswD+OSL+pFr2dwAuRcQL1S/N0Yj4m5b18XkAs4MuOV5VVtq2tCQ6gKcA/DVacgwTfXwaDRzDQZzZdwF4PyI+jIgbAH4KYO8A+vH/SkS8CuDSpxbvBXC4enwYix+MgenQx1aIiHMR8Ub1eAbA7ZLorTmGiT42YhDBfj+A00t+PoN21qEPAL8i+TrJA4PuTAdbq7r+ADAJYOsgO5MgS47306dKorfyGHZTtl3xBbrOdkfEnwH4EoBvVH+itlYs/j/WxjxqrZLj/bJMSfT/05Zj2G3ZdmUQwX4WwANLft5eLWuViDhbfT8P4GW0s1T21O0qv9X38wPuzx9oU8nx5Uqio2XHsJdl2wcR7K8BeIjkZ0muBvBVAEcG0I+OSI5UF0hAcgTAF9HOUtlHAOyvHu8H8MoA+7KstpQc71QSHS06hj0v2x4Rff8CsAeLV+Q/APC3g+iD6N8fA/hd9fV2G/oI4CUs/gl3E4vXOZ4BcC+A4wDeA/AfADa3sI//AuAtAG9iMbC2Dahvu7H4J/qbAE5UX3vadAwTfWzkGHq4rFkhfIHOrBAOdrNCONjNCuFgNyuEg92sEA52s0I42M0K8b+SmWD6Hm3QcAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(arr, cmap='gray')\n",
    "plt.plot(p[0], p[1], 'rs')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.268308Z",
     "start_time": "2020-08-06T18:56:09.265573Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "class OpenCVCheckerRefiner(CheckerRefiner):\n",
    "    def __init__(self, hw_min, hw_max, cutoff_it, cutoff_norm):\n",
    "        super().__init__(hw_min, hw_max, cutoff_it, cutoff_norm)\n",
    "\n",
    "    def proc_arr(self, arr): return grad_array(arr)\n",
    "\n",
    "    def refine_point(self, arrs, p, W): return checker_opencv(*arrs, W)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.279600Z",
     "start_time": "2020-08-06T18:56:09.269183Z"
    }
   },
   "outputs": [],
   "source": [
    "refiner = OpenCVCheckerRefiner(hw_min=3, hw_max=10, cutoff_it=20, cutoff_norm=1e-3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.475625Z",
     "start_time": "2020-08-06T18:56:09.280621Z"
    },
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "file_img = Path('data/dot_vision_checker/SERIAL_16276942_DATETIME_2019-06-07-00:38:48-109736_CAM_3_FRAMEID_0_COUNTER_2.png')\n",
    "img = File16bitImg(file_img)\n",
    "ps_c_p = torch.FloatTensor([[858.7428864 , 633.38643518]])\n",
    "bs_c_p = torch.FloatTensor([[[878.5844935448612, 614.1667623245968],\n",
    "                             [838.8941064355354, 613.1104586040717],\n",
    "                             [838.9778791064288, 652.5319092059557],\n",
    "                             [878.5655638183706, 653.6357474223086]]])\n",
    "ps_c_p_refined = refiner(img.array_gs(torch.float), ps_c_p, bs_c_p)\n",
    "assert_allclose(ps_c_p_refined, torch.FloatTensor([[860.234, 635.174]]), atol=1e-3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.705409Z",
     "start_time": "2020-08-06T18:56:09.477052Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(620.0, 650.0)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQsAAAD8CAYAAABgtYFHAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAGLFJREFUeJzt3X+s3XV9x/Hnq73t7Q+UbqIBaSdlVogQivWOyZIxHQsiYeAGy1jYCCVLxUB1m4mDEI3TkP2hm6mykdQyFNfpSFXGHBDY5hb+sNUCReWXq4D0IkiBFVKgv+x7f5zvldPLPZ/v536+53vOqXs9koZzvj/f53sO7/s95/35oYjAzKzOnGEHYGaHBycLM8viZGFmWZwszCyLk4WZZXGyMLMsWclC0hJJmyQ9LOkhSadL+oSkJyVtq/6d07X91ZK2S3pE0nvbC9/MBkU57SwkfQm4OyI2SJoPLAL+DNgdEZ+Ztu3bga8ApwFvBv4deFtE/KzPsZvZANXeWUg6EjgDuAEgIvZFxK7ELucDX42IvRHxGLCdTuIws8PYWMY2y4GdwI2SVgL3AB+u1l0p6RJgK/CRiPhf4Fhgc9f+k9WyQ0haA6ypnr5zzhz/fNIvkvq+bv78+clzjo+P91w3b968nuvmzp1bFA9A6jNTt2+p1HGbtIZu67jdnnjiiWcj4o2l++ckizFgFbA2IrZIWgdcBVwHfAqI6r9/A1yWe+KIWA+sB5g7d24sXLiwe13P/Uo/BHX7lb5ZbcRat2/dcVP/A6bWLViwoOe6Y499Tb4/xPHHH99z3dKlS3uue93rXtdzXSoBASxatKjnutTrTL1nqf0gnaD279+f3DcllYxTx53NH9nLL7/8x7MKavq5MraZBCYjYkv1fBOwKiJ+GhE/i4iDwBd49avGk8Cyrv2XVsvM7DBWmywi4mlgh6QTqkVnAg9KOqZrs98DflA9vhW4SNK4pOXACuA7fYzZzIYg52sIwFpgY1UJeRRYDXxO0ql0voY8DnwAICIekHQz8CBwALgipxLSj++YqWMcPHiw9fP385hNvt6Ufm1q6xY7de1T6+pe54EDB4riGRvr/bGve52p319+9rPeH/PUfpB+LamvKKnXsnfv3uQ5ZysrWUTENmBi2uI/SWx/LXBteVhmNmpcgjCzLE4WZpbFycLMsjhZmFkWJwszy5JbOm1ddwktVcIrLbU1aU7e1qDGg2jiO13q+qVKf6kSHdS3fOwldQ1S8UC6xWlpPHWvMxVvqsVpk3J3at9UqbfutcyW7yzMLIuThZllcbIwsyxOFmaWxcnCzLI4WZhZlpEpnXaXh0p7RjYpN7ZxzrZGTmpLk3Omypypcm3q2taVP/ft21e0b1uvM3XOJr2e2+hdW8J3FmaWxcnCzLI4WZhZFicLM8viZGFmWZwszCzLYVc6bWsOjzbmI2lr3pAm+5aWgZvEU9oDtK70nDrunj17eq6rmzAppbQc2aR0WjoQdb8H7PWdhZllcbIwsyxOFmaWxcnCzLI4WZhZFicLM8viZGFmWbKKxpKWABuAk+lMhHxZRHy7WvcR4DPAGyPiWUnvBv4FeKza/esR8cm6c/RjNOvSUcHr9m2ry/wwuqGXdhevey2pfVPrUl2+60b3TsU0jAmMm4wgX6p05O8SuS1M1gF3RMSF1UzqiwAkLQPOAp6Ytv3dEXFu/8I0s2GrTYWSjgTOAG4AiIh9EbGrWv1Z4KN07jbM7BdYzn3TcmAncKOk+yRtkLRY0vnAkxFx/wz7nC7pfkm3SzqprxGb2VDkfA0ZA1YBayNii6R1wCfo3G2cNcP29wJviYjdks4BbgFWTN9I0hpgTfW4KHgzG5ycO4tJYDIitlTPN9FJHsuB+yU9DiwF7pV0dES8GBG7ASLiNmCepKOmHzQi1kfERERMOFmYjb7aZBERTwM7JJ1QLToTuDci3hQRx0XEcXQSyqqIeFrS0ar+75d0WnWO59oJ38wGJbcashbYWFVCHgVWJ7a9EPigpAPAK8BFMcu6aGmpMqVJF/Vh3Pk06S6e2jfVxbq0rFq3b+mI2HUfm9RxU/GWdvmGdLyp0cYXLlyYPG5KW0MHzFZWsoiIbcBEYv1xXY+vA65rGpiZjRa34DSzLE4WZpbFycLMsjhZmFkWJwszyzIyo3t3l7raKBUdbg2/2urBmCo3NpmcNxVvat9Uz8i6kbTreqX2UtezNKV0kuI6pb1vU9eoHz25u/nOwsyyOFmYWRYnCzPL4mRhZlmcLMwsi5OFmWVxsjCzLCPTzqK7JtxGW4q6mnMb3eKbaPJaUvuWdglv0kU9NWt5aRsCSL+WVFuKum7oKaWfhboZzVPvWenM7aX79eI7CzPL4mRhZlmcLMwsi5OFmWVxsjCzLE4WZpZlZEqnuSWpVNmrSZfcNrqwNzlmv7sX5xy3ycTSKaVlwbrSX+mI4k3K0qUjitddv1RMe/bs6blufHy857rSLvy9+M7CzLI4WZhZFicLM8viZGFmWZwszCyLk4WZZRmZ0ml3aSlVRkqVp9oqnZaes0nv0CavpTSmJhMGp9bX9bjsJdWrFMpHFE8dt+6cqcmPS3vXQvln7MUXX+y5rslkzDPJurOQtETSJkkPS3pI0uld6z4iKSQdVT2XpM9J2i7pe5JW9TViMxuK3DuLdcAdEXGhpPnAIgBJy4CzgCe6tn0fsKL69+vA9dV/zewwVntnIelI4AzgBoCI2BcRu6rVnwU+CnTfJ50P3BQdm4Elko7pa9RmNnA5X0OWAzuBGyXdJ2mDpMWSzgeejIj7p21/LLCj6/lktewQktZI2ippa1tNm82sf3K+howBq4C1EbFF0jrgE3TuNs4qPXFErAfWA8ydO9fZwmzE5dxZTAKTEbGler6JTvJYDtwv6XFgKXCvpKOBJ4FlXfsvrZaZ2WGs9s4iIp6WtEPSCRHxCHAmcG9EnDm1TZUwJiLiWUm3AldK+iqdHzZfiIin6s7T9qC4TY7fxkTNdfuWrqtTWm6sO2eqNJi6fqlSZZNJiFPnbNJTM7U+VSKuO25pr97UwMSpSadL5FZD1gIbq0rIo8DqxLa3AecA24GXa7Y1s8NEVrKIiG3ARGL9cV2PA7iiaWBmNlrc3NvMsjhZmFkWJwszy+JkYWZZnCzMLMtIdFGXdEi9PVVzLm1/UNekvHRk5ibtIdrqop56LaWTJjcZKTrV/iDV5qGu7UGqHUFpl++6c7bVHih13tLPZr/bWfjOwsyyOFmYWRYnCzPL4mRhZlmcLMwsi5OFmWUZidIpHFrSS5Xw2hppu3TE5ybd10vLo3UjUKfKnKnJhpt0fU+Nep1al5r0ty6e0velyXtWOgJ66hpA+n1pYxTzEr6zMLMsThZmlsXJwsyyOFmYWRYnCzPL4mRhZllGonQ6Z86cQyZxLe3Z12TS5DZG0x7GRM1NjpsaKbrunC+//HLPdakSXmq/utG9S9/v0tIypEugdfumpHrmln4WmpTCZzxXX49mZr+wnCzMLIuThZllcbIwsyxOFmaWxcnCzLKMROl03rx5HHPMMVnblvYObXvi5ZnU9forjalu8NxUCTS1b26vyXvecw/7x187GOzYnjFW3rXyNct37dpVdM7SyYKhvcGHU0rLtVA+QXRqvybXbyZZn1ZJSyRtkvSwpIcknS7pU5K+J2mbpDslvbna9t2SXqiWb5P08b5GbEM3U6IAOLCgfNZzG325dxbrgDsi4sJqJvVFwAMR8TEASR8CPg5cXm1/d0Sc2/dozWxoapOFpCOBM4BLASJiHzC9GdtioLy5opmNvJyvIcuBncCNku6TtEHSYgBJ10raAVxM585iyumS7pd0u6STZjqopDWStkraWtes18yGLydZjAGrgOsj4h3AS8BVABFxTUQsAzYCV1bb3wu8JSJWAp8HbpnpoBGxPiImImKiSZt6MxuMnGQxCUxGxJbq+SY6yaPbRuACgIh4MSJ2V49vA+ZJOqpP8doImLd35mrL2B4n/V9kte9uRDwtaYekEyLiEeBM4EFJKyLif6rNzgceBpB0NPDTiAhJp9FJSM+1FL8NwTu/9c7XLOv3vJo2enL/FKwFNlaVkEeB1cAGSScAB4Ef82ol5ELgg5IOAK8AF0VNX+0FCxZw4oknZgWSakOQOk1dm4c2JqZNxQrp2nuTSZNL21K88sorPdc988wzyXM+//zzRcdtMklxSun1azKUQUrdcduYfLvf7UmykkVEbAMmpi2+oMe21wHXNQvLzEaNm3ubWRYnCzPL4mRhZllc6zIA1t18M0tmmEpw1/g4a847bwgR2ajxnYUBzJgoAJYkBpK1/19G4s5iehf10lGb+z0R7JRUWWv+/Pk919VNhpvSZLTx1GTDJY444gieffbZ5Dapa5QaubrJCOilpdXUOeuGDUiVI9saXb6NsmoJ31mYWRYnCzPL4mRhZlmcLAyAFxctmnH5CwsWDDgSG1Uj8QOnDd8nPvCB5FiZZr6zMLMsI3FnMTY2xhve8IafPy8tiZWWn6B8pO3UORf1uLWfkhohLFVqq+sOvnjx4p7rSsuqbQ1Q1KTXaZOyay91PTVL420j1rpzunRqZkPhZGFmWZwszCyLk4WZZXGyMLMsThZmlmUkSqeSDinNlQ6QmypFpnqHNtm3dNJaSA/omzpuXUkstW+qBJp6nU0mlk6VI5sMnptSWuJs63XWaTJAc7/368V3FmaWxcnCzLI4WZhZFicLM8viZGFmWZwszCyLk4WZZclqZyFpCbABOBkI4DLgHDqzpx8EngEujYifqFMwXletf7lafm/dObrr26XtD8bHx3uuq6s5l3bBTu3XpM6dqvenXiek2xGk2n60UeuvO25K6r2GZqNp91LXVqLf3b6ntNGFfVjtLNYBd0TEicBK4CHg0xFxSkScCnwT+Hi17fuAFdW/NcD1fY3YzIai9s+ppCOBM4BLASJiHzB9QozFdO44oHO3cVN00tpmSUskHRMRT/UtajMbuJw7i+XATuBGSfdJ2iBpMYCkayXtAC7m1TuLY4EdXftPVssOIWmNpK2Str700kuNXoSZtS8nWYwBq4DrI+IdwEvAVQARcU1ELAM2AlfO5sQRsT4iJiJiIjUMnJmNhpxkMQlMRsSW6vkmOsmj20bggurxk8CyrnVLq2VmdhirTRYR8TSwQ9IJ1aIzgQclreja7Hzg4erxrcAl6ngX8IJ/rzA7/OXWC9cCGyXNBx4FVgMbqgRyEPgxcHm17W10yqbb6ZROV9cdfM6cOSxcuPDnz0u7LafKgnWl0dTkvalSZGpd3TlTZbpUia6uG3Xqa11p1/e6kmFp1+7Sa1Cn9DPUVhf1utdSWrZuq5Q7k6xkERHbgIlpiy+YYVOqKsgVzcIys1HjFpxmlsXJwsyyOFmYWRYnCzPL4mRhZllGYnRvOLQ8VFrCKx0VHNIjW+/bN70rzKtSPWTrSqepkliqxFn3WlKTH5eO/F1Xoivt6dpkFPPS0bSblBtLJ+1ua6TyQfKdhZllcbIwsyxOFmaWxcnCzLI4WZhZFicLM8viZGFmWUamnUV3G4nSbuhNuounztndfX66VNf2Om3NLp66RkcccUTPdc8//3zPdcPoLl6nje7Zde0oPLq3mVkNJwszy+JkYWZZnCzMLIuThZllcbIwsywjUTqNCPbv3//z56nu4qWjINeVkVJdzVMzpqVKsnUjRadKp01GmU6V/1Ld11OvpW6S4tIR0JsMK1CqyWjZpZ+xJqN7l8bTb76zMLMsThZmlsXJwsyyOFmYWRYnCzPL4mRhZlmySqeSlgAbgJOBAC4Dfh/4XWAf8CNgdUTsknQc8BDwSLX75oi4fPoxu82ZM+eQyXxLy1OpUmRdr9NUaTBVyk2ds0npr7REDOlyZOq4pSN0NznnqGmrFFl33DYmRh5Wr9N1wB0RcSKwkk4yuAs4OSJOAX4IXN21/Y8i4tTqXzJRmNnhoTZZSDoSOAO4ASAi9kXEroi4MyKm/hxvBpa2F6aZDVvOncVyYCdwo6T7JG2QtHjaNpcBt3fvU23735J+s1/Bmtnw5CSLMWAVcH1EvAN4CbhqaqWka4ADwMZq0VPAr1Tb/gXwT5JeP/2gktZI2ipp6+7duxu+DDNrW06ymAQmI2JL9XwTneSBpEuBc4GLo/o1JSL2RsRz1eN76Pz4+bbpB42I9RExERETqaHezGw01CaLiHga2CHphGrRmcCDks4GPgqcFxEvT20v6Y2S5laPjwdWAI/2PXIzG6jcXqdrgY2S5tP5H3818F1gHLirKt9MlUjPAD4paT9wELg8InqPBMtre52mynCpdaUTDUO612mqpFhayq2Tirf7Ws2kycDFbUhdo9T7WXf92igbNilxNlFaHm2rZ+5Msj45EbENmJi2+K09tv0a8LVmYZnZqHELTjPL4mRhZlmcLMwsi5OFmWVxsjCzLE4WZpZlJEb3hkPrxQsWLOi5Xar2XteNOiVV5963b1/Pdan2GXXa6l6cqr2/8sorRcdsawLj1PvZ1kjbTUbhLm3X0GS09tJr3+82Ib6zMLMsThZmlsXJwsyyOFmYWRYnCzPL4mRhZllGonQ6Z84cFi5c+PPnqVJlqot1aYmubt/x8fGe65qUFEtLeHUlsVT39tQ5U6XnutJfk9JgL3XXdhijXpe+L8Po+u7SqZkNhZOFmWVxsjCzLE4WZpbFycLMsjhZmFmWkSidTh/dO9XrNFUWbFK+S5UNSyc/bjLycuq11I3unYq3dHT0OqnXWlpurHs/S0ddb6KtUcNL923r/ZzxXH09mpn9wnKyMLMsThZmlsXJwsyyOFmYWRYnCzPLMhKlUzi0zJMq/ZUO8FpXRiotczbpHZoq15b2vAWYP39+z3V79uzpua5JqbeNkmJdT+FBlg2bqvsslH52S9eVyLqzkLRE0iZJD0t6SNLpkj5dPf+epG9IWtK1/dWStkt6RNJ7+xqxmQ1F7teQdcAdEXEisBJ4CLgLODkiTgF+CFwNIOntwEXAScDZwN9LKh+j38xGQm2ykHQkcAZwA0BE7IuIXRFxZ0RMNafcDCytHp8PfDUi9kbEY8B24LT+h25mg5Tzm8VyYCdwo6SVwD3AhyPipa5tLgP+uXp8LJ3kMWWyWnYISWuANdXTvWvWrPnBLGNv21HAs8MOoovjSRu1eGD0Yjqhyc45yWIMWAWsjYgtktYBVwEfA5B0DXAA2DibE0fEemB9dYytETExm/3bNmoxOZ60UYsHRi8mSVub7J/zm8UkMBkRW6rnm+gkDyRdCpwLXByv/vT6JLCsa/+l1TIzO4zVJouIeBrYIWnqFuZM4EFJZwMfBc6LiJe7drkVuEjSuKTlwArgO32O28wGLLedxVpgo6T5wKPAauC7wDhwV1Uj3hwRl0fEA5JuBh6k8/XkiohIF8yrryMjZtRicjxpoxYPjF5MjeLRqDVeMbPR5ObeZpbFycLMsgwsWUj6c0kPSPqBpK9IWtC17nOSdnc9v1TSTknbqn9/Ooh4JH1R0mNd5z212lZVjNur5u2rhhzPuyW90LX84wOKR5KulfTDqtn/h6pth3V9esXT+vVJxHR313l/IumWatthXaNe8cz+GkVE6//oNMp6DFhYPb8ZuLR6PAF8Gdjdtf2lwHWDjgf4InDhDNufA9wOCHgXsGXI8bwb+OYQrs9q4CZgTrX8TUO+Pr3iafX61H2mu7b5GnDJMK9RIp5ZX6NBfg0ZAxZKGgMWAT9Rp8/Ip+mUYAftNfEktj0fuCk6NgNLJB0zxHgGYaZ4Pgh8MiIOAkTEM9W2w7o+veIZlJ7vmaTXA78N3FItGupnaIZ4Zm0gySIingQ+AzwBPAW8EBF3AlcCt0bEUzPsdkF1u7ZJ0rIZ1rcRD8C11Xk/K2lqRuRjgR1dh5ixCfsA4wE4XdL9km6XdFK/YqmJ51eBP5S0tTrvimqXYV2fXvFAi9enJqYp7wf+IyJerJ4P8zM0Uzwwy2s0kGQh6ZfoZNblwJuBxZIuAf4A+PwMu/wrcFx0erTeBXxpAPH8MZ2esycCvwb8MvCX/TxvH+O5F3hLRKykc/1uGVA848Ce6DRh/gLwD/08bx/jafX61MQ05Y+Ar/T7vH2MZ9bXaFBfQ34HeCwidkbEfuDrwF8BbwW2S3ocWCRpO0BEPBcRe6t9NwDvHEA8vxERT1W3iXuBG3m1t2zbTdhnFU9EvBgRu6vHtwHzJB3Vdjx0/hp+vdrmG8Ap1eOhXJ9e8Qzg+qRiojrXacC/dW0/rGs0Yzwl12hQyeIJ4F2SFkkSnSbjfxsRR0fEcRFxHPByRLwVYNp3ufPojJ/RdjwPTZ23WvZ+YKon7K3AJdUv2u+ic4s301engcQj6ehqGZJOo/M+Ptd2PHT++ryn2ua36IxjAkO6Pr3iGcD1ScUEcCGdHw+7hyUb1jWaMZ6ia1T66+ts/9G5k3iYzgf+y8D4tPXd1ZC/Bh4A7ge+BZw4iHiA/wS+Xy37R+CIalsBfwf8qFo/MeR4ruy6Ppvp3IUMIp4ldP46fR/4NrByyNenVzytX5/UZxr4L+DsadsO5Rol4pn1NXJzbzPL4hacZpbFycLMsjhZmFkWJwszy+JkYWZZnCzMLIuThZll+T+Hy90HtnTsQwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(img.array_gs(torch.float), cmap='gray')\n",
    "plt.plot(ps_c_p[:,0], ps_c_p[:,1], 'rs')\n",
    "plt.plot(ps_c_p_refined[:,0], ps_c_p_refined[:,1], 'gs')\n",
    "plt.gca().set_xlim(845, 875)\n",
    "plt.gca().set_ylim(620, 650)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ellipse"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-07-31T01:42:07.667986Z",
     "start_time": "2020-07-31T01:42:07.663616Z"
    }
   },
   "source": [
    "Note that torch doesnt have a `polygon2mask` function, so just convert to numpy and then back to tensor for now"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.710261Z",
     "start_time": "2020-08-06T18:56:09.706433Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "class EllipseRefiner(CPRefiner):\n",
    "    def __init__(self, cutoff_it, cutoff_norm):\n",
    "        super().__init__(cutoff_it, cutoff_norm)\n",
    "\n",
    "    def it_preproc(self, p, b):\n",
    "        bb = ps_bb(b)\n",
    "        bb = stackify((bb[0].floor(), bb[1].ceil()))\n",
    "        W = p.new_tensor(skimage.draw.polygon2mask(*torch2np((tuple(bb_sz(bb).long()), b-bb[0]))))\n",
    "        self.bb, self.W = bb-p.round(), W\n",
    "\n",
    "    def get_bb(self, p, b):    return self.bb+p.round()\n",
    "        \n",
    "    def get_W(self, p, b, bb): return self.W"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-05-17T21:08:31.844283Z",
     "start_time": "2020-05-17T21:08:31.835913Z"
    }
   },
   "source": [
    "### dual conic method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.719024Z",
     "start_time": "2020-08-06T18:56:09.711295Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "@numpyify\n",
    "def fit_conic(arr_dx, arr_dy, W=None):   \n",
    "    one = arr_dx.new_tensor(1)\n",
    "    \n",
    "    # Condition array points\n",
    "    ps_cond, T = condition(array_ps(arr_dx))\n",
    "        \n",
    "    # Form homogeneous coordinates of lines\n",
    "    ls = grid2ps(arr_dx, arr_dy)\n",
    "    ls = torch.cat([ls, -(ls*ps_cond).sum(dim=1, keepdim=True)], dim=1)\n",
    "    \n",
    "    # Form linear system \n",
    "    A = stackify((ls[:,0]**2, ls[:,0]*ls[:,1], ls[:,1]**2, ls[:,0]*ls[:,2], ls[:,1]*ls[:,2]), dim=1)\n",
    "    b = -ls[:, 2]**2\n",
    "        \n",
    "    # Get weighted least squares estimate\n",
    "    aq_inv = wlstsq(A, b, W)\n",
    "    \n",
    "    # Get conic matrix\n",
    "    Aq_inv = stackify(((  aq_inv[0], aq_inv[1]/2, aq_inv[3]/2),\n",
    "                       (aq_inv[1]/2,   aq_inv[2], aq_inv[4]/2),\n",
    "                       (aq_inv[3]/2, aq_inv[4]/2,         one)))\n",
    "    Aq = torch.inverse(Aq_inv)\n",
    "\n",
    "    # Rescale conic matrix to take conditioning into account\n",
    "    return T.T@Aq@T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:09.727590Z",
     "start_time": "2020-08-06T18:56:09.719964Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "@numpyify\n",
    "def ellipse_dualconic(arr_dx, arr_dy, W=None):\n",
    "    Aq = fit_conic(arr_dx, arr_dy, W)\n",
    "    return conic2ellipse(Aq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.050124Z",
     "start_time": "2020-08-06T18:56:09.729717Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXMAAACECAYAAABiWsuHAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAVP1JREFUeJztvV2IZWt63/df9bmrqrvOmTkzSIfpIZKJbnQjCULi4CCEg8FRTHwTjO3EzIVhbiZgE4dEShDkJpDcOGNQCAw44IDBTrBBQhiGZOJc5EaRjU2CLWQrIZYlja05mj7T3fVdtVcuuv+rfutfz7v2ru69q85p7Qc2e+/18X49X//ned/1rq7ve21oQxva0IY+37T12A3Y0IY2tKENvTttjPmGNrShDb0HtDHmG9rQhjb0HtDGmG9oQxva0HtAG2O+oQ1taEPvAW2M+YY2tKENvQf0Tsa867o/3nXdb3Rd95td1/3cqhq1ocelDV/fX9rw9v2l7m3XmXddty3pn0j6Y5J+W9KvSfozfd//49U1b0MPTRu+vr+04e37TTvvcO+/Luk3+77/fyWp67q/IelPSmoKxkcffdQ/e/ZMdCBbW6+Dg77v1fe9tra21HXdcOxN2cP1eY7X8BzPd1031GOaz+ejdvh+1pV1s76KeG22rzqfZVf1VP3iOffr5uZmKLPVZp/7rd/6LX3yySetyu/N18PDw/74+HjEu+wv25Lt4THem+fm87kkaXt7eyjb19zc3Aznfe/29vaoHVU9FT+n2u3rs69ZlttVlZ3tsSxub28P/PQx1re1tXVn/Pj9/PlznZyctIXqnrw9Ojrqv/CFL5S8qcanolZbkyq5bVGlW5U9mKqvurdq0yJKe5TnFtmPymZU/3/3d3/3k77vvzzVlncx5l+R9M/x/7cl/RtTN3z1q1/Vt7/9bV1dXWk+n2tra0v7+/uSNBybzWba29vTfD7X9fW1pFvltUDP53Pd3NyMO7Kzo52dHd3c3Oj6+nq4Zj6fa3d3V3t7e5I03H95eamrqytJt0Zxb29vZPS3trYGg0CDYqPB6+gMrKD8bG1taWdnZ6SkqdgkX7O9vT2M0cXFha6vr4d+7ezs6PDwUJJ0cnKiy8tLbW9va2dnp1nWzs6Ofvqnf3qKTffm6/Hxsb72ta9pf39fOzs7A59SmG10zVeP99XV1XBMes2Pvb097ezs6PLycujz5eWlJOno6Eiz2Wzg29XVlT799FNdXl4O/JvNZnry5MnQFkm6vr4e6vGYeCzJv2y7r3Eb+r7XbDbT7u7uHZ6Tr2lg+r7X7u6udnZ21Pf90K9Xr17p6upKx8fHevr0qa6vr3V+fq75fD7oxcHBgQ4PDzWfz3VxcaG+74extn784i/+4hSb7s3bL3zhC/rGN74x9I0OLEFYOlz/9n2+fmdnZzRebP/u7u5I31qO3+WaL/7tc7x+e3t75Px9TV7fAhP5O4ntZd+3tra0u7s79LuSEdaVY0hHL0m/8Au/8M9afDK9izFfirqu+7qkr0vSs2fPmp7KBqjrurLD0q1wsOMWaOkWnVp4XObOzo52d3eHMnx9GmoKqgfYBp8K7vsTRbFPbgsRI7+pjBYuC3TXdYNDIkON3G5ubnRzczMIqftsY21jUSnYqrZvIF+Pj4/v8CqVvTpHfr4pc2Qw6LDtbP2bwn59fT0oD49fXV0NjoEKk0bCcmCycff5KkpqRXCsJx28v+1QXLdlU9LgLKroKw0GHVL24W2JfP3ggw+GdlvWKv2lrm5vb5cRL3nedd3AW/+3XPuY9cfEyIS6laDObeJ3gi+WzciHbc72LxizURkEcdkut3eZSOW+9C7G/HckfRX/n705NqK+778l6VuS9JM/+ZN9IlF67qmw29emoNhoG+VQcCyAVhga152dnUHwqgGm0TRZYdKYuww6FRpiOpDKmLhvRBOu12WwX2yjfzNy2d7eHglTouMFBv3efP3444/7FD6i4RRq9pnf7IPHkcbCqC7REBWSSOv6+nrgcRpzX5fpkkRPU2Ex+2oetRxXOhL3Yz6fjwxgGoFK6RP12zgu4agX8pZ8/cpXvtK7fXRMrjcpIxL3h5TIPR2zpCF6SyfhexIoVW2iHPB8hY5b30mps9n3alyob7Q1BIxZBsvP9PAUvYsx/zVJP9Z13Y/qtUD8aUl/dtFNRpGS7hhxn6ewkyo0ZKXlgKUiGqVVYRVTNhXyYp2JFCkcLYNFZtCAuA1uOw2U67ES2VEwXEzHxjYZkd5XWN/QW/G1pew5Lqx/kdPO85QHj7sVmwbR/LeiV0h7mQgl0bV0l+/pQBaND/lGo2I9qAAL29M6viQyf2veEjhM9TXRdCVraVRN7gOjWddZyQPbw98ex6pNCYrcRgMF3lPxNcvKyIR9clTvekgEeFOOoAIILXprY973/XXXdf+RpG9L2pb0P/R9/48W3DMollMGrYlJ5rpaaRdf77JZD9Ga847SXfRtQy9phEDS8GY6IA2OnVSiZZdhI911t2G/88HMc9OQMwVjg+66E9W6HqOa7e3toR4LRApyg0f35quVgZEQ+UCj6HYm35K3mct2BMb+M3LiZKeP8Zt8T0dPY5ptIl/NRyLIRH45LuSNr6nSZ3ZGabAoT5WBq8auRW/DW48L5abVz5x74FhnVMTUYqZKKN/Mtad8cRwq+Z6SP3/74/SmidEBI0SXm6CAdTPNQlvGsTI/q3FNh7Esf98pZ973/d+R9Hfuec8dxSJNIZEKmVQK0Aq3qjqkW2ObVKFKllUpcaJflk0D0qqzYmQiuaoP2Z/WOC5Db8NXGl3X2UIaeV9rHFo8ldq8oZPO8papd+oc+1jJmOtZVgFpzLOdyecsz07mvnRf3mYaKNvA76o/SSkji65JOV4UGbAMfy8zllP2IY12nvPxnM+r2rRILhKc3YfWPgGaxA6lN5PGiPTy8nKEZOnNmUPzAGZ6guXnwPj6ra0tHRwc3Jl0qZAjGWZvLmm0Osd1EkUn8uCKCuaHWZ/z+qb5fLyKYXd3d9RGphs8UZgocFkP/zbEFUdEX2y/VKPVqcminEtZZIQ5/rw+eUA5XKTMbisNG1dDWCb5yQk/lmE0ajnwKg6uzvB4cEVERoPLRFqroEx7SHUaLSMigopMZyRSn4pqKlDgMrMNdKa0EUnmhc9dX1/fkR1GT25Ttsttof76GHWCbeG4UAar+YWHypm/M1VKkwbBVKEt3kODRcQwhaLNeDqQluBk2FihB/aH11tR2U4b8mzXFAKiwfa9TCG4zDQ490Xm9yUaMPax1Y+pcpKmophly2ih7KmxzvNVmi4VniuyMgyvqGXoKzly2x6Cn1U7K8NJakXIBDWWTRrLyilI41UxrTmJFjixjvi8y8rUDdtNMDllZxLpsw8JCAhEyOep9lfR67L04MY8ByK9auaCu268TM/X2ZMy955erFIqDxRzWilYRHSJAlim1z0nYry+vh7l4asxqMK2VAT22wLKMn2Ny2HePSd2W5NRq6IqhMyxpXJmhOZPzh0QPU2tDa6MXPKsmlidckAcs0RozoNyFRXPe/6Ck7CZJuPqoypKoVxzfiAjkCqyWRWl3KRhzPFMveEcGcc/ddPAxGVWixpo6BI8+OMybB98rcthnyrQxzqq1Wt5H1euEeQRyFVOq0qtMfq7D4AxPagxT4YQTfOanFjKyTAf67rbh0uqcKUylB5sph/M+KkwSrodbH/bkKZQzee3k6puX4sx6elNLcNDp8XUkB/AoSK4jHwYY13UQpGVs2R/OfbuW/UwRguhVcdYJtvEdAvrtEGsjAXTDFTw1n2sk4Y50T7TcbyeiNbG3mkZ/09EX0WV70o5hlk+ZczXVw42U17sfxpil5m8Tl20oUyd4cNnrIeUa+YrZ5GouhqHFkBxfZ5UbbUrx4MgocpOLKIHNeY5KGmYUslbQiLdRVQtBec1/s3cZIUoeX9rOWJe53a2Ui+J3NgHnsvH8pPZHqtqgsj/iYLYrnXnzVs8qJxVEseOx1huCndVNv+3ooVWG8kXlmk+VGNd/a7+V45sKs3WkkvPpeR4rctJVytH2OaqjRWKlsY5boKf6iGepCo9kQ7Ssi+NlxFXeeqUj4zi2M48PtW+dLI25tUYtexWHl+Wtw+eZskQxKjVT76ZpjrpyUcPVAuJ08szVM+ycgUG/7td+eBEenYrGZdPZfmul+fo1OzBq2glFYN9TsdoNMn+cM36OqiFxD3O1Zhlv9hn9oURVE6IMmpJSmV1SoSU/KGRrNB0oq+qvmyjdDe1lOUTXVeGnuPVksFV09bW6+02MiXpc6YqKuLvlIcq0mpFcS6Lqagce0Y/+ZBYgrHkD4+xLZygzn4lD+1AnCJ026mD5CvbXFFGc8vSoyJzH8tzPM57q3IqqsKmZdHhomMVUubvyrHweJaRx1rjkkTF9/1Tk4XrQm5Z/jL1VIiOypS0iM85RlX5b5NTbhlOl1dFR8sg9KptLacx1S7/vs+Kh7ehKhpdxshUzmYquqrKrcBKK0Lm/dSDVgRR3ctyF8nyFN9XQZ/p1Sxd1432z5DGD8nkpECigJZi+RqjL+eMuScHc55GzqwrvbPbtb+/fwdZOF9vJfQeK33f38mjS7eG157a7eOSNBpzl0fBNFWTSD7O8SQqzjmIVZPr8Ji4zZXiuo9UFCJUzzWkIk6hVZZR5Vuzfl8v1bsT8tvHuVLF93lSMmWSSsiHRdJYsB+5nK/ircegMjjr4K35wTmDKnKo+mWi/rmfOZmccweJzPnJsUsnbT5Z9zlxTF4T+GQElX3NPlb2pyVrbHfFo8pYs0+f2Zy5dHeJYStso2GYMuKkRKkVMud1LcHJcDfDNd+TTsbX5SSVNH4qjJOvFHYb8uxnhRj8ybW0eY3P+dpVbciU5D5kjjUV0Ne2eNlCRS2euU7zKoU/y0k0vag+31OVW+U8W0Z2Sm7Jz7w+jUA6Mx5ftUF3H+24SMkP3lOVw7Iyes30B8thPQnMsg7aiGqn1dSHKvWVzp16XjmsZdFzBT6nzr1NJPkoSxM5kNUqixw0opY0RmQQ9zlxPj69fuahEpFXH17XQsNTeXsqaiIaGvG8r3I2/k6jUSGlXJfONq+LKgRSGXRTRh3phKmQFZ/YtzR0lcGmsmdZlDPeT4PAe9j+llGbQtIZedFQWS5YPnfJ5DVpgFZJjGhZfkveSZTfjC4qcJbOMRE5P3YwPEa+cTfVqq7UqwRrOQcljVM3tAGJ5NOh5FhkH7Me8/W+81uPYsw9AFWYk4OcRqkV3nB1CnfX88oODnIqWKX0VPCuu91rmQiwMgLpLNgu1i9ptJVtPtHK9pgqtJ5KT4WnQfFe6usy5pXzq9re6lOiJgp4Fa20xjTb42tdFuWC9dEx8zq2o4oqW46B/eu6brQLYH4y5UMnz2VqGb1SN9ZBiV7Jq3xi2X2dKovXVMYy042+LhG5f3NCm0bc19OoJ/+oY6kTLWdTOexqktM8y2dg6IBazo1lc5yXoUd5AnQKbfm8j/k/76lQSA4InUY1IHm/r69QVOt8ywlkzisFnugrDUAVurMe5nRzlYG9OWf02Y51onIiEFMK+dR1aVTzGvfB1/s8c+WVLLFMPjm7qN0sJ8tm3f5dGZuk5AnrTtnNfhCEVO1cB6XTIVn+eLzSQX9Xcm1K/rMMj32ObWUDWKfv4zbIaUArHlV9aKFqjpHPJaBr3Zd9r2yXy1+WHmVvFi7N8+9ELtLdySci7cqb8k1DLLMyvolqfT09p1/0II2NIT2/UdeUgNEZtQy4j/PJt1R2Gw+3i5O9VYqA6JNPxK2DaFhT4ZKvFcp2H6r9dSSNxobGOfevYb3ZhmrnPSNxpvD4pCW/Sa3teH2OSp38Jj+JHrlnD7d/MHH+JgFCZXBXQeQDjXNGJIkw7WD5yf5IYwftdCVl2ACF/MwJzQrV8sloAiBfW61H5/UsLx9mYh9TTknp1BPUtNA371kmH296FGPeUm4Kewp+CmwKlylD0UQPFaXyJwLkNVlWdU2lgFPIj4xdRNW4LLovlfCzRBUKTWfG88tEF4kE0/BV1BrLqTYtM+FKeec1FSqn8W+1MeVtyuGskqp+5e+WQePvROtSPXmffU2w1ELVeT4dQ9UH8rFyDFM05USXjYSnrvvMI3OSBys34JHublwk3d2dzDsW5nsv6ZF5ff5P9JZIIcNrIzgiMUYb+bi476HAExXkXuqZ866Ehe1mqqoaX16/ztUsLQRRhaz+Ld2i6Xwqt1KiRCytcUnlqJQzc74ZMfI+t8/RnzR+qXQiK8paTr7zHCcViXazX5XDcDsrJLkqMsKteJifBEPkRQVyXBaJOptGOe/nfIYpc9IJfDiXxqgonW3KIFF4tS95trcCg4z+mGI2HyvAcd9o69Fy5tUxhmZppExVmJbb4ma5acT5Ow06j/PeTLHQ8JtBuYY2jQL7yTIYLmZonvdn3dVY5gqRPLdqxa+ihGx7hSj9zTHIa6fqJLUMRl47hbRb6DwdRKYTWm3xtYwOKVt2sFl/K587Resy6ES4pr7vR0YtjTnTYa1+EPjk5L3P+zv5mpFMXkdwYGI6pAI/1XxEy3FNjVXl1HwPd2dN8Fo5j8+8MbehToMljSdGabCSLAiJTphL47X+9n1pQFrCIt3uPU4GVUa7MiRECpW3JUpIA115eCOLfCx9iuG5PGtdlMqex3O8qUSt9vo+IlErQoWIkiojm7vpZf08VjlM8ynLNmonsq/mUShvVVrFspBLVVtjRBC0jBO8LyVqrfrh+lNXKj4Teeey5DTOrKdyKlVbWmNAvvlBKFJGAJQR7oVDXWddVQRhec2J8q673bwt+5Z6eh+ePviuiZx0lO4qoSeprq+vm28qn89vdyXkU5QZYlVIq6VYppw5p4BKGoWdXEmRZbNfTNe4H/7s7++P2uLz7k+FCunI8slE103FqCKEVZKdEgXTfWmt5/b4pROtngLOdBb5Wi2Rc5+zv+lYPX5pHKvtSdke88d1cJvTRP5pZD1GlfOjXtCZcdLN7WZ710Vp2EiMoP0/9S8nHHNSOCMy98n9ZSrT7WE0k2M8tX8/21sBK1OF6KW7YKW1HXPysEoVJ1DLMV40d9KiR0+zpBeiJ04Fqe5rHauOV8rW8voUNl7v36yLAlB5aF/DcCvr9nUtI1AhyEUh37qNeEWMkEitflV85P8W//NctoHlpJNtpSRSNljGVNibfZhCiCmjec8iJa7asA5EzramnJIqY5YfGnE6YoKzqt48n5FUpcMcR85JVPrZcp6kjAwSYKRsMj28DK8q3rG++6TbHu3lFMwlphGvQkZeQ4/nc2QMkVLmK/1CCbahCrl83Pm/ylsmWq+MVG5pS8GucorcDZKTblNKkuNDo0BUtC5KdJEOmuNfUSX8lAXOTWQEQpTPMipD6aWc29vbwzj7k/nNbH9rYtKTYpX8Vk7BPOMyydb12QaCga7rRksZ15VmybGZcqg+z/cLVJO0LCcNtDSecM4XNST/bViJyltG/+rqapCDjNByV1Qfl25TrS4rr/FvR3m5GIP8zRRwyin56nuXpUdB5lOenoZoCpW1kDmvTwNDT14pEetOQ5LInN+V8SICSaYxjCTqyDb62qm1rNmm1v91UxXCsi2L2lOh22WIZbccFo1dpjjyukV1V3K4KLqojrX4uajfVV51mfvelpIfLWNejcOicV6m3qpO6lXlPBPQJP/p+JIfBAWM2NNepXMnX6o5lynZmpKFZXTH9ChLE42Q6MWy4fyfebYKfZFJGZbZ0/OanITlgF9dXeny8rJkuu8jGqjK4oQPhdr/nfc3UvRLmD0X4DEgKqGBZ/9yfHNsTMxrrpKMjhLtkBc0pNUqBLad5WbfeK5C5nkv25Z5TtbPCK1aXUFFTUOT4+w+ElXm/czLp9xWZfl31kfQsWq+dl03vFC6NTeR40P+57xRi8i7TC9UfK30IJ0O57oYaVdRjCNw8oXRbD5Y5m9fx4iLToDtzQ8jlirtd19DLj3SfuYevOp1ZmnMeZ90l8mJgKkUVMbMwdPwZL1Mr7DtJBvjKWPedd2dJxS55e3V1ZWurq60t7d3R0BdZxURMOxNZ9NC8YlQVk1pCNNYs04ua8uoaFlymTm5XF3TSn2krNGQVEpWOcqsu+W4eR/TZsnTSuYX5WA9DqsmtzWd8tT1bm+ld+xDpdPsTysSyDZU15HnlRGvdI33MnXJNpufrdTnFKBIPq1jEvvBkTn3DalyYalYPre1tTV6Cw8dg8ug4pjSyLeUQdIoB2mqJiLsfZ1/r5CBjzM3Sq9tRc5ldkQ6WTavqZxbheKr8tZBGRFVyDvPZZ66uq6lJC2Fr5B3RiRWKis7V0OxvmqcU6l9zLzM6ISomoaiNTZJBApcQcPz0t2wf1Xk/ky9x5b6LI0nGVtpytYkdEtOK97l8j6eu7y8HJbwOtqdcniuz5GIt4owWs9Ik+9jdX+mHJ4BTLW6xVRt03Cfh/weHJlvbW1pb2/vjpekIGR6g4aKu9rRmHPwq/AkjbF0F636U00s5dpUX5eCmc7IwkGj7j7s7OwM5+1oMpROgyGNHz6woFYTPvl/XQY9+5xtqfgpaRjTdES8LiMvlpl1SWMkmOkVtpNP3nIPn5bTdb1ul6NK9o/L2HKspxyE28EVT2yP0xzp/Nj3dTzZ674RqVbLME3ppLnKa0oGqxwzzxFQ0Yj76W9+LFNnZ2e6vr7WxcXFkDL1J/Pb5IlfjL63t6f9/f3hN/XcPEh58bHsh/dSsnGmTHEJ8sXFxajfu7u72t3dXZpfjzYBSqpQKL/T4CaSy3Am60ojwTArQyAiufn89p2cFJw0CNluG3obKKMaIjMqKctIlG7hqwydVD9A0hKq/L1qqsquHDaNcxqMKSLfXE72l2XmhkzpvKmQLp8Kx34lIjd4sNywjckL972qk+NGw1XpRAVQlpkcfxeqotmK8hqPVxU9VTLuuniNgY+kkT4mMrdBvLm5GRD52dmZLi8vdXl5qYuLixFIqyInG9v5/HZrEPPZjjJXrTHvXYEOE2V8yvZ5m2qO1bLjL30G9jOv9jqg4kgaDCoRFjuZxoF1mej9u+52p0Ma54uLi1Fo5mP+b2Gxp6cxdx25k9xsNhsQ+MHBweBxnTaysBjpGRm4L60HWCxAXGLntthJ+KErTjb7dVqrpDReLJ+CzDy02+W5k2WIRtiGrxXac06jMua5OyNzoX61oQ27f/vc9va2Dg4ORq9AdGhPmaKC5pPH0t35BdaXRjDLlW6jPvJhlWT5Z6Qx5XSZ8vH/1vXU5+qJauqAdGsDaMwvLi6Gcb+6utLFxYXOzs50dXWlly9f6vLyUufn5zo/Px/Jg3QbZfizv7+v7e1tzWYz7e3taTab6ejoaNiJ1Xy/uroaTYRaX6v5FfadSy0TMFiOuCOs+86l1IvowZ8ArZByi6oc8H3q4f/qmgrJ0ZBbQCwwNzc3uri4GBl4IrF0TmSWdPt0a9UvblfLPhC5VBFIhvSVg5vK5a2KFhmSFmpnCqai7DPH2+VyPBhOJzLnN6MwG13Wk23nKgZ/bKhdB9MkdGDSeHUPDfJ9DDD7uM4IK+tkaoJG15QGm21LI5fymHpI4nVMR1RRs9Mp5+fnurq60unpqa6urgZjTn2lzvlT9XF3d1fz+fzOU9qMLChjdGJpf+icfX+ulrK8cEvtRfpBerSHhlqpACIXKtr29vbggenhmO8kwmZdrJuI3AJib+88m9H31dWVTk5OdH19PXh6h2/X19cjTy+NX0pt5X/y5In29/c1m810eHio3d1dHR8fD6jKaNDCkA+hUID8UAKXAdLIZMplNpvdGQejx8egnDicIuaI/b8VwpqvieSo7LyGUZXHdG9vb4haqnDaxntvb2/49vWOEsxHLzklb/2fY0EH7rZzEpbflAtfz+WC63DUdHJTlCmEdLSmdKhuN+eMTETtacgdJZ+fn+v6+lqnp6c6PT3V5eXloKc25tbjdOCulwZ0a2tLl5eXA39fvXqlvb09nZycaHd3V0dHR9rf3x/GxKhdGu9hn/lx9jXnatiGlA/awmXo0Yy5yUykt0skyQGX7s7yEr3Rc5uYzyJjJQ2TJE5LEH2fn5/r1atXury81PPnzwcv73SMc3GmzKc6LJvNZprNZrq8vNRsNtP+/v5ImNxWC0emjXwt96GhwWHOjorktetUjFxyuUrKFEv+bqVhqijKBpA0FWFR8NlXG1caBZ+Tbvf2oUOsjE2mwzw5ZeU13zgHQiXd3d3V1dXVCAVWy3KJECujTkrUtmq+Ws9yTiCv4blKBsjjComnPtPoWl6pt0x50pifn5/rxYsXurm5GYw5Zb5Cu3TGXdcNPNrZ2dH5+fnwEhij9L7vNZvNRnKXoMqpTy544Ic2ji/AYRTncj/Txjy9dXU+c+aSRh5QGoeazEUa4aZn4+DT0xth29MbjZ+dnens7Eyffvqprq+v9f3vf38w8GnMXYcVnrm0q6sr7e/v6+DgQBcXF0PItr+/r/39/QHhua9GzszDZa6Z48I+Z+heLbV8CKKCVMT20NGm0lflEqmT5zbcvteojAiOTp5RXxUCm5gWyfH18ZTV3d3d0TIzLo/LMZJulddOiOVmTp/y1hqnVRHRv4lGt+W4E0hR93wso0jXR/5XKRqCLhvs09NTnZycDCg80XLOyVSRTM5fGATM53Odnp4O8yVG3V7hYkPsvHem0bL8Cm37f4IX1rMMPcrSRCtVpTxpvKTxADB04T3VU6W5xJConYy3h3/16pUuLi706aef6gc/+IFOT0/1/Pnz4ZiRua+3MeeDAK5/Nptpe3t7SK0cHh7q6dOn2tvb08uXL7W/v68PP/xQT5480eHh4Shk29nZ0f7+/mjNa4ZgNhhOKTFfawVKlHif/Nt9iQao4qvbId1dX84ycvUCoxfeL91NpTDiYv705cuXgzHN9BQnEX0s29mKCExW8NlsNnJMTtkZMBDp0nkxGmM6puu6UTrO92b0wLFdNZkHHp/WevZE3IwwKjTOyISGjg7NIIV56evr64GvL1680MXFhV68eKFXr16NnKb10BEyly6zza6Lr5rk3JntVaa5vLzaddi+JNBiPzkJTkDJ8XLbrUuOrpehR1maSKryp5WS5zl68BRoK2xew4HzgFrRPIlyfn4+hG3MxZ2eng7nmY6R7qJjoz0z0Qxxnmxvb0/X19dD/m1nZ0eXl5eDgFxfXw/Hq7RUNWbuVyvyWbfSm5Ytn0qdx7KsTLtlaiEdtflKQ2pjKo1XsBAB52Qy66AMVWmdLM9987nkZYbp2eeUYZ9jxLNuXrKOKs2TxBTZVFlJHMdW2k26feze3/4wCrLOMQXmyNfzIwYbNKoEBVw5Q2DCa33OzpUy4faQ/zkOlJNqnJimWZbPD76ahTPKu7u7g6Eziq6Qna+vOkavxskmCkUiN6Nxp0yMyL/3ve/p7OxM3/3ud/X9739/yMFdXV3p008/vTPxSSEykRFbW1uazWba3d0dljrt7e3p9PRU+/v7w2z7xcXFaNng7u7usCRpPp8PYSKNdRpnhsNuA1do8Np1GAEaxSo8p5GisFvBmJKQ7qK1KoRNfnLC6+LiYpi8fvXq1eBErdgOmckfR1O+huk7y4/HlYbWqM/LUJ0+cztyWWs+mej+Zx7dfeaYMCpdlLJcBVXRCr/JVx9jZEHepuPyfdQnH2OZHn/riyOtk5OTYWy5YyGXjnqOyssPW8bcyxeN+rkEmfp2fn4uSaNlquwXDX3agoxGUpd9zH3JDMUietSlie5ohbxbHcj8mj+ZU0zvXoWBVqgKkb948WI0K85VLOksKiKzmM/1hKhzfjYwNt4Os90+CktGGll/Nb4W2nUjcwov25NUtdGUiK7Fa5aVqJyKaJRlvlHWuHIoPzlGGdExFPeHaN+To+aZHbTvb4ETGr/qQ5RsANMay1UQeUrDQ8NdEVMspExBSHcnsNnPROVMfxCVu2w70q2tLR0cHAzrxf18x2w2GxlIlum0pee5bFAtU3RIPsYy2H6OnfvBVW6VQ6Pdqz7L0IPnzLnUqeu6O8uGKATM+fq7MgY5uZCTYVR2C8D5+blOTk50cnKi58+f6+zsTN/73vd0cnKiTz75RM+fP9fl5eWA7mzUM7zmN9vqftgJWACNsh2R2PAYRRhVsLzLy8shd+fvHEdThocpJPfx9Pch8yZDTVPFNz4Zy+tIRPJWFs4h+LwRsCevjcz7vh8enrJR9YS05zJ8jHMU6UjcBudRvQLK5c7n8yFnbmNm9M80nlddbG1tjfbVNrDxRBpRKn9TRzwm60ToTA+xLTRcWbcBiHR3zXWVUsmUQzoN9p9yc3h4OOShPd4HBwfa2dnRkydPRnzNnDnnkBh5WY78AJJz9DTY5l3XdcO8GZE6UbXHg8crsMC+G3g4dbs0r5a+cgVkgSVZoHMWn0aTqJtIMxGJy+HEQqUYnrw8OzvTycmJfvCDH+js7Ey///u/r5OTE3366ad6+fLlYBA4gz6FNt0Wt7HrbjeatzG3IHJN83w+18HBgQ4PD3V0dCTpNnQ3siQa4xvFczxtcBge3zf39jbEMJzG3MbHbeP1/nCpoHnuNqcxz3NURvOL8xxpMJzO4pN+Tve1FI5OKtvtvtDQS+M19V6WyBUPvobpinxwjPVXY+expa6sktwmok5+V231fTlpXfUnkXim1nwt9dnlm4993w+R0P7+vo6OjoY14Xt7e6MoqSX/lDE7aUfqTtNZvhg5WzclaTabDcaXupdpFo8NF3JkZCJpZLOWpUefAE3kyIElQ7NTiZwWlZ8hstGS927wSpWzs7Mhj83JkEQXrfQADQdRFxGGc3F83Pjk5ESSdHZ2Nsxi28AkcqnCr0QCSVSaVVOGlC1+VCisRYvCeI6plYqoyryTbh8MssJ7EyWudqCyUSarupkG4yqqy8vLwUFzElvSEP4bAUq3O4j6m9FT5qQ9fg9JjGpowCnfOcmXcplzJ1WE5m/qWYICGjvLN5f52iF7ZZjnL8zfrCuJKSsaUtsLR21cYtpa2UM+Mlde8XJRtLxSY9513Vcl/Y+SfkhSL+lbfd//la7rvijpb0r6EUn/n6Q/1ff986VrHtcxCLQHLw1ECgmNtNMPnHjht3Q7MZErV05OTvTixYsRIv/BD34wLHXig0EUKv9vkc+5XS7Limpk4Xbu7+/r/Pxch4eHA7JuhW1GmFU9rQdfsr2r5mumwKrQOlFJZfhzrCn0VHTz1GkqO2KnzhjR2DEeHBzo6OhIh4eHQ5rFRjYdTEvxHfn5t+XPzw3c3Nzo8PBw9HTo1taWjo6OBuNnHtmB57yRHYbrJ4LjGFVoc5V8vbm50cnJyZ0UHR0yI59sSyV7VZSRuXTbAUdcdJoeLz6kc3h4OERYR0dHo0nsCjiknWC7vIjC/OTqM4M8t4dRb/LF4CEnPnN8eI4bt1HGlqVlXgx5Lekv9X3/45L+sKRvdF3345J+TtJ3+r7/MUnfefN/IVWIu5W6MOVgVddnqJlIkV43lyQS1eXj/tnGqTZX7ef1DKeyDa7bv5mjy9B2Ud0Vci/uWylfp1B/y1AmsX2LZCPHkpNi1XI1ojT+ngq/p1BxxZtc5pZpAX/SSZNfy7aHYxS0Ur5yHCteLJO+WxZdtnQtI2MaS/PR6ZQ09K2xrY4z3WPAZH65Lv7OlErV59b4pF3KSMXX3CeFthCZ933/XUnfffP7Zdd1vy7pK5L+pKSfeXPZX5P0v0v6z6bK8gQDJ7CYP+S3dHfSzspSDQK9ue9znpuG05NjL1680Keffqrnz5/rk08+0dnZmV68eDFMUOXELJ1Da3ATZVAYmfe3sT45ORnKms1murq60vHxsebzuZ48eaKnT58OqNMz9dkGRig5bpwkLRDSyvjquvjYtJE0EXWFwlvCXCG67Du3ODUyPz091dnZ2ShfyfmIJ0+eDPnyatldGteqTWwHlxienZ3p5uZGR0dHOjs7k6Qh7DdS50Mq+/v76rputAY+xyRz8E4VZo4V96+Mr6yXyDINa8uQVc7G19NIp9GugEw6Rq9O2d7e1tHR0TDx6f2IWimQqs0sn5EzHYKzBV71kik8f0vj7TNaRD1JfeX9XMyxiO6VM++67kck/ZSkX5X0Q28ER5L+hV6HddU9X5f0dUl69uzZKIRMY5zKw2+nKlLYiVar/BwnRGnUma/mk502tClsbseUF50Sal5fpQm8wU+uRWbayWPnNlAoadwTFVQTaqvk64cffjjwqJr0YpuSUgam2mlKNMyIxr+N1jjp6U2ymCdfhCorQ+42MJ3ldjiv6iWoJtdHhM6US2ulAykjHMt3NcH8pu0/onfg6/Hx8QhMsY5KB6rzUXZ5Xa6QqYAay6CDSdR838n+dNCUY85BbW1tDWCKji33YXE/Ul+TiLgr3SEqXzayWdqYd133RNLfkvQX+75/EYzpu64ra+z7/luSviVJP/ETP9FLt4/yVnnACnVWBt+dTITvPBvL8TF+OOHpDx/Vn1pa1xKyHPgK0VPQvATJOXz/ppF3Dq+lHOSDy7eCJwpiP1bJ12fPnvUUvgoxJh+zLURGeT2NNz823HTQ3LGQD/McHBwMeVArYEVTIXHrHp83WnOksL+/P4o8M2z3agw+ps88L/nGh6uoA55oLVJZK+Grc/2th8DctwpMVECIKzp8DZd9ZkTCvjHSJDL3g0HpeO5LCR7dPuo2UTrvqfLw0t2J/Mw8mCg/GSksS0sZ867rdvVaMP563/d/+83hf9l13cd933+367qPJf3espVasNl4Gp2WYNhD0oMzTPH9mbOkMc8JM4fl3tS+tYKligQqg1Uxj22jAeZDS34KlMbp4uJitFtbGmLpFiU4+nDbMsXSMLQr5SujGROjg8oZ8p50Tl3XjXjB8JMPkKQht4Ezcjs4ONDBwcEovVKhtynFyWjHv1PG+r4f5Ovw8HA4z4eSaNTn89vlbHTCJvI0UR4nC7mZ1Kr46onb3EaW591W7qHCMaNue+w4/jzvvvic8+IeF/LAE9tGzHbQLR5WaZVKhyv+0lltbd2uMScvGOlVzpjlpc3yp0rL3Ged+UKz373u2V+V9Ot93/9lnPplSV978/trkn5pibJGm2C1FJ8elt8ZeqRRz3DV1zDNwm9OWrUMZoUyKqruy+tTcNgfTt5xEq9qc07+VQ9XSeO30xRCuzK+WgnZTtfVmhxaNJ4mjhWFnvxP5J7ywCc7Wyscqj6x3mxP61rypHpSmPKc8pqOjB/p7pa3aRzfHFsZX11nVS+phYTZP6aTOG5Zbjo8O2ROOOZk5KLJ5NQL2hjeV82X5Fi3+pT3+n+uOc9P1s+ymZJbhpa56o9I+nOS/u+u6/7hm2P/uaT/WtL/1HXdn5f0zyT9qUUFeTKKu90RreYkaKIhI7BE3UNndm5f3uDrc2MePurt/HiFxtMwV6kOGhq2M9FnogL31UbcL73Y2trS2dmZdnd3B7TuB06k8fa4Xs7m8qzQXPbnibVENG/aszK+2phzzbWftnN9GXHlGKcSVOOfjpd8pUx5TKj0foCkUtRsQ6LJiueZ8qCTYZRV5UatoNxGlZNqLYNSpbCo/Kvm63w+Hy0iIFUOKMdx6rrUMxo+GzKPDYELDV5rV9EpqvjakocEQBzn1HPyzO1x+1rOgW2u+JpPrS6iZVaz/B+SWqX920vVAqqQVut8HquuI2UY2kLH/BDVJS2DHJe9JwXDx2iYKLSJqqtyU5haAu3rIh+7Ur5WY7hISZIWCW3LwVYRiMtLZJSTefepc9F1lVxV91eoLPWiOlb1r7h/pXxtofJlDUw6pgQ6JOuv+2+gx1QDx+xtn2y+7/W+h/aFWxZk25g/b+XSUw4zlUP5WJYebddEIkjp7ksAjFrdUZPRTAqZES5z5tIYDTFHz8mzTA1UxHZNGVdTi4lEcEToHAMiXbbNntrhpctLtOg8rL26x5wKsGoykmDf0/i6jzmpZ/Tl80k0Wh4vRh9E5hxX5qmrflfjQFTOa9JoVPM6blfOTyR643hRNruuK2XDdbh/pvumju5L8/nrXQRzbxGPU+Vk+D+jR0mj/pksp9Qxgxr3nVtYZFqkomo8KtlaRORZbkvdAgOJ0H2suobOP8/blq0Mma+SbMgsFDkY+Z+P1Xpw8mk5ClKFahORcfCYj16EwKi0VYhY9aGKMvzf6ZKcA6ABoOFzealUVvacJ6Cis+1TCvC2REWcKr8VbSXqyetZPo0txyjTGdVniiiPyyB2trGSx6mIiv2ucskcn3QWbK/5zD6vkuxA2EZOuPuajEQyQiFvfSx508p5p6Ej2iVVxyo0LN3fqFcOfZmxXmat+yK7s0iGSI/6DlAKR2UQ7QWZW/d5zuS7THt4GmqjHd5nj5/IKFEU28vjleHJPlboLimFOENIbsPZMvLZLraBXj1XAq2aWkiEik3lNJ98LQ00j6WMZLlVTpp5+tznh+PDdlblZ//uMxbp2NLYJqBo8cT3MW9dGa5lAMl9yTx1+a6bslYhdUaTUnuFV+47VDlB24B0eLQBCXpYj/uR/5cZq4zWUw4roJDtqlA3nVZGYVn3fRzPo+xnbrLBJhFRm5He7CeX39lgeWmWJxyIeruuG22Owz2Q89FfDnzFcCJ9tyfP5++WsfD/RLIM0TxTn8jHRsCUqMl1e0c3PyyTbVwlUfFdT4VOK2OeE7p9f7uZUsvhuyw67TT+uQKoKsN1tJAteb7sOJCviTw5NlPzIlVfsy10zutw0tRPpnukGnX6vNOc3E+FkaGkcgOsRPSuv+u60Ysg/LFdaJHtgOvI9F9+J7nvrSjC/UkQ4P5zm2TaA0Z+vJ/66z56/JahR3uhMxWQlDO8eYzl+DsN4hQSo/HmhEVLWVsGnW3I49nG/M1jKSBsG41TGuoslyiUkUM1Lusw6O4/nY8RZRqpjLAqVM5yq/7yP8eNwu8xySeHW1HTfdB3RVTQNOZJVFi2MftdIXufI2pbl5Nm2eRXto3fbG+eo0GtUCeRqqNm85YrWzh2Oe9EokFvyVILqefYp05WiP/m5vZFJdbZjBzIW6bKclzogJahR3uh883NzbBrnSkROVMoZJ50m+PM1Qq+Ntdnu34vU/M2qEatNJaVAadA0zG4vqqvlaOqUFs+SMINg7jhvqQ7/TGzuVUnXzPHeu3pV634RsgeR/OF/HI/OanF6MLHzNcp51iNm3eg5H4WfKqWm5dJd3cqTMNTOZcpY5+RgZdFtiYn+74fbSvhvbJZvzTes6iV822lGd6VrIPpfNkuGrRMQ1Cv2FbpNZ/8Dl0TnbJ1wE95+kOdthHnFhwEi1knHUXL0baMuvvD9rH/vpdPbNu4W/58DaMP20OOKetddo259Ej7mSdCTKSUoUgOeqLyTL1MIdBEv4milm1/hTpIlSGq+pFGoNW+RZMpFfJsIcJ1kA2x28Nv8naKX1WZLYPr/3zIhI9ZV8iXe4AwGsv/VdvZpik+JIqrEHW2L0FHhv8sZ50IvEXscxUVmFpoPK9ZlBZKnWA6lGV7riuXCaacLaIpPWGfKG+pj5Rjpngru3IfPt7nCdAHz5lz+RgFJHOfPk5vyjDP/6mQXIJHxjtc434Y3tv69PR0eEuIN0bi29yrUMzf1R4NUw6oEgxHJ9w/xJED35SSb8PJtwlxjHPihePmR6BXTYw62AYKOPnh8csQ3tfznO/LXPj29vawJ4ff73h+fj5C/l6y6n13zs/PB9TneukMGe2kgWgRkRv3gvELEoi+csvj3Pa45eBs2DKdxvFcB6UumYg+8/rqOPlKQ+5xcz004FWkQ1Tr8TDPM5XRSq/w2DLO2W3xtgHewZRywz2dco6PCDvTQXQKqR+tqK5FD47MbQArA8R8nFTnM1vejpMFVEgzywKyt7en+Xx+x0h6EtUIrhUNVCiyQiNVyJTlMRdIQckP0QBTUK1oggLP+ph7XzUleqwMkvvsULTai8LEsJnzHF13O6FtJe/7frRplYmTSDSYCQRMrYirJZPZfxoeO2IbQveVUUKFzFvjRkeZUVjL4K6CaNASIFS60GqH+5YrXWgsczkto+ZM4zBdkyt96KDdVn6bquMVz5m/5yvq3A7ro8uqljynPWAU24pkFs3nJT3KBChDUOZUc+KO/ym4VXgjjZ8g83++FV3SkF/1Dnqz2Uyz2Uw3NzfDLnecha+QR6ZlMppooRMqhtvmNhweHg5vwPGbaqqllt5Eym/ISWGhI2s5lFUrPfueq25a1yfqyCiiKt9OgH21g3PO3A6RbfKmZd7KwdexLZSzKbRGIorM3RD5ImHL+c3N7bsl+WpCaTz/kkaT48B+VWO0SjKipONtRStE3jxH5865nsphVk6J+ky5YXRcGdTKMbaMetVvpjqta7YZOUfT9/3oocR0zpSxtFFuD52P6b7zII9mzDPM9ScHQBqjtDRQFDSjB/83WmMI61TFwcHBsLPd0dGRuq7Ty5cvB6GwkrWig+o7HQ+NAoWVxtyG/OnTpzo6OtLx8bGePHmig4ODwaA7UvC+J34RMYXW9WSImuO3LnJ/p9A2x4H8rJyB/+fa9cqYSxomiW1EaUD8ur6zszPt7++r7/th4p1pIcpWCyGbiCaJ1vx6OqfMGE3c3Lx+IcnLly/16tWr4SUoNprkkceBy+NaeeD7Irhlqeu6YW8gyxCj1pT7fMiP7akm7as5qyw7c8/WG9/Pcm1MqxVMblf1YZvZRkfGfPGFddJt8GS7+UzHTR6mEWf7WuBvGX0iPYoxl8YvRbbHMlm5UrBJ6bGTORkm0ZDaSDCPznf9cacy5+UqxNZSninUYSFxyicNgY8xvVIpais0q5CIf68DlbP8TJO12uVrFk1StSKbdIx93w8GnFvL2vAw7OWj/1K9XnoRcnNbcoIu02JM9TE/zu16iWQTkebYVO1ad74825aT3OSr+U9dacltC4FXAIj/Pc7+n4DGSwO5IV8LCVftMt8sR/zmS9b5QGMC0oz2KtBHJL6MvC1DD7400cvT6HEpIMxZt0JNTnZKt3sp5wJ7l0sFd72z2UxHR0e6urrShx9+qL29PZ2cnIyYYHTu9Ea2J5GJj6UAuo0WiqOjI+3u7urp06eazWb64IMP9MEHH+jo6GhA6EdHR8NLaivk4rJTmYmeqvFbB9lY2bASMbs9zGUS1dqg54Sny2qlj6jY5uf29raePHkyGM3z83Ntbb3ef/ri4kJnZ2dD2oDpGMtENQFvSiPEJZF+OfTR0dHw35OflrmzszNdXl4Oryv0Xvp2RJZP1+PxSINUtWmdEZf77PFII+7+0fmY90ylpUOgMef4Z8TmD5fqHh4ejuTeKbSLiwudnJzo5uZmeIUfXzTD8lMnXJ5TY+bjzs6Onjx5MiB0voSd6JuTt3t7e4NTSafG/4wKJd3rAaGKHg2ZJ6om0dNX9zP8qFZGVLli57Q8gDmZYWS+v78/rHiYz+dDmqM1e1/1jf9TgBmWOy3ANe/Vi2mrENRGcAoJV//XheQqZFVdw5Cyui4jn/xPQ07n4P/cDImoO1E5N7bKCKGF4tKJMhQ3IudvIkdOwhKV25hPRVhVm/LYOqjlzKTpt8aTzykPNOLpGKp+pcwzTUnAYOIqIZ6rVnaxDjpp6mj1mZJ1Oh/XOwUOVsnDR8uZS7cTlOnhpTpd0DJEDHl5HcNYh0ce3IODA0m3L5k+ODjQ1dWVTk5OBoY613p9fa3z8/PRHuxsY4WYyFTmuff29vTkyRPt7u7qS1/6kp4+faqPPvpIX/rSl3RwcKDj4+NhEnQ2m5XCTIPFCTSi1ZZBuu+kyjJUKSMnbF0nox6OVY6by8hjzP1b8egcrq+v9eTJE/X96z3F0wCZj9JttCZp9G5Hlsf6Ofbmp3lKZE7Ef3Nzo9PTU11fX+vk5ESXl5c6PT3V6enpgCbdFveP0WM6s0TBdELrWqXE8SOf2caWQZoyVK2okiuQUo7dFjtvT0x6hdrl5eUAzE5OTkbO0zLZSo2YZ57gdASwvf36hdGMrL3cNY26AYXb437SDqUtsxNhP5lPv4++PniaJQ1d7r/SQiAcfIc19HoUDg8qEbyZZcXhbLSRuF8MYc++vT1+ejBDK64ZnTLmFrzDw8PhiU6nWb74xS/qgw8+0PHx8TAZ6llzT+Ay9CQqoUInanG7kqr5h1WQJ4HIr/zwmqRE6zRinNimY5A0RE6WCa9MshLRAFoJJQ3LUC0HTolkvSQaklzayiiP6Ru+AtDvBuXLw30dPxXia6U1FhnUdyXKMoEM9Y/I1m2lbricLLeaJ8tP1Q6mM/xtI8oI3OPul4RwtUmmOxwNewLdixO8WIILCyoDneje9WUmIUFKynw6nc/0BKg0ZmQyzedJRAa5TrVlRBPV+DqiJkmDkd3e3tbTp0+1s7MzPABgRb+8vFTXdbq4uBgtO5oy5lyx4plwI/IPP/xw+D4+PtbTp09HLxzmQ0M05nw8Ol9aTaLAcIyZw1wH0RkT4VbpEl6ffUjE6f8ZARBVpTF3v/3YeBpIn7NBtdFldEclZRol0ylclSHd7q3PPfONxHMPfdct3a7CSIfF8c2xNl/z/KoodY3/czK0RZWDdnk+TmNIufVDhpeXl8O40yhXq2EMGjzuNOaUIfKN7xPl1hqMxswfOgl/e9lrrmTJbAPHNZ1w5bzuE3E9appFUukt83p7RK4TNhNSuDgAGbZ3XTeawGBO1eX6abK9vb1h6aJfA+ZUi4WCy7VIFA5PhDF18tFHH2k2m+nLX/7yMPF5fHysvb29YY25l7dRgbiLotMIHMsUDioF5xTWRR6LbFfLgaSBJk/cdiMiP1nHezzOfX/7hvu+7wfnacV/+fLl6LFvvsjEDpqT5O4LjYQV3MsNicZ9zv208fATp063eFLWu+lZbmwE7ADcB+nuSq00DlPzJqugyggyLZHommmEjJxJGVlI40jI553qpE77iW2+otAps6OjoyEKkzQ4UBthtpErV8xXpsns0G2cnTK7uLgY+OlvO2uCvXwQzOPJ73SQVaZhWXoUZC7dncTJcPFtDE+GKzb46SQScWXu9PDwUJeXl0OoxryZvfCUMWc+zytSjMxtsP2wECc+E+klml1WaddttBcR2zo1uVhdz2NTZefYMMVmo7y3tzdy4JWzd5umDCP5UfHGZRBp23jYkPhYa2VFS+4zGqnGY51URb75m+1stY854iw/bUDe6+tsIJ2+sDFn6pH10GjbYbNPXGRg8JWbg9EYcxLddiAn1WkX2K+MYipnVo3zfejBjbkNZK6vNSIjYpFul705fJoK7VLJ7GGNlFyWwzDmQD3hcnV1NXjqi4sLHR4e6urqaljG6JwnPT0V0+V5rxBPinkppNMre3t7Oj4+HtD60dHR6KEXqb1igEjFRBROdORziXjXSUQ1pkrZFzmr7JPL5LLFHHcbb+dTnUbhnIzD6VQs5jXZVobarsftsoL3fT+8lJvIzKjbaNwIjimoqk7/zsf+OS6Z3lk1UfeYOsx0VaZKfY0NZMvw+2MUnWCGZTElZoctaQBHno9yaqRaIZTEtjHt6zbZTnlHS09en52dDROs/vaciMuQbifUCRjdf6L2KrLPdi1Dj47MM6xoobjqeiumjTwVw9c4/M68KwXUAuDrvEJha2tr5HEvLy+HkJphvz98YMSTJ157zAlQp1KIzPNhkyoPWnn6Cg0R5SWSeghjbpoyVj4v1ei9JcRT6NjG1rJhBeKWwJQVptmIzhM0VMaM17rtnHhjftUpARv+ap6jRRVir9DxOnjrPrLfdGS8JtvoNnGMfY2BWSsaYT3+T2RsEOWVXFzXzftbS3uTeMx8Io9zQ7R872y1HNJyStlhX5l2SV3mvZWct+hRNtoimkpFyxyTdDdES8+eubZcY+ocpweOdZjhfFjp+PhYu7u7w5LF6+vrO8icQsw3GSUyd+6W+6h/8MEHw6qWqbXlST7PazJnyVAzlYTGbh2UwmcFrPL1RizSeN6EEZuRN8fa9xmd0zAybHWqy7JBp0wnwrwpw+zKCTEKy1DaRtyRYD6sYoTuc6w70SvRWjUBloBkXcg8Ix5pvGyO48F7KAeZJqpSShkh0ZESMfPj40b/1mFGTlyEkA7CRKNqwGag5u0WjL6Nyo3U+ayL6yFvPAZ0+Gm7fD/b4vZ4/JalB98ClzP1XFbkiajsLD0qBSQNNsvnva7LA20jwGVrFgRue7u/vz9atzqbzYZ1yjn5mGiLy5kyPN/d3R0mO5knd3RQoU46P5fn+pmScV/4OHuWsS7FbyGJNLhsR6Jc8yEpDYaNC50/76NseXztVHI1FBFnIjm3N9visL9Ke0i3D67QUDiKY4oxn+6tKOtvGYJ1UNd1ow3JOM500JkyqtIq1WR39pnXpc4zsmGE7N1OuXsmU6ipV0mu07LhPLifGHZ6hZPXNvjkv+2KyySx79WKlWocqkh2ET3Km4bYWKLaZTx5Iqb07h5YekIPGO9hTtffRsc5EeL0i5cr+jFxk405BdBpGgqTjbHROlMruQwqxyoNho/z0X0ajxbdJ8R/GyI/kjf5m0T+sBxp/K7LHAvel1EdlZjhMh1d8ruV4yUI8XhzMo73ELXx/kTzPsf200F5XCo58Ll1RVmmdMSVMTIxzZGoPUEa+5I6baNqvtmYpzzZdnh1Cffb8cIFr/1nVOr7OdZOh9kR23hfX1/r7OxsZMwpR8krEvvKxRLJx/ssP5yiB0+zON8svR7Ms7OzIXRmKDeFjlL4jbg9KPbSXm7IEMsCwtDH37lH8c3NzbA8kExMRtpLEznnE2U09BYut40Clzk+99NjZMWmcSC5HS67Mq7rIiJsIiyi9SqHmu338Ywm8hjDa64kSF7weQFPRPpJUPPckVIVXThqI8/oOAhIsp0edxuLyqBbH5h3TXDg/6QpHVkF0fhQVzhGHmfyaWrSNiMoflOfONmZsuu6bLQlDevQnXbhLprmMZd7ui3+8BkSt9/8ODk5GaLy5CP7mOjcdoI2g/s9+R7KlcfjbSKuR50AbQlgdm4KzbXKo5BPKUF1nMbQjKXAJlI0Y+n9bUi4pzUNUDXBVn18TSuPXjk6f68bhbdoCkHyOHOfvJZOYVGfW6GqdPcVgdybh3nVVu7adVXGjFHQIoPaCptbIKXVx4oewkmznmXkdapd1p0pIoCpyjXZkNrRumw6GkbimcLjxDSNuXWaS0xp5FtpE7a/xfMcl/vYrUX04MacXoqoicf4kmeSB1m63avc+StpjGY8SH5YJInePQffE5XcP8PtbDHP54hELFBMpZCqlStpwJh6Yf0MFTkuFAbP+jOF8FCUaIWGkgLMMZPG64OdD02EKo2jNxptlsPUllEcl8F1XTdsrZB88LeVWrr7/AIdj/vCvvq6VGwje5ebkYd5R+TaUm6fay1jfVeq5kDcB7a9Ojdl7NMQZpSRxrJatui5CKNnOm5Pfvt5Dtsa1kGD7eWitkNE3Vxqyt9sJ/Ux73P7DSJc9iLnP+XEK3qUCVB3lMbUx1MpfZ909yUVVYc5sGmkKyGjkWEelY8NV8i4QpoMk02cjOHypDQKLDdzpLm6Io1ihQQ4pmnk1kmJUnOcSRV6yfuznBbim0KLGRFRhvj0puvmJJ+kkbxSCZnySgVdBoGmMcsU05TCu46HQOYeK7eN34vuaxn0amx83nLLvpGXnMuSNFqlxFSPJ76tA1xy6GuIvmmwK50iIq/0zZTXUbb8zTmARRHbsvQo+5lTobi8x9dI9aQThYH3pQD4Y8Tl8kj04hW6sYfOndtogDPXl44oDYlzdpliYXqGIWI6Jf+3Acpwj8rAsXHf1qX0XdeN5gWsHNJtjrwyCKYKpVXnzAuiI8oS50VSVnKMcpVROhsaVtaTcwJTip0Ik30xCkwwQkfutrMuUob86+BvppH8qSbyKh0mVQ45Zd5jy4jZxPHnsbyOjtRzXdSxysakQ26NbTqX6nxGLTTktAGsh5SAbFl68DQLB5QTAmSmNE6pZLjJTtJLp0FgCicnqaj8ro8hklMU3KK37/s7RsDlGwnkJAuF1d+cyGRbaPB8jUNJEvPwRN9pRBj1ON20DiRHIZVuJ5Yy8qmUPJFaoi5eb2NO/lN+qGA0MtUyMmns0M3PRMMeWyLFVLIKYbXQHM+7zJxHocF3+yvDyPLXMT+SxpM6wnSEiXqa/GDEwmvcF/KA+pZLICsD7GPUZ5/LlWZT0VvWwfw5204Z8jn2kbKejos6yHp4LR3JZ9qY26iQqRw4IneTlbOlLJUCS3cNagomB5n3kLlEfkbfmaqwYNIAVYaslV8koq+MXu4p4X67XexTjkHVh3URhT0jBp9fVjhb6aB0DFPXkGhMLAfZzpQ7l+X7MmJIQyXpjlFh+VWOmJEU+54GNB1g5VTWQTQ2/k+ijuWnuk6qI+2MWlqIOK9z2ZwjqXSIRpnOvOJrK8qpxjpXtVRl+XdlZ1rHaX+WpQc15vP5XOfn53cW9ksabTtLlCfdRd1SvRomGZiGl4xKISS5TVTwRAckpxgYAXDvj/x2GURkOTlq4fMkjtvjuh2is62VQaJTqQzDKsjI3GPJiUaPP6OJqfrpoCohp/OzoWmhrGo8PN6ZRpFun0T19R5bR0lsE2Wpks00uJzgT0Ul7zxX4wl+R4k0/BnNrMtB9/3t1rM0mKYqXVhRy5BnJMI5qsoZsP/8GNEn8jUlyrfuMPW5CCB4PKhDRP6MkDk3xvsk3QFtLee0THuSHhyZs2NEMVMe0ecXdW4RWnM5PJYGg9+J4vNelsGyMmWSdTHk53VsS5VDTkOyDCp/KMr6M68sLZ7E9DW8tqqH1/lYhWSqCErSHQTVIjqVZZWdv6vPovvczuzfojqnyn8XqnSg0ptlZa/Sl6nxaZU7FYlPjXPyyPMhrTLz/iq3XvF9GaqurWRhWXqUnLnRBr1YtYTNRK9M1M3BqCYyc8D5yRxlVScVhTPmLXS2tbWlw8PD0tAzCmk5lHR07Fv23/dWDwiZqpRSGopVkPP6RN/ObxOhpmEnKsl5hSqaqDYySqdWKXSmVIjM6WArhERUXo1trnyQxgg/kXslY/5U6Uf3nfJNuXQ5b4PkFtH29usXZDvSYvmcrKTMZh/93zxgGsT9J89sXFMOUk+yDM7TePxSJxIUmciflhyZj7k0keNRva6R/ZfuRn8VP92m++rqozzOX+WXknE+5u9EdNlJht6V8TaZKRRGtoHlmTgJVtXjUMt7RPgel8NQkFQxPoUnFT6pFfqzP3Q260LvTgnlXEPVxzQCbmv1yfFgv/i71a+WIXA72DbWZXmsxozncxKO/ZPqCGuZuYMWaFgH+p5qw/7+/iidZSKvMwKryiGoWTbKqeQhy/S1Ld1ImcmIJ+9LveV1BIfcOoB9S+fN9iawZF2sm225z6P+D27MvWc4Z8UTXeUDEM5r+Xobx0oppgQl0UCF4pLZnDQhsqZw8LVlJv/m8kdHH1zWRWPNNhIhJLOpFMw3tvq1bgPQ9/3oMeW+74cx4V4ofESbhpBGjpO9Uw7ZlAiOZVZO3OeY4++69sZZ6QQYlme5STxOvkwZ+sohVyDASP4+yn5fct1eBpvn3AePQRpCynKFyKsyWznsCuRJY0PLa3JehPdl2WxLykzyw3KQc1x5PcuiHakmu1NWTfn8wyJ6lL1ZuFyotV6cDLIweJG/lZ7eMBnJgZ1iIr9Nvsahm6ThibJ8IrBCkSwjjblTTDc3N6P9aNhXC3WVJrESc801jSY3CKuQxzoM+3w+H96zyYk8t0u6dTbcqD+NICeCPUY0EuwDJxNp8KrIhBEO20d5YJiehj5l0YgpDVqldJUxZx9MLIupmypnzaVs+VzDqsk8SQPMujgByL4yuslVJu4ny6qcZ1LlFOmIaUdoB9I+WJc5pnRKHGfzMCdaSZkN4H3sq9tJ3c9+EqDwukX0KGmWKnTi96JwrTUp1Zo0zPyzNFYyrtH2uby+qnMKPWdZKSjZpxyHlmD4dwoax6VyZv7QCK2S2EcqJetyKiYdbI5F1f8cm+qeCs1l23jdfcYhkZl/t2RuyrBX8t3qrz9VKoO8Xgc54uIukFkXdWBKrxNYtaKZrCOdcyUj6TwqGUq+8WGsyhCT0tlU/K14xAhsiketeqvxmKKljXnXdduS/p6k3+n7/k90Xfejkv6GpI8k/X1Jf67v+8sFZWh3d3dAXQwT6flbqRN7KntWe/2tra3Ry55psLquG+oj6vIA58M/RLld1w37ORBZEGFLtxvkmzwhKN0yKoXHbaMBlm4n1BgVWBB9rZWLnpuRCqMfzg9sb2/feZ3WKvjqPvNxaSpACniicSobnS13EOQEcpW+aBmVlCnKDMsiuq5AAud4jM4rR+S2tFZpteTb5xIo5OQtx8r8bJT5zny9ubnRixcvRm3MMXcbsq+JrBOZpiNkn1lfPrtRfapIKyN80hRI4jnqVMtZZTlE+dV4UC45HhURtC1D99mo4y9I+nX8/28k/bd93/+rkp5L+vP3KKtJU4hpGTRXXd8qh/e3UAVD+UqI0ojwPlMLkeX1LeRRofQK5S/Tb7brDa2Nr9WYmhahjSnHzvtbKKl17m0px79CUy1UmmXkvVPIrDV+i85pBXw1Mq92C2Qb6CAznUJ+MX1RGTGmayrZnzLQVdtbZVX1JFV6XunmsnQfPXjbepZC5l3XPZP070r6ryT9x93rmv6opD/75pK/Jum/lPTfT5Uzn891eno6PCBEY0TkNJvNhslOU2XMyNQM79NDcmKN9ebElxFPIuZKEYnuKpTB9hCtGY3mbohd141QjttXvVUpEaDfZOMoge0kGuCSqlXxVbq7lawpc/9bW1tlyJ5ILHOg5AHz5UTrLDMnnH0sqULtiZ5yDH0fJ7lNrq/rxg+4GJkx4qqAQI4LEW71QEzl9FbFVxtz5or5qSIk3st2Ul9zsjLvqRA6+e4xoY6kTLWcjlTLQZZZAZFqGSozC5T/lhHOx/9pHzI/zlz7MrRsmuWbkv5TSU/f/P9I0qd939va/rakr1Q3dl33dUlfl6Rnz56NDDk75/8O29J7L4OMqlUGnLiq7qXy5eQGGZlrh32+UiauyuCxnDShoXI7K4XNui0EKWjeNtX7ySTRsL2hb2oFfH369OkdZFb1MR23x5DUQn/sP4875TF1jvzPtkwhbJ7P+9i/Sh4LAzvc7zZR5pPSYFcGrmrPG/qmVsTXXLtN8JH5+mo+KtvncWk5Vo51GjZ/T8lMK5LLhRJVhEeDXkWx1p/cZIyGvJIV/k9QReSf+7ZId1f2TdFCY9513Z+Q9Ht93//9rut+ZumS31Df99+S9C1J+qmf+qkexyXVD7aYOZk3fNOekcHNMCoNNwWHTKShz9xYtWkXEVa2PZGZ78l+xbjcMSjMO5NygpZjQYOQjshGnVEGxmJlfP3hH/7h3uNDXmVbq/Gggklj504EbgVIvpBniapSDnw+VxKk86lWxqQTJ1Gm0iCwXPYn76+iD/OXwID3O2pj/1bJ148//rjn1hutVGBuLVvJqnnYmhzPcXA5jEr8P/WNUXyi/sp2VOi5Qt3pGDIC6Pu78zqsO/vHfnI8U2fY1tT9KVoGmf8RSf9e13U/K2km6VjSX5H0Ydd1O2+8/TNJv7NUjdHg3CUtjSWVko/JOx3hTeVT2YiaK6RoY7GzszPsrmiFo5JYgLw0kQrrsrjDIZWRHphKbmK7fT7XMbvf0q2x4DptC1e+uqvrbneN9L0hNCvlq41gtVqGwpvps6ScAPZSy/n8didLKiq3Ak6FItLhBDR5V6XpPJae4HbZVG6+6Nj9yrQDc8kulzzL9lZ6QOfs+60zfX/3lXWr5OvW1usnmv27cjo2aBzTNPg+V4GXrM+UIIxtSGfr9whQj1pL/xJFsx3UV7bDlPrK85n6qwAgZTeXv2b/7fiWNeTSEhOgfd//fN/3z/q+/xFJf1rS/9b3/X8g6e9K+vffXPY1Sb+0dK2L65RUe7YWKuJ970JVaDVFFQpdlrLsCtlM1Zf3TY1X0a+V8/W+/W/R2/Cxuqdqz6KJs1X1YarcZep4W1leB18r5EiiE8t2T+nr21Cid1OCpPuUlSCLE6SLKKO3Ze5x3VP/34a6+wjNm7DtP+lfL3X6Q3q91OmLkv6BpP+w7/uLBfd/T9KJpE/etsGPTF/S+9H2f6Xv+y/7xIavG75WtOHro1K2fcTbiu5lzFdBXdf9vb7v/7UHrXRFtGn745W/Ttq0/fHKXyf9QWv7fdaZb2hDG9rQhj6jtDHmG9rQhjb0HtBjGPNvPUKdq6JN2x+v/HXSpu2PV/466Q9U2x88Z76hDW1oQxtaPW3SLBva0IY29B7Qgxnzruv+eNd1v9F13W92XfdzD1Xv21DXdV/tuu7vdl33j7uu+0dd1/2FN8e/2HXd/9J13T998/2Fx25ri7qu2+667h90Xfcrb/7/aNd1v/pm/P9m13V7i8q4R12fC95u+Hrvuj4XfJU+/7xdBV8fxJh3r7fj/O8k/TuSflzSn+m67scfou63pGtJf6nv+x+X9IclfeNNe39O0nf6vv8xSd958/+zSg+1y+Xnibcbvi5JnzO+Sp9/3r47X/MR1XV8JP2bkr6N/z8v6ecfou4Vtf+XJP0xSb8h6eM3xz6W9BuP3bZGe5/pteD+UUm/IqnT6wcQdip+/EHl7Yav7ydfP2+8XRVfHyrN8hVJ/xz/m7u2fdao67ofkfRTkn5V0g/1ff/dN6f+haQfeqx2LaBv6vWueX62eOld896CPpe83fB1IX0u+Sp9Lnn7Ta2Ar5sJ0Anquu6JpL8l6S/2ff+C5/rXLvMztxSow655j92Wzypt+Pr+0ueNt6vk60O9A/R3JH0V/99ql8WHpK7rdvVaKP563/d/+83hf9l13cd933+367qPJf3e47WwSWvb5bJBnyvebvi6NH2u+Cp9bnm7ul0u19rMW/o1ST/2ZoZ2T693c/vlB6r73tS93sLsr0r69b7v/zJO/bJe7zgnrXinyFVR//C7XH5ueLvh673oc8NX6fPL25Xy9QGT/D8r6Z9I+n8k/RePPemwoK3/ll6HY/+XpH/45vOzep3L+o6kfyrpf5X0xcdu64J+/IykX3nz+w9J+j8l/aak/1nS/h803m74+n7y9X3h7bvydfME6IY2tKENvQe0mQDd0IY2tKH3gDbGfEMb2tCG3gPaGPMNbWhDG3oPaGPMN7ShDW3oPaCNMd/Qhja0ofeANsZ8Qxva0IbeA9oY8w1taEMbeg9oY8w3tKENbeg9oP8f0vFZ29GoceYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "file_img = Path('data/dot_vision_circle/SERIAL_19061245_DATETIME_2020-05-17-20:48:11-005188_CAM_1_FRAMEID_0_COUNTER_1.png')\n",
    "bb = torch.FloatTensor([[520, 975], [561, 1016]])\n",
    "arr = bb_array(File16bitImg(file_img).array_gs(torch.float), bb)\n",
    "arr_dx, arr_dy = grad_array(arr)\n",
    "fig, axs = plt.subplots(1,3)\n",
    "for a, ax in zip([arr, arr_dx, arr_dy], axs): ax.imshow(a, cmap='gray')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.169672Z",
     "start_time": "2020-08-06T18:56:10.051600Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fa81749e828>]"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPsAAAD6CAYAAABnLjEDAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAKu1JREFUeJztnXl0XfV17z9bs2R5kDzItuTZso0NWDbGYErADCGEgIGVocBLX9omIVlt30teu94LaVZfx6yVvtUm7VrNSh5taEiaQlISJtctBsJgDLaxsbHxPI+SJVmWbCRb4+/9oWuefO/+yb+rK92rq7M/a2lJ+umc8xvO2Tr3fM8exDmHYRgjn5xMD8AwjPRgxm4YEcGM3TAighm7YUQEM3bDiAhm7IYREVIydhG5W0T2isgBEXlssAZlGMbgIwN9zy4iucA+4OPACeBd4GHn3C7fPuPHj3dVVVWXtWn95+To/4N6enoS2nJzcxPatGOKSFCbb3+tTeu7u7t7wPv60I6pjd23bhq+uceTqh+G1k/o+obu69s2lFTObTL7d3V1JbQlc85CbOXYsWM0Njaqi5EX3FMiy4EDzrlDACLyNHA/4DX2qqoq1q5de1lbe3t7wnYlJSXq/hcvXkxoGzNmTEJbR0dHQlteXuJUtTaAzs7OhDZtobVxtrW1BY2ntLRU7Vvjww8/TGjTLhLfumkUFBQktGn/TLW10AzL989LG6dmNNoaaWPUDMbXj4Z2HrV5a+fn/Pnz6jG1sY8bNy6hrampKaGtsLBQPaaGNs74/W+55Rbv/ql8jK8Ejvf5/USszTCMYciQC3Qi8qiIbBaRzdp/NsMw0kMqxn4SmNbn96pY22U45x53zi1zzi0rLy9PoTvDMFIhlWf2d4FqEZlFr5E/BDzS3w4ikvBspT1D+YSQ4uLioG1DxRpNA/DtP2rUqIQ2TW/QnjO1Z+mioiK1b+0ZOVTs0cjPz1fbNW1Be37U9tfGo62FD+0ZV5u3pqn4rg1tnNocR48eHbSdNh7t+gN9PbRrSBujtuah5xYSn+P7E1QHbOzOuS4R+QPgJSAXeMI5t3OgxzMMY2hJ5c6Oc24NsGaQxmIYxhBiHnSGERHM2A0jIqT0MT5ZnHMJ4oMmaPmcXTQ0IUUTxDRhxydeafjEvHg0xwetraWlRd1fE3a09dBEJW0+PmcT7ZiacKYR6pji6ydU0NLG43Mc0hyPtG2186gJb5rg6BNVtWurtbU1oU1bi9Bxg77G8ee3P3Ha7uyGERHM2A0jIpixG0ZEMGM3jIhgxm4YESGtajwkqoWaWuwLl9RUT01J1ZRQLRT2woULaj/aGwKtb03Z1RTXZEJCQ92Hy8rKEto0hd6nkocou76+Nddhbc19xwxdX+1Ni6+f0LcL2ry1a0jbrrm5We1bu7Y0l1dtPKEh2qCvUTJvruzObhgRwYzdMCKCGbthRAQzdsOICBl3l9VcPDVhBsJdKrXtNPHK50qquUpqQkqooKWNJ5m4bG082rxTTWIZmg/t3LlzCW0+oUjr5+zZswltmuin4RMctXULTWaqCbVaPz53We38asfUhEBN5PXFs2v9J5MQ1O7shhERzNgNIyKYsRtGREjpmV1EjgDngW6gyzm3bDAGZRjG4DMYAt1tzrnGkA1FJEFI0UQlnxCSjMgWjyZk+ITA0Pjm0IolmtjjE2E0QSs0Ft8nXmlonmihop+25r75aLHZ2rprx9TaUo2b1wo9hAqTvmst1LNN21/r2ye6JSNma9jHeMOICKkauwPWisgWEXl0MAZkGMbQkOrH+JudcydFZBLwsojscc692XeD2D+BR6G31pthGJkhpTu7c+5k7Hs98Cy9xR7jt7GKMIYxDBjwnV1ERgE5zrnzsZ/vAv6iv31ycnISxDdN3PCJPaGhiKGlnZMpl6uhCSaaAKRVBPXNUQv/DA0J1cbj60fzCNTEHk2Y1MQnX7UUTVTVBFjtnGmCoa/yTKgXmyaIhYYg+7weQ6sGa9slUxFJm492bfhI5WN8BfBsbKHygH91zv1nCsczDGMISaX80yFg8SCOxTCMIcRevRlGRDBjN4yIkPYcdPFCjCYg+UQHTbjQhB1NcNFCCUM9p0AXUrRjaqGa2nx8QpMm7IQKjqHrA7pYpIlsoeWZffn8tPmEejNqAqovFFabjyYE+sSvEHyCbmios7YWmqjqqz6Uanlnu7MbRkQwYzeMiGDGbhgRwYzdMCKCGbthRIS0q/HxymUy9dlDXTdDkyn6VP/QGGHNDVZTpbVa7L462qHuk5pKnkxFmFRi9rW+fVVMtPlo22qKdjI15EPXLXSOoZVsQI/Z1865ppxr4/b1E5pc1Yfd2Q0jIpixG0ZEMGM3jIhgxm4YESHjJZtDkzaCLrJpYlpojG8yroa+UsHxhApaySSH1LbV5h0a2+/bX1uP0ComPtdjn1twSD8avnOrCWKaK2poLH0y16UmBIZeg5oQ6OsnpEpNf4Kd3dkNIyKYsRtGRDBjN4yIcEVjF5EnRKReRD7o01YuIi+LyP7Y97KhHaZhGKkSItD9GPgH4Cd92h4DXnXOfUdEHov9/o0rHUgr2RyaTBHC45s1MU0TN3yeepoXnNaP5tGkiVzJiDBa36FVZpJJoBlapjgVr0Uf2rbaeLT5+AQozYstdEzaedQETK0P3/7adR1a4cZXejuZ5JIaV7w6Ynngm+Ka7weejP38JPBASqMwDGPIGegze4Vzrjb2cx29mWYNwxjGpCzQud7PmPpnUnorwojIZhHZfObMmVS7MwxjgAzU2E+LyBSA2Pd634Z9K8KMHz9+gN0ZhpEqA/WgewH4AvCd2PfnQ3eMFySSEZU0wSW0ykZoG4SX0dWEKq1NE+18QpMmVGnbagJSMiV9NS+40GSXoecB9HXTznlo2Ksv4aQvSWM8oUkbNTHMd85CyzOHrmWq16WPkFdvTwHvAPNF5ISIfJFeI/+4iOwH7oz9bhjGMOaKd3bn3MOeP90xyGMxDGMIMQ86w4gIZuyGERHSGuIqIglihhaaqOV2g3BPstByxj4PK01wCa2MEhpy6wv9DK0koqGtj1ZC2rdtaFUWTUT0eX2Fei5qa5SMl+CHH36Y0DZ27NiENi1Pnya8addLMtV1tDmG9uMTAkPKWvcXOm13dsOICGbshhERzNgNIyKYsRtGREi7QBcvxGjCji83XGjRAE3g0Pb1lRkOLfkcms8stIABhIfDamvk82LT0Pb3rUc8mrjoE+g0NCFQy0EXmnMQdM+6ULFTO99am0/80sYZmqcvVHQOHafloDMMw4zdMKKCGbthRAQzdsOICGbshhER0qrGawknQ2OJfWjqo6Z6hiZyhPAKIaFVXTRVOJk4fk391vrWXI99brnattoanTt3Lmg733xCz+WYMWMS2srKEpMW+yrHaPPUFHqtLfQNiu98a28itOSUoaWhfS7O2tuW+HUzNd4wDDN2w4gKZuyGEREGWhHmz0TkpIhsi33dM7TDNIY9PT3ktbd/9JWTRNEIIz0MtCIMwPecc3+TbIchCSd97rKhrrWaS6QmzDQ1xde+6EUTlbQ4aE2w0eKqNXzupVo1G03g0/rRhB1NiAOor09MCHz8+PGEtubm5o9+Lr1wgWmNjVScOsXkc+eY0NrKhNZWRrW3U9LZSU6cqHUxL4/WwkKaS0poLC2lYfRojpeXUzd1Kg1jxkAfMWny5MkJfU+YMCGhbebMmep8QhN9hsapa2uuiYignzOtH01c1K5/TZgE3VU4/pz3VwUnJAfdmyIy80rbGSOPgo4Olhw8yIKTJ5l/8iQVff5xNBcX01BayqHx4/mwsJC2ggIu5uXhYgac393NqI4OSjs7KW9tZVZjI8uOHCE39g/hw8JC9lVWsreykp0zZoBi7Mbgksqrtz8Qkf8KbAb+yDl3dpDGZGSQ3O5uFu7bx+Jdu5h3+DD5XV1cyM9n/9SprL/qKo5NnMje4mJaldebGn3vuHnd3Uw9e5a5584xp66O+SdPsvTQIVi3jpOVlexctIgdixfT6slUZKTGQI39B8Bf0lsJ5i+BvwV+V9tQRB4FHgWoqqoaYHfGUFN69ix3v/Ya1+3YQWlbGy2lpby7eDEbKis5WFFBT5/HDu1jawhdubkcmzCB+qoq3l64EJxjUnMzNYcPs/zYMe5au5Y7Xn2VPQsWsOmGG7igfIw3Bs6AjN05d/rSzyLyj8DqfrZ9HHgcYMmSJeEZ7Y20MKa+nuvWrmXepk3Q08Pu6mreralh/8yZuJycy57ZBx0R6svKWFtWxvZ77mFCQwNLN29m8bZtLNq5k1NvvMGWu+/m+FVXXfZ8bwyMARm7iEzpU9jxQeCD/rbvs1+CkBJaCcSHJnSFxr37vLtCEy9qIpnmOaXFNvtESK1vzYtNO2ZdXV1C2+HDh9V+GvfsYdW2bazcs4eenBxemzePn0ycSH1xMbS0wPvvA6jGrt3Zk8lBoCWCPHbsGABrxo+ncOVK7jp2jM8ePcp93/8+u6dM4Rc33MDx8eM5ffp0wr4AM2bMSGibNm1aQpsmymrCm3YefR6XmvCmrVFoLL2vuo12rSdTEeaKxh6rCLMSmCAiJ4A/BVaKSA29H+OPAF8J7tHILD09LH3vPW5/6SWKOjt5Y948nq+p4VxJCfWnTmV6dAC05+by4qxZbF2+nFv27OG+rVv5k+ee480FC3j7/vtpD9QLjMsZaEWYHw3BWIwhpqypifuef56ZR46we/JkfrZiBafGjcv0sLx05eby60WLeGfuXFZt3codu3ax5Ac/4MVPfYr91dWZHl7WkdZAGCNzLNqxg3tffBEnwgurVvFceXnWPAdfKCzk5zfeyKY5c/jShg18/qmneOeGG3j5jjvoDsypb5ixj3hyurr45Jo1XLdlC8enTeNXn/kMLePGwZEjmR5a0hyeOJH/+6Uvcdcrr7Bi40amHz/O05/7HOc8zi7G5aTV2Ht6ehJCETWBweeppAkXoR5RWpvPiy00waPmlRf6WsonwmgCnXbM2trahLZ9+/Zd9ntxWxv3//SnXNvUxL/NmcNP582jJ7ZNQ0ND0DG1cYaG9vrQvPdGjx6d0KaNsbm5me3TpnF9Tg7/bfNmfueHP+Svb7qJA8q6aQlBtXOutYXuC+Elm7UwXG07LUQbwiodWYhrBBl77hxf/elPuaq5mb9ZvJgnFyygJ4m3HMOddysr+dbKlXTm5PAXb7zBvIMHMz2kYc/IOfvGR5SfPcuj//IvjG5t5Vs33MDrI9SZ6fjYsXzztts4MXo0n3/mGRbu3ZvpIQ1rzNhHGONaWvjSz35GQUcH//TII+wqL8/0kIaUc0VF/PnHPsapyZN5+NlnWbB/f6aHNGwxYx9BjGpt5XefeorCzk5+9PDDnIpIcElbQQFPPPTQRwY/K+agY1xOWgW6nJycBK81TVDQwklBF680MSO0hK7P+yi04okm2oUKWr5QRK1vTaiK94zL6+zkN3/yE0afP8+frljBnjNn4MwZjnhU98bGxoQ2bd21NdLOmW8+2rbaWmrrFtp2qf9v1dTw7bfe4pGf/5zHbr2VfYqHpCboauPRwoV9FXe061LrJ7T8tY+QY5pAN9JxjlWrVzOvuZnvXncde8aPz/SIMsL5wkL+YsUKOnNz+eaGDRR5/jFEFTP2EcANmzaxeMcO/nXBAjZOmZLp4WSUxpIS/s/y5Uxsa+PBZ5+FJHzHRzpm7FlORV0dd73yCnurq/nFvHmZHs6wYPf48fz46quZv28f17/7bqaHM2wwY89i8rq6+MyvfkVbcTHP3X9/1ri/poM1s2ezr7qau15+mYmK5hFFzNizmI+tW8ekxkaeW7WKNiUkM9KI8PyqVXTm53Pfiy8i9nE+/b7x/amFV9pGUyM1V1JNXQ2NUfftr6mmoTXJte18bxy0+agJJHbv5ub169k4dy7rS0uhro6TJ08mbKYp+RBeIUd745BMDHVobXrtbYmWL8DnqqvFn+8aN46nli7li+vXM/2VV3hz3jzGKVF+mtvz+CRETu261Oatbae1+WrQh1Yq8mF39izlsxs20J6Xx7+tWJHpoQxr3po7lz0VFXxmyxaKPa/OooIZexZStXcvVx8/zpqlSznvqX1mxBDh6eXLGd3ezic/CEqoNGIxY882nGPFCy9wprSU1xcuzPRosoKj48ezYdYs7tq5kxJP0cQoEFIRZpqIvCYiu0Rkp4h8LdZeLiIvi8j+2Hc9s70xqFTt3UvFsWOsWbKELkvcEMxzNTXkd3ezeN26TA8lY4RcLV305oV/T0RGA1tE5GXgt4FXnXPfEZHHgMeAb/R3oO7u7gQBSosb9rleakKXJnBoooW2nc9VMTReW3Pd1MaoJYz0uV6ePZuYfv9Un9xwn/j3f+fcqFH8vKCAzjhXWE2MS8bFM7Q6j7a+PtFOE940tHEmc860GPm+lVVOA+9UVLD0zTdZc801dPRxpZ0+fXrCvto58+VZ0K7XVERizSYgrLR0f+LpFe/szrla59x7sZ/PA7uBSuB+4MnYZk8CD1zpWEZqTGxspPrIEd5ZtozOwDcBxv/nmVmzKG5vZ0lEn92TemaPlYFaAmwEKvqkk64DKjz7PCoim0Vks6+2mhHGsvffpzsnhy2LF2d6KFnJ3nHjqJ04kWXbt2d6KBkh2NhFpBT4JfB159xln3Fc72cH9fODc+5x59wy59yy8hEeWz2U5HR3s3THDnZXV/Oh8l7YCECEzYsXU1VXx2RP/vmRTJCxi0g+vYb+M+fcr2LNp0VkSuzvU4DEhyZj0Jh99CijLlxg69VXZ3ooWc22hQvpEeGaCGa1CSkSIfTmid/tnPtunz+9AHwB+E7s+/MDGYAWj+5L2qgJFKECUqgABLqQElplxucZF49vjlp7U1MTc7dv52JeHpvGjaOzqYkzZ84E7evzxtII9ZbTBCmfN6K2RtoxQ0su+zzGNOFOK1d9RoSDlZVctWcPL1x/PaCLcZMmTQoatw/t2vKJpfH4Ek5q5zJ+jfoTl0Pu7L8B/BZwu4hsi33dQ6+Rf1xE9gN3xn43hgLnuPboUXZOn06nvW5LmR1z5jClqYkJQ1nHbhgSUhHmLcDngHvH4A7H0Jjc3My4tjZ2j9DEkelmT6wu3NwTJ2gcxhVxBhvzoMsC5sfes++dOjXDIxkZ1I8bx7mSEuaeOJHpoaQVM/YsYNbp0zSXlNBglU8Gh9hz++xhUsgyXaT1ATAvL+8yrybQwxh9HnSh3nKakKIJSFpSQd/+mvChleoNrejiE2u0basaGjhaXk5Hn2OHJrZMRlTyrXs82pr7hCFt3UNDNZOp4qMl6tQEukti3MExY1iyfz+dDQ3qNaitr6/Et+bxpo1d219bH19SzZCwWUs4mcXkdnQwpaWFYxFNIjlUHI/ddKYp7skjFTP2Yc6YM2fIdY5TY8dmeigjikulqiuU124jFTP2Yc6Y2Pv0RqXwoTFwzo4aRbcIE5SP8CMVM/ZhzuhYPMGZ0tIMj2Rk0ZOTw9mSEsZHyNjTXrI53sNMEy18Yo/mURWa303bTss95utHE0dCq4ZoIoxP7Invpygm2NU7R3cfLzFNCEwmH5lGaIhraJUYCA//DA0rTiY/m9Z337bzRUWUXLyonjPteL7QU83jTds/NN9cqecfuybexo/J58kIdmcf9hS2tXExL4/uEVRuebjQlp/PqAjlpbMraJhTcPEiF8xFdki4mJ9PURKxA9mOGfswR5yjx4o/DAk9Il4/8JGIGXsW4MzYh4wolY4wYx/m9OTmkhsoXhnJkdvTE6lPTWl9GHTOJSikWhyzT/XU0NRizY1Vi3fWXCxBV5s1JVSrQhLqautz1Y0f04W8PEq6uhLUe03ZTSbOXNs2lTjzZAh1g03m7YI2Tm3ufdtKOju5UFCgXi+aQu9zJw7NGRBaOcZ3bWjjjB9TSgknjczSXlREYVeX3d2HgJKODto8r0BHImbsw5y2mOfcWE8KZWPgjL1wgQ+T+BSZ7ZixD3POxQI2ouTplQ4KuroYe/EiDRHyTEylIsyficjJuFRVxiBzydgnmrEPKpd84qMUc5BKRRiA7znn/ia0MxFJECk011Zf5Q1NfNDcTkOTHPpiyjWxSBMNNQFJE2tCxTBIFP06qqroyslhVmsr7/dx79XEnlCXU1//oeJiqDjoO2Z/Lp190c5DMmW2tWvjksg1K1beumnSJCYq89FKO/uESe060Mapxc2PVv7ZJJMbIN5++rsGQnLQ1QK1sZ/Pi8ilijBGGujJy6O2rIxKJZusMXCmNTbSnZNDbXk5EzM9mDSRSkUYgD8Qke0i8oSvsGPfijBa+mPjypyYMIHpDQ2QRNYZo3+mNzRwqqyMrgiV0UqlIswPgDlADb13/r/V9utbEWa8ZVsZEAenTGH0xYtMjlBWlaEkt7ub2XV1HJwyJdNDSSsDrgjjnDvtnOt2zvUA/wgsH7phRpt9sayy8yKWIHGomFFfT2FXF/sqo/U0OuCKMCIypU9hxweBK5bGdM4lJNPTRAefeKV5vGniiNaWqtijoQkumieXNh9f1Q9NnDwxfTpnx45lcW0t22++GdAFJG1frdoJhItKIVVIwL9mof1oaNv5vOo0ATU+uSn05jCoqaujR4STc+cyqrhYzWugHc8XZ56MyByPJmz6hEBNfIsXdPtb2xA1/lJFmB0isi3W9sfAwyJSQ28swRHgKwHHMgaCCLuqq1m+dSsF7e10eP5RGGFcs28fh6qqaFPcT0cyqVSEWTP4wzF87Jo/n9/YvJn5Bw+yY+HCTA8na5l05gwVTU28XVOT6aGkHfOgyxKOVFXRPHo0S3fsyPRQspplH3xAtwg75s/P9FDSjhl7luBycnjvmmuoPnSIsRFKfzyY5PX0sGznTnbPmcP5CNa4T3u+o3iBRfNy8gkUmvilCRyh5XJDhTjfmLR+tO00Yees5zWaJrxdEpoO3HILt73zDrfs3s2WiYmuII2NjQltWmgvpJawMlRgg/AqNaGiqi9Rp7ZuU+Nq491w5Aij29rYd8stVFRUfNSuCXTa+vgqtUyYMCFof02U1a4hX4i3dm3FnwurCDNCaCkvZ/+iRdS88w6FKcaURw7nuG/XLs5MnMihCH6EBzP2rGPTrbdSfOECdx46lOmhZBVX19Ux8+xZ3l25EiKaqTeas85iamfM4OjcuTywZ0+kMqOmgjjH57Zto7GkhN1Ll2Z6OBnDjD0LefOTn2Rsezv37duX6aFkBTccPcrspiZ+uXgx3RFOy53WmYtIgsCiiR4+AUgTODThLTRHmk/M0DyiNGFI86DTRCVNhNHy14EuXpWXl1/2e0d5OVurq7l/3z52L19OU6zooybQ+cJ4tbn7BKiQfX1ej6G55bRzrglVExVh0tc+depU8js6eOTFFzk1cSKHbrqJhcp2WphpMjnxtLGHhiBr12oy4cI+wVLD7uxZyvMf+xhOhE+//rpFw/XDXevXU37uHM/deScuos/ql4j27LOY5tGjWXPjjSw6coSle/dmejjDkqraWm7esoWN117LkaqqTA8n45ixZzHramo4NGUKn33tNcrM0eYyirq6eGT1as6VlrLm1lszPZxhgRl7FtOTk8PPPvEJBPj8Sy9ZuulLOMd/37+f8pYWnrr3Xi5EKINsf5ixZzlnxo7lF7ffzpxTp/i9AwcyPZxhwYMnT3J3XR2vrlhhH9/7kHY1Pl651FxgfXXTNZVcQ1MzQxMn+tCUWE0J1VRcrcKHLzZaG/ukSZMS2mbPnv3Rz82zZ/NeezsPvvYaF6urWb9o0Ud/873Z0No1F97QBJo+BVlDWyPtTYumsPvU+IWxSMB5x4/zewcO8P6sWRz6/OeZEzdPbX8t9ly7XnzKtzZ2zU1Z2047pi9ppKbGx6v5VhEmAqy/7z52Tp/Ob77xBtdG1Ltuen09X/6P/+B0WRk/ufPOyHrK+bDVGCG43Fx+9IlPcHTSJH5n7VrmnTiR6SGllclNTfzeiy/SWlTE91etoj1CZZ1CMWMfQXQUFPDDe++lYexYvrp6NVcdPZrpIaWFuefP87Vnn6U7J4d/WLWKlgiGr4YQUhGmSEQ2icj7sYowfx5rnyUiG0XkgIj8XETsX+kwoLWoiL9/4AFOl5XxlTVruPX06UwPaUi5urmZ723bRmdeHn/34IM0KqGuRi8hAl07cLtz7sNYltm3ROQ/gD+ktyLM0yLyQ+CL9KaX7pd4oUsTyZIppayJPZrAoble+pIxauKVJkBp4qI2H03c84kwY2Our1eiSlGZ+/azZt487nv8cf7kgw+oKS1l7XXXQZ+/a/Hfhw8fTmjTRFHNxTOZajSaq7AmWGpi2syZMz/6efnu3Ty0fTtNY8aw5utfp7y8nL6OxZMnT07YX0tCqbVpwqQvzlwTYLVrI9TN2Ceqhrpx+7jind31cskJPD/25YDbgWdi7U8CDwT3agw5HSUlPPf7v8+m+fO5b8MGfuellyjy+MlnG7nd3Ty4bh2ff+UVDk6dync/+1nOx8UPGIkEvXoTkVxgCzAX+D5wEGh2zl36F38CT0koEXkUeBRg2rRpqY7XSILu/Hx+eued1JaXc++GDUxraOCf77qL432ytGQbE1pa+MJ//icz6ut549prefbmm+mJUFWXVAgydudcN1AjIuOAZ4EFoR045x4HHgdYunSpRWykGxFeue46Dk6dym+/9BJ/9MwzvLp0Kf88bRqdWWQkOc5x7+HD/NbevXTn5vJP99zD9jlzMj2srCIppxrnXLOIvAasAMaJSF7s7l4FnByKARqDw+EpU/jOQw/x4Pr13LVlC4v27uXJa6/lvcmTL3uWH47Mb2riyzt3Mq+5mQ9mzuQXK1fSHKFSy4NFSEWYiUBnzNCLgY8Dfw28BnwGeBr4AvB8SIfxgpomOiRT9UOLwdZED0308wkhoYkXQ8VBTXzyeTpp66EJWppXnU9A6svma6/l9J493Pz00/zx229zcOZMXr79dk5OnaoKYk1NTQltra2tCW0+DzpNQNLix7U6gHN7erj9179m0a5dnC8t5Zef/jQNd9zBVBGmxm2rrYevIkw8ofkPfOdMW/dQwVK7LjVxD3TxN5nknyF39inAk7Hn9hzgF8651SKyC3haRP4K2EpviSgjCzi+YAH/8JWvcP2WLdy2bh1ffeIJ9s+ezepFi3rrymX4Tl9VX8+dmzdTc+AAnfn5vL5yJW+vWEFnYSEVw/xTyHAmpCLMdnrLNMe3H8KKOWYt3bm5bFi+nK2LF7N8yxZu2riR/3HoELVlZaxbuJBNc+fSmsbySIUdHdTs389NH3zAzNOnuVBQwFs338yGG2+kzRNHYCRHdBNyGQC0Fxay7qab2HD99czcuJGP7drF59av59Nvv83eykrWV1TwfkUFDUNgcOPa2ri6ro5lx4+zuLaW/O5u6srK+NUtt7DpqqsojlhJ5aHGjN0AoDM/n3cWLOCdBQuobGxk2cGDLD14kC/HfOxPjxrFngkTOFRWxp7iYk6MGkVLQUHYR37nGNXZSWVrK7NbWpjf1sb8+noqY05NTcXFrL/mGrZWV3Okj2AYrbKLQ09ajd05l5AAURO5NAEHdGEotOJIMpU3tCSN2jG1Ni20URP8fGG8mqCljUebjyaSVVdXq/1oglhzc/NHPx+67TYOOUfe/v3MPHiQmQcPsuzYMW7t42/fUVDAubIyLhYX01ZYSHNZGWs/9Sk+sXo1Y1paKLp4kaILFxjb3ExRHyHqYnExp2fM4K3bbuPEvHk0VFZSUFRECdC3ZKUWEjrFc7fXREzt/PrEr3i0c+u7XrRrWDu/2nnUxq1570FY2G1/4rLd2Q0/IpypqOBMRQVbbroJnKP0/Hmmnj7NuDNnGNvUxJjmZora2hjT0kJB7GKeXFtLQXs7F4uLaS4r49jMmTSXlXG2rIzTU6bAjBkZFwGjiBm7EY4IH44ZwyHFNbXvK6Qnv/zlfg8z1gw9I1iIq2FEBDN2w4gIaf0Yn5OTkyC6aN5LvnLGmkASUsYWdPHKV7JZO2ZoyGOx8m5a8/LTxCfQx66NR2uLrxwDfs82LZRWE0A1wVHr21d5Rlsjbe6aUKWtha+SjkZopRetLbREN+jXRqhQm0xosLZ/Mrn/7M5uGBHBjN0wIoIZu2FEBDN2w4gIaRXoenp6EkL6tBA/LWk/6CWSfcUW4tE8onyeSpoIFBrWqYlK2hx9nk5auyYgaZ52mqikhXmCHkqrrbvWpo3RN5/QMtChZa2TCaUN9aTURLJkQmFDBT5NxNQEXe1a840pmRBXu7MbRkQwYzeMiGDGbhgRIZUiET8WkcMisi32VTPkozUMY8CkUiQC4H86557pZ1/DMIYJIWmpHKAViUgaEUlQTUPdSyG8xLKmfmv9+JTdUNVUc4nU1FFNEfeVnw4tDR2qxmtjBF3Z1dZDc3dNpoqPNs7QNw6hpZ0hXHnXzk8y16CG5nqsrXtoDgLfHEPe6qRcsllEckVkG1APvOyc2xj707dFZLuIfE9E9BEahjEsCDJ251y3c66G3vzwy0XkauCb9BaLuB4oB76h7Ssij4rIZhHZ3NjYODijNgwjaZJS451zzfTmi7/bOVcbqwPXDvwznkyzzrnHnXPLnHPLJkyYkPKADcMYGCFq/MRY2Sf6FInYIyJTYm1Cb1HHD4ZumIZhpEoqRSJ+HasWI8A24KtXOlB3d3dCmWRNKPIJcQ0NDQltmrusJgAl47KquU+GVpQJSQoI/uSFmmAT6parkUzcvCbuaH1r2/nizDUhUhPEQl1wfS6ryaxxyDG1WHhf36H5CpJxbdXQ5hPv8t2fQJdKkYjbQwZoGMbwwDzoDCMimLEbRkQwYzeMiJD2hJPxQoombviqpYR6RGnCTmgiR9CrhoQKb5qQp8XN+4QUbUwhwgwkJwpp22pj19Zc69uXcFLzLtPWMjRm34e2nqHVgkITjPYnfsWjjT1U0PVVrdHGFL9tfxVh7M5uGBHBjN0wIoIZu2FEBDN2w4gIaRfo4r2nWlpaErZLJlxSE4Y0ry9N9Esm/FMjNFRTE2F8wpnmiRYqPiWTCFJLaqh5I4Z6DvrWTBt7aOiqJg5qSUdBvza06yA0saW2vr4EpZqAql2X2lqGVjkaDOzObhgRwYzdMCKCGbthRAQzdsOICGkV6EQkQZDQhCKf51Sop5Ov79B+NBFH27a5uTloO+1458+fV/vWhC5tjloIpiYK+XKpaWKc5rmliXHaMbWwVR+aMKqJacmUM9ZEv1DhTJtPaE47ICFs29e3tkba9eITjrVt44/ZXxit3dkNIyKYsRtGRDBjN4yIEGzssXTSW0Vkdez3WSKyUUQOiMjPRSQ8RMkwjLSTzJ39a8DuPr//NfA959xc4CzwxcEcmGEYg0uQGi8iVcCngG8DfxjLKHs78EhskyeBPwN+0N9xuru7OXv27GVtmjrqixvW1FmtfriW5DAZZVdTv7U2TZ3Vxq6pq75EkKGuktocQ11BQVevtXOhrZvWj09B1tDmrq2btp0v1ltrD3Wb1vrW3HJ950Z706Otkfb2RqsW5HNx1trjz1nKFWGAvwP+F3DpzI8Hmp1zl2Z/AqgMPJZhGBkgJG/8vUC9c27LQDqwijCGMTwI+Rj/G8AqEbkHKALGAH8PjBORvNjdvQo4qe3snHsceBxgyZIlAyoIaRhG6lzxzu6c+6Zzrso5NxN4CPi1c+6/0FsG6jOxzb4APD9kozQMI2VScZf9BvC0iPwVsBX40UAOogkKPjdWza1RE000cUQTpJJJxqgJUFrf2r7adtp4IFz80oQ3bT4+t1xNvNJcaDVhUjs/vnh2rZ9Qt+dkKqiEVtLRREhN+NLOo0/kDU3eqYnJGporNOi20l+CyXiSMnbn3OvA67GfD+Ep5mgYxvDDPOgMIyKYsRtGRDBjN4yIIMlUuUi5M5EG4Gjs1wnASHnxPpLmAjaf4U5/85nhnJuo/SGtxn5ZxyKbnXPLMtL5IDOS5gI2n+HOQOdjH+MNIyKYsRtGRMiksT+ewb4Hm5E0F7D5DHcGNJ+MPbMbhpFe7GO8YUSEtBu7iNwtIntj6aweS3f/qSIiT4hIvYh80KetXEReFpH9se+JGQmGKSIyTUReE5FdIrJTRL4Wa8+6OYlIkYhsEpH3Y3P581h7VqdQG6yUcGk1dhHJBb4PfBJYCDwsIgvTOYZB4MfA3XFtjwGvOueqgVdjv2cLXcAfOecWAjcCvx87J9k4p3bgdufcYqAGuFtEbiT7U6gNSkq4dN/ZlwMHnHOHnHMdwNPA/WkeQ0o4594EmuKa76c3NRex7w+kc0yp4Jyrdc69F/v5PL0XVSVZOCfXy6V8UvmxL0dvCrVnYu1ZMZdL9EkJ90+x3y+lhEt6Puk29krgeJ/fR0o6qwrnXG3s5zqgIpODGSgiMhNYAmwkS+cU+8i7DagHXgYOkt0p1P6OQUoJZwLdION6X29k3SsOESkFfgl83Tl3WT2jbJqTc67bOVdDb/ak5cCCzI5o4KSaEi6etNZ6ozd11bQ+v3vTWWUZp0VkinOuVkSm0HtXyRpEJJ9eQ/+Zc+5XseasnpNzrllEXgNWEJhCbRiSUkq4eNJ9Z38XqI6piQX0prl6Ic1jGApeoDc1F2RZiq7YM+CPgN3Oue/2+VPWzUlEJorIuNjPxcDH6dUgsjKF2qCnhHPOpfULuAfYR++z1LfS3f8gjP8poBbopPd56Yv0Pke9CuwHXgHKMz3OJOZzM70f0bcD22Jf92TjnIBr6U2Rth34APjfsfbZwCbgAPBvQGGmxzqAua0EVqcyH/OgM4yIYAKdYUQEM3bDiAhm7IYREczYDSMimLEbRkQwYzeMiGDGbhgRwYzdMCLC/wNMmkP+zB2KlgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "e = ellipse_dualconic(arr_dx, arr_dy)\n",
    "ps_e = sample_ellipse(e, 100)\n",
    "plt.imshow(arr, cmap='gray')\n",
    "plt.plot(ps_e[:,0], ps_e[:,1],'-r')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.173800Z",
     "start_time": "2020-08-06T18:56:10.170689Z"
    }
   },
   "outputs": [],
   "source": [
    "# export\n",
    "class DualConicEllipseRefiner(EllipseRefiner):    \n",
    "    def __init__(self, cutoff_it, cutoff_norm):\n",
    "        super().__init__(cutoff_it, cutoff_norm)\n",
    "        \n",
    "    def proc_arr(self, arr): return grad_array(arr)\n",
    "    \n",
    "    def refine_point(self, arrs, p, W): return ellipse_dualconic(*arrs, W)[:2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.182320Z",
     "start_time": "2020-08-06T18:56:10.174773Z"
    }
   },
   "outputs": [],
   "source": [
    "refiner = DualConicEllipseRefiner(cutoff_it=20, cutoff_norm=1e-3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.375674Z",
     "start_time": "2020-08-06T18:56:10.183267Z"
    }
   },
   "outputs": [],
   "source": [
    "file_img = Path('data/dot_vision_circle/SERIAL_19061245_DATETIME_2020-05-17-20:48:11-005188_CAM_1_FRAMEID_0_COUNTER_1.png')\n",
    "img = File16bitImg(file_img)\n",
    "ps_c_p = torch.FloatTensor([[ 89.74539133, 996.88665043]])\n",
    "bs_c_p = torch.FloatTensor([[[ 89.5465538485512, 1029.648805759194],\n",
    "                             [ 57.6381616678298,  996.914281777708],\n",
    "                             [ 89.9441082645408,  964.144356548655],\n",
    "                             [121.8893325612110,  996.858987492743]]])\n",
    "ps_c_p_refined = refiner(img.array_gs(torch.float), ps_c_p, bs_c_p)\n",
    "assert_allclose(ps_c_p_refined, torch.FloatTensor([[ 95.89870615, 997.13709847]]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.650228Z",
     "start_time": "2020-08-06T18:56:10.376850Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(982.0, 1012.0)"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAREAAAD4CAYAAADLqNJwAAAASXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjIuMSsyODU3LmdmODBkNTZlNjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvoazhzAAAAAlwSFlzAAALEwAACxMBAJqcGAAAGqRJREFUeJzt3X2wXdV53/Hv78oyQhJIQq+UFwc3IgGTgVgqJc1AGZOxMeMpxLE9uM0Yxx5rYuMpdseZ4qlbJ249Y1J3GvJSu3LA4GTAIYEE4mIMJQmk04AtU2wLY0CYOEgBCYPe0BsCPf3jbNlHl7ufpbvW1b1H499nRnPP3Wu/rLPOvo/22c9aaysiMDOrNTbTFTCzo5uDiJk1cRAxsyYOImbWxEHEzJq8ZqYrkFmwYEGsWLHisNaV1Fs2ExmorD4HDhyYxpr8WG0bZWXZPkvbZsbG+v9/a2m/bL8tat9nqf1K5X0mU5/HH3/8hxGxtOpAjHgQWbFiBWvXrv3R79nJM2vWrN6y/fv395a1nFTZB/ya1/Q37a5du6qPmb3PV155Jd32ta99bW9Z1kYvv/xyVX2g/g9+zpw5vWV79uxJt80+l2OPPba3LGu/0h9z1n7ZH/QxxxyT7jc7j2rrM/69vOlNb/pB1UE6/jpjZk0cRMysyUh/nYHD/26XrZd9ZSldpmaX8pmXXnqpt6z0nmq/mpUufffu3Vt1zNo2gLzts/1mbdRyT2Tfvn1V+y21wbx586q2Lb2XrL61X22n+h6hr0TMrImDiJk1cRAxsyYOImbWxEHEzJo4iJhZk5FO8Uo6pJfljh07qvaTpRlb0l21vV1ruzJDnror9VjNZOnf2bNnVx+zNr25c+fO3rJSL9lMVt+sPqX0eW3bl1LHWdtn6f5MqZfsZBX/CiRdL2mLpPVDy06QdI+kJ7qfi7rlPyvp7yTtk/Sxcfu5WNJjkjZIunpK34WZzZjD+a/0BuDiccuuBu6NiJXAvd3vAC8A/xb47PDKkmYBfwC8FTgTeLekM+urbWajohhEIuJ+BsFh2KXAjd3rG4HLunW3RMQ3gPHXWecCGyLi+xHxEvDlbh9mdpSrvbG6PCKe6V4/CywvrH8S8PTQ7xu7Za8iaY2kdZLWbdu2rbJ6ZjZdmrMzMbgzOWWd8SNibUSsjojVCxcunKrdmtkRUhtENks6EaD7uaWw/ibglKHfT+6WmdlRrjbFewdwBfCZ7ufthfW/AayUdBqD4HE58K8P50DDKbAs3VU70rSUusvSaFlar2Wyo9oRmNnI4dK2WXqzJfWZpY4zWftl5wHkn1k2MVN2zFIKt2XbTJbOrh25Xpsa7lMMIpJuBi4ElkjaCHySQfC4RdL7gR8A7+rWXQGsA44HDkj6CHBmROyQ9GHga8As4PqIeGRK34mZzYhiEImId/cUXTTBus8y+Koy0X7uBO6cVO3MbOS527uZNXEQMbMmDiJm1sRBxMyajPQo3gMHDhyStszSXdnI2KyslBbNUoLZJLpZyrmUFs1Sqi2jeGvTkNn7LLVfNmK0NnVcSlHWTuScpY5bPrOsbO7cuel+a59nk/2tTPUDvHwlYmZNHETMrImDiJk1cRAxsyYOImbWxEHEzJo4iJhZk5HuJyLpkNx9ljPP+itkOf7SzNfZfrN8e3bMXbt2pccslfcpDZHfs2dP1TFrH7xdOmamtk8QwJw5c6r2m73PUn+OWrt3707Ls/Oz9r20PBB9Ir4SMbMmDiJm1sRBxMyaOIiYWRMHETNr4iBiZk1GPsU7nEbNUnu1s2KXhpXXPvA7S22WUpRZejh7qHlpKoAXX3yxtyxLCWbvpTTLfu0Q+WOPPbaqDGD+/Pm9ZVnbZyny0qz1tdMIlNLy2fmXtX3299Ay+/xEfCViZk0cRMysiYOImTVxEDGzJg4iZtbEQcTMmox0ivfAgQOHpNayUbPZKMssFVaarTwrz9JzWeq4lOLduXNnb1mWasxSuABbt26tKsuOWRpxnM0Un422zWbZL42oXbZsWW/Z0qVLe8sWLFjQW1aa7T1Lm2Yp59L5l7VRJkufl9LKk1W8EpF0vaQtktYPLTtB0j2Snuh+LuqWS9LvStog6duS3ji0zSuSHu7+3TGl78LMZszhfJ25Abh43LKrgXsjYiVwb/c7wFuBld2/NcDnhrbZExHndP/+VVOtzWxkFINIRNwPvDBu8aXAjd3rG4HLhpZ/KQYeABZKOnGK6mpmI6j2xuryiHime/0ssLx7fRLw9NB6G7tlAHMkrZP0gKTL+nYsaU233rrt27dXVs/MpkvzjdWICEmHM8DkdRGxSdLrgb+S9J2IeHKC/a0F1gKcfvrpdQNXzGza1F6JbD74NaX7uaVbvgk4ZWi9k7tlRMTBn98H/gb4+cpjm9kIqQ0idwBXdK+vAG4fWv6eLktzHrA9Ip6RtEjSMQCSlgC/CHy3od5mNiKKX2ck3QxcCCyRtBH4JPAZ4BZJ7wd+ALyrW/1O4BJgA7Ab+LVu+RnA/5R0gEHg+kxEFIOIpEOGNGe5+Np+GVl/hNIxs/4TWX2yvhOl8qw/x7PPPpvud/Pmzb1lL7ww/t75j2X9T0rvJevbk30u2SznWX8OgG3btvWWZdMarFixords4cKF6TGzfiTZeZK1D+T9PbJjHnfccb1lpf5Ek1UMIhHx7p6iiyZYN4ArJ1j+f4Gfm3TtzGzkudu7mTVxEDGzJg4iZtbEQcTMmjiImFmTkZ4KAMrD5g/KhlRn+yjN9p6lMLP0XDb9QOkhzlm6dePGjb1lmzZtSvf7/PPP95Zls8i3PBy6drb3bLh6VlfI07Glofe1svpmM6+XUryHe/6Plw0ZKT3EfrJ8JWJmTRxEzKyJg4iZNXEQMbMmDiJm1sRBxMyajHSKNyIOGUWbpQSzGcCzlGrpgd3Zw6OzVGOWOi6lKLMU75YtW6rKIE/7ZaNbW5RSmH2yka+lkcPZeZKlq7M0bWm0dzaiNisrPRB90aJFvWW1D+2e9tnezcwyDiJm1sRBxMyaOIiYWRMHETNr4iBiZk1GOsU7XpZmy1KCWborK4N81Ge2bZZWziZbhvo0bulhX7VtlI0kLaVwa0dXZ/stpeWzh4xn+81S688991x6zGxi5Cx9Xpp0OktnZ6NxszR3aeT6ZPlKxMyaOIiYWRMHETNr4iBiZk0cRMysiYOImTVxEDGzJofzQO/rgbcBWyLirG7ZCcCfAD8F/D3wrojYqkHi/1oGD/XeDbw3Ih7qtrkC+ES32/8SETceTgWH891ZX4bafgXZMG3Ih2pnZVnfgNKw++yB1NnDmEvDymv7gtQOrYf8c8nqU+oLkql9uHbWT6T0QO9seodsOonSfrN2yPp7ZH1ISp/ZZB3OlcgNwMXjll0N3BsRK4F7u98B3gqs7P6tAT4HPwo6nwT+OXAu8ElJ/RMlmNlRoxhEIuJ+YHyIvhQ4eCVxI3DZ0PIvxcADwEJJJwJvAe6JiBciYitwD68OTGZ2FKq9J7I8Ip7pXj8LLO9enwQ8PbTexm5Z3/JXkbRG0jpJ60rduM1s5jXfWI3Bl7b6L7Cv3t/aiFgdEatL4wrMbObVBpHN3dcUup8HR4VtAk4ZWu/kblnfcjM7ytUGkTuAK7rXVwC3Dy1/jwbOA7Z3X3u+BrxZ0qLuhuqbu2VmdpQ7nBTvzcCFwBJJGxlkWT4D3CLp/cAPgHd1q9/JIL27gUGK99cAIuIFSf8Z+Ea33qcioj+f1hkbGzskPZbNbp2lRbM0Y+kBz1mKLds2S4tmaVqAnTt39pZl6bnStAa1qe5su9JUAFl9s22z9svOg9K22XtpeQh7Vp7VN5u2AGDx4sW9ZVn7ZWWl9pusYhCJiHf3FF00wboBXNmzn+uB6ydVOzMbee6xamZNHETMrImDiJk1cRAxsyYOImbWZKRne4+IQ0YqZmm02nRhaebrUtq0T236F+pHDpfqmo3UrR3hWxptm6WOa0cOZ9uVyrP2axmtnJ1H2cjhbJb40n7nz5/fW5bVt/Rw8snylYiZNXEQMbMmDiJm1sRBxMyaOIiYWRMHETNrMvIp3uFUaZZOrE3FliZqzlJsWRq3tgzqR5qWUp+1kx/XjooF+OF7f8iBua/efmz3GMu+tKx3u5a0claeTWBcar9Mdv5l7VdKt2afd3ZuZiN1p3oUr69E7IiaKIBky+3oM9JXInao//XQQyye4H+f58bGOPfUU2egRma+EjmqTBRAAJYWesCaHUkOImbWxEHEzJo4iNgRNbZ74lOsb7kdfXxj1Y6opTcu7S+sz6jaCDmqgsjs2bN7y7L+HtkUAqV+IllOfc6cOb1lWX+E0jH7+is8P3v2xNmZWbMYGxsr9p/I+ivU9rMpzfbeMmt77TEztf1sSnXNPu+Wvim1/U+yY061oyqI/KR726pVbNmypbyi2TTyF1Mza+IgYmZNHETMrImDiJk1cRAxsyZN2RlJVwEfYJDx/0JE/I6ks4HPA/OBvwf+TUTskPRTwKPAY93mD0TEr2f7j4hDUnJZird2+HfpQc1ZOjGrT1Y2/JDyiWTDw7Oy0hQDmex9HokpGKA8632t2mHwWeq9ZYb07JilNsjOldoU+b59+6q261MdRCSdxSCAnAu8BNwl6SvAHwIfi4j7JL0P+A3gP3abPRkR57RV2cxGScvXmTOAByNid0S8DNwHvB04Hbi/W+ce4Ffaqmhmo6wliKwHzpe0WNJc4BLgFOAR4NJunXd2yw46TdL/k3SfpPMbjm1mI6I6iETEo8A1wN3AXcDDwCvA+4APSfomcByDrzoAzwCnRsTPA/8OuEnS8eP3K2mNpHWS1m3fvr22emY2TZqyMxFxXUSsiogLgK3A4xHxvYh4c0SsAm4GnuzW3RcRz3evv9ktP32Cfa6NiNURsXrBggUt1TOzadAURCQt636eyuB+yE1Dy8aATzDI1CBpqaRZ3evXAyuB77cc38xmXusAvFslLQb2A1dGxDZJV0m6siu/Dfhi9/oC4FOS9gMHgF+PiBeynY+NjR3ywOOdO3dWVbI0arZUh5qyLK1cShdmV2A7duzoLSulq2sfdF37sPTStlnq+EgdM2v7bFT28ce/6pv3YZe3pIezzyzrRtDyEPvJagoiEfGqm6MRcS1w7QTLbwVubTmemY0e91g1syYOImbWxEHEzJo4iJhZEwcRM2sy0nOsRgR79uz50e9ZSitLhWVlpfTbrl27qvY7b9683rITTjghPWbWU3fr1q29ZXv37k33m5XXTuJcmjS5dgRrbWoY6kdBL168uLes9JktXdo/q312Lgx3YZisbNR27QTjNXwlYmZNHETMrImDiJk1cRAxsyYOImbWxEHEzJo4iJhZk5HuJwKH9hfI+gdkZVk+vdTnIOubkg2pzrYr9TnIpjwY7jczXjacH/I+JpnSQ6czWV+QrI0ypakdshnSlyxZ0luWfS5ZPxDIp37I+mWUpjXI+rVkfWmyNqpt9956TOnezOwnjoOImTVxEDGzJg4iZtbEQcTMmjiImFmTkU/xDqcXsyHp2dDnLBVWSl9m5bVD5LOh4QDLly/vLasdsg952m/btm29ZVmKvJSizMprp2+YO3duesxs5vUTTzyxt2zZsmW9ZaUh+1laufah8JCf19nnnaWGS10BJstXImbWxEHEzJo4iJhZEwcRM2viIGJmTRxEzKzJyKd4h0fZZmmrLH3ZMnN4tm0ppdqnNNt2lk7Mjll6L1lqOUvxvvjii71lWfoX8lRtVt9sVGw2KzvkKeDaFO+iRYvSY071yNiDamdmz0aYl9Lyk9V0JSLpKknrJT0i6SPdsrMl/Z2k70j6S0nHD63/cUkbJD0m6S2NdTezEVAdRCSdBXwAOBc4G3ibpJ8G/hC4OiJ+Dvhz4De69c8ELgfeAFwM/A9J/T1pzOyo0HIlcgbwYETsjoiXgfuAtwOnA/d369wD/Er3+lLgyxGxLyKeAjYwCEBmdhRrCSLrgfMlLZY0F7gEOAV4hEHAAHhntwzgJODpoe03dssOIWmNpHWS1mVPgjOz0VAdRCLiUeAa4G7gLuBh4BXgfcCHJH0TOA7I77y9er9rI2J1RKxesGBBbfXMbJo03ViNiOsiYlVEXABsBR6PiO9FxJsjYhVwM/Bkt/omfnxVAnByt8zMjmKt2Zll3c9TGdwPuWlo2RjwCeDz3ep3AJdLOkbSacBK4OstxzezmdfaT+RWSYuB/cCVEbGtS/te2ZXfBnwRICIekXQL8F3g5W79tKOFpHQo9LBshvQs117qW5Hl27O6Zf0GSnn6bCh71lemNEQ+63+S9YPYvXt3b9nevXvTY2ay95K1X9aHBPL+Hllfmax9SjPMZ+dY1lemtN+sfbM+TJnavid9moJIRJw/wbJrgWt71v808OmWY5rZaHG3dzNr4iBiZk0cRMysiYOImTVxEDGzJiM9FUBEHJIOzdJ+WeouS+OWUpS1D/RuSetlqeOsrJS6y46bPcw6S/EeKdn7nD9/frptNvN61ga100lAnnbOzr/aNC3UdyMoTd8wWb4SMbMmDiJm1sRBxMyaOIiYWRMHETNr4iBiZk1GOsUL5QduH44sxVaapbt2RveW2b+ztF9WVkrdLVy4sKo+WSq71D61n1/2mZVS5Nkxax+8XXofWYo3S7eWzpOsTtmDuY9UWnnC/U3p3szsJ46DiJk1cRAxsyYOImbWxEHEzJo4iJhZk5FO8Uo6JLVWm9KqnVAZypMC98kegl06ZvY+s/Rcqa7ZKN99+/b1lrWkyLP6ZqnPbBLn0uTdWXl2zKx9SmnRLL2epaSzzxry9s3Szll9p/rh474SMbMmDiJm1sRBxMyaOIiYWRMHETNr4iBiZk0cRMysSVM/EUlXAR8ABHwhIn5H0jnA54E5DB7c/aGI+LqkC4Hbgae6zW+LiE+VjjGc786GnWe5+Gy70sO1sz4HLcPVM7VTAZT6MmR9ErI+JlmfjdID0bNts5n2s/YrHTMrz/pWZOdC1o8G8oewZ+ffkZrtPWuD2ukt+lSf6ZLOYhBAzgVeAu6S9BXgt4HfioivSrqk+/3CbrO/jYi3tVXZzEZJy5XIGcCDEbEbQNJ9wNuBAA6G5QXAPzbV0MxGWss9kfXA+ZIWS5oLXAKcAnwE+K+SngY+C3x8aJtfkPQtSV+V9IaJdippjaR1ktZt3769oXpmNh2qg0hEPApcA9wN3AU8DLwCfBD4aEScAnwUuK7b5CHgdRFxNvB7wF/07HdtRKyOiNULFiyorZ6ZTZOm7ExEXBcRqyLiAmAr8DhwBXBbt8qfMrhnQkTsiIgXu9d3ArMlLWk5vpnNvKYgImlZ9/NUBvdDbmJwD+Rfdqu8CXiiW2eFulvjks7tjv18y/HNbOa1TgVwq6TFwH7gyojYJukDwLWSXgPsBdZ0674D+KCkl4E9wOVRytVxaNqtNj2XbVd6CHY203ntg5GztGdJy4zkWfmePXt6y7IZ0rP2gfppGGofag71KdXscyl9Zln6PKtvKcXb8nnXHnOymoJIRJw/wbL/A6yaYPnvA7/fcjwzGz3usWpmTRxEzKyJg4iZNXEQMbMmDiJm1mSkZ3uHQ9NRWToxG4VamxqGPA1ZOzq4lJrLUnC1I1QhHxk7b968qv2WRoTWplSztGjpfWbp1uyYtZ91qU5Z+ryUOq59GHg26tizvZvZSHEQMbMmDiJm1sRBxMyaOIiYWRMHETNrMtIp3og4JLVWO2lyadRnpnZEaJZOLY2izN5LlvIrDYquTbdmo5Xnz5+fHrM0yrdPNolzlo6GPN2ftUHtA+NLsvqUuhhk9c3OzeyYtZ9JH1+JmFkTBxEza+IgYmZN1PJd70iTtBN4bKbrMc4S4IczXYkhrk9u1OoDo1enn4mI42o3Hukbq8BjEbF6pisxTNK6UaqT65MbtfrA6NVJ0rqW7f11xsyaOIiYWZNRDyJrZ7oCExi1Ork+uVGrD4xenZrqM9I3Vs1s9I36lYiZjTgHETNrMlJBRNJHJT0iab2kmyXNkXSapAclbZD0J5LqHx83NfW5QdJTkh7u/p0zjfW5qqvLI5I+0i07QdI9kp7ofi6a4fr8pqRNQ+1zyRGuw/WStkhaP7RswjbRwO9259K3Jb1xhutzoaTtQ231n6a6Pkmd3tl9bgckrR63/se7NnpM0luKB4iIkfgHnAQ8BRzb/X4L8N7u5+Xdss8DH5zh+twAvGMG2ucsYD0wl0H/nv8N/DTw28DV3TpXA9fMcH1+E/jYNLbLBcAbgfVDyyZsE+AS4KuAgPOAB2e4PhcCX5mhNjoD+Bngb4DVQ8vPBL4FHAOcBjwJzMr2P1JXIgxOxmO75/jOBZ5h8FDwP+vKbwQum8H6/OM0Hnu8Mxic9Lsj4mXgPgYPUb+UQbvA9LZPX32mVUTcD7wwbnFfm1wKfCkGHgAWSjpxBuszLSaqU0Q8GhET9Qa/FPhyROyLiKeADcC52f5HJohExCbgs8A/MAge24FvAtu6kxRgI4MrhBmpT0Tc3RV/ursc/u+S+sdcT631wPmSFkuay+B/1VOA5RHxTLfOs8DyGa4PwIe79rl+Or9eDelrk5OAp4fWm67zKfuMfkHStyR9VdIbpqEuJZNuo5EJIt3JdimDS6h/AswDLh6l+kj6VeDjwM8C/ww4Afj301GfiHgUuAa4G7gLeBh4Zdw6AUxLzj6pz+eAfwqcwyD4/rfpqE+f6WyTwzGuPg8Br4uIs4HfA/5ipurVYmSCCPBLwFMR8VxE7AduA36RwSXnwTE+JwObZrA+/yIinukuh/cBX6RwqTeVIuK6iFgVERcAW4HHgc0HL8m7n1tmsj4RsTkiXomIA8AXmMb2GdLXJpv48dUSTN/5NGF9ImJHRLzYvb4TmC1pyTTUJzPpNhqlIPIPwHmS5mrwJKCLgO8Cfw28o1vnCuD2GazPo0Mngxh8t13fv4upJWlZ9/NUBvcfbgLuYNAuML3tM2F9xt1j+GWmsX2G9LXJHcB7uizNeQy+oj4z0Q6moz6SVnTnEZLOZfD3+Pw01CdzB3C5pGMknQasBL6ebnGk7wxP8i7ybwHfY3Di/RGDO8Sv797EBuBPgWNmuD5/BXynW/bHwPxprM/fMgis3wIu6pYtBu4FnmCQITlhhuvzR137fLs7IU88wnW4mcHXpv0Mvr+/v69NGGRl/oBBxuE7DGUlZqg+HwYe6drvAQZXutPVRr/cvd4HbAa+NrT+f+ja6DHgraX9u9u7mTUZpa8zZnYUchAxsyYOImbWxEHEzJo4iJhZEwcRM2viIGJmTf4/7RvS9+VTAXoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(img.array_gs(torch.float), cmap='gray')\n",
    "plt.plot(ps_c_p[:,0], ps_c_p[:,1], 'rs')\n",
    "plt.plot(ps_c_p_refined[:,0], ps_c_p_refined[:,1], 'gs')\n",
    "plt.gca().set_xlim(80,  110)\n",
    "plt.gca().set_ylim(982, 1012)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-05-03T20:36:21.604772Z",
     "start_time": "2020-05-03T20:36:21.600609Z"
    }
   },
   "source": [
    "# Build"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-08-06T18:56:10.851923Z",
     "start_time": "2020-08-06T18:56:10.651429Z"
    }
   },
   "outputs": [
    {
     "data": {
      "application/javascript": [
       "IPython.notebook.save_checkpoint();"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Converted control_refine.ipynb.\n"
     ]
    }
   ],
   "source": [
    "build_notebook()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.9"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}