Positive Semi-definite Matrix.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Based on this tweet: https://twitter.com/gabrielpeyre/status/1350376375312650240"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.special import factorial\n",
    "from celluloid import Camera\n",
    "from IPython.display import HTML\n",
    "from IPython.display import Video\n",
    "from matplotlib.animation import FFMpegFileWriter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "grid = np.indices((1, 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(2, 1, 3)"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "grid.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_psd(kernel_fn, n=2, **kwargs):\n",
    "    grid = np.indices((n, n))\n",
    "    A = kernel_fn(grid[0], grid[1], **kwargs)\n",
    "    return A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gaussian_kernel(x, y, s=1):\n",
    "    return np.exp(-((x-y)/s)**2)\n",
    "\n",
    "def laplace_kernel(x, y):\n",
    "    return 2*(x==y) - ((x-1)==y) - (x==(y-1))\n",
    "\n",
    "def hilbert_kernel(x, y):\n",
    "    x = x + 1\n",
    "    y = y + 1\n",
    "    return 1 / (x+y-1)\n",
    "\n",
    "def pascal_kernel(x, y):\n",
    "    return factorial(x+y) / (factorial(x)*factorial(y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_figure(n, width=4):\n",
    "    fig = plt.figure(constrained_layout=True, figsize=(n*width, 3*width/2))\n",
    "    gs = fig.add_gridspec(nrows=3, ncols=n)\n",
    "#     im_list = []\n",
    "#     line_list = []\n",
    "    for i in range(n):\n",
    "        ax = fig.add_subplot(gs[:2, i])\n",
    "        #im = ax.imshow(np.zeros((2, 2)))\n",
    "        #im_list.append(im)\n",
    "        ax.axis(\"off\")\n",
    "        \n",
    "        ax = fig.add_subplot(gs[2, i])\n",
    "        #line, = ax.plot(np.zeros(1), color=\"k\", marker=\"o\")\n",
    "        ax.axis(\"off\")\n",
    "        #line_list.append(line)\n",
    "    return fig#, im_list, line_list\n",
    "\n",
    "def make_plot(\n",
    "    A_list, width=4, fig=None, \n",
    "#     im_list=None, line_list=None\n",
    "):\n",
    "    n = len(A_list)\n",
    "    \n",
    "    if fig is None:\n",
    "        fig = generate_figure(n, width)\n",
    "    \n",
    "    for i, A in enumerate(A_list):        \n",
    "        ax = fig.get_axes()[i*2]\n",
    "        ax.clear()\n",
    "        im = ax.imshow(A)\n",
    "        ax.axis(\"off\")\n",
    "        #im_list[i].set_data(A)\n",
    "        \n",
    "        ax = fig.get_axes()[i*2 + 1]\n",
    "        ax.clear()\n",
    "        eigvals = np.sort((np.linalg.eigvals(A)))\n",
    "        #line_list[i].set_data(np.arange(len(eigvals)), eigvals)\n",
    "        ax.plot(eigvals, color=\"r\", marker=\"o\")\n",
    "        ax.axis(\"off\")\n",
    "        \n",
    "    #fig.tight_layout()\n",
    "    return fig"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABIgAAAG4CAYAAAAubtzJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO39e9xld13f/b++SQgkcvZAjXKQCgKKIoiSqUU0MlFrrbaW1lSrdzmLqKBVb7QeK/XQ4rmCWE8t+Lvxtt7VX1VuCij1BB7gp3iqoICIoEECSoRkMt/fH+tKScJMyJCZta6Z9Xw+HvuxZva+sj+fneuada393t/DmHMGAAAAwH6dt3UDAAAAAGxLQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcuuLkHj7/xPqtucfbCvz1/zXI99XcetV6xF99lvVrV3V729lXrXfCqN6xa7/hVb1213jx27YrFzu2dBV9w/CfGqXz92uehyy950JrluvJxl65aD6iXP/Opt/g8tPY56L4/9sQ1y3XXV6738i5+44q/S6sLr3rnqvXOe/u69cY7rlm1Xtes9/2bK9aq6tixVcv9/Jt/4JSuhR553j89ty8OgdWd7D2ZEUQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzl1wcw++8G/PX6uPqi676LpV6z3jgc9brdZTe9Rqtare1F1WrXe3Llm13gWvWrVcx69662q15rFrV6u1FJzr1jtFl1/yoFXrPf8Nr1i13uUr/tO58nGXrlcMzhH3/bEnrlrvrq9c95x88RvX+51z4VXvXK1W1XlvX7feeMc1q9brmnWvF+aa9Y4dW69W1XXrvgcBOKyMIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB27oKbe/Cpv/Ootfqo6hkPfN6q9S676LrVaq392p7aut+7N3WXVevdrUtWrXfBq9ardfyqt65XrJrHrl213qm68nGXrlrv8nV/tHr+G16xWq21X9va3zs4E+76yrlqvYvfuO45+cKr3rlarfPevl6tqvGOa1at1zXrfu/myvU6dmy9Wtetd41eNVeuB3BYGUEEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdu+BmH33xXVZqY/HUHrVqvWc88Hmr1brsoutWq1XrvrZa/3v3ptb92bxbl6xW64JXrVaqquNXvXXdgofclY+7dNV6l6/3o9Xz3/CK9Yq17mur9b937MPFb7x21XoXXvXOVeud9/b16o13XLNaraquWfd7N1eu17Fj69a7br1r2blirWrV1wZwQmNs3UFlBBEAAADA7gmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcxfc3IN3e9nb1+qjqjd1l1XrPbVHrVbrGQ983mq1qi676LpV6639+tb83tW6P5t365LValVd8KpVy3ETVz7u0tVqXb7uj1bPf8MrVq239utb83vHdi686p2r1jvv7evWG++4Zr1i11y7Xq1qrlyvY8fWrXfdutd6c816a7+2OVetB+eMMbbu4MwZ646lGecdjv+XRhABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAO3fBzT74qjes1UdVd+uSVeu9qbusVuupPWq1WlXPeODzVq132UXXrVpv7de35vdvzZ/LWv/fHdu58nGXrlrv8pV/tJ7/hlesWm/t17f294/FeW9/56r1xjuuWbVe11y7Wqm5Yq2qjh1bt951614LzZXrrfn65pyr1apW/97BGTPGyvXWHW8yzj9/xVorv7YLL1y13skYQQQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADs3AU39+Dxq966Vh9VXfCqVct1ty5Zrdabustqtaqe2qNWrfeMBz5v1XqXXXTdqvXWfH1rf+/W/tlkP6583KWr1rt8vVN6Vc9/wytWrbfm61v7e3eYjXdcs27Ba65dtdxcs96xY+vVqrpu3WuFuXK91V/fnOsVW/u1HV/xtbEvY6xcb93xH+P889etd7vbrljrdqvVqupOt1+33kkYQQQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADs3AU39+A8du1afVR1/Kq3rlrvgletV+tuXbJesepN3WXVek/tUavWe8YDn7dqvcsuum61Wmu/trW/d3CmXPm4S1etd/m6p/We/4ZXrFZr7de29vfulFyz7rXQXLlex46tV+u69X6XVs2V663++uZctd6ar28eX/m1zePr1oMzZJx//rr1bnfbVeudd/v3Wa3W8fe982q1qq6+1x1XrXcyRhABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnLrjZR+dcqY2DcseuXbXe8aveulqtC161Wqmq7tYlq9Z7U3dZtd5Te9Sq9Z7xwOetVuuyi65brVat+9oW37hyPTgzrnzcpavWu3zF0/rz3/CK9Yq17mur6pm3/EvnNetem3Ts2Lr1rlvvd85csVa16murmitfN6/++o6v+Prm8fVq1ervediRse54jHH+yvVud7tV6x1/3zuvVuuv73en1WpV/fmRsWq9kzGCCAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANi5C7Zu4EbmXLfcsWtXq3X8qreuVqvqgletWq67dcmq9d7UXVat99QetVqtZzzweavVqrrsoutWrQe8d6583KWr1bp83VN6z3/DK9YteCqOHVu33nXrnpPnmvXWfm0rX1eu/vqOr/z65vEVa6382tiPMdYtd97K9S68cNV63en2q5a7+l53XK3Wnx9Z93v36n/+zFXr1Zed8F4jiAAAAAB2TkAEAAAAsHMCIgAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ0TEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAIAAADYOQERAAAAwM4JiAAAAAB2TkAEAAAAsHMCIgAAAICdG3POrXsAAAAAYENGEAEAAADsnIAIAAAAYOcERAAAAAA7JyACAAAA2DkBEQAAAMDOCYgAAAAAdk5ABAAAALBzAiIAAACAnRMQAQAAAOycgAgAAABg5wREAAAAADsnIAJO3RhXNMZrGuP4wfGKrVsCdsZ5CNiScxCwtTNwHhpzztPRGrAXy4nn2dXFN7j36uqxzfncbZoCdsV5CNiScxCwtTN0HhIQAadmjNdU9zzBI69tznut2wywK2OcVz24ekF15xN8hfMQcOa5FgK2dobOQwIi4NSMcbwaJ3hkNqdpq8DpNcbdqqPVpxwc3+9mvtp5CDjzXAsBWztD5yEnMOBU/elJ7n/dql0A56YxbtMYD2+MpzfGb1ZvrH6s+uTq56rPrV5/kv/aeQhYw8nONc5BwFrOyHlIQAScqpec4L6rq6et3Qhwjhjjno3xuMb4r9Wbq1+svqJ6e/XV1UOqD2zOf9mcz6m+suW8c0POQ8BanlZdc5P7nIOANT2tuul0sFt9HjLFDLjlxviQ6pXVH1TvW92jJaV+mkUZgVtsjIuqh7dMG7u8uv/BI6+rfv7g9qLmfOvNPMcV1dNzHgK2MMavVw+qzs85CFjbGHeu3lJdVd2p03QeEhABt8wYo3p+9bDqAc15sikeADe2nD8+rCUM+pTqEdXtqne2jBa6PhT6g1yYAIfdGLeprqz+P835+K3bAXZojMtbrp0ua84Xna6nveB0PRFwzvu86pHVk4RDwLt59xE931j9VUsg9Cm9a6eNP6ye1XJR85LmvOlUMYDD7mOrO1b/79aNALt1pDpevex0PqkRRMB7NsYHVL/fMrXs7zfn8Y07Ag6TJRx6dnXxCR79m+p/tARCz2/O16zYGcDpN8Y3VF9TvV9zvmXrdoAdGuMFLeegjz6dT2sEEXBLfGd1++qxwiHgRsY4v+UccaJw6E3VPZrzpou5ApzNjlYvEw4Bm1iuvT6u+s+n+6ntYgbcvDH+QfU51Tc35+9t3Q5wSIxxx8b40upV1fuf5Ks+QDgEnFPGuEvLFDPTy4CtfHh1h+pXTvcTC4iAkxvjDtX3V79bfcvG3QCHwRj3aoxnVK+vvuPg+Jcn+erXrdUWwEo+seU9lIAI2MqRg+Ovnu4nFhABN+fp1QdXjzEKAHZsjNEYRxrjJ6pXV0+ufqZ6aHP+/epLq5suNn119bR1GwU4445Wb+s0LwwLcAoubZnG/yen+4mtQQSc2BhHqidV39Ocv7Z1O8AGlq2c/0n1lJYpFVdV31597412M5zzuY1RN97F7GnN+dy1WwY4Y8YY1eXVi5vz2q3bAXbrSPUrnYEdx+xiBry7MW5bvbx6n+ojmvOvN+4IWNMYd64e2zJS6O7VH7UsRP2jzfn2LVsD2MwYH9pyPnxSc/7HrdsBdmjZXfpN1Vc057ef7qc3ggg4ka+q7l99mnAIdmSMv1t9SfWvWgLiF7eMJPzvdjAE6OjB0fpDwFYednA87QtUl4AIuKkxHlB9dfXc5vy5rdsBzrBlysTDW6aRfUZ1rPrx6jua8xVbtgZwyBxtWfPj1Vs3AuzWkera6jfPxJObYga8yxjnVb9U3be6f3OebGci4Gw3xoXVo1qCoQdXb27ZtfA/Nuefb9kawKGzrMn25pYP0J6wdTvATo3xkurC5nzYe/za94IRRMANPbFlVfzPEw7BOWqM960e3zJ17JLq9w/+/p+b82+3bA3gEPu46g6ZXgZsZQmqH9rygd4ZISACFmPcvfqW6vnVczbuBri1xriiG+8q9j3Vh1afX11UvaB6dPX/Wl8I4D06Wh2vXrR1I8BuPai6XWdo/aESEAF1/Rok31+dVz3hTGyZCKxoCYeeXV18cM89q3/fMmf9x6rvbM5XbtQdwNnoaPXS5rxq60aA3TpycPzVM1VAQARU/bPqH1RPac7XbNwLcOs9vXeFQzf0puZ8zNrNAJzVxrhLy7SOb9q6FWDXLq1e15x/dqYKnHemnhg4SyzrkXx39estU1CAs9myE+E9T/LoB63ZCsA54pNa3jdZfwjY0pHO4PSyEhAB9R+qu1SPac7rtm4GeC+N8X6N8b3Vb7esk3Eir1uxI4BzxdHqbdXLtm4E2KkxPri6e2dwelkJiGDfxnhky4K139qcv711O8B7YYwLG+Mp1auqJ1TPatmR8OqbfOXV1dNW7g7g7Las03h59aLmPLZ1O8BuXXpwNIIIOAPGeJ+WN5J/WP3bjbsBTtUYozH+UfXK6hnVr1Uf2ZxPas4fqB5bvbaaB8fHNudzN+sX4Oz0oS3Tdk0vA7Z0pPrb6v93JotYpBr26xuqD6ke3pzv2LoZ4BSM8ZHVd7Ssi/H71ac158/d6GuWMEggBHDrHD04CoiALR2pfr05rz2TRYwggj0a46HVU6pnNuf/3Lod4BYa426N8QPVy6sHVU+uPurdwiEATpej1R8356u3bgTYqTEuqj66Mzy9rAREsD9j3Kb6weqN1Vdt3A1wS4xxu8b4yuqPqv+j+q7qQ5vze8/0J0kAu7VcM31SRg8B23pIdZtWCIhMMYP9+fLqI6vPbM63bt0McDOWxVH/SfVtLVNCf6b68ub8X5v2BbAPD6tun4AI2NaRg+OvnelCRhDBnoxx3+rrqv+7Of/b1u0AN2OMh1S/WP1E9TfVJzfnZwiHAFZztLquevHWjQC7dmn1R835l2e6kIAI9mKM86ofaFn9/skbdwOczBiXNMaPVL9e3a96fPXRzfnCTfsC2J9HVi9tzqu2bgTYqWU0+ZFWmF5WAiLYk0dXn9AyPeWNWzcD3MQYFzfGv2lZZ+hzqm+v7tOcP9Cc123bHMDOjHHX6qHVC7ZuBdi1e1cfUP3qGsWsQQR7MMYlLW82X1T90MbdADe0fDL0OdW3VHevfrL6iub84037Ati3T2r5MN36Q8CWrl9/yAgi4LT53uq21eObc27dDHBgjIe1/MJ/TvWX1Sc052cLhwA2d7R6W/WyrRsBdu3SlnPR761RTEAE57ox/nH1WdXXNeertm4HdmuMKxrjNY1xvDFe3xi/1DJc+J4tW9c/tDlfsm2TAByM7DxavbA5j23dDrBrR1rWQltluQEBEZzLxrhz9X3Vy6tnbNwN7NcYV1TPbgmDRvVB1d+rfqq6b3P+SHMe37BDAN7lPi3na9PLgO2McYfqga00vaysQQTnum+r3r/6Bz4Bg009vbr4BPc/uDn/Zu1mALhZRw+OAiJgSx/bMqhntYDICCI4V43xiOqx1TOa87c27gb2a4y/0/JJ9IncY81WALhFjlavth4csLEj1axeulZBARGci8a4qPqB6tXV12/bDOzUGKMxvqBlUcGTLQ7/uvUaAuA9GuPC6hMzegjY3pHqd5vzrWsVFBDBuenftMyff1xzXr11M7A7Y9yren71w9Urq39d3fTf4tXV01btC4D35GHV7RMQAVsa47yW89Fq08tKQATnnjEeVH1F9UPN+aKt24FdGeP8xvjillDo0uoLq0c0539omfL52pbRRK+tHtucz92sVwBO5JHVddWLt24E2LX7VXdu2fF2NWPOk416B84ayw5JT29Zz+Ta6m+rezfnX23aF+zJGA+ofrAlGPq56gnNaQoZwNlkjJdWx5rz723dCrBjYzymZQfcD2vO/7VWWSOI4Gz37ttnX1jdrvqULduC3Rjjwsb4murl1X2rz23ZOVA4BHA2GeOu1UMzvQzY3qXVm6s/WrOogAjOfifaPvu2B/cDZ9IYH1P9evVN1X+tHtCcz8nwXICz0WUtH7YJiICtHal+de1rSgERnP1Otk227bPhTBnj4sb49pZtR9+v+kfN+TnN+RcbdwbAe+9o9daW4B9gG8toxvu18gLVJSCCc8HJprGY3gJnwhiPqH67+vKW6Z0PaM6f3rQnAG6dMUZLQPTC5jy2dTvArj3s4CggAk7Zd5zgPttnw+k2xp0a41m9a2ebT2zOJzTnW7dsC4DT4r4to69NLwO2dqRlN8XfWLuwgAjOZmOcV31Wy65lr8/22XBmjPEPq9+rHlN9e/WRzfkLm/YEwOl09OAoIAK2dqR6RXO+fe3CF6xdEDitnlR9QvWvmvOHt24GzjljfED1XdU/r36nZa2h1T/NAeCMO1q9qjn/ZOtGgB0b44LqY6sf2qK8EURwthrjQ6tvrX62+pFtm4FzzBijMf5Fy6ihf1z9m+pjhEMA56AxLqw+MaOHgO09sHqf6le3KG4EEZyNlqllP1xdUz3OltpwGo1x9+qZ1ae1/HJ+THP+3rZNAXAGPazlDdkLtm4E2L0jB8fVF6guI4jgbPXF1cdXX9Kcf7Z1M3BOGOO8xvjCllFDj6i+pPr7wiGAc97RlgVhX/yevhDgDLu0ekMb7Ug9DDyAs8wY961eUb2w+gyjh+C9NMYV1dNbdq358+pt1f1aPkF+vHUoAHZijJdV1zTnx2/dCrBzY/xx9VvN+dlblDeCCM4mY5zfMrXsnS1vYIVD8N5YwqFnV/esRnVJSzj0rOpy4RDATozxvtXHZP0hYGtj/J3qQ9poellZgwjONl/SMi/185rzDVs3A2exp1cXn+D+TxG8AuzKZS0fFAiIgK1denDcLCAyggjOFmN8WPXN1U9Xz9m4Gzh7LduH3vMkj95jzVYA2NzR6qrKLpXA1o60bEL08q0aEBDB2WCZWvYj1dWZWgbvvTHu381/KrPJgoAAbGCM0RIQvbA5j23dDrB7R6rfaM53btWAgAjODk9t2YL1yc35xq2bgbPOGOc3xpe3fCJz7+q7WwLXG7q6etrarQGwmQ+r7p7pZcDWxrht9ZA2nF5WAiI4/JYRD99U/VT14xt3A2efMe5TvaT69urnqg9vzi+pHlu9tpoHx8c253M36xOAtR09OL5g0y4A6qOr21a/umUTtrmHw2xZK+WXq7/b8qb2TRt3BGePMc6rvqj6lpad/55cPccUTQCqGuP/W31Yc95n61aAnRvjqdV/qC5pzj/fqg27mMHh9uXVx1b/XDgEp2CMD6l+qHpEy6ihxzbnn23aEwCHxxgXtvyO+NGNOwGoZQezP9kyHCpTzODwGuPDq2+ofrJ63sbdwNlhjNEYj69+u2Ue92OqfyAcAuAmLq3eJ+sPAVtbFsw/0sbTy8oIIjiclqllP1K9rfpCU2LgFhjj7tUPtqwp8T+qRzenXckAOJGj1XXVi7duBNi9e1SXtPEC1SUggsPqK6qPqR7VnH+xdTNwqC2funx+9V3V+dUXVs8UrAJwM45Wv9qcb9u6EWD3Lj04bh4QmWIGh80YD6y+vnpec/7Ext3A4TbGB1Y/Xf1w9YrqI5vz+4VDAJzUGO/XMg3Z9DLgMDhSvb36na0bERDBYTLGbVqmlr2letK2zcAhtqw19DnV71afXD2l+sTm/ONtGwPgLHBZNRIQAYfDkeplzXls60YERHC4fFX14OqJzXnl1s3AoTTGB1Q/UT23+sPqQc35nc15fNvGADhLHK2uqn5j60aAnRvjfaoHdQiml5WACA6PMT6q+trqx5vzv27dDhxKY/yT6pXVP2wJVD++Of9w26YAOGss69YtmxnMed3W7QC79zEta2huvoNZWaQaDocxLmyZWvbm6snbNgOH0Bh3rb63+pzqt6pPas5XbtsUAGeh+1UfnOllwOFw5OD4a5t2cUBABIfD01qGFn5mc75562bgUBnj06tnV+9XfV3175rz2m2bAuAsdfTg+IJNuwBYHKn+4LC8BzTFDLY2xkdXX109pzn/29btwKExxp0a44ern6n+svrY5vxG4RAAt8Ijqz9qztds3Qiwc8uU10s7JNPLSkAE23rX1LIrqy/ethnY0BhXNMZrGuP4wfHftaw19HnVN1cf05wv37ZJAM5qY9y2+sRMLwMOh/tU79shWaC6TDGDrX1N9ZHVZzTnX23dDGxijCtappBdfHDPPVsWoP6z6tLm/PWtWgPgnHJpy+8aARFwGFy//tChCYiMIIKtjPGQlrWHfqw5f2brdmBDT+9d4dANXSccAuA0Olodq35h4z4Aagmtr6r+YOtGrjfmnIDdnL4AABHwSURBVFv3APuzDHH+jequ1Uc051s27gi2M8bxapzgkdmcPsgA4PQY4zeqq5vz4Vu3AtAYv1O9vjk/detWrufCG7bxtdVHVI8TDrFrY3xmdfwkj75uzVYAOIeN8f7VgzO9DDgMxrhT9eEdoullJSCC9Y3x0Oorqx9pzv++dTuwiTHu1hjPq36qen31jpt8xdUtUzAB4HS4rGW0qoAIOAw+ruWcdGh2MCsBEaxrjNu17Fr2xuop2zYDGxhjNMbnV79ffUZLCHSf6tHVa6t5cHxscz53sz4BONccrd5S/ebWjQC0LFB9vHrZ1o3ckF3MYF1fXz2g+rTmvGrjXmBdY9yrelbLRfovV49pzusX5XvuwQ0ATq8xRsvvnhc253VbtwPQEhD9TnO+betGbsgIIljLGB9X/evqPzXnz23dDqxmjPMa48nVK1t+GX5R9fAbhEMAcCbdv/qgTC8DDoMxzm+ZYnaoppeVEUSwjndNLXtD9WXbNgMrGuP+1Q+2BEM/Xz2hOV+7bVMA7MwjD44v2LQLgMUDqjt2yBaoLiOIYC3fWN2venRzvnXrZuCMG+M2jfHV1Stafvb/ZcvUSuEQAGs7Wv2v5nzN1o0AtHxwWocwIDKCCM60MS6tvrx6dnMa2sy5b4yHVD9UfWT1vOqLm/NN2zYFwC6NcdvqES2/lwAOg0urv6z+eOtGbsoIIjgTxriiMV7TGMerl1RvbgmJ4Nw1xkWN8a0tuzG8f/VZzfnPhEMAbOhIdXHWHwIOjyPVrzTn3LqRmxIQwek2xhXVs6t7VqNlpN4dqk/fsi04o8b4hOq3q6+o/lP1gOb8f7ZtCgA6Wh2rfmHjPgBqjPev7tMhnF5WAiI4E57e8knVDd324H44t4xxp8Z4ZsuF93nVZc35uOa8atvGAKBaAqJfac6/3roRgOphB8dDt4NZCYjgTLjHKd4PZ6cxPr363eqx1X+oHticL9q2KQA4sHxS/+BMLwMOjyMtoxp/Y+tGTkRABKfflSe5/3WrdgFnyhjv3xjPrX6mekt1aXN+eXNevXFnAHBDn3xwFBABh8WR6rea82+3buREBERwOo3xBdX7Vsdv8sjV1dNW7wdOpzHGwRpbv199dvV11UOa82XbNgYAJ3S05YOM39q6EYDGuE310A7p9LISEMHpM8aXVT9cvbB6TPXaah4cH9ucz92wOzg1N9yJbzl+UcuIoedUr6o+ujm/sTmv2bZRADiBMUZLQPQ/mvO6rdsBqD6quqhDukB1LbsrAbfGcgHy76qvrH6i+rzmfGdLWARnn3ftxHf9Yuv3rL6nemf1lOp7XGwDcMjdv7ok08uAw+PIwVFABOekMS6onlk9unpW9SRvnDkHnGgnvqorm/M7124GAN4LRw+OL9i0C4B3ubR6fXO+futGTsYUM3hvjXG76nkt4dC/rZ4oHOIccbId9y5ZtQsAeO8drf6wOV+7dSMAB450iEcPlYAI3jtj3LH62eqzqi9tzn/TnHPjruDWGePCg7W0TvazbCc+AA6/MW5bPSLTy4DDYowPavkQVkAE55Qx3r96UfXwlvWGvmvjjuDWWXYn+8zqd6t/X/129Y6bfJWd+AA4W/y9loVgBUTAYXHpwfHQ7mBWAiI4NWPco/ql6sOrf9Sc/2XjjuDWGeNBLTvv/VR1bfWpzfnRLVMn7cQHwNnoaMvvtF/YuA+A6x1p+QD2FVs3cnOGWTFwC43xgJZPom5ffXpz/tLGHcF7b4y/U31TSxD0V9XXVc9qzmOb9gUAt9YYv1W9rTkfsXUrAFWN8WvVNc358K1buTlGEMEtMcbHVf+zOr96uHCIs9YYt2uM/7P6o+oLqu+o7tOc3yccAuCsN8YHVB+d3cuAw2LZ3OjBHfLpZWWbe3jPxnhky/SbN1VHm/PVG3cEp26MUX129W3Vvar/Vv3r5vyjLdsCgNPskw+O1h8CDouHVLfpkC9QXUYQwc0b459W/716dfXxwiHOSmN8TPWS6nnV26rLmvMzhUMAnIOOtkyd/q2tGwE4cOTgeOhHEAmI4GTGeEL1f1UvrT6hOf98447g1IzxQY3xo9WvV/etHlc9uDlftG1jAHCajXFFY7ym+vzqdtU/27YhgP/t0urVzfkXWzfynphiBje1TMX56pYFfP979ajmvHrbpuAUjHFx9eXVV7ac57+1enpzvm3TvgDgTBjjiurZ1cUH91xcPbsxsgMnsKnlveWRzpJpr0YQwQ2NcV7Lor3fVP2X6rOEQ5w1xjivMT63+sPqG6qfre7fnF8lHALgHPb03hUOXe/ig/sBtvQh1d06C9YfKgERvMsYt6l+tPqS6ruqz2/Oa7dtCm6hMS5tmdf8n6u/aNlt75825x9v2xgAnHH3OMX7AdZy6cHx0K8/VAIiWCxTcn6q+tzqa6qnNOfxbZuCW2CMezbGj7d8KnH3lq3rH9qc/3PTvgBgDcvo77ee5NHXrdkKwAkcqf6meuXWjdwSAiIY484tc0I/rXpic35zc86Nu4Ibu37xzTGOHxz/j8b4t9UfVJ/ZMi3yvs35o8JNAHZhjDtWP1ndubruJo9eXT1t9Z4AbuxI9WvNedNz1KE0vA9m18b4wOrnq/tX/6I5f2LjjuDdvfvim1WzGtVzq/+zOX1KCsB+jPFh1f9T3af6suovW9YcukfLyKGnWaAa2NQYd6iuqr65Ob9263ZuCQER+zXG320ZOXS3lsWoX7BxR3Biy7a99zzBI29szg9cuRsA2NYY/7BlM5FrWnabffHGHQG8uzE+qXph9anN+fNbt3NLmGLGPo3xUdUvtwxJ/iThEIfWGLfvxOFQLeEmAOzDslvn11U/Xb2qeohwCDjEjhwcf23TLk6BgIh9uPH6LX/esqDvtdXHN+fLNu4O3t0Yd2+Mb6tefzNfZVoZAPuwrDf0U9XXt+zY+fGmVwOH3JHq95rzqq0buaUu2LoBOOPeff2Wv9Oyfsu3Nefvb9YXnMgYH1s9tfrsg3t+svrd6iu78RpEFt8EYB/GuF/LekMfWn1J9T02FAEOtWWHxYe1XMufNYwgYg+e3o3fWNeyuO+XbdALvLsxLmiMz26MX65eWn1q9R3VvZvznzXnN1aPrV7bEm6+tnqsxTcBOOeN8Y+ql1V3rS5rzu8WDgFngQ+r7tIyc+WsYZFqzl1jjOrvV794kq+YzSkkZTtj3Kl6dPXFLesM/XH1XdUPN+dfb9kaAGxq+fT966qvrX6j+sfN+afbNgVwC43x6OoHq/s35x9s3c4tZYoZ554x7lJ9XvWElu3rj3fi0XLmrbONMT6kJRR6dHWH6n9WX1r9THNet2VrALC55QOU/1J9evWj1ROb82+3bQrglFxa/VX1v7Zu5FQIiDg3LKOFPrZ6fPXPq4tahiP/q5YpOd+X9VvY0vIzeqR6SvVZLcHl/1V9R3P+5patAcChMcb9W9Ybunf15Or7TCkDzkJHql9tzuNbN3IqBESc3ca4Q3VFy2ihB1Vvb9nZ4lnN+Vs3+LprWtYiukfLyKGnWb+FVYxxm5YFp59SPbR6S/Vt1fc2559t2RoAHCpjfGbLddzVLesNvWTjjgBO3Rh3bZnJ8pytWzlV1iDi7DTGR7WEQp9b3b767er7q+c259u2bA2q66c6Prbl088Prv6o+s7qR5vz7Vu2BgCHyrLe0DdUX1P9est6Q6/ftimA99IYn1r9bPVJzfnirds5FUYQcfYY4+LqUS3TyB5WvaNlis4zq5cafsyhMMZ9Wrbg/YLqfaoXVU+sfvZsG2IKAGfcGHdu+ZT906ofrr6wOd+xbVMAt8qR6rqWwPusIiDi8Fvmoj+++vzqztUftEzX+bHm/KstW2PHxriiG09b/M/VR1b/sDpWPbf6zuZ8xWY9AsBhNsYDWtYb+pDqSdX3+8APOAccqX67Of9m60ZOlSlmHE5j3LZlId8nVJ9QXVv9ZMtooZe4eGBTSzj07G688HnVX7dMI/uPzfnG1fsCgLPFGP+4ZYeyt1ef3Zy/tHFHALfeGBe0rDn6o835RVu3c6pOtPU3nHljXNEYr2mM4wfHKw7uv3djfEv1p9WPt4zO+Krqg5vzc5rzF4VDbGaM0RgfXn1v7x4OVb2lOb9WOAQAJzHG+Y3xb1s++Pu96mOEQ8A55CNa1sj9la0beW8YQcT6Tjz64p3VH7ZM0bmu+unqWdULrNvCppa1ET65+pTq8pYFp09mNqfgHQBOZNnA4TnVp1b/qXpSc75z26YATqMxvrD6vurezfknW7dzqgRErG+M11T3PMEj11XfWP0n23+zmWUnlQe3BEKf0rIg+vnVW6sXVM+vvr76oBP8169tznut0icAnE3G+IiW9YbuUX1x9SyjwoFzyjIQ4gdaNqp5bfW05nzutk2dGgERZ8Yy9/KDq79b3fvgdv2fH3KS/8roC7Yxxt2qoy2B0NHq/Q4e+Y3q5w9uL23OYwdff6JRcFdXjz3bfgkAwBlx480crqzu2LIux2c35y9v2RrAaXeOvD8QELF49x2Z3nPaOcYdu3Hwc8Mg6J7deJe8a6vXVK+u/l51hxM8o9EXrGOM21SX9q5pYw8+eOQvWkYIPb9leuNf3MxznPq/GQDYgxO/UTpePbk5/+M2TQGcZmPcrrqkZWDET/auD5lv6Kx6j3v6RmucbNHhM0W901tr+SV+z2ocHJ/dGP+iMe7RGI9ojH/VGN/cGD/eGC9tjCtbpty8vPq/q2+r/ml1l5ZRF99WPab6pIPnu6g579ucn9qyM9nVN+ni6uppZ+w1cu57T/9mxrhnYzyuMf5r9ebqF6uvaNk95atbRrZ9YHP+y+Z8zs2GQ1VzPrc579Wc5x0chUOwd+fytYJ66t34uS5ujHs1xsc2xqcfXCd+VWN8R2M8p2V9oZtu5nBey+9dzpSz+WdKPfUOW70x7tAY92uMT26ML2iMr2mM72+Mn2mMlzfGX1Z/2zIA4hc7cThUy4fJZ43TM4Jo7eFU6r2n//6C6qLqdgfHm/75pn//1uqut6Cz61rmUr66+uOD27v+POdbb8FzXN+j0RecPif/N/PdLT/rl1f3P7j/T6ufaxkl9MJT+rkFOJmz7VpBPfVuWG+M86v3rT6gutvB8US36x97n5NU+puW0bj3PsnjlhM4Uw7bz5R66h3WestO2XdtGfXzwS3rip7oz3c8wTNfWb2++rOD4+tv8PcfaRlNdFNn1Qii0xUQvaYTLzp8ZfVFt77Au/neTpzQXdmy6N31xk0ef2///oyWX5o39eaWUS/nHdzGrfjzDf/+xE78A3l19QudPOi5/s8XnOC/fW89vncFQX/6v9dggcPk5OegWnbI+8XetZbQH1gUEzjtbv5a6ElnoOL3dfJrIfVOb70TXcve9JrxRG7ua777JPXeXP3rTny9eNPjqTz2pdWdTlDv2uqqg15O1O911V+2hD5vOjie/DbnMkL85P8ezqo3SmcV56Bzud7a76effIO/39x57FQe+86T1Htzy8jCmzufvafjie77kk58zjt2cLvdTe4/Xv15Jw5+rv/zG5rzHSd/xdYgusGzjOPdsl+UnNhs+aG8/nbbm/na32wZynb97R2n+OcTPfbLnXjrbr/EOTuc/Bw0q9v/7wtWgDPFtRBnr2d18gDoLc15/JSf8Rx5o3RWcQ6CW+rf9+4h0BtPy0CIc2CWzJkeQfSG6pNvfYF38z868fCtN7SseXNDN32B783fX9KJt7R+ffVxLaHODUOeW/7nE30D1v7UxS9xznY+qQS2dvPXQo88AxVf0MmvhdQ7vfVOdi17Sy6iT/Y1LzpJvT+rPr53v3a86fHmHrvp18zqT1r79+Q58EbprOIcdC7XW/v99GU3ue/mznW39LFfOEm9P6uOdOvOcTc+zjm9N7gV5py3/lZXzHr7rHmD29tnXXFanl+9M1vvXTVfM+v4wfHM1XJzO923Lf7NuLm5ud3wdq5fK6inntvhvp3rP1PqqXeY651Dt9P9TVgvYFDPzc3thjf/Ztzc3La+nevXCuqp53a4b+f6z5R66h3meufI7fRMMQMAAADgrGWbSQAAAICdExABAAAA7JyACAAAAGDnBEQAAAAAOycgAgAAANg5AREAAADAzgmIAAAAAHZOQAQAAACwcwIiAAAAgJ37/wOoFIonYl9K8wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x432 with 8 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 10\n",
    "gaussian = generate_psd(gaussian_kernel, n=n, s=3)\n",
    "laplace = generate_psd(laplace_kernel, n=n)\n",
    "hilbert = generate_psd(hilbert_kernel, n=n)\n",
    "pascal = generate_psd(pascal_kernel, n=n)\n",
    "fig = make_plot([gaussian, laplace, hilbert, pascal], fig=None)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 7.84 s\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABIgAAAG4CAYAAAAubtzJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeZhdVZk37F9lICSMATESG4J0mKGZQRRlCJBGlO6mFaEdoFGxRUWBtl9RW7pbxTlOSIu+H9raDAIa5UU0IjIIKKASQBCZEkCIIQGSGMhUqfX9caiYoiqBrKTqVLLv+7qeq5J99nPW2rvOWWefp/beqyNJCQAAAACNNaTdHQAAAACgvRSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgYi2O+uss1JKyUEHHdTurgCraNy4cSml5Jvf/OaAtmvcAPqyqmPSCSeckFJKTjjhhB7Lp02blmnTpr2gdQGa6KCDDkopJWeddVa7u8IapEAEAABrib6KVwCwJgxrdwfgnHPOycUXX5yHH3643V0BABpk8uTJ+dWvfpUZM2a0uysA0HYKRLTdE088kSeeeKLd3QAAGmbevHmZN29eu7sBAIOCS8zWAqeeemruuuuuLFiwIH/84x/zla98JRtvvHGvU4w33njj/Ou//muuvvrqPPLII1m0aFEef/zx/PCHP8z+++/f63mf7zr9a665JqWUXsvf+ta35sYbb8zjjz+eBQsW5OGHH85PfvKTHHvssT3W22233XLhhRdm2rRpWbhwYR5//PH85je/yRe+8IUMG/aX2uSK7iXyd3/3d/nOd76TP/zhD5k/f37+/Oc/59e//nXe+973pqOjo1e/vvnNb6aUknHjxuXkk0/OHXfckQULFuRPf/pTzjvvvGy88cYr39HAGrHddtvlk5/8ZG699dY8/vjjWbhwYaZPn57zzjsvL33pS3utv/w17C9/+ctz1VVXZc6cOZk3b15+8pOfZO+9937Bba/quJEkI0eOzL/927/l1ltvzbx58/LnP/85d999d770pS/lxS9+ca91P/jBD+a2225b9vw33XRTjjvuuFXbSUC/GzduXC666KLMmjUrCxYsyK233pqjjjqqxzq19xV6zWtekxtvvDHz58/Pk08+mUsvvTTjx4/vc91VGTeWHw/33XffXHHFFXniiSeW9bGUkm222SbbbLNNSinLYqDvAwes2PLfsXbYYYdMnjw5TzzxRObPn59f/OIXOfzww3usv6rf4ZLkwAMPzOWXX55HHnkkCxcuzIwZM/LLX/4yH/3oR3ut+0KPc1b1+I11kzOIBrmvfvWrOeWUU/Loo4/m61//ehYvXpyjjz46++23X4YPH54lS5YsW3ennXbKJz7xiVx//fX50Y9+lKeeeipbb711jj766Bx55JF53etelylTpqxWfz7xiU/kQx/6UB588MFccsklmTt3brbccsvsu+++ecMb3pBLLrkkSas4dPPNN6eUkssvvzzTpk3LxhtvnPHjx+eUU07JRz7ykXR2dq60rU996lPp6urKzTffnEcffTSbbLJJDj300Hz5y1/Ovvvum7e+9a195n3mM5/JxIkT8//+3//LT3/60xxyyCE5+eSTM378+EyYMGG1th94fsccc0z+5V/+Jddcc01uuummLF68OLvsskve/va353Wve1322WefPPbYY73y9t9//5x55pn52c9+lq9+9asZP358jjnmmLz61a/OEUcckRtuuOF5217VcWPTTTfNNddckz322CP33HNPzj///CxevDh//dd/nZNOOinf//738/jjjydJNtlkk/z85z/PXnvtld/85jc5//zzM2TIkEycODEXXXRRdtlll/z7v//7mtmJwGoZN25cbrnlljz44IP5zne+k8022yxvfOMb88Mf/jCHHXZYrr322urnPuaYY3LkkUdm8uTJufbaa7PHHnvk9a9/fQ455JC84hWvyL333rts3dpx44ADDsiZZ56ZG264Ieeff35e9KIX5d57781//Md/5P3vf3+S5Itf/OKy9adOnVq9PUD/eNnLXpZf/vKX+d3vfpfzzjsvW265Zd74xjfmxz/+cf7pn/5p2femVf0ON3HixPzoRz/KvHnzcvnll+fRRx/NZpttlp122imnnHJK/uu//mvZuqtynFN7/Ma6p4jBGQceeGAppZR77rmnbLLJJsuWDx8+vFx33XWllFKmTZu2bPnGG29cNt98817P89KXvrQ8+uij5e677+6xfNy4caWUUr75zW/22f4111xTSusUomUxe/bs8sgjj5SRI0f2Wn/5tj/3uc+VUko5+uije6236aablo6OjmX/P+uss0oppRx00EE91tt222175XZ0dJRvfetbpZRS9ttvvx6PffOb3yyllPLQQw+VrbbaatnyoUOHLttf++67b9t/r0KsS9HXODJ27Niy3nrr9Vr38MMPL52dneXcc8/tsfyggw4q3d797nf3eOzoo48upZRy77339su4ccEFF5RSSjn33HN7PH+SsuGGG5aNN9542f+7x5gPfOADPdYbMWJE+fGPf1yWLl1adt9997b/ToRocnSPSaWU8tGPfrTHY0cccUQppZQf/ehHy5adcMIJpZRSTjjhhB7rTps2rccx1vLrllLKUUcd1eOxU089tZRSys9+9rMey1d13Fh+PDz55JP73Ma++iaEGDyx/Dj0mc98psdje++9d1m8eHF58skny0YbbVSSVf8Od9lll5VSSvmbv/mbXjnPfZ5VOc6pPX4766yz2r7PxRqNtndArCC+8Y1vlFJKectb3tLrsVe84hW9CkQriy996UullNKjcFJbIHrwwQf7HDyWj+4C0eGHH/68fVvRF70VxZ577llKKeXf//3feyzvPgh729ve1ivnxBNP7PPLpxBi9eL5xpHnxu23314eeOCBHsu6DzCeWwTqju6x6NWvfvWyZWti3Nhiiy1KZ2dnefTRR8uoUaNWmr/ZZpuVJUuWlFtuuaXPx//mb/6mlFLKpz/96bb/ToRocnSPSdOmTStDhgzp9fj06dPLrFmzlv2/pkD03CJQkjJkyJBy3333lVJK2XrrrUtSN250j4e//e1vV7iNCkRCDO7oHoeeeuqpsuGGG/Z6vPs7y1vf+tbnfa6+vsN1F4i22267leauynHO88XKjt8UiNatcInZILbnnnsmSZ+XVfzqV7/qcXlZt1e84hV53/velwMOOCAvfvGLM2LEiB6Pv/SlL80jjzxS3acLLrhg2T2RLr300lx33XX55S9/2esGj9/97nfzvve9Lz/4wQ9y2WWX5Wc/+1luvPHGPPjggy+4rc022ywf+MAH8prXvCbbbrttNtxww17b0pdf//rXvZZ1b/Po0aNfcPtAvTe96U058cQTs/vuu2f06NE97ju2aNGiPnN+8Ytf9Hnfs2uvvTYHH3xw9txzz1x//fUrbXdVxo199903Q4cOzfXXX59nnnlmpc+77777ZtiwYcvuDfJcw4cPT9I6TRxov6lTp6arq6vX8kceeSQHHHDAaj33dddd12tZV1dXbrjhhowfPz577rlnHn744dUaN2655ZbV6iPQfr/97W8zf/78XsuvvfbanHjiidlzzz3z7W9/O8mqfYe74IIL8o//+I+5+eab893vfjfXXHNNbrzxxjz66KM9clblOKdbzfEb6xYFokFsk002SZLMnDmz12NdXV29Zv76+7//+1x22WVZuHBhrrrqqjzwwAN5+umn09XVlYMPPjgHH3xwr8FmVZ122ml54IEHctJJJ+XMM8/MmWeemSVLluTKK6/MGWeckQceeCBJcuutt+ZVr3pVPvzhD+f1r3/9svt+3HPPPfnP//zPXHzxxc+77bfeemu23Xbb3Hzzzfn2t7+dJ598Mp2dndl0003z/ve/f4XbMmfOnF7Luu93NHTo0NXZfOAFmDRpUk477bQ89thjmTJlSh599NEsWLAgSXLiiSdmm2226TOvr7EuSf70pz8l+cuYuCKrOm5suummSdLrgKovm2++eZJkv/32y3777bfC9Z5bkALao69jgaR1PLC6xwIvdKxanXGj+7mAtdcLHStW9Tvc5MmTc9RRR+WMM87ISSedlH/5l39J0vojefe9HJNVO85J6o/fWLcoEA1i3WfljBkzpsdsZUkyZMiQbL755j3e8B/72MeyePHi7LPPPrnnnnt6rL/lllvm4IMP7rGs+y9ry1eGl9c9qDw358tf/nK+/OUvZ4sttsiBBx6Y4447Lscee2x22WWX7LLLLlm8eHGS1llOr3vd67Leeutl7733zt/+7d/mve9977IZRa6++uoVbvvb3/72bLvttvmP//iP/Od//mePx17+8pcvu0EjMLhsscUWOfXUU3PnnXfmFa94Ra+/nB1//PErzB0zZkyfy1/ykpckSebOnbvStld13Oj+AvlCZubobnvSpEk544wznnd9YN31Qseq1Rk3+jqbEli7vNCxYlW/wyXJlVdemSuvvDKjRo3K/vvvn9e+9rV517velSuuuCJ77rlnfv/736/Scc7qHL+xbjHN/SB22223JWlNY/hcL3/5y5edmtxt/Pjxufvuu3sNLB0dHX0+x1NPPZUk2WqrrXo9ttFGG2X77bdfaf9mzZqVyZMn541vfGOuvvrqjB8/Prvuumuv9RYvXpxf/vKXOeuss3LqqacmaU1FvTLdU8V+73vf6/XYQQcdtNJcoH223XbbDB06ND/96U97HVy89KUvzbbbbrvC3AMPPLDPqei7D4y6x8QVWdVx45ZbbsnSpUvz6le/OqNGjVrpc3ev+6pXvWql6wHrvr7GkyFDhiw71uoeq/pr3Fi6dKkzomEtsNdee/V5huBzj2tW9Tvc8p555plcc801OeOMM3L22WdnxIgROfLII5Os2nHO6hy/sW5RIBrEuq9J/fCHP5yNN9542fLhw4fn7LPP7rX+9OnTs91222XLLbfssfyss87KLrvs0mv9+fPn5/e//31e+cpX9rj+fciQIZk0aVKvgWS99dbLoYce2ut5hg0bls022yxJll3feuCBB/boc7fuSvrzXQc7ffr0JOlVMd9jjz1y5plnrjQXaJ/u9+6BBx6YIUP+8hGzwQYb5Bvf+Eavwvbytt9++5xyyik9lh199NE5+OCDc9999+UXv/jFC2r7hY4bs2fPzsUXX5yxY8fmc5/7XK/i1AYbbLBsHJs1a1YuuOCC7LvvvvnIRz7S55ezbbfd1unX0AATJkzIUUcd1WPZe97znowfPz4///nP8/DDDyfpv3HjiSeeyBZbbJH111+/ehuA/rfpppvmox/9aI9le++9d970pjdlzpw5mTx5cpJV/w536KGH9vn+f+73rFU5zlmd4zfWLS4xG8Suv/76nHfeeXnnO9+Zu+66K9/73veyZMmSvO51r8vcuXPz6KOP9rgB4xe+8IWcd955ue2225at+8pXvjI777xzLr/88hx99NG92vjsZz+b888/PzfeeGMuvfTSLFy4MIccckiGDx+eqVOnZo899li27siRI3P11Vdn2rRpufnmm/PQQw9l/fXXz+GHH56dd945P/zhD5dVvs8444wcccQRufbaa/Pggw9m/vz52WWXXXLkkUfmySefzNe//vWVbvu3v/3tfOADH8gXv/jFHHLIIbnvvvuy3Xbb5bWvfW2+//3v57jjjltDexlYk2bOnJmLLrooxx9/fKZOnZqf/vSn2WSTTXL44Ydn4cKFue2225bdgP+5fvzjH+fzn/98jjzyyNx+++0ZP358jjnmmCxYsCBve9vbnveSi5px4z3veU923XXXvOtd78rBBx+cKVOmZPHixXnZy16WiRMn5uijj152Q9r3vOc92W677fKxj30sb3nLW3LDDTdk5syZGTt2bHbaaafst99+Oe6445YdZAHrpssvvzyTJ0/O5MmTc//992f33XfPUUcdlSeeeKJXkbs/xo2rr746++23X37yk5/k+uuvz6JFi3L77bfniiuuWMNbCqyO6667Lm9/+9uz//7758Ybb8yWW26ZN77xjRkyZEje+c535s9//nOSVf8O9/nPfz7bbLNNrr322kyfPj2LFy/O3nvvnQkTJmT69Ok97vX6Qo9zVuf4jXVP26dSEyuOjo6O8v73v7/8/ve/LwsXLiyPPvpoOeecc8rGG29c5s2bV2677bYe659wwgnltttuK/Pnzy+zZs0q3//+98uuu+660imhTzrppPK73/2uLFy4sMyYMaN87WtfK5tttlmvae6HDRtWPvCBD5Qrr7yyPPTQQ2XBggXl8ccfL7/85S/LO9/5zjJ8+PBl6x5++OHl/PPPL3fddVeZM2dOmT9/frnnnnvKl770pWXTv3bHivq20047lR/+8Idl5syZZf78+eXXv/51edvb3rbCabW7p4wcN25cr200DaMQ/RN9vR9HjhxZPv7xj5f77ruvLFiwoDz88MPlnHPO6XNcSXq+P1/+8peXq666qsydO7fMmzevTJkypeyzzz692l1T40aSMmrUqPKhD32o3H777eXpp58u8+bNK3fddVf5whe+ULbYYose6w4fPry8+93vLjfeeGOZM2dOWbhwYXnooYfKz372s/K+972vbLbZZm3/nQjR5FjZez1JrzGoZpr7E044oRx11FHlpptuKvPnzy9PPfVUueyyy1Y45fSqjBsv5Hhl1KhR5dxzzy2PPPJIWbJkyUq3Vwgx8LH8OLTjjjuWH/zgB+XJJ58sTz/9dLnhhhvKEUcc0StnVb7DveENbygXXnhhuffee8uf//znMnfu3HLnnXeWj3/84+VFL3pRr+d+occ5q3P81u59LtZotL0DoiLGjx9fSinlwgsvbHtfhBBidcIBhhBCCCHWlXi+QrUQgzncg2iQGzNmTK/rRUeOHJkvfvGLSbLs2lUAAACAWu5BNMi9//3vz/HHH59rr702M2bMyEte8pJMmDAhW221Va688spceuml7e4iAAAAsJZTIBrkrrrqquy+++454ogjstlmm6WzszP33ntvvvzlLy87iwgAAABgdXSkda0ZAAAAAA3lHkQAAAAADVd9idnSGeOr8q5eMLQq7/Q7j63KyzWjq9LG3PJ0Vd6w+x+ryuuaM7cqr3QuqcpLceIY/etn5bJ+ff7aMWji2D2q8maffEBVHtAeU887o9/bqB2Hxl/wrqq80XdXpWWDmZ1VeevNrTvGGPr04qq8jkV17dXmZenSurzOVc8rnXW/g+o+1h7n1bbXVddeqe1nV1dd3gC399MF/1vX3gt0WMfr+/X5gbVbzfcxZxABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDDatNvHrB0Kq8CSOXVuVN2u2SqrzTc2xV3syMrsobk7FVecPur0pL15y5VXmlc0ldg6XU5cEaNnHsHlV5Ux6bWtleVVpmn3xAXSIw6I2/4F1VeaPvrmtvg5mdVXnrza37zB/69OKqvI5Fde3V5mVp3bFlOuvySmfF76G2j7XHXbXtddW1V2r72dVVlzfQ7QE0hDOIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpuWG3i6XceW5U3abdLqvImjFw6oO2dnrrtm5nRVXljMrYqb9j9VWnpmjO3Kq90LqlIKlVtwcrMPvmAqryJdW+1THls6oC2V7t9wMAZfXdd3gYzO6vy1ptb8RmcZOjTi6vyOhbVtVebl6V1x3rprMsrnXW/h6p+1h4L1e6Trrr2Sm0/u7rq8ga6vUrV+wVgLeMMIgAAAICGUyACAAAAaDgFIgAAAICGUyACAAAAaDgFIgAAAICGUyACAAAAaDgFIgAAAICGUyACAAAAaDgFIgAAAICGUyACAAAAaDgFIgAAAICGUyACAAAAaDgFIgAAAICGG1adec3oqrTTc2xV3qTdLqnKmzBy6YC2V7t9M1O3P8dkbFXesPur0tI1Z+4q55TOJXWNlVKXBysx++QDqvIm1r3VMuWxqQPaXu32Aatug5mdVXnrza37XBz69OKqvI5Fde3V5mVp3bFXOuvySmfd76G6nzXHJ7VtddUdC5XaY6iurrq8gW6vUv1+cUwKNIMziAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAablht4phbnq7Km5nRVXmn59iqvEm7XVKVN2Hk0gFtr3b7avfnmIytyht2/6rndM2ZW9VW6VxSlZdS6vJgJWaffEBV3sS6t1qmPDZ1QNur3T5osvXm1n1ODX16cVVex6K69mrzsrTuWCiddXmls7Ouvdp+1h4v1LTXVddWqe1jV1dd3kC3V6l+vzhGBFgZZxABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANJwCEQAAAEDDKRABAAAANNyw6sT7H6vKG5OxVXkzM7oq7/QcW5U3abdLqvImjFw6oO3Vbl/t/qz5/Q27v6qpdM2ZW5VXOpfUNVhKXR6sxOyTD6jKm1g3VGbKY1MHtL3a7YN1wdCnF1fldSyq+5yqzcvSumOTdNbllc7OuvZq+1n7+V3bXteqt1dq+9jVVZc30O1Vqt8vA3zMVgZ2vwC0izOIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpuWG1i15y5dQ3eX9femIytypuZ0VV5p+fYqrxJu11SlTdh5NIBba92+2r2Z+3vrva1UvvaLJ1L6hospS4PVmL2yQdU5U2se7tlymNTB7S92u2DwaRjUd3nRm1eltYdK6SzLq90dta1V9vP2s/T2va66torNf3s6qpqq3qf1LZXqWqfJNW/g2plYPcLwNrGGUQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADTesNrF0LqnK65oztypv2P1VaRmTsVV5MzO6Ku/0HFuVN2m3S6ryJoxcOqDt1Wxf7b6s/d3VvlZqX5u174WUUpcHKzH75AOq8ibWvd0y5bGpA9pe7fZBf+hYVDn+L6377E5nXV7p7Kxrr7aftZ9vte111bVXavvZ1VXT2MC1tRrq98kAH9OUgd0vAE3hDCIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhhtWnVlKXVrnkqq8rjlzq/KG3V+VljEZW5U3M6Or8k7PsVV5k3a7pCpvwsilA9Ze7bbV7sva313ta6X2tVn7Xqh978HKzD75gKq8iXVvt0x5bOqAtVe7bfC8ltZ9lqazLq90dta1V9vP2s+b2va6Ko8ta/vZ1VWXV9NebVuV6vfJAB9jlIHdL7XKQO8XgDZxBhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADTcsAFvsZS6tM4lVXldc+ZW5Q27vyotYzK2Km9mRlflnZ5jq/Im7XZJVd6EkUsHrK3abavdl7W/u9rXSu1rs/a9UPveg5WZffIBVXkT695umfLY1AFrq3bbaJDOVf9MTJLS2VnX3tK69qrH/9r2uiqP9Wr72dVVlzfQ7VWo3ycD/JlfBm6frI4y0PsFYC3jDCIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhhvW7g68YKXUpXUuqcrrmjO3Km/Y/VVpGZOxVXkzM7oq7/QcW5U3abdLVjlnwsilA9ZWUr9ttfuy9ndX+1qpfW3Wvhdq33uwMrNPPqAqb2LF223KY1MHrK2kfttY+5TOzrrEpXWfi9XjcW17XZXHXrX97Oqqyxvo9irU75MB/gwuA7dPVkexXwD6hTOIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpuWLs70O9KqUvrXFKV1zVnblXesPur0jImY6vyZmZ0Vd7pOXaVcybtdklVWxNGLq3Kq22vZtuS+n1Z+7urfa3UvjZr3wupe+vBSs0++YBVzplY91bLlMemVuXVtlezbbTZ0rrPqdpjk+r2uiqPhWr72dVVlzfQ7VWq2i+Vv4NqZWD3Sa1ivwAMKs4gAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGi4Ye3uwKBVSl1a55KqvK45c6vyht1flZYxGVuVNzOjVznn9Bxb1dak3S6pypswcumAtle7fTX7Mqn/3dW+Vmpfm1lclwZr2uyTD6jKm1j3VsuUx6YOaHu128caUHmskKV1n1Ppqjw2qe1nV1dd3kC3V6l+v1Tm1SgDu09qlYHcJ60GB7Y9gO5STakAACAASURBVDWho6Myr/K8norDDWcQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADScAhEAAABAwykQAQAAADTcsHZ3YJ1TSl1a55KqvK45c6vyht1flZYxGbvKOTMzuqqt03NsVd6k3S6pypswcumAtle7fbX7s+Z3l9S/VvJ4ZR4MErNPPqAqb2LdWy1THps6oO3Vbh/LWVr3uZGuymOFymOMdHXV5Q10e5Xq90tlXo0ysPukVhnIfdJqcGDbA1heR8fANjd0aF1ibV7FYYoziAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAaToEIAAAAoOEUiAAAAAAabli7O8CzSqlL61xSldc1Z25V3rD7Vz1nTMZWtTUzo6vyTs+xVXmTdrukKm/CyKUD2l7t9tXuz9rfXx6vS4O13eyTD6jKm1j5Vpvy2NQBba92+9ZJXZWf3ZWf+enqqssb6PYq1e+XyrxaZWD3S41in6xZta9NYGB0dFTm1Z0v0zF0aF3e8LryS8fI9avysnDVU5xBBAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBww9rdAVZTKXVpnUuq8rrmzF3lnGH3VzWVMRlblTczo6vyTs+xVXmTdrukKm/CyKUD2l7t9tXuz9xUlwZNNfvkA6ryJtYNlZny2NQBba92+wazUvkZnK6u2gYHtr1K9fulMq9WGdj9UqPYJ2tW7WsTGBgdHZV5dee9dAwdWpe33vC6vFGjqvKy2SZ1ebNXPcUZRAAAAAANp0AEAAAA0HAKRAAAAAANp0AEAAAA0HAKRAAAAAANp0AEAAAA0HAKRAAAAAANp0AEAAAA0HAKRAAAAAANp0AEAAAA0HAKRAAAAAANp0AEAAAA0HAKRAAAAAANN6zdHaBNSqlL61yyyjldc+ZWtTXs/qq0jMnYqryZGV2Vd3qOrcqbtNslVXkTRi4d0PZqty9frEsDVs3skw+oyptYN1RmymNTB7S92u0bEF1ddXmVn8HV7VUq1f2szKtVBna/1Cj2yZpV+9oE1kkdQ4fW5a03vC5vww2q8vKiuu+b88dvUtfePaue4gwiAAAAgIZTIAIAAABoOAUiAAAAgIZTIAIAAABoOAUiAAAAgIZTIAIAAABoOAUiAAAAgIZTIAIAAABoOAUiAAAAgIZTIAIAAABoOAUiAAAAgIZTIAIAAABoOAUiAAAAgIYb1u4OsJYpZdVTOpdUNdU1Z25V3rD7q9IyJmOr8mZmdFXe6Tm2Km/SbpdU5U0YuXRA2zuiKgsYKLNPPqAqb2LdUJkpj00d0PYGRMVnYpKkq2vN9uN5lOp+VubVKgO7X2qVgdwva8k+qVb72gQGt46OurShQ+vyhteVNTpGjarKy4vqvv/9eYdNq/IeO7Buf+YHq57iDCIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhlMgAgAAAGg4BSIAAACAhhvW7g7QAKXUpXUuqcrrmjO3Km/Y/VVpGZOxVXkzM7oq7/QcW5U3abdLqvImjFxalQesm2affEBV3sS6oTJTHptalTd0y7r2VklX1wA08hel8vM0XZV5tcrA7pdaZSD3y1qyT6rVvjaBdVNH5XkoQ4fWNTdy/br2NtukKm3++Lq8xw7sqMp74LivVeUNPW3Vc5xBBAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwCkQAAAAADadABAAAANBwHUlKuzsBAAAAQPs4gwgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIAAAAABpOgQgAAACg4RSIgFV2fJJpSZY++/P49nYHaCDjENBOxiCg3fpjHOpIUtbA8wANcXySbyTZYLllTyd5R5KL2tIjoGmMQ0A7GYOAduuvcUiBCFgl05Js08fy6UleNqA9AdYFxyc5O8nWSR5O8qH0PrAZkWS7JDs9G/+WngdE3abHOAT0P8dCQLv11zikQASskqXp+9rUriRDB7gvwNqtr79+LUjy7STz0ioG7ZjWgU73+NKV1sFLRx/PZxwCBoJjIaDd+mscUiACVsmfkozpY/n0+KsZ8MINS+uMoS1X8PjCJPcmuSfJ75f7eW+Su+Ov90D7OIMIaLf+GoeGrUYu0DB7JtkkrYr18pXpp9O6LARgRTqS7JZkwrPx6iQbrWDdrrTOKupaweMfSt/X3RuHgIHwoST/N8mo5ZYZg4CB9KEk30qy3nLL1sQ4ZBYz4AXZOsmPkjye5L1pVae7nv3ppoxAXzNpjE/yziTfTTIzye1JJj27/DtpjSd9eTgrLg4lrfHmHTEOAe1xUVr3TkuMQUB7XJTkpiSdWbPjkEvMgOe1SZIbk7w0ySvTurwDoFtf9xLqyl/+CvXHJFc/Gz9P8uhK8swEBKwN3pxWoXvHJH9oc1+AZvpZkvWTHLgGn9MlZsBKrZdkclozCE2M4hDQ0+gkX0nvWcWGJHkiyQFJ7ltBbncR6PlmMQMYbHZMsiTJA+3uCNBY45Lcuoaf0yVmwEqdn+SQJP+c5Nr2dgUYJNZP8oYkP0jrxvWbr2C90VlxcajbRfnLLGUvi+IQsHbYMa3iUGe7OwI0UkeSrdL649qapEAErNDHk7wprb/oX9jmvgDtNSTJYUm+mdb9hC5Jsm9aZw89toKcNX3QAjBYuLQMaKcXJxmR5KE1/LwuMQP69I4kH07y9SSfbHNfgIFzfHpe8vWNJC9KclxaU9LPTXJZkgvSOquwK8lvYlYxoDmGpnXp/RXt7gjQWOOe/bmm/xinQAT0cmSSc5NcmeSUNvcFGDjPvWn0Nkk+kdZ9Nq5Iqyj0oyQLn5PnXkJAk7wsrXs03tPujgCNtfWzPxWIgH61V1qXjtye5Ni0pqwGmuEL6X2z6SSZkeSY58m9KApCQDPs8OxPBSKgXfqrQOQeRMAy49I6S+CJJK9N6xIRYN02Iq2b0N+eZMwK1vmrgesOwKC347M/3YMIaJet07rsf+4afl4FIiBJsmlal5SNTOsSsz+1tztAP3txkrPSurnh+UlKktkrWNfNpgH+Yse0btb/VLs7AjTWuPTP8ZkCETTY8UmmpXUZ2Ywk2yf5+yS/b2engH61a5L/L62Div9IckuSQ5PskeTU9D5z0M2mAXraMS4vA9qr+56Pa5oCETRU981ot0lrIFg/rULR2Db2CVhzli8AT0vy6SRXJbkzyRvTKhLtkOToJNc8m3NRWjMYTk9rdrLpz/7fvYUA/kKBCGi3rbPmp7hPko60zioHGmZaWsWh55qe1uwcwNrrubORdXsiyWeTfD0ujQCosXlal+OeluSLbe4L0Eyj0jrD+4Np/QFwTTKLGTTU1qu4HFh7fCp9z0Y2P2v+QAKgSbpvUO0MIqBd+msGs8QlZtBYKxpQ3IwW1l4jkrw3yVYreHxFywF4YRSIgHbrLhD1xyVmCkTQUFf1sczNaGHtNCytewXdl+TLSRauYD0FYIDVs2OSBTGeAu3jDCJgjXpJkjckuStuRgtrsyFJ3pzWX7K/nuSPac1I9raYjQygP+yY5N60jp0A2mFcks60ZqFe09yDCBroy2nNWvYPaZ1xAKxdOpIck+S/kuyc5LdJXpPkx89Z7+z8ZRrUD0UBGGB17ZjkN+3uBNBoWyd5NK2Zatc0ZxBBw7wurbOHPhbFIVgbPHe6+s+k9eXksrSmIf3HJPukd3HoorRmJBz67E/FIYDVMyKt8dT9h4B26q8p7hNnEEGjbJTk3CR3pjXVNTC4PXe6+m2SfCDJn9K6tOyiuMwBYKCMT6vorkAEtNO4JDf203MrEEGDnJ1kbFpnHCxpc1+A53d2+p6uflGSCwa4LwBNZwYzoN2GJPmr9N+N8l1iBg3x8iSnJDknyS1t7gvw/LZI6y9EfTFdPcDA6y4Q3dvWXgBN9pIkw9N/l5gpEEEDDE/rMpU/Jvlwm/sCrNzwJKdl5fcIM70ywMDbIa3x95l2dwRorP6c4j5RIIJG+D9Jdk3rDKL5be4LsGJ/m+SOJJPSurb8X2O6eoDBYse4vAxor+6zyxWIgCo7JPlIku8m+VGb+wL0bbskV6Q1E1lHWlPWH5VWoegdSaandTPq6c/+34xkAANPgQhot/4+g8hNqmEd1pHk62mdCv2+NvcF6G3jJP+e5NQkC5OckeQr6XkT+YuiIATQbmPTmg1WgQhop62TPJn+uyrEGUSwDnt7klendZnKzDb3BZru+CTTkixN60yg89K60enpSb6d1llEk2KGQYDByAxmwGAwLv17L0pnEME66iVJPpPkmiTnt7kv0HTHp3Wj+O4p68clOTnJH9K6nOy3beoXAC+MAhEwGGyd/pvBLHEGEayzvpJk/bS+hALtdXb+Uhxa3ogoDgGsDXZMMi/JjHZ3BGi0/i4QOYMI1kFHJ3l9WjMd3d/mvkDTdeQvM04819YrWA7A4OIG1UC7bZRkdPr3EjNnEME6ZqMkX01rquzPtrkv0HS7JLkhrSJRX/rzAx6ANUeBCGi3/p7BLFEggnXOJ9OaaeMdSTrb3BdoqvWTfCLJbUm2T/LfSZ5+zjpPp3WWHwCD24ZJtooCEdBeCkTAC7L87EinJJmS5Ja29gia67Akv0ur+HNBWn91PiWtou30JF3P/nxHTF8PsDbY/tmfCkRAO3UXiNykGlih7tmRtknrDd2R5KBnlwMDZ4sk30lyVVrF2kOS/HOSJ559/KIkL0sy9NmfikMAawczmAGDwbgki5P8qR/bUCCCtVxfsyONenY5MDD+Oa0vDscm+a8kf5Pk2nZ2CIA1Zse0Ltt/oN0dARpt6yR/TFL6sQ0FIljLrWgWJLMjQf9Y/pLOPya5K8n5z/7cI8lZSRa1rXcArGk7Jnkwrb/cA7RLf09xnygQwVpvRTcpMzsSrHnPvaTzpUl2enbZQUl+37aeAdBfzGAGDAbj0v/f8RSIYC13WR/LzI4E/aOvSzo7khye/j3dF4D2GJJkuygQAe01NK0/TCoQASs0MsnrkzyS1umGZkeC/jMqrb/c9MUlnQDrpnFJ1k/yh3Z3BGi0sWkVifq7QDSsn58f6EcfSutSl4OSXN/ersA67dVJ/r+0zhbqi0s6AdZNZjADBoOBmOI+cQYRrLW2S/KBtKbVVhyC/rFBkq8kuS6t4tDH0rqEc3ku6QRYdykQAYNB91nsLjED+nROkoVpFYmANe+QJHcmOSXJF9Oauv6jaV3COT0u6QRogh2TzEryZLs7AjRa9xlEj/RzO4OuQHTaaaellJLjjjtutZ/rf/7nfzJz5syMGjVqtZ9rr732SiklJ5100mo/F6yu1yc5IslHksxsc19gXbNRkv9O8vMkS9K6vOy0JM88+/hFSV6W1nXgL4viEMC6zAxmwGCwdVrF6meeb8XVNOgKRPvss0+S5De/+c1qPc/ee++dN7/5zfnUpz6VZ55Z+W7cYIMN8sgjj6SUksmTJ/e5zm9/+9tMnjw5H//4x7PBBs+dwwYGzoZJvpDktrS+xAJrzuFJfpfk5CSfS7JHkhvb2iMA2kmBCBgMBmKK+25lMMVWW21Vdthhh9V+nilTppSnnnqqrL/++s+77mc/+9lSSilLly4tDz300ArX23fffUsppZx55plt30+iufGZpJSkvHwQ9EWItTmOT8q0pCxNykNJuTqt99bd8f4SQgiRsllanwunD4K+CCGaHXcm5fsD01b7N3ZNx3bbbVeWLl1azjvvvOddd6eddiqLFi0qV155ZbnttttKKaVsvvnmK1z/7rvvLtOnTy9Dhgxp+3aK5sXOSVmclG8Mgr4IsTbH8UmZn9aBf3d0JeUHSRkxCPonhBCi/XFAWp8PrxkEfRFCNDvmJeWLA9DOgFxidswxx+TKK6/MrFmzsmjRotx7770588wzM2RIz+YPPvjglFLy6U9/usfyESNG5P/8n/+Tu+66KwsWLMj06dOX5c+bNy+33357j/VPOumkDBkyJN/97neft2/nnHNOOjo6cvrpp+eOO+5Ikuy5554rXP/iiy/OuHHjcthhh73QzYc15twk85J8sN0dgbXc2WnNULa8jiS7J1k08N0BYBAygxkwGGyS1j0y+3uK+6Sf70E0ZMiQXHjhhfne976X8ePH59JLL825556brq6unH322Tn//PN7rL/XXnslSW677bZly0aNGpVrrrlm2b2EvvSlL+Waa67JRz/60XzrW9/KRhtt1GP9JDnssMPS2dmZX/3qVyvt3/HHH59DDz005557bu65557ceeedPfrRlxtvbN2N4vDDD3/hOwLWgDclOSjJmUmeaHNfYG239SouB6B5dkzrjwbT29wPoNnGPftzrb8H0Ve+8pVSSilnn312GTp06LLlw4YNKzfccEMppZSddtpp2fL//d//LaWUsv322y9b9p3vfKeUUspHPvKRHs996KGHlm6nnnrqsuWjRo0qS5YsKXfcccdK+7bhhhuWP/7xj2XWrFll0003LUnKEUccUUop5eKLL15h3sYbb1xKKeXmm29u+2lmojmxSVJmJOVXSekYBP0RYm2NDZPytfS8tGz5mDYI+iiEEGJwxA+Tcscg6IcQotnx2rSOU/cdmPb654n322+/snTp0jJ58uQ+H3/HO95RSinlxBNPXLbs7rvvLvPmzSsdHR0lSdl///1LKWWFzzF9+vRSSikHHnjgsmXbbbddKaWUKVOmrLR/n//850sppZxyyinLlm255ZallFL+8Ic/rDT3mWeeKTNmzGj7C0U0J76U1o109xoEfRFibY2DkvJgWu+ly9P7HkTz07o3Ubv7KYQQYnDEH5JyySDohxCi2XFKWseqLx6Atoaln7z3ve/NkCFD8swzz+Sss87q9fiuu+6aJOno6EjSupRs++23z0033ZRSSpLk3e9+d5LkE5/4RJ9tPPHEE9lqq60yderUZcs233zzJMlTTz21wr7tvPPOee9735vf/e53+drXvrZs+YwZMzJ79uyMHz8+G220Uf785z/3mf/kk09mzJgxK3x+WJP2SPLutKa0/22b+wJro5FJPpnkfUnuS/KqJDclOT6texFtndYpux9KclGb+gjA4LJekm2TPP8dTQH617gkC5PMGoC2+q1AdMQRRyRJ/umf/mml6z38cOtKut133z1Dhw7Nb3/7l6/ARxxxRGbPnp1f//rXfeZuueWWue+++zJ//vxlyxYsWJAkWX/99VfY5le/+tUMHz48p512Wrq6uno8duedd+aQQw7JHnvskV/84hd95o8cOXJZO9CfOtK6MfUTST7S5r7A2ujlSf4nyfZJvpzWPbyeefaxi6IgBEDf/jqtL0puUA2029ZJHknrFJ/+1i83qR4xYkRe/OIX57rrrktHR8dK4+qrr07S+wbVI0aMyJgxY/LII4/02cYOO+yQLbfcskdBKUkef/zxJH85k+i53vSmN+Xggw9Oklx11VUppfSIQw45pEd/nqujoyObbrrpsnagP/1zkgOSfCDJnDb3BdYm66V11tANSUYkOTStM4ieWVkSADxrh2d/KhAB7dZ9tvtA6JcziLovG3vRi170gnO6CzLdBZ/Ozs4sWbIko0eP7nP9f/u3f+uxfrcZM2bk8ccfzw477NArZ6ONNspnP/vZLF68ON/5zneWXcq2vG222SaHHXbYCqe632GHHTJkyJAel7VBf9gsyaeT/CLJt9vcF1ib7JXWWUO7JvlGkjOS9H3BMAD0rXuK+z+0tRcArQLRTwewvX65udHUqVNLKaX8wz/8Q5+Pv/KVryxDhgxZ9v/bbrutPPPMMz1mO/vd735XSinl0EMP7XmTplNOWTaD2XMfS1IuvfTSUkopf/3Xf91j+Re+8IVSSimf+tSnVtjvvffeu5RSVjgL2oknnlhKKeXd7353229WJdbt+FpSliRl10HQFyEGaxyf1sxjS5MyPSmXJWVxUh5Nyt8Ogv4JIYRYO+NbSXlkEPRDCNHsGJ7Wce5ZA9dm/zzx4YcfXhYvXlxKKeWnP/1p+dznPlcmTZpUvvvd75b777+/PPTQQ8vWXW+99cqiRYt6TR3/pje9qZRSyoIFC8q3vvWt8qlPfarccMMNZebMmeXuu+8uS5cuXTZF/fJx3HHH9ZqhbNdddy2LFy8uDz30UBk1atQK+z1ixIiyZMmSsmTJkrL++uv3evzCCy8sS5YsKX/1V3/V9heLWPdi+S+7XUm5YhD0SYjBGsen90xkJSm/SMroQdA/IYQQa2/8KilXDYJ+CCGaHdukdXz7zwPXZv89+T777FMuvfTSMmPGjLJ48eIya9ascscdd5Svfe1rPc782WuvvUoppfz3f/93r+d4z3veU+6///6yaNGi8tBDD5VzzjmnbL311mX27NnlN7/5TZ/tDh8+vMyYMaP86le/WrbsuuuuK6WU8vd///fP2++77767lFLKfvvt12P5xhtvXJ555pkyefLktr9QxLoXfX3ZfTqm3RZiRTEtvYtD5dnl7e6bEEKItTvmJOUrg6AfQohmx6vTOr49dODabP9Gr2q8+c1vLqWUcsYZZ6xwnQ9+8IOllFL22GOPNdbue97znlJKKQceeGDb94FY92JafNkVYlViafp+zywdBH0TQgix9sZL0vo8efcg6IsQotnxlrTGo/ED1F7Hs/8YdIYMGZLNN988s2bN6rF8woQJ+cEPfpAnn3wyO++8c55++uk+80eMGJE//OEPueOOO3L00Uevdn/WX3/9PPDAA7npppvyhje8YbWfD55rafqeVrArydAB7gsMZh1J3pvkC+n7PTM9ycsGskMArFMOTnJNksOSXN3ergAN9+EkH08yMsnCAWivX2YxWxN23nnn3HrrrZkyZUoeeOCBDB8+PHvssUde9apXZdasWfm7v/u7FRaHkmTRokV5y1vekkMOOSSjRo3KM8+s3uTG22yzTb7+9a/nW9/61mo9D6zII0nG9bF8oKY0hLXBNknOT3JIktvSmoZ41HKPP53kQwPfLQDWId0zmJniHmi3rZPMzMAUh5LWH2IH5RlE22+/fT75yU9m//33z+abb56urq5MmzYtV1xxRSZNmpTHH3+83V2ENeriJG98zrKnk7wjyUUD3x0YdN6R5PNpfWi9P8k3kxyf5Oy0PjwfTqs45P0CwOr4YpK3Jdmo3R0BGu/HSTZLsv8AtTdoC0TQJHskuSXJr5NsGV92YXkvTfJ/k/xtWqf6nxRn1gHQf36SZPMk+7a7I0Dj3Z3kriQDdZObQXuJGTTFiCTfTjIryWuTPNne7sCg8uYkX0kyPMm7k/x3/FUDgP61Y5JftLsTAGmdOPDjAWyvr/t7AgPoY0l2S+tUZsUhmur4JNPSuln7tCTvTDI5yXeS/C7J7knOjeIQAP1rVFr3hHT/IaDdNkuyQZKHBrBNZxBBG70qyRlJvpbW6czQRMcn+UZaH4BJ60bU/51kSVrvjy+mNZsfAPS37Z/9qUAEtFv3BEYDeWsFBSJok42S/E+S/7+9+w+SpCwPOP5dDu44jiAxdyByt7cnGE7KRCktkiiXEKmgqJhCI7AIBAgk0YQkBARyEmJChQBaSBK8yB/8NOYkwVSilqix8gNKIknFIJhgiGZWKX6ECCh4P4F988czUzf0de/e7s722zP9/VS9Nbsz+2y/d9u70/3028/zP8BFmeci5XQlu5JDPWPEbZfX1j8dSVKL2cFMUlOMdx9NEEktcC3xS7+B6FYmtdV4xfOH1DoLSZLgCOJ252/mnoik1usdI9d5i5k1iKQMTgTOBa4G/jnzXKScVgHbKl6zU5kkqW7riVp4O3JPRFLrjQNbgSdr3KYJIqlmK4l6K/cBH8w7FSmrk4m2nUuBnYXXtgAba5+RJKnt1uPtZZKaYS31XzA1QSTV7AbgQOAMogiv1DYHA58CbidqcL0WOAuYIopRTwHnAZuzzE6S1FZjxC1mJogkNcE49SeIrEEk1egM4J3A+4nW3VLbvAf4Y6Io9cVELa4XgP/EhJAkKa9xYDkmiCQ1wzjw2Zq3aYJIqska4E+Bu7Azk9rnEGL13InAPcA5wH9lnZEkSS/W62Dm+5Ok3JYRx8/eYiaNoDHgFuIX7iziNhqpLX6RqDV0HHAB0bnPg29JUtPY4l5SU6zuPpogkkbQ+cCbiJPjTua5SItpktjHXwAeBr5KJEcfAF4DXIcJUklSM60nugV9N/dEJLVejhb34C1m0qJbD1wFfAa4MfNcpMU0SXToW9H9fDVwKJEgOgdIeaYlSdIesYOZpKZY2310BZE0AvpXUdxHdCs7L+uMpMV3JbuSQz1jwLGYHJIkNZ8JIklNMU6sun+k5u2aIJIGrLeKYoL4BVsGLCVuMZNG1VJ2XekoGq94XpKkpjgQeBkmiCQ1wzjwOLCz5u2aIJIGrGwVxb7d56VRdAyxUm6s4vW6l8ZKkjRXR3QfTRBJaoJx6q8/BCaIpIGrWi3hKgqNmh8mVsvdTSRBrwa2FL5mC7Cx5nlJkjRXdjCT1CRryXOR1QSRNGBV94m6ikKjZBJ4EDgLuAZ4NXApUWtrirhneqr7+eYcE5QkaQ7WE7dy2G1WUhOMk+f8l4mEBwAAC2pJREFU0S5m0gDtB+wgCvL2327jKgqNinXAJuAtwL3A8cD9fa9vxoSQJGn4rAf+m2gwIkk5rSJW53uLmTTE9gbuIE6gr8VVFBpu/Z34OsDpwMXA14E3AL/efby/6htIkjREjsDbyyQ1Q64W9+AKImkgxoAbgROAc7sfX5R1RtL89Trx9YqtTwC3ElcU/ho4H3g0y8wkSRq8vYHDifc4ScqtV7vWBJE0pK4CzgQuI5JD0jAr68S3F/AE8K76pyNJ0qJ6BbAPriCS1Aw5E0TeYiYt0AXErTfXA3+YeS7SIFR13FtZ6ywkSVp8k0Q3ToimC5MZ5yJJEMfizwJPZ9i2K4ikBTiNqDf0V8BvZp6LNAjHAM8DS0tesxOfJGmUFG+pPqT7OVg7UlI+uVrcgyuIpHk7HrgF+HuigO901tlIC7OGOBi+G3gG2F543U58kqRRU3ZL9Yru85KUS64W92CCSJqX1wOfAv4DOAnYmXc60rwtBy4n6i68A/g94k3pHOzEJ0kabVW3VFc9L0l1GCdPi3vwFjNpzl4JfI4o2HsCsdpCGkYnAx8i3oQ+SdTSerj72mZMCEmSRtt3iE6dZc9LUg7LgYNwBZE0FF4GfAFIwJuBx/NOR5rVJNABXug+TgKvBf4JuB14EtjQff7hiu8hSdIouqPkOW+plpTTmu5jrgSRK4ikPXQAcCewCjgW+GbW2UizKxbfnABuBZYQiaHzgJuwfpYkqX2WAG8DHiGaM6whTsg24gpaSfnkbHEPJoikGU0ShQrHgR3APsBbgX/LOSlpD5UV39yHuC3ylcD3a5+RJEnNcA7wKuDngU9nnosk9fQSRNYgkhqmuPpiOZEkWpltRtLcVBXZ3B+TQ5Kk9loB/AHRudPkkKQmWUuUhng00/atQSRVKFt9sQxbn6r5xoBTiDeXMhbflCS12YVEXcn3556IJBWME8mh5zNt3wSRVMHWpxpG7wDuI7qSPQZsL7xu8U1JUpsdTCSG/hK4N/NcJKkoZ4t7MEEklfqVGV5z9YWa6HjgX4C/BfYlbpFcR9RYmCIKUU8RhaktvilJaqsPEivCvVgiqYnWkvd80wSR1OcAYuXFx4AHgK2F1119oabZQLSs/wJwEJEQOpLYj6eJZNA6olvLOkwOSZLaaz1wLvBnwLcyz0WSisbY1VExFxNEUtfrgK8C7wIuBY4iDiKmcPWF8psEOkRdoQ5wGfB54C7gcOB9wI8CN1Nde0iSpDa7irjYd0XuiUhSiYOBpeRNEI0BKeP2pUY4H/gw8L/AqcA9eacjvUixo17PM8DvE1dCt9U9KUmShsgG4qLK7xCJIklqmqOJ2mhvAz6XaQ4miNRqBwI3AScRbU7PBp7KOiNpdx1gouT57xD3KUuSpJl9BXg5sdq22MBBkprg3UQB/R8Dvp5pDntn2q6U3U8QdVpeDlwAXJd3OtJulhIr2qqSQKtrnIskScPq3cRx31mYHJLUXL1u2dYgkhZZsX7LJ4C7ieVzx2BySM2ykqgx9G3gVuC5iq+zo54kSTPbB/gj4H7g45nnIkkzGQe+R5SRyMUVRBp5xfotE91xL/Bm4PtZZiXt7kjgt4AziFb1dwIfIRJGxRpEdtSTJGl2vwocRhzzTWeeiyTNJHeLezBBpBa4kt2L+0JUiTc5pLpNEvvkOPEGsJGoe3UBcfC6jVg1dB3wjUJsMc6OepIkVXsJcDnwd8AXM89FkmbTO87PySLVGmn7AluJHb1oGlhS73TUcmXdyKaJe30fBT4K3AA8Wf/UJEkaOVcClwCvA+7LPBdJms13gduBX8s4B1cQaSStAt7XHWXJIcifnVX7lK1m24t4M5igutaQJEmam9XEbdt/jskhSc23AvgRogZpTiaINFJeBfw2cDqxeujTwNe6z1m/RbkcDpxGdTeyl2JySJKkQbqCuEj4u7knIkl7oAkdzMAEkUbEccCFwAlEDZebiRouD3VffxDrt6heBwGnEMnKo4lbybYDy0u+NvcbgSRJo+THgTOBD+N7rKTh0JQEkW3uNVT629VPAZuIZcNfAo4iWoOvIW4te6gvbjOwjqg5tA6TQ1q4/n2x0/18fyIhdCfwCPAnRHvdi4g/+r9ErF7r52o2SZIG6xrgaeLioCQNg16CyFvMpD1ULPC7FngvkWU9m0j67MgzNbVMcV+cAG4jkkXLiITR1cAniNVrPb3EpKvZJEkarP4uoXsBH8dutZKGx1rgeeCxzPOwi5kabwL4aeB64IdKXp8iVgVJdXmYKH5Z9Axxm+M99U5HkqRWK+sSugU4Dy/CSBoOtwEbyH9ea4JIWfRf5SmuolhPJIR6Y033+YTt6pXHKuBN3XEccFjF17kvSpJUvw5xQbFoivwnW5K0J/6x+3hsxjmANYiUQe8qzwSxA04AtwBfAZ4gbsm5gTgZ/zJRT+jVVN+PmbuQl4ZbWS2hA4ATgY8A9xP75SeBk4EHgKcqvpf7oiRJ9VgCvIHoVlbVJXS84nlJapJJ4I3E4oje+UhOyeGYhNSB9EL3cXKAcSsgHQnpLZB+GdL3IKWSsRPSzZDOhnRYxbZ+UIj5wRzm6nAUR9k+9TyxPydIWyF9EdIlkF4PackMce6LDofD4XAsfMx0bLkK0umQ/gLSk8T773OQtlF+bNlpwL/H4XA4ZhoNPK8Y3D+sw+IkGIxb3Lj57pSnQ9pSiNsB6U5IfwPp39n15t0/pkueS905L9b/i2P0x57sG6sh/Ryk34C0ieoDyqch/QykpQvcnsPhaNcYhvd844xrclzZMelWSHdAupddF3Aeg3QTpF+A9JKKuDZeuBmGn7FxxjU9rq5t7kssiniMxiW4F/5N5vtH2bi5xZ0JaT9I+0M6ANKBkF4KaSWkgyC9l90TNlshXQrpaEhvhPSzkI6H9HZIJ0E6BdL/Ub5TPgNpM6TPQrqLSPh8C9ITkLZXxCTiF+JrkD4D6XpIF0M6FdJPQToU0lRFXGeW/xeHo2qU/c5sg3Q7pNsg/SukZwuvP8XCkpUOh8PRP5pyrGCccU2N2wLpfEhHQHoNpJ8kjkvfCumdkN5D9THpNKQvQ/oApKMgjVVss0N7L9w04WdsnHHDHjfobW6DdCOkK4ik9uchPUD5IoriyHU+MpAi1R3KC8NNA1tniNuP8iJIOeK2zRC3fIa47d2Py4onL6uIS0QLu7G+uLGKr80pAQ8BzxLdmYqPlzC/otGT2GlCg1X1NwiiLtA3iNpWD/Z9/MQMcVNY1FLS3IzCsZBx7YybzzFwAnay61h2sY9pbQIxO/8GGWfcwuNmiy3+3aPv470oPzeGOPd/HHgEeLTweA1wcEnMFHnOR/YexDepKgA3BnxshrgLGxS3aYa4i2aIu77v82Km7eIZvueHeHGqjr6PL6d850rd7zldMTZVxE0DbweeI3bq4uOXgENL4r5NdBSrcirlb0SzFertJYGquphJc1X1N2ia6sKVEPtdWbJy44DmJak9RuFYyLh2xs3nGBjg2u5j8fIzwGVUH8ueRlxg3dY3ep//A7C6JM4mELPzb5Bxxi08brbY67ofl53Hf6AibppYODJd8fpzNO98ZMHLkDoVy6I6xg1FXI6lew7HIEeH+e374JJ0h8MxmNFhON7zjTOuyXEeW85/DMvP2DjjmhyXa5sNOx9Z+DcZlnsKjRv8TtmwndnR0uEBpcPhyD2G5T3fOOOaHNeL7eCx5VzHsPyMjTOuyXG5ttmwMZhvVHeCwbjBxjkcwz7c9x0OR+4xLO/5xhnX5DjH/Mew/IyNM67Jcbm22ZQxkCLVkiRJkiRJGl5Na5wlSZIkSZKkmpkgkiRJkiRJajkTRJIkSZIkSS1ngkiSJEmSJKnlTBBJkiRJkiS1nAkiSZIkSZKkljNBJEmSJEmS1HImiCRJkiRJklru/wFi4EawCpFEhQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x432 with 8 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%%time\n",
    "with plt.style.context(\"dark_background\"):\n",
    "    fig = generate_figure(4, width=4)\n",
    "    ax = fig.get_axes()[1]\n",
    "    ax.text(0, np.mean(ax.get_ylim()), \"$eig(A)$\", fontsize=20)\n",
    "    for i, title in enumerate([\n",
    "        \"gaussian\", \"laplace\",\n",
    "        \"hilbert\", \"pascal\"\n",
    "    ]):\n",
    "        fig.get_axes()[i*2].set_title(title, fontsize=20)\n",
    "\n",
    "    video_file = \"psd_multiple.mp4\"\n",
    "    moviewriter = FFMpegFileWriter(fps=10)\n",
    "    with moviewriter.saving(fig, video_file, dpi=100):\n",
    "        for n in range(2, 20, 1):\n",
    "            gaussian = generate_psd(gaussian_kernel, n=n, s=3)\n",
    "            laplace = generate_psd(laplace_kernel, n=n)\n",
    "            hilbert = generate_psd(hilbert_kernel, n=n)\n",
    "            pascal = generate_psd(pascal_kernel, n=n)\n",
    "            make_plot([gaussian, laplace, hilbert, pascal], fig=fig)\n",
    "\n",
    "            ax = fig.get_axes()[1]\n",
    "            ax.text(0, np.mean(ax.get_ylim()), \"$eig(A)$\", fontsize=20)\n",
    "            for i, title in enumerate([\n",
    "                \"gaussian\", \"laplace\",\n",
    "                \"hilbert\", \"pascal\"\n",
    "            ]):\n",
    "                fig.get_axes()[i*2].set_title(title, fontsize=20)\n",
    "\n",
    "            moviewriter.grab_frame()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<video controls width=\"100%\">\n",
       "  <source src=\"psd_multiple.mp4\" type=\"video/mp4\">\n",
       "</video>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "HTML(f\"\"\"<video controls width=\"100%\">\n",
    "  <source src=\"{video_file}\" type=\"video/mp4\">\n",
    "</video>\"\"\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}