{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Final Processing of GAMA-12 Blind source catalogue\n",
    "\n",
    "![HELP LOGO](https://avatars1.githubusercontent.com/u/7880370?s=100&v=4>)\n",
    "\n",
    "\n",
    "The final processing stage requires:\n",
    "1. Quick validation of blind catalogues and Bayesian Pvalue maps\n",
    "2. Skewness level\n",
    "3. Adding flag to catalogue\n",
    "4. Merging MF catalogue with XID+ flux densities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import seaborn as sns\n",
    "from astropy.table import Table,hstack\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import pylab as plt\n",
    "\n",
    "from astropy import units as u\n",
    "from astropy.table import Column\n",
    "\n",
    "import herschelhelp_internal\n",
    "from herschelhelp_internal.utils import gen_help_id\n",
    "import numpy.core.defchararray as np_f\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Read tables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "cat=Table.read('./data/dmu22_XID+SPIRE_GAMA-12_BLIND.fits')\n",
    "cat['RA'].unit=u.deg\n",
    "cat['Dec'].unit=u.deg"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<i>Table length=10</i>\n",
       "<table id=\"table4509288936\" class=\"table-striped table-bordered table-condensed\">\n",
       "<thead><tr><th>HELP_ID</th><th>RA</th><th>Dec</th><th>F_SPIRE_250</th><th>FErr_SPIRE_250_u</th><th>FErr_SPIRE_250_l</th><th>F_SPIRE_350</th><th>FErr_SPIRE_350_u</th><th>FErr_SPIRE_350_l</th><th>F_SPIRE_500</th><th>FErr_SPIRE_500_u</th><th>FErr_SPIRE_500_l</th><th>Bkg_SPIRE_250</th><th>Bkg_SPIRE_350</th><th>Bkg_SPIRE_500</th><th>Sig_conf_SPIRE_250</th><th>Sig_conf_SPIRE_350</th><th>Sig_conf_SPIRE_500</th><th>Rhat_SPIRE_250</th><th>Rhat_SPIRE_350</th><th>Rhat_SPIRE_500</th><th>n_eff_SPIRE_250</th><th>n_eff_SPIRE_500</th><th>n_eff_SPIRE_350</th><th>Pval_res_250</th><th>Pval_res_350</th><th>Pval_res_500</th></tr></thead>\n",
       "<thead><tr><th></th><th>deg</th><th>deg</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy</th><th>mJy/Beam</th><th>mJy/Beam</th><th>mJy/Beam</th><th>mJy/Beam</th><th>mJy/Beam</th><th>mJy/Beam</th><th></th><th></th><th></th><th></th><th></th><th></th><th></th><th></th><th></th></tr></thead>\n",
       "<thead><tr><th>bytes27</th><th>float64</th><th>float64</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th><th>float32</th></tr></thead>\n",
       "<tr><td>5728</td><td>177.33249773350468</td><td>-2.8626750655047006</td><td>46.99437</td><td>51.617813</td><td>42.285637</td><td>24.519966</td><td>29.18767</td><td>19.781033</td><td>7.3959165</td><td>12.345512</td><td>2.9401817</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.9984494</td><td>0.9996332</td><td>1.0003422</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>6193</td><td>177.47688850399337</td><td>-2.753637517003424</td><td>51.88781</td><td>55.96602</td><td>47.426086</td><td>27.64024</td><td>31.665827</td><td>23.390295</td><td>1.2530155</td><td>3.3373992</td><td>0.31067345</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>1.0021628</td><td>0.9996473</td><td>0.99841005</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>9576</td><td>177.61792998658532</td><td>-2.755769776516541</td><td>50.22379</td><td>53.389</td><td>46.41272</td><td>60.492035</td><td>64.15296</td><td>56.363983</td><td>53.880367</td><td>57.85861</td><td>49.16425</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.99980855</td><td>1.0003723</td><td>0.9990191</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>11951</td><td>177.43870031367155</td><td>-2.798325812505323</td><td>35.613663</td><td>40.001526</td><td>31.563438</td><td>26.074862</td><td>30.359695</td><td>21.792027</td><td>7.168383</td><td>11.612774</td><td>3.176016</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.9984573</td><td>0.9988439</td><td>1.0031193</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>13253</td><td>177.48710389475912</td><td>-2.775530044816809</td><td>41.376884</td><td>45.664272</td><td>37.065506</td><td>18.41389</td><td>22.501331</td><td>14.166896</td><td>6.5003223</td><td>11.309651</td><td>2.5763736</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.9988781</td><td>0.9983494</td><td>0.9987865</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>16636</td><td>177.51418531436926</td><td>-2.8010771166096355</td><td>30.843939</td><td>35.30931</td><td>26.058865</td><td>19.496435</td><td>23.797674</td><td>15.104127</td><td>1.905052</td><td>4.6192346</td><td>0.55968374</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>1.0027441</td><td>0.9991835</td><td>0.9985236</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>18022</td><td>177.49956774378163</td><td>-2.746520568003721</td><td>27.84859</td><td>32.46442</td><td>23.372057</td><td>9.111092</td><td>13.182294</td><td>4.733516</td><td>1.9789332</td><td>4.283366</td><td>0.59092546</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.9985813</td><td>1.0005714</td><td>0.99898714</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>18715</td><td>177.5833549025818</td><td>-2.7575925842354785</td><td>29.730644</td><td>34.355747</td><td>25.205248</td><td>23.45012</td><td>27.798838</td><td>19.390047</td><td>2.2101705</td><td>5.0710087</td><td>0.601719</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.99868774</td><td>0.9981934</td><td>0.9991434</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>25862</td><td>177.59221428983025</td><td>-2.74489278772619</td><td>35.422165</td><td>39.50032</td><td>31.021103</td><td>44.73535</td><td>48.559315</td><td>40.788452</td><td>28.819529</td><td>33.691803</td><td>23.399258</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>1.0001178</td><td>0.99950063</td><td>0.99873555</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "<tr><td>27920</td><td>177.48763655281635</td><td>-2.8043057209669406</td><td>13.002316</td><td>17.77816</td><td>8.310316</td><td>2.8451154</td><td>5.900669</td><td>0.96971977</td><td>1.7667637</td><td>3.9922001</td><td>0.4775168</td><td>-0.0062363395</td><td>-0.007360109</td><td>-0.008557404</td><td>0.0018146746</td><td>0.0023701238</td><td>0.0038019842</td><td>0.9983479</td><td>1.0002941</td><td>1.0012697</td><td>2000.0</td><td>2000.0</td><td>2000.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<Table length=10>\n",
       "          HELP_ID                   RA         ... Pval_res_350 Pval_res_500\n",
       "                                   deg         ...                          \n",
       "          bytes27                float64       ...   float32      float32   \n",
       "--------------------------- ------------------ ... ------------ ------------\n",
       "5728                        177.33249773350468 ...          0.0          0.0\n",
       "6193                        177.47688850399337 ...          0.0          0.0\n",
       "9576                        177.61792998658532 ...          0.0          0.0\n",
       "11951                       177.43870031367155 ...          0.0          0.0\n",
       "13253                       177.48710389475912 ...          0.0          0.0\n",
       "16636                       177.51418531436926 ...          0.0          0.0\n",
       "18022                       177.49956774378163 ...          0.0          0.0\n",
       "18715                        177.5833549025818 ...          0.0          0.0\n",
       "25862                       177.59221428983025 ...          0.0          0.0\n",
       "27920                       177.48763655281635 ...          0.0          0.0"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cat[0:10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Look at Symmetry of PDFs to determine depth level of catalogue"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/scipy/stats/stats.py:1713: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.\n",
      "  return np.add.reduce(sorted[indexer] * weights, axis=axis) / sumval\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa8AAAGoCAYAAADxbmq5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XmcY2WV//HPuUlq641uuoGGbjZB9kUsWYRBHDdwFERxBFHBGQd1VFxmfImjo/PDWdxxw3FAEZxBYERQVFDBARFl68YGmkV2pGmWXuju6q6uqiT3/P64N6mb5CaV6q6kKtXf9+sVO3nuTXIqUjl1n3vueczdERER6STBZAcgIiIyXkpeIiLScZS8RESk4yh5iYhIx1HyEhGRjqPkJSIiHUfJS0REOo6Sl4iIdBwlLxER6TjZyQ5ggqldiIh0OpvsADqBjrxERKTjTLcjL5EaP7z9z3W3vf2IXdsYiYhMFB15iYhIx1HyEhGRjqNpQ5k2Gk0Pisj0oiMvERHpOEpeIiLScZS8RESk4yh5iYhIx1HBhmzT6hV56PovkalNR14iItJxlLxERKTjKHmJiEjHUfISEZGOo+QlIiIdR8lLREQ6jpKXiIh0HCUvERHpOLpIWTqOuseLiI68RESk4yh5iYhIx1HyEhGRjqNzXiIp1LBXZGrTkZeIiHQcJS8REek4Sl4iItJxlLxERKTjKHmJiEjHUfISEZGOo+QlIiIdR9d5yZSlHoYiUo+OvEREpOPoyEtkHNR5Q2Rq0JGXiIh0HCUvERHpOEpeIiLScZS8RESk4yh5iYhIx1HyEhGRjqPkJSIiHUfXecmkUycNERkvJS+RCaCLl0XaS9OGIiLScZS8RESk4yh5iYhIx9E5L2kbFWaIyERR8hJpIRVyiLSGpg1FRKTj6MhLJpymB0Wk1ZS8RCaBphNFto6Sl2wxHWGJyGQxd5/sGCbStPph2kmJaOrTUdk2wyY7gE4wrZKXmf0S2BdYPdmxJMxnasUDUy+mqRYPTL2YFM/YplpMWxrPanc/fqKDmW6mVfICMLMl7t4/2XGUTLV4YOrFNNXigakXk+IZ21SLaarFM92oVF5ERDqOkpeIiHSc6Zi8LpjsAKpMtXhg6sU01eKBqReT4hnbVItpqsUzrUy7c14iIjL9TccjLxERmeaUvEREpOMoeYmISMdR8hIRkY6j5CUiIh1nWiWv448/3on6G+qmm266deqtadP0O68p0yp5rV49ldqaiYi01rb8nTetkpeIiGwblLxERKTjKHmJiEjH0UrKIuOQz+dZsWIFQ0NDkx2KdLienh4WLVpELpeb7FA6kpKXyDisWLGCWbNmsfvuu2OmBW9ly7g7a9asYcWKFeyxxx6THU5H0rShyDgMDQ2x/fbbK3HJVjEztt9+ex3BbwUlL5FxUuKSiaD/jraOkpeIiHQcJS8R6Qjuztlnn81ee+3FwQcfzF133ZW638jICGeddRYvfvGL2Xffffnxj3/c5kgj//Ef/8Fee+3FPvvsw69+9atJiWE6U8GGiFAoFMhmJ/7roFgskslkJuS1rrvuOh5++GEefvhhbr/9dt7//vdz++231+z3b//2b+ywww489NBDhGHI2rVrJ+T9x+P+++/n8ssv57777mPlypW8+tWv5qGHHpqwz0J05CXSUZ544gn23XdfzjjjDA4++GBOOeUUBgcHAVi6dCmveMUreOlLX8rrXvc6nnnmGQAuvPBCXvayl3HIIYfwlre8pbz/mWeeycc+9jFe+cpX8olPfILf/va3HHrooRx66KG85CUvYWBgAHfn4x//OAceeCAHHXQQV1xxBQA33XQTxx13HKeccgr77rsvp59+OqVV2XfffXfOPfdcjjnmGH70ox9N2M/+05/+lHe9612YGUceeSTr1q0r/4xJF110EZ/85CcBCIKA+fPnA3DNNdfwmc98pmb/m266iWOPPZaTTz6Z/fffn/e9732EYbjVsZ566ql0d3ezxx57sNdee3HHHXds1WtKJSUvkQ7zpz/9ibPOOot77rmH2bNn8+1vf5t8Ps+HPvQhrrzySpYuXcrf/M3f8KlPfQqAN7/5zdx5553cfffd7Lfffnzve98rv9ZDDz3EDTfcwFe+8hW+/OUvc/7557Ns2TJ+97vf0dvby1VXXcWyZcu4++67ueGGG/j4xz9eThh//OMf+drXvsb999/PY489xu9///vy6/b09HDLLbdw6qmnVsR+6aWXlhNk8nbKKaeM+XM//fTTLF68uPx40aJFPP300xX7rFu3DoB//ud/5rDDDuOtb30rzz33HAAnnngi5557bupr33HHHXzlK1/h3nvv5dFHH+Wqq66q2eejH/1oauyf//zntyhW2TqaNhTpMIsXL+boo48G4B3veAff+MY3OP7441m+fDmvec1rgGi6buHChQAsX76cT3/606xbt46NGzfyute9rvxab33rW8tTWUcffTQf+9jHOP3003nzm9/MokWLuOWWWzjttNPIZDLsuOOOvOIVr+DOO+9k9uzZHH744SxatAiAQw89lCeeeIJjjjkGgLe97W2psZ9++umcfvrpW/Rzl47skqor9gqFAitWrODoo4/mq1/9Kl/96lf5x3/8R/77v/+74Wsffvjh7LnnngCcdtpp3HLLLTUJ9bzzzpvQWGXrKHmJdJjqL0Ezw9054IADuPXWW2v2P/PMM/nJT37CIYccwsUXX8xNN91U3jZjxozy/XPOOYe/+qu/4tprr+XII4/khhtuSP0SLunu7i7fz2QyFAqF1NdNuvTSS/nSl75UM77XXntx5ZVXVoydf/75XHjhhQBce+21LFq0iKeeeqq8fcWKFey8884Vz9l+++3p6+vj5JNPBqLknDzSrCftM6320Y9+lBtvvLFm/NRTT+Wcc86pGGsmVtk6mjYU6TB//vOfy0nqsssu45hjjmGfffZh1apV5fF8Ps99990HwMDAAAsXLiSfz3PppZfWfd1HH32Ugw46iE984hP09/fz4IMPcuyxx3LFFVdQLBZZtWoVN998M4cffvgWx3766aezbNmymlt14gL4wAc+UN6+8847c+KJJ/KDH/wAd+e2225jzpw55aPLEjPjjW98YzlB/+Y3v2H//fcH4Oqrry6fC6t2xx138PjjjxOGIVdccUX5CDLpvPPOS429OnFBNEV5+eWXMzw8zOOPP87DDz+8VZ+b1NKRl0iH2W+//bjkkkt473vfy95778373/9+urq6uPLKKzn77LNZv349hUKBj3zkIxxwwAF87nOf44gjjmC33XbjoIMOYmBgIPV1v/a1r3HjjTeSyWTYf//9OeGEE+jq6uLWW2/lkEMOwcz44he/yE477cSDDz7Y5p8aXv/613Pttdey11570dfXx/e///3ytkMPPZRly5YB8IUvfIF3vvOdfOQjH2HBggXl/R599FFmz56d+tpHHXUU55xzDvfee2+5eGNrHHDAAfz1X/81+++/P9lslvPPP1+VhhPMGk0LdJr+/n5fsmTJZIch09gDDzzAfvvtN2nv/8QTT/CGN7yB5cuXT1oMneod73gH5513HgsWLKgYv+mmm/jyl7/Mz3/+87bHVOe/p6ZPjk3T77ymfn4deYnINuF//ud/JjsEmUBKXiIdZPfdd9dR1wQ77rjjOO644yY7DBknFWyIjNN0mmqXyaP/jrZOy5KXmS02sxvN7AEzu8/MPpyyz+lmdk98+4OZHZLY9oSZ3Wtmy8xs2k3qSmfq6elhzZo1+uKRrVJaz6unp2eyQ+lYrZw2LAD/4O53mdksYKmZXe/u9yf2eRx4hbu/YGYnABcARyS2v9LdV7cwRpFxWbRoEStWrGDVqlWTHYp0uNJKyrJlWpa83P0Z4Jn4/oCZPQDsAtyf2OcPiafcBuj/SZnScrmcVr4VmQLacs7LzHYHXgLUtoAe9bfAdYnHDvzazJaa2VkNXvssM1tiZkv017CITHf6zou0PHmZ2Uzgx8BH3H1DnX1eSZS8PpEYPtrdDwNOAD5gZsemPdfdL3D3fnfvr75+Q0RkutF3XqSlycvMckSJ61J3r23THO1zMPBd4CR3X1Mad/eV8b/PA1cD6q0iIiJAa6sNDfge8IC7f7XOPrsCVwHvdPeHEuMz4iIPzGwG8FpAF7eIiAjQ2mrDo4F3Avea2bJ47J+AXQHc/TvAZ4DtgW/HXZwL7t4P7AhcHY9lgR+6+y9bGKuIiHSQVlYb3sIYParc/T3Ae1LGHwMOqX2GiIiIOmyIiEgHUvISEZGOo+QlIiIdR8lLREQ6jpKXiIh0HCUvERHpOEpeIiLScZS8RESk4yh5iYhIx1HyEhGRjqPkJSIiHUfJS0REOo6Sl4iIdBwlLxER6ThKXiIi0nGUvEREpOMoeYmISMdR8hIRkY6j5CUiIh1HyUtERDqOkpeIiHQcJS8REek4Sl4iItJxlLxERKTjKHmJiEjHUfISEZGOo+QlIiIdR8lLREQ6jpKXiIh0nJYlLzNbbGY3mtkDZnafmX04ZR8zs2+Y2SNmdo+ZHZbYdoaZPRzfzmhVnCIi0nmyLXztAvAP7n6Xmc0ClprZ9e5+f2KfE4C949sRwH8CR5jZPOCzQD/g8XOvcfcXWhiviIh0iJYdebn7M+5+V3x/AHgA2KVqt5OAH3jkNmA7M1sIvA643t3XxgnreuD4VsUqIiKdpS3nvMxsd+AlwO1Vm3YBnko8XhGP1RsXERFpffIys5nAj4GPuPuG6s0pT/EG42mvf5aZLTGzJatWrdq6YEVEpjh950VamrzMLEeUuC5196tSdlkBLE48XgSsbDBew90vcPd+d+9fsGDBxAQuIjJF6Tsv0spqQwO+Bzzg7l+ts9s1wLviqsMjgfXu/gzwK+C1ZjbXzOYCr43HREREWlpteDTwTuBeM1sWj/0TsCuAu38HuBZ4PfAIMAi8O9621sw+B9wZP+9cd1/bwlhFRKSDtCx5ufstpJ+7Su7jwAfqbLsIuKgFoYmISIdThw0REek4Sl4iItJxlLxERKTjNJW8zOwYM3t3fH+Bme3R2rBERETqGzN5mdlngU8An4yHcsD/tDIoERGRRpo58joZOBHYBODuK4FZrQxKRESkkWaS10hc0u4AZjajtSGJiIg01kzy+l8z+y+iju9/B9wAXNjasEREROob8yJld/+ymb0G2ADsA3zG3a9veWQiIiJ1NNVhI05WSlgiIjIl1E1eZjZA+jIkRtTZaXbLohIREWmgbvJyd1UUiojIlNToyGu2u28ws3lp29XlXUREJkujc14/BN4ALKV2dWMH9mxhXCIiInU1mjZ8Q/yvWkGJiMiU0kx7qN80MyYiItIujc559QB9wHwzm8votOFsYOc2xCYiIpKq0Tmv9wIfIUpUSxlNXhuA81scl4iISF2Nznl9Hfi6mX3I3b/ZxphEREQaaqY91DfN7OXA7sn93f0HLYxLRESkrjGTl5n9N/AiYBlQjIcdUPISEZFJ0Uxvw35g/3hZFBERkUnXzJIoy4GdWh2IiIhIs5o58poP3G9mdwDDpUF3P7FlUYmIiDTQTPL6l1YHISIiMh7NVBv+1sx2A/Z29xvMrA/ItD40ERGRdM20h/o74Ergv+KhXYCftDIoERGRRpop2PgAcDRRZw3c/WFgh1YGJSIi0kgzyWvY3UdKD8wsS/oKyyIiIm3RTPL6rZn9E9BrZq8BfgT8rLVhiYiI1NdM8joHWAXcS9Ss91rg02M9ycwuMrPnzWx5ne0fN7Nl8W25mRVLqzab2RNmdm+8bUnzP46IiGwLmimV7wUucvcLAcwsE48NjvG8i4FvUaeNlLt/CfhS/JpvBD7q7msTu7zS3Vc3EZ+IiGxjmjny+g1RsirpBW4Y60nufjOwdqz9YqcBlzW5r4iIbOOaSV497r6x9CC+3zdRAcTXjR0P/Dgx7MCvzWypmZ01xvPPMrMlZrZk1apVExWWiMiUpO+8SDPJa5OZHVZ6YGYvBTZPYAxvBH5fNWV4tLsfBpwAfMDMjq33ZHe/wN373b1/wYIFExiWiMjUo++8SDPnvD4M/MjMVsaPFwJvm8AYTqVqytDdV8b/Pm9mVwOHAzdP4HuKiEgHa5i8zCwAuoB9gX0AAx509/xEvLmZzQFeAbwjMTYDCNx9IL7/WuDciXg/ERGZHhomL3cPzewr7n4U0dIoTTOzy4DjgPlmtgL4LJCLX/c78W4nA792902Jp+4IXG1mpfh+6O6/HM97i4jI9NbMtOGvzewtwFXjWZDS3U9rYp+LiUrqk2OPAYc0+z4iIrLtaSZ5fQyYARTNbDPR1KG7++yWRiYiIlJHM0uizGpHICIiIs1qZkkUM7N3mNk/x48Xm9nhrQ9NREQkXTPXeX0bOAp4e/x4I3B+yyISEREZQzPnvI5w98PM7I8A7v6CmXW1OC4REZG6mjnyysfNeB3AzBYAYUujEhERaaCZ5PUN4GpgRzP7N+AW4N9bGpWIiEgDzVQbXmpmS4FXxUNvcvcHWhuWiIhIfc2c84Koi3xp6rB3jH1FRERaqplS+c8AlwDzgPnA981szJWURUREWqWZI6/TgJe4+xCAmX0euAv411YGJluu1MQrag8pIjL9NFOw8QTQk3jcDTzakmhki7knbsS3xJiIyHTSzJHXMHCfmV1P9J34GuAWM/sGgLuf3cL4ZJy86r4OvkRkOmomeV0d30puak0osjXqHVwpgYnIdNRMqfwl7QhERESkWc2Wyss0kzwPpsIOEek0Sl7TRDL/eMpYeZunTDG6EpiIdBYlr2milHzcK5NWMimlJi7iMSUwEekgDZOXmS0CTgX+AtgZ2AwsB34BXOfuatA7xVQkMSUjEZmm6iYvM/s+sAvwc+ALwPNE13u9GDge+JSZnePuN7cjUBkfJS4Rmc4aHXl9xd2Xp4wvB66K1/TatTVhiYiI1Fc3edVJXMntI8AjEx6RtIwabYjIdDFmwYaZHQ38C7BbvL8B7u57tjY0mQilc19jtYjy8v9oylFEpr5mqg2/B3wUWAoUWxuOtEK9KsOa/Sj9ZaIEJiJTWzPJa727X9fySKRlxjNdqHZSItIJGlUbHhbfvdHMvgRcRdSkFwB3v6vFsYmIiKRqWG1Y9bg/cd+Bv5z4cERERMbWqNrwlQBmtqe7P5bcZmYq1uggpWnAZqYPNWUoIp2gmcUor0wZ+9FEByKtYdZc8YWVbk3uLyIymRqd89oXOACYY2ZvTmyaTeXKyvWefxHwBuB5dz8wZftxwE+Bx+Ohq9z93Hjb8cDXgQzwXXf/fFM/jdRlBoxRdaikJSKdotE5r32Iks92wBsT4wPA3zXx2hcD3wJ+0GCf37n7G5IDZpYBzidasXkFcKeZXePu9zfxniIisg1odM7rp8BPzewod791vC/s7jeb2e5bENPhwCOl82xmdjlwEqDkJSIiQINzXmZ2spnNc/dbzWyBmV1iZvea2RVxt/mJcJSZ3W1m15nZAfHYLsBTiX1WxGMiIiJA44KNf3P3tfH9bwHLgBOA64DvT8B73wXs5u6HAN8EfhKPp66hWO9FzOwsM1tiZktWrVo1AWFtm5y4E4caIIpMafrOizRKXpnE/b3c/Tx3X+HuFwMLtvaN3X2Du2+M718L5MxsPtGR1uLErouAlQ1e5wJ373f3/gULtjqsac0MAqtfDl/KW0pgIlOXvvMijZLXTWZ2rpn1xvffBGBmrwTWb+0bm9lOZlF9m5kdHseyBrgT2NvM9oiXXTkVuGZr309GNaoqVN4SkU7QqNrwg8CngD/Fjz9qZpuAnwHvHOuFzewy4DhgvpmtAD4L5ADc/TvAKcD7zaxAtELzqe7uQMHMPgj8iujo7yJ3v28LfjZpIK6cFxHpSOZNzBGZ2Rwg6+5rWh/Sluvv7/clS5ZMdhgdIazzf7u7Y0AQ6KIvkUnS9C/fNP3Oa+rnb6bDBu6+Ppm44guYpUPV+3ul9IeMA8XQaeYPGxGRydBU8krx6wmNQtrCPTriqk5J7umJKnQI6x2iiYhMokbtob5RbxNR1w3ZBqhllIhMRY0KNt4N/AOJNbwSTmtNOCIiImNrlLzuBJa7+x+qN5jZv7QsIhERkTE0Sl6nAENpG9x9j9aEIyIiMra6BRvuvtbdB5NjZnZY60OSdhqrojBqGaXKQxGZWsZbbfjdlkQhbVHZHqq5hOSMXhOmBCYiU8V4k5dqz6aBLakgVMW8iEwl401e/68lUUjb6a8QEelkjQo2ysxsF2A3YK2ZHQvRYpOtDEwmj3v60Vm9cRGRdhszeZnZF4C3Ea1kXIyHHVDy6nCBxQUZ8ePq+0Zlsiq1jQoMTFlMRCZRM0debwL2cfe0i5WlQ40mnyhjFVPOaXn8P9V5KnQwdzXvFZFJ08w5r8eIlzKR6cfMtmgqUAdeIjKZGvU2/CbRH9+DwDIz+w2JVlHufnbrw5N2UCGhiHSaRtOGpUVillK7krG+70REZNLUTV7ufgmAmX3Y3b+e3GZmH251YCIiIvU0c87rjJSxMyc4DpkkniwxrLtP+ljUNqo1cYmINNLonNdpwNuBPcwsOW04C1iT/izpFMlVk53asvnyfok7FUUaZjhROX0pgamIQ0TapdE5rz8AzwDzga8kxgeAe1oZlLRHdcsns/SjLKCcqKL9rGZciUtE2qnROa8ngSfN7E3ALkTfUyvd/bl2BSftZzSeRdTFySIyFTSaNjwU+A4wB3g6Hl5kZuuAv3f3u9oQn4iISI1G04YXA+9199uTg2Z2JPB94JAWxiVtEFh6t/h6R1+lzhpmtVOHae2kRERapVG14YzqxAXg7rcBM1oXkrRD1FnDEut7lcbjG+md50vre6Wt7ZVc+0tEpJUaHXldZ2a/AH4APBWPLQbeBfyy1YFJe5hZaiJqVLzR8PVQ93kRab1GBRtnm9kJwElEBRsGrADOd/dr2xSftMGWJqpGryci0koNu8q7+3XAdW2KRUREpCl1z3mZ2cGJ+zkz+7SZXWNm/25mfe0JT0REpFajgo2LE/c/D+xFdLFyL1EJvUxz9WYSS105ovZQXrvNJ3YaUkSkWqNpw+SZi1cBL3P3vJndDNzd2rCkncyMAG+YeLzqfjEujc9YlMSi81xW3q62USLSSo2S1xwzO5no6Kzb3fMA7u5mNubf1WZ2EfAG4Hl3PzBl++nAJ+KHG4H3u/vd8bYniNpQFYGCu/c3/yPJljAzDAjdCevsk9b3sOCQTSSu5DZd9yUirdIoef0WODG+f5uZ7ejuz5nZTsDqJl77YuBbRKX2aR4HXuHuL8RVjRcARyS2v9Ldm3kfmUA2Vn8oEZEpoFGp/LvrjD9LNI3YkLvfbGa7N9j+h8TD24BFY72miIgINK42PKbRE81stpnVTAduob+lsiTfgV+b2VIzO2uMOM4ysyVmtmTVqlUTFM62LRtUTgKWCjSqlc6PjYRQCCuLN9yd0J1iWFvUISJbTt95kUbThm8xsy8SddNYCqwCeoiqDl8J7Ab8w9YGYGavJEpeyWR5tLuvNLMdgOvN7EF3vznt+e5+AdGUI/39/fqW3EqlnoWZICrgKIQNElficdGjWy7wmrZSoYPh8TkwnQQT2Rr6zos0mjb8qJnNBU4B3gosBDYDDwD/5e63bO2bx9eSfRc4wd3LC1y6+8r43+fN7GrgcCA1eUlrmBkeesNy+dRxtYYSkTZotCTKUcBt7n4hcOFEv7GZ7QpcBbzT3R9KjM8AAncfiO+/Fjh3ot9fxrYlxRv18paOukRkIjWaNjwDON/MHiKaOvxlXKzRFDO7DDgOmG9mK4DPAjkAd/8O8Blge+Db8ZdaqSR+R+DqeCwL/NDd1QhYRETKGk0bvg/AzPYFTgAuNrM5wI1Eyez37l5s8PzTGr2xu78HeE/K+GNorbBtRqNpxnrbxjsuItNPo/ZQALj7g+5+nrsfD/wlcAvRObCatb5keiktMFk7HtUf1lYROsW4ZVRl5eHoGmDV46V/ky+VfNzMePVricj017CrfElcuLEzUcHGL7UkyvRV8eVvRgYndAgZXYCymEgiQdVJsSJRosoGpReL0l8YJ7BMeRyiGsTSvapENcZ49YO0cR2FiUxfjQo25gAfAE4Duhgtld/RzG4Dvu3uN7YlSmmryvxlZAwsdIZT+h6GXtsGyoF8CLmU5ZiLIdHqzVWZZbxVjWM9R3lLtgVrN41MdgiTptGR15VErZ3+wt3XJTeY2UuBd5rZnu7+vVYGKFPDlixYWS+JKLmIyNZqVLDxmgbblhJduCwiItJ2zZ7z2oWoo0Z5/3odL6TzVV/eVTrX1RVEXeRDr9xW9GgqMDkd6O4MFaNOHbnAKsYLYXSuLJupHC+9bvXr1BuviLlq3ImqkZLTk6WuIFZ+zujPkJz+HH2tynNvo/unj4tI+4yZvMzsC8DbgPuJzsdD9Lur5DUNlb+I42/nMJEkMoERxF/0I8UoaZWUCzIs2l7aVAihGDq5TPRFX1pupQgUi04upd41jBcWq06ioYOV1w6rlNY/sfw6VUu2lPas7iBSWs/MUiY2S9tq3rf8HCUxkXZq5sjrTcA+7j7c6mBk6jCDMKz9ti4tXFlM+SIHUsedKIEFQe23ezF0Minjpeeljo3zeq6614U13F9n5kSmsjGv8wIeI+6MIduWeglii9o8NeobNSWNLzAddclk+eHtf57sECZFo1L5bxL9cToILDOz3wDloy93P7v14YmIiNRqNG24JP53KXBN1Tb1MZApp/UTfeN/B7WsEmmNRqXylwCY2Yfd/evJbWb24VYHJlOXx9WFKafEyoUT1VOLYehYUDXu4IwWTYxVSVjaVrrmrHo8fsnacTPcnaBq3OMTaGn7R4Ueta9f+rstrbqxcjz5mSiBiUy0Zs55nZEyduYExyFTkJlF5enx41JlXQjkMha1gCpt86ivYcFL7aNGexyG7oyEMFyM7pfGi+5szsOmvFc8J3TYXHAGa8adzQVnYMTJh5XjG0dC1g6FbC541bizerDIppGw4r035UNWby4yMFIZ00gRBoZDBvOV+4cerRidD0ufw2ifxtBHF+1MjjuJysaU7iQisuUanfM6DXg7sIeZJacNZwFr0p8l042ZlSsPi155hJHNGJnAGRwJyyXqfA1kAAAgAElEQVTwEO1T8Ogvo+QXeOgwVHAycel66cjNHQbzTjaIjugKiRcr7W8p48NER1NDiRLHzQVnqOB0ZYyhwuhR3eaCM1ws0pONxkvvPVx0RjY7M3JGMXENWyGMklhv3OMqedyVD4liqvqsimHp2jIqNiavCRORidHonNcfgGeA+cBXEuMDwD2tDEqmnnrtocysInEl1RsvevoXeSHue5i2f9pZ1qLDSMrcZenILW18MF877kTXraVVURbCuJlw9XPGKL9XopJ2qldx+PYjdm1zJO3T6JzXk8CTZvYmYBei38uV7v5cu4ITERFJ02ja8FDgO8Ac4Ol4eJGZrQP+3t3vakN8MoVkjYqpQ4jO8XRljELoVW2jomk0K7eNisZDdzbnnYwZPdnK9lAbR6JjtVndQbm4wt15YXORkaIzf0aWbDA6vmawwLqhkJ1nZelKnIBbtanAU+vz7LV9F7O6M+Xx9UNFHl07wu7b5ZjXN/qf/lA+5NENeRbMyDK/L1OOqRA6qwcL9OUC5vZmKmLdlA8JzOjLVf4M+fhwszvbXMsqEdkyjaYNLwbe6+4Vi06a2ZHA99Fqx9uUUmeNUpVhPpGsskG0bEroMFzwyv6HcQFHxpx80RksUN4wVHRm5gKK7qwfKpb7Cw7mQ7bryRC688zGAoW4Kdn64RF2mJGhK2M89kKezYWQ0OGFoSI7zswypzvgvueGWT1YpOjw7MYCi+fk2HNuFw+vHWHlhjyhw3MbC+wwI8v+O3Tx/MYiz2wsEDqsGiwyuztgn+27GCk6azaHOLBhKGTNYJFdZkfJczBfOpcW/QyzugIyFhV0lBTyTneGiqIWqGxxpSQmsuUaJa8Z1YkLwN1vM7MZLYxJpqjSl22A15z/SraNSiuh35T3ioILiPZbu7lQecRGdNT2/MZCVDmY2N8dVm4osDEf1vQ8XLkhz10DxfJrlMb/vC7Pk+vyFaX9pQQ2mA/JBaPn7EKHdUMhj76Qpy83usxmSFTYsXqwWHEkV3rOpnxIT7Z28bJ8WCo2qe2TqLQlsnUaJa/rzOwXRGt6PRWPLQbeBfyy1YHJ1GVmeEoFhZlVNPJNGm+ZeEj6lfD1eiqWEmP15mQiqx4PLL2oJBtY6ntn6/RgtMT/JjVzHYqIbJlGBRtnm9kJwElEBRsGrADOd/dr2xSfiIhsoa3teziVqxUbdpV39+uA69oUi8g2T904RJpTd2bDzA5O3M+Z2afN7Boz+3cz62tPeDJV1ft+TbsmqtH+kL4WV7JTRcU4UWeO6m2hj3bESCoWixTDkGKxWLl/GBWQ5IuVE4fuzlAhpJAyPzlSDCnWua4sLd5iqcqwarz0sP54zVuISJVG0/IXJ+5/HtiL6GLlXqISetmG9WSr20NFt96skaxpcHeKYVSBmOQeJY+1m4usH46SQqml03Ah5KE1Izy8doShQpSUimGUVG58fJBrHhjguY0FCmE0PlwI+dk9q/jWTY/x4LMD5Ith9PqFAr+/5yG++N8/4477HyVfKMbvG3Lvihf43M+W8+vlzzASv0chdJ4ZKHDxH9fxq0c3MlwICcMopk0jRe5YsZmlKzdH46UWV2FUyv/8pmLUIirRUipfhIGRZIur+BICh6FidRsttZMSGQ9L++sWwMz+6O4vie8vA17m7nmLSqfudveDU584ifr7+33JkiVj7ygTpnTdVvV/RaE764ZC8sXahSsHR4oMF72m20UYOoP5kJUDhYrXyxgMF0KWrhwqX0cFMLcbZmSK/PrB1WwcHj2yWjgryz6zQ/5wz0O8MLCpPD5/zixefvhLufupF3h2/VB5fHZvjrcesSfPbCzy3KbR1+nOGCftOwsH1g2NjgcG/Tv3MqMrYKTygI65vQE9GcOrjjW7MtCVMcKq8YyVyulrj021OvM2q+n/1/fc72D/14t/3rJAJumcV1M/f6NzXnPM7GSio7Nud88DuLtbtMysCEGdysPAoguX06oDhwqe2r7phaEiqweLNePPDBR47IWRmlL7P68b4eGVayhUTeWtXLeZR5Y/UHPosnr9AL+8d2VN6fqGzXnufHqQXLby12G46DyxboQ5PbXl8WsGi+RS5kiHC05XYDVJpxhCMahNRqHrPJfIlmiUvH4LnBjfv83MdnT358xsJ2B160MTEZHJ1MpVmrf2qK7uOS93f3fV7bl4/Fl3f9VWvatMK10Zq2moWzr/1ZWpukDXnVwQTclVj+eLpY7zlePPbyqwbnOx5hqyDWueY9PKh/FCvmJ85PknGPjjdRQH11eMFwZWs/aG7zCy6omK8XBkiMeX3Mj651ZUjruz/JlNPPnCUMW4u/PcxjxPrRuuKLpwd4byIeuHiikFJc5IfC6u4rWgfK6sWumCbRGp1bBUvpqZXeDuZ7UqGOlMuQBygZEPo2mz0vdtXy6gLxdNv20Yjir4CqHTnQ3oyji9oTEwEq2dtXZzkWII3dkAd2e46KwbiookXtgctY7aMBwytzdDUBjioT/exqpnV0T7rnqS7p33Idu3Hetvvpj1d98AYZHNf/oDMw96FT17H8HgPdczsOyX4CEDd1zNnJedxOxXnElhzVMMPX4XGz3k6fvuYMc99uNFR76WfNDN6oFhwHng+c3sNGuQY/ecTW8uw2A+ZP1QdG7u8XV5Dt2pl5ndQflnt7yzYSRkQV+Wnlzc3xAIQyiEId3ZoKLzRkjUWipjtR05SglM579EKo0reQH9LYlCOlrpyzYXOAWjXFlYGu/ORAluKFGgYWbkMhCGIas2FSvOmpkZawbzXP/IpuicUGLbylVreeK26zAPKYblxk4MPb6MtTd9n8CLeGEkHi+yafn/se7W/yUTBHgxPkILiwws/SnFkSG65i7Ei3GBSBjy/GP3k5+1M7N2O2D0fd1ZuWGEZU9v4qCFM8pZpBgn1AdXD3HADj2jjXiJjqbWbi6yMJupOcc2UgjpzdVOekQ9IOt/zjo3JjJqvB1snh/PzmZ2kZk9b2bL62w3M/uGmT1iZveY2WGJbWeY2cPx7YxxximTwMxSCzTMjJE6fZ02F9LbQK0fCqM1xKrGRwY3EgSWSFyR4tAAhEWK+ZHK8fwwFPPRv8nxkWEyfXMIi4XK8WKRnnk7pVRPwvYzcqnZoy8XpDbZzWXSf7YtTUBKXCKjxpW83P34cb7+xUCj55wA7B3fzgL+E8DM5gGfBY4ADgc+a2Zzx/neso3Rd7vItmPMaUMzezHwcWC35P7u/pdjPdfdbzaz3RvschLwA4/OVt9mZtuZ2ULgOOB6d18bx3A9URK8bKz3FBGR6a+Zc14/IuqocSFQexHO1tmF0Y71EDX+3aXBeA0zO4voqI1dd526TSS3BaVzMmkVclkzRlKuCOvKREdM1eO9OUt/ne4eisXaXvCZrm7CMCQIAsLElGIQBIQYQSYgTLSICoKAcHiQXK6rouO8ASMb15OdORescmJiYDiqeAyqz2EVPXW8+rq0ktGfK31xFHdPnYbUOS+Byu+8+Tulfi12hGQZ/paUzTczbVhw9/909zvcfWnpNu53Spf2q1hvuaPUkybufoG797t7/4IFCyYoLBmPZCuj3qyRUovAvL4Mc3qCiv9jHdi+L8viObmoyi4eD92Z35fjyMV99MZtqCx+o+zMuSx86avpnjmHIJMtj2MZZr30jWTn7IhlRv8my+S6yc5fRNDdV/7mtyBDdvYCikMbCEc2l4O3TIZM3xwGXljN4LrVeBj1w8hYVAU5MOI8u7FQ7m8YxKtEj4Tw7ECh3OKqxICNI2HNeGDEY9TsX8rLaSX11Z+1bJuS33mztps32eFMmkaNeefF555+ZmZ/b2YLS2Px+ERYQbRGWMkiYGWDcZmiSn35zIzubEBfrvLvDzNjdk+WnWfnoi/vuMcfZmzXm2HfBd3kMjCUD1m3OWSk6CyYkeU1e81kh74MRQ8Zjhvjds+ax8KXncCsXfYiLOQpDm6gOLyZTO8seg98Nb17HwWZLEHPDDzXTZDJEcxaQHa7hWRm70DPnv10v+hwglwPhc0DFDatJZg5j65d9qd77yPJ9M5kcN1qXljxCL0ZZ4dZ3bx4hxnM6M6yfijkkbVRUciMroCdZ+foywUMFpwn1uUZKTpdGZjTE9CTCwjdGBhx8mHcqiZj5DJBvPZZKRE5GYjL5KPikNDTE5jylkik0bThUiqPgj6e2ObAnhPw/tcAHzSzy4mKM9a7+zNm9ivg3xNFGq8FPjkB7yctkPaFGpiR1jgqExi5jDFS1dIpExg92YBn8pUz04EZO8zM8vSGyguRzYzeeQtZ//BS8LBiPDt3IZm+2RBWvVZXD1277F87JRcW6d31QKyrr2q4wA4zAmbP6q4c9+joaW5v5a+PE13TtmNK9WExhGyXpU4HZtB0oMh4NVqMcg8AM+tx94oWA2bW08yLm9llRMUX881sBVEFYS5+/e8A1wKvBx4BBoF3x9vWmtnngDvjlzq3VLwhIiLSzDmvPzQ5VsPdT3P3he6ec/dF7v49d/9OnLjwyAfc/UXufpC7L0k89yJ33yu+fb+5H0emkmxATdsogN5cQG+2dsOsroCdZmZqxuf2Zjhgx+6aE6Fz+ro46MADyWQqn9M9Yxa7vfxEsj0zKsYzvbPY/sC/INs3u2I8yOaYu2AhXd2VR1hmMLsnV7H0S8mMriD1Z+vKWOqUn1nUNT+tDVSIzmOJjFfdI6+4Ae8uQK+ZvYTR6cPZgBajlLLSfxi1lYTRlnxI+SLlbADzejPM7c2waSTkuY0F3KEvZ8yem2O37XKsHyqydOUQw0Vn51lZZu3QRTGEY3efyZX3reO5jQUO2KGH3bbbHQ8Xc9TRL+e6X/yCJ598kl3372fRgS/DgH1e/zcsv/LrrPzjb1jw0uNZfML7CLI5dj7yjTx353Wsuvdm5u1xAHu9+jSyXb24Bbyw+jmeW7mS+XNmctTBL6avpwszY7jgrB8OmdkV0L9LLzPj5OXASDH6uRbNzjGjK6AQRp9JNuNkzOjJGl2ZqPNIsejkMvG0qo1WWhaBIG4DFVj9acTyZ63KQ2mzSVoepa5G63mdAZxJ1BLqTkZ/bzYAl7j7Ve0IcDy0ntfkKi2oWDvuFMJoZeHoC9fK45vzIRuGipX9/Hx0Xa+08ftXDce9/ka3DQ2PcM+zQ4RYRbVhcXiINZvzOAGWSxxZFfN0BU4umyHIdiWCDdmhL2CHGdmaI7r5fRn2mNtFpuqQqzdr7Dw7VxNTNojW9zIql0gxoqPP6KOoWvcrPsqrPjdWeqSEtU2YMut5JbUxeW3del7ufglwiZm9xd1/PGFhybRV7xovi9f8qvlCjtf8ShvPh6W1wqr3j5oAV3eb8iBLkM3VvH/Q1Q2FbO1vQyZHV1dKWycLmNeXI5OyVtf2fZmaxAVRIqq+xgsgG1hqQ10bIxOlFXU02F1kmzTmOa9k4jKz/2ttOCJTmbKHyFTR6JzXPdVDwItL4+5+cCsDExERqafRdV5PEJ3f+ldgM1Hy+h3wxtaHJdNNWgsoIJ5uq91Sb2mQbJD+OtkgvZ1UYKNrYVW3gfL4faouOSMfQm/aeDEkuiqrUqm6sHq6L1o8M6XNU0r8lTukV2OoSENkVKNzXiea2cnABcCX3f0aM8u7+5PtC086Tb0klQviaruqjbO7ArIWrYuV3LTjzAyzugMefyFfXtPLgMVzulg4K8dtTw2SLzqFeA2sXed1ccSiXn58/wbWDRUZKUbvuXhOF28/uI9rHtzAUxsKcQcMY/u+DH/14pkse3aYx9aOkA+j/ftyAYft3MvmkZCnB6L3zlh0EfW83ix9WWMwXrCslBS7M5QrDyt/5vRMUy8xG/E5w/hO9aKUpe1KYDIeU61KcKI0bMzr7leb2a+Bz5nZe4CuRvuLlMu/E5WHBlhgdBEdjYwUKfcMtMCYnckwoytgbbxi8qzugExgzOyGeX1ZHn9hhE0jIQtn5eiOL7o6eVaOe5/bzNMbChy6sIft+6L/lP/xmPnc+udBbn1qM3+5xwxePD+qMPzYy+ez7Nkhfv6nAV6+ay+H7dxLYMaL5nXz5LoRbnhskH227+Jlu/SWizIWbZfjT6uG2K4ny74LusnG4725qOoxlzF2nJktj5dWPA4CmB3/DEmBQVemdrz0WVTnuuqjOa+6oyQm27Ixu8q7+ybgY2Z2CHBU60OS6aCUxKoFZuSClGnCwJjbm6k5MssGxm7bdTGY95r9D96pl30XVI4HZhy1ax+H7NRbFY/xkoW97Dmv9u+v3bbr4oxDa8dnd2d4+a4zaioJM4Gx08wsPVUdiM2gJ2v05WrbQJWa+6aVwGcbXNdVjxKXbOsaVhua2bFmtk/8cBYw08z+qvVhiYiI1Neo2vBrRKsYZ+NGua8CrgM+ambHufvH6z1XZEuUuktUF0pE021W7tKRHO/NGpsLVa2YgJldARtHwprx7XoCNgyHNe8xuztgMB/WrMHVkzWKYe25ulLrq7RY6/1sIjJxGk0bvgY4EOgFngZ2cfdBM/s88Ecqu8yLNK30Be+MXtRsBrm4kqHoyXZSRnfG6MtBvuhsGA4pOvRljRlxO4pCCKsHC+TDqM3UnJ7oP+tiCM8O5BnMO705Y8GMbDnhPLU+z7qhYtQdY1aOXCaqVnx2Y4E1m4vkAthxZrbc4mpgOGTDcEjGoouVe7KVra9KCbM0nly+JBd30q9OYAH1izfqUXsokUij5OXu7mZW+lu09LsY0lxDX5FU5fM+5XZRleMZnJ6sJRJb9MXflYlaNIVh9LzS/rnA2XFmNuopaJTPUQUZ2Hl2jnzRK8cNdpuTY8eZUeeNclsng4WzsszrzZRbWZVSzuzugFldAUGQ2D+OqTsTn89KjFv883XlKsdLSsud1e1hqPZQMk7TtaqwnkbJ6xdm9jugB/gu8L9mdhvwCuDmdgQn05uZpV4fZRbXnVvldVLl/YPKoxgzI3DIBLXXVQVxmXu1IDAyYc1w+bqz6iINM6vbMDcTWGp7qKDcHmp8iUjtoUTG1ug6r0+Y2VHRXb/NzF4EnEyUyK5sV4AyvdX7oo63jmO08RYRmV4aFWyYu99aeuzujwJfTtlHKxGJiEhbNTp3daOZfcjMKiZSzazLzP7SzC4BzmhteCLttyXTc+34G05/JoqManTO63jgb4DLzGwPYB1R5WEA/Bo4z92XtT5EkbGZRf8xV1XNY2Z0Z6KuHtXf/X05GC7UlsH35YyhgteUzZf6JFafKnOvrJqs3OaJcat4TjL2Zuicl8ioRue8hoBvA982sxwwH9js7uvaFZxsu+q1mQrKrZi8fI1VqR4jKtyIkk6pJDYqsgjoMacQOvkw6miRDaLxXBCNDeadIDmeiUrzN42EmEVVhaXzc6E7hWL0nrnMaLFGslyktFpyWVyAEmAVlwpU56O0c4CGEpdItTHbQwG4ex54psWxiNSo12bKzFKvkTKz1M7zZtG1Vtmq1lRmRleGcvPfpFzG6M1ZzXhgRld2tIy+Ot5ckJ6EMqQcmUFqRSIoaUnztrUyedD1WiIi0oGUvESoX2S/JQc+9QorVG8hMnGamjYU6TTlNkpVY6VpuOqehF3xGpMj4eg2A7rjdk8jRS8XdgQ2+ldfcrpxtKDD8JDyRc3l97XoXFypq0fpOY3iF5F0Sl4y7SS6T9UtiAjwypWVS+2e4gKOsNw7cLQNVDGMijWS56cyVtln0BLFG0WHuP1izbpcoys8qw2UyJbQtKFMW2bJpFKZEOq3YIoKNKo3l3sWphVWlI6wqrYFVQmtvD/piSv5eiLSmJKXTHvjXuhxy95li54lsrW2xUpDUPISEZEOpOQlMiFUSyjSTkpess2qV+mXDdInATN1xiG9PD70Uuuo5NKU0b0wHq/uiVhaoFN9DEUaa2nyMrPjzexPZvaImZ2Tsv08M1sW3x4ys3WJbcXEtmtaGadsm8prdFWNB2bkgtpVjs2MbMbIBLXjacUapU4fpSRWykilPUtl9mkJDJTARBppWam8mWWA84HXACuAO83sGne/v7SPu380sf+HgJckXmKzux/aqvhEoJR4oFh14ZeZkTUIi57aHsotfaLQoCa5QZSI0hJlssy+YhyVgIg00sojr8OBR9z9MXcfAS4HTmqw/2nAZS2MR0REpolWJq9dgKcSj1fEYzXMbDdgD+D/EsM9ZrbEzG4zszfVexMzOyveb8mqVasmIm6RlhjvLKBmDSVN8jtvYN3ayQ5n0rQyeaXNetT7fTwVuNLdi4mxXd29H3g78DUze1HaE939Anfvd/f+BQsWbF3Ess0KrLKAwz06J5WxaAmV2vNi0bbkuFn6lGFpW6k9VHXLquj9aos3IHG+TCSW/M6btd28yQ5n0rQyea0AFiceLwJW1tn3VKqmDN19ZfzvY8BNVJ4PE5lQpaKL5FpbyfGMVf6ylLpqlJJYJi7wKO1fOo9VSopWyl7J16b+4pUVY9T2YhTZ1rUyed0J7G1me5hZF1GCqqkaNLN9gLnArYmxuWbWHd+fDxwN3F/9XJGJZmapiaLUNqp2PE5iVLd7Gq1AnIi1ulS8IVKpZdWG7l4wsw8CvyJah+8id7/PzM4Flrh7KZGdBlzulX9y7gf8l5mVFsT9fLJKUWTKadArUUQmXku7yrv7tcC1VWOfqXr8LynP+wNwUCtjExGRzqUOGyITIm1SMb0IQ0S2npKXSJV6FYPVHTeSanKUe9TmifpFGOPJa0qBIpW0GKVIlcAMC6LFKpOrKgeBEXi0onIxkU0qFpqsbvUUJ7AAj4ouzCpXea7qsDERxR0i2wIdeYmkMDMygSUWjhwdzwZW/sVptiCjlNMCKisHq0vyqwWmxCWSRkdeIg3USxyWOIJq+rVavL/ItkRHXiJtovZQMtHmzeia7BAmjZKXSAPVbaNKMga5Om2j0rrHl47g1B5KZGJo2lCkATPD3QmsskVTKRlliRJSubDDRssxjNGlTaLh0fZQ1efSSty95txXqTKx3uKZItsiJS+RMZSTScohkMUZqjoJWbzBnNQTZ2oPJbJ1NG0ospXqJxaVCoq0ipKXiIh0HCUvERHpOEpeIk1KK5hoVyVgqWhDlYciESUvkSaVFqscLW+vXFiyZv/ELTmWofHFz/WU3ktJTETVhiLjUloluVD0muu1SqXxpcc1z6V+uXvTbaYYf6WiyHSk5CUygerllOojMBHZOpo2FBGRjqPkJTJOHl+UnDYF2OgIq1HRRb32UI1iENmWadpQpEkVBRqJtlGldbmi8dLO9dZWHr2Tdt4qrT1USWlU57tElLxEmladjJJto6q7QJlN7NGRijREKmnaUGQCpCUW5RqR1lHyEhGRjqPkJSIiHUfJS2SKqVd1qM4aIqOUvESalLZCMtTvjmH19mfs4ovq0nlH7aFEklRtKDIOpYQUevV4vEpyVVYp7e/VpfRbSO2hRCJKXiJtoGQjMrE0bSgiIh1HR14iWyDRSKPtSm2mktOH1cuz6EhPpjslL5GtkFwGZcx945ZSTb/2GBmolMTqjSuJyXTW0mlDMzvezP5kZo+Y2Tkp2880s1Vmtiy+vSex7Qwzezi+ndHKOEXGy2z0Vjlu5Vv1WOn+2K9tTe3X8DW26tkiU1/LjrzMLAOcD7wGWAHcaWbXuPv9Vbte4e4frHruPOCzQD/RH5JL4+e+0Kp4RbZUvaOvrU1AW0tHXTKdtfLI63DgEXd/zN1HgMuBk5p87uuA6919bZywrgeOb1GcIiLSYVqZvHYBnko8XhGPVXuLmd1jZlea2eJxPhczO8vMlpjZklWrVk1E3CIiU5a+8yKtTF5pkxbVsys/A3Z394OBG4BLxvHcaND9Anfvd/f+BQsWbHGwItNJvWIO6Xz6zou0MnmtABYnHi8CViZ3cPc17j4cP7wQeGmzzxWZKhq2gUrdv/7JqLRtRv3WVHVfB53zkumtlcnrTmBvM9vDzLqAU4FrkjuY2cLEwxOBB+L7vwJea2ZzzWwu8Np4TGRKMhtNMOVkY6PjtfvXVhTWS1yl4XpJsnr/0nuLTGctqzZ094KZfZAo6WSAi9z9PjM7F1ji7tcAZ5vZiUABWAucGT93rZl9jigBApzr7mtbFavIRKmXNLakIjHt6KnRCs062pJtiY3nosmprr+/35csWTLZYYjUSHbAaFa9ZFTdFHis/aXjNP3/4jT9zmvq51dvQ5E2qfcbOd5802j/afS3qEhDag8l0kbJnoil+5aYU6wZT3uNOvuLbEuUvETaIJmI0voOls5ljZW4tnR/kelGyUukzRoeUU3A64hsC3TOS0REOo6Sl4iIdBwlLxER6ThKXiIi0nGUvEREpOMoeYmISMeZVu2hzGwVsAlYPdmxJMxnasUDUy+mqRYPTL2YFM/YplpMWxrPandvavFdM/tls/tON9MqeQGY2RJ375/sOEqmWjww9WKaavHA1ItJ8YxtqsU01eKZbjRtKCIiHUfJS0REOs50TF4XTHYAVaZaPDD1Yppq8cDUi0nxjG2qxTTV4plWpt05LxERmf6m45GXiIhMc0peIiLScTo+eZnZPDO73swejv+dW2e/opkti2/XtCCO483sT2b2iJmdk7K928yuiLffbma7T3QM44znTDNblfhM3tPieC4ys+fNbHmd7WZm34jjvcfMDmtlPE3GdJyZrU98Rp9pcTyLzexGM3vAzO4zsw+n7NO2z6nJeNr2GZlZj5ndYWZ3x/H8v5R92v171kxMbf1d22a4e0ffgC8C58T3zwG+UGe/jS2MIQM8CuwJdAF3A/tX7fP3wHfi+6cCV0xyPGcC32rj/0/HAocBy+tsfz1wHdH6ikcCt0+BmI4Dft7Gz2ghcFh8fxbwUMr/b237nJqMp22fUfwzz4zv54DbgSOr9mnb79k4Ymrr79q2cuv4Iy/gJOCS+P4lwJsmIYbDgUfc/TF3HwEuj+NKSsZ5JfAqs5YtJ9hMPG3l7jcDaxvschLwA4/cBmxnZotSUccAAAXISURBVAsnOaa2cvdn3P2u+P4A8ACwS9VubfucmoynbeKfeWP8MBffqivO2vl71mxM0gLTIXnt6O7PQPTLBuxQZ78eM1tiZreZ2UQnuF2ApxKPV1D7S17ex90LwHpg+wmOYzzxALwlnnq60swWtyiWZjUbc7sdFU8JXWdmB7TrTePprpcQ/SWfNCmfU4N4oI2fkZllzGwZ8DxwvbvX/Xza8HvWbEwwtX7XpoWOSF5mdoOZLU+5jedoYlePWrW8Hfiamb1oIkNMGav+66uZfSZKM+/1M2B3dz8YuIHRv1YnSzs/n2bdBezm7ocA3wR+0o43NbOZwI+Bj7j7hurNKU9p6ec0Rjxt/YzcvejuhwKLgMPN7MDqcNOeNskxTbXftWmhI5KXu7/a3Q9Muf0UeK40bRL/+3yd11gZ//sYcBPRX5ETZQWQ/GtqEbCy3j5mlgXm0LopqzHjcfc17j4cP7wQeGmLYmlWM59hW7n7htKUkLtfC+TMbH4r39PMckSJ4lJ3vypll7Z+TmPFMxmfUfxe64h+j6ub0rbz96ypmKbg79q00BHJawzXAGfE988Aflq9g5nNNbPu+P584Gjg/gmM4U5gbzPbw8y6iE4UV1c0JuM8Bfg/d2/VX4RjxlN1nuREovMZk+ka4F1xNd2RwPrSdPBkMbOdSudLzOxwot+XNS18PwO+Bzzg7l+ts1vbPqdm4mnnZ2RmC8xsu/h+L/Bq4MGq3dr5e9ZUTFPwd216mOyKka29Ec1n/wZ4OP53XjzeD3w3vv9y4F6iqrt7gb9tQRyvJ6rGehT4VDx2LnBifL8H+BHwCHAHsGeLP5ex4vkP4L74M7kR2LfF8VwGPAPkif46/lvgfcD74u0GnB/Hey/Q34b/dsaK6YOJz+g24OUtjucYoimue4Bl8e31k/U5NRlP2z4j4GDgj3E8y4HPpPx33e7fs2Ziauvv2rZyU3soERHpONNh2lBERLYxSl4iItJxlLxERKTjKHmJiEjHUfISEZGOo+QlIiIdR8lLOo5VLm+zrN6yF2bWZ2aXmtm9cTuxW+JWR8nXWG5mPzKzvnh8Y/zv7ma2Od7nfjP7QdxtIm0ZkGVm9uo6MdRdVsTM/sXMnk68xusT2z5p0bIefzKz103UZycyXWQnOwCRLbDZo15yY/kw8Jy7HwRgZvsQXZBc8RpmdinRhbfVXSQedfdDzSwDXA/8NXBpvO137v6GJmIowP9v7+5dowjiMI5/n0oj2NjFTmIEY5NKi5BCrWyUgIWNf4E2YqH2Nr5haSGkFotoYRHQQkkhxCgSEbGQYMDCIrVY6GMxe2Rd7i7JgXJ7eT5wxe7O7s5Vw84Mv4ertt9J2g+8lfTcdqfCy33bd+s3SJqiVEU5BhwEXkg6YvvXNt4XsSvkyytG2TjwrXNg+7M3a8zVLQGHez2kGjSWGaB6uweLFTkHPLL90/YapVrE8Z2+O2KUZfCKNhqrTbU96dNuHrgm6bWkm5Immw2q4q1nKKWWupK0FzgBLNZOzzamDbdMKVD3WJHLVVTGvDZTwIc1HiZiaGTwijb6YXu6+s31amT7PSVN+g5wAHgj6Wh1eUwlg2kFWKcUoG2aqNpsAOu2V2vXlmp9mLb9pV+He8SKPAAmgGlKjcV7nebd/k6/50fsNlnzipHmEtexACxI+k0pLPuJ7a2bdda8xoGXks7abqYFbKlXrIjt77U2D4Fn1eHQxcNEDJt8ecXIkjTTmYqromGmgK87fY5L5Mh14MYAfegZK9KIypijVCWHEutxQdIeSYeAScqaW0RUMnjFKJsAXkn6QImtWKF8AQ3iKbBP0mx13FzzOt/jvhngInCqy5b429U2/lXgJHAFwPZH4DElc24RuJSdhhF/SyRKRES0Tr68IiKidbJhI1qvqkBxq3F6rd9OxH/Qh06id9Np2xv/qx8Ru0WmDSMionUybRgREa2TwSsiIlong1dERLROBq+IiGidP6IwSWU4PNBAAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 3 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "skew=(cat['FErr_SPIRE_250_u']-cat['F_SPIRE_250'])/(cat['F_SPIRE_250']-cat['FErr_SPIRE_250_l'])\n",
    "skew.name='(84th-50th)/(50th-16th) percentile'\n",
    "g=sns.jointplot(x=np.log10(cat['F_SPIRE_250']),y=skew, kind='hex')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For 250 $\\mathrm{\\mu m}$ depth is ~ 6mJy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/scipy/stats/stats.py:1713: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.\n",
      "  return np.add.reduce(sorted[indexer] * weights, axis=axis) / sumval\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAGoCAYAAAATsnHAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3XmcZGV59//PdWrpZfYNGGbAAUE2WZRhABFxiT9BCW4YRRLFvJ7oYzS45PGFSfyZPGaPmrhrNCrojyBKjBIFoxgQUUEGHHAEZJNlYIRhhll7q6pz/f44p6pPVZ2qru7p6jrd/X2/Xk1X3+dU1V3FdF193/d1rtvcHRERkawJet0BERGRNApQIiKSSQpQIiKSSQpQIiKSSQpQIiKSSQpQIiKSSQpQIiKSSQpQIiKSSQpQIiKSSfled2AKVPpCROYC63UHsk4jKBERyaTZOIISSfXvtzyS2v7GUw+d4Z6IyHRQgJI5T4FLZHbSFJ+IiGSSApSIiGSSApSIiGSSApSIiGSSApSIiGSSApSIiGSS0sxl1mmVNi4ic4sClMxbuj5KJNs0xSciIpmkACUiIpmkACUiIpmkACUiIpmkACUiIpmkACUiIpmkACUiIpmk66Aks3RBrsj8phGUiIhkkgKUiIhkkqb4RBq0m1pUGSSRmaMRlIiIZJIClIiIZJIClIiIZJIClIiIZJKSJKTndL2TiKTRCEpERDJJAUpERDJJAUpERDJJAUpERDJJAUpERDJJAUpERDLJ3L3XfZisWddhiczldHLV6JMpsF53IOs0ghIRkUxSgBIRkUxSgBIRkUxSqSOZVnN5nUlEZpZGUCIikkkKUCIikkkKUCIikklag5Ip0VqTiHSbApTINGgVsHUBr8jUaYpPREQySSMoaUtTeSLSKwpQAigQiUj2aIpPREQySdXM5xmNlLJByROCqplPSFN8c5QCkYjMdgpQs4QCjojMN5riyxAFIdHU37yiKb4JaATVRQo4IiJTN+tGUGb2PWBlr/vRYCXwVK870YL6NjXq29Sob517yt3P7nUnsmzWBagsMrON7r6+1/1Io75Njfo2NeqbTCddByUiIpmkACUiIpmkADU9Pt/rDrShvk2N+jY16ptMG61BiYhIJmkEJSIimaQAJSIimaQAJSIimaQAJSIimaQAJSIimTTrAtTZZ5/tRAVj9aUvfelrNn91ZI5+5nVk1gWop57KUiktEZHums+febMuQImIyPygACUiIpmkACUiIpmkDQtl3iqVSmzZsoWRkZFed0XmsP7+ftauXUuhUOh1V2YdBSiZt7Zs2cKiRYtYt24dZtp9W6afu7N9+3a2bNnCYYcd1uvuzDqa4pN5a2RkhBUrVig4SdeYGStWrNAofYoUoGReU3CSbtO/salTgBIRkUxSgBKRGefuXHzxxRxxxBGccMIJ3H777U3n7Nmzh5NOOqn2tXLlSt797nfXjn/961/n2GOP5bjjjuONb3zjTHa/5rbbbuP444/niCOO4OKLL0b7600vJUmIzDHlcpl8fvp/tSuVCrlcbloe69prr+W+++7jvvvu45ZbbuHtb387t9xyS905ixYtYtOmTbWfTz75ZF7zmtcAcN999/H3f//3/OQnP2HZsmU8+eST09KvyXr729/O5z//eU477TRe/vKX873vfY9zzjmnJ32ZizSCEumRhx56iKOPPpo3v/nNnHDCCZx//vkMDQ0B0V/mZ511FieffDIve9nL2Lp1KwBf+MIXOOWUUzjxxBN57WtfWzv/oosu4r3vfS8vetGLuOSSS/jRj35UG3k85znPYc+ePbg773vf+3j2s5/N8ccfz5VXXgnADTfcwAtf+ELOP/98jj76aC688MLaSGDdunV86EMf4vnPfz7f+MY3pu21f/vb3+ZNb3oTZsZpp53Gzp07a68xzX333ceTTz7JmWeeWXsf3vGOd7Bs2TIADjjggNq5J510UupjrFu3jksuuYQNGzawYcMG7r///v16DVu3bmX37t2cfvrpmBlvetOb+Na3vrVfjyn1FKBEeujXv/41b33rW7nzzjtZvHgxn/nMZyiVSvzJn/wJV111Fbfddht/+Id/yF/8xV8A8JrXvIZbb72VO+64g2OOOYYvfvGLtce69957ue666/joRz/KRz7yET796U+zadMmfvzjHzMwMMA3v/lNNm3axB133MF1113H+973vlpQ+MUvfsHHPvYx7rrrLh588EF+8pOf1B63v7+fm266iTe84Q11fb/88svrpuCqX+eff/6Er/uxxx7jkEMOqf28du1aHnvssZbnX3HFFbz+9a+vJRzce++93HvvvZxxxhmcdtppfO9736udmxx1NVq8eDE///nPeec731k3XVh1/fXXp76m5z3veamvYe3atR2/Bpk8TfGJ9NAhhxzCGWecAcDv//7v84lPfIKzzz6bzZs389KXvhSIptZWr14NwObNm/nABz7Azp072bt3Ly972ctqj/W6172uNgV3xhln8N73vpcLL7yQ17zmNaxdu5abbrqJCy64gFwux4EHHshZZ53FrbfeyuLFi9mwYUPtw/akk07ioYce4vnPfz4Ar3/961P7fuGFF3LhhRdO6XWnrdW0y3b72te+xle/+tXaz+Vymfvuu48bbriBLVu2cOaZZ7J582aWLl3a9nkvuOCC2vf3vOc9Tcdf9KIXtQ1w+/MaZPIUoER6qPEDzcxwd4477jh+9rOfNZ1/0UUX8a1vfYsTTzyRSy+9lBtuuKF2bMGCBbXb73//+3nFK17BNddcw2mnncZ1113XdgG/r6+vdjuXy1Eul1MfN+nyyy/nwx/+cFP7EUccwVVXXVXX9ulPf5ovfOELAFxzzTWsXbuWRx99tHZ8y5YtHHzwwanPc8cdd1Aulzn55JNrbWvXruW0006jUChw2GGHcdRRR3HfffdxyimntHyNUP9+pwWT66+/PjVwDQ4O8tOf/rSube3atWzZsqWj1yBToyk+kR565JFHaoHoiiuu4PnPfz5HHXUU27Ztq7WXSiV+9atfAVFm2+rVqymVSlx++eUtH/eBBx7g+OOP55JLLmH9+vXcc889vOAFL+DKK6+kUqmwbds2brzxRjZs2DDlvl944YVs2rSp6asxOAG84x3vqB0/+OCDOe+88/jKV76Cu3PzzTezZMmS2iix0RVXXFEb+VS96lWv4vrrrwei7SjuvfdeDj/8cACOPvroln2urrtdeeWVnH766U3HqyOoxq/G4ASwevVqFi1axM0334y785WvfIVXvvKVLZ9bJk8jKJEeOuaYY7jssst429vexpFHHsnb3/52isUiV111FRdffDG7du2iXC7z7ne/m+OOO46//uu/5tRTT+UZz3gGxx9/PHv27El93I997GNcf/315HI5jj32WM455xyKxSI/+9nPOPHEEzEz/umf/omDDjqIe+65Z4ZfNbz85S/nmmuu4YgjjmBwcJAvf/nLtWMnnXRS3TTb17/+da655pq6+7/sZS/j+9//Psceeyy5XI4Pf/jDrFixgqeeeqrtSHF0dJRTTz2VMAy54oor9vt1fPazn+Wiiy5ieHiYc845Rxl808xmW97++vXrfePGjb3uhswBd999N8ccc0zPnv+hhx7i3HPPZfPmzT3rw1zzne98hwcffJCLL7646di6devYuHEjK1eunPF+tfi31tGC1Rz9zOvotWsEJSJzxrnnntvrLsg0UoAS6ZF169Zp9DSDHnrooV53QSZJSRIyr822KW6ZffRvbOoUoGTe6u/vZ/v27foAka6p7gfV39/f667MSprik3mreh3Ltm3bet0VmcOqO+rK5ClAybxVvchTRLJJU3wiIpJJClAiIpJJClAiIpJJClAiIpJJXQtQZnaImV1vZneb2a/M7F0p57zQzHaZ2ab464Pd6o+IiMwu3cziKwN/6u63m9ki4DYz+4G739Vw3o/dXfVJRESkTtdGUO6+1d1vj2/vAe4G1nTr+UREZG6ZkTUoM1sHPAe4JeXw6WZ2h5lda2bHtbj/W81so5lt1EWVIjLX6TMv0vUAZWYLgf8A3u3uuxsO3w48w91PBD4JfCvtMdz98+6+3t3Xr1q1qrsdFhHpMX3mRboaoMysQBScLnf3bzYed/fd7r43vn0NUDCzmd+sRUREMqebWXwGfBG4293/ucU5B8XnYWYb4v5s71afRERk9uhmFt8ZwB8AvzSz6v7Nfw4cCuDunwPOB95uZmVgGHiDq7S0iIjQxQDl7jcxwba+7v4p4FPd6oOIiMxeqiQhIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZpAAlIiKZ1FGAMrPnm9lb4turzOyw7nZLRETmuwkDlJn9JXAJ8GdxUwH4/7rZKRERkU5GUK8GzgP2Abj748CibnZKRESkkwA15u4OOICZLehul0RERDoLUF83s38FlprZHwHXAV/obrdERGS+y090grt/xMxeCuwGjgI+6O4/6HrPRERkXpswQAHEAUlBSUREZkzLAGVme4jXnRoPAe7ui7vWqznEffxNNMCsl70REZk9WgYod1em3n5KBieIbrsrUImIdKLdCGqxu+82s+Vpx919R/e6Nfs1BicREZmcdmtQ/w6cC9xG9Fmb/JvfgcO72C8REZnn2k3xnRt/V1kjERGZcZ2UOvphJ23zlWseT0SkK9qtQfUDg8BKM1vG+BTfYuDgGehb5lWDU/V7MvFBcUtEZP+0W4N6G/BuomB0G+MBajfw6S73K9OSo6ZkCnknoyll8ImIdKbdGtTHgY+b2Z+4+ydnsE+ZlxaHOhkxKTiJyGTt2DfW6y70TCeljj5pZs8D1iXPd/evdLFfc46Ck4jI5EwYoMzsq8AzgU1AJW52QAFKRES6ppNafOuBY+MtN2SSqm+bE6VMWmIYVb2YV6MrEZFmnWy3sRk4qNsdmU0s8ZXWVm1vjOmhQxg6oTthotKEx8f0J4CIyLhORlArgbvM7OfAaLXR3c/rWq8yzhprajS2AR6mRxuP/6MRk4hIe50EqL/qdidmu9RgY+hiKBGR/dBJFt+PzOwZwJHufp2ZDQK57ndtdtBISESkOzopdfRHwFXAv8ZNa4BvdbNTIiIinSRJvAM4g6iCBO5+H3DARHcys0PM7Hozu9vMfmVm70o5x8zsE2Z2v5ndaWbPnewLyKyJpvdaZEQoUUJEJNLJGtSou49V06PNLE9nqytl4E/d/XYzWwTcZmY/cPe7EuecAxwZf50KfDb+PivVBRezttHG4ztEb6vV2pIlkzR9KCLzWScjqB+Z2Z8DA2b2UuAbwH9NdCd33+rut8e39wB3E00PJr0S+IpHbgaWmtnqSb2CjHES9fk6iDCNMSx5fxGR+ayTAPV+YBvwS6ICstcAH5jMk5jZOuA5wC0Nh9YAjyZ+3kJzEMPM3mpmG81s47Zt2ybz1D3XSZASEUlKfubt2Tl/Ny/vJEANAF9y99e5+/nAl+K2jpjZQuA/gHe7++7Gwyl3aRpAuPvn3X29u69ftWpVp08tIjIrJT/zFi1d3uvu9EwnAeqH1AekAeC6Th7czApEwelyd/9myilbgEMSP68FHu/ksbPGfepTc2lVpBwlTIjI/NZJgOp3973VH+LbgxPdyaK5rS8Cd7v7P7c47WrgTXE232nALnff2kGfMqVVcHL31OCT/hjN56oEkojMZ51k8e0zs+dWEx7M7GRguIP7nQH8AfBLM9sUt/05cCiAu3+OaD3r5cD9wBDwlsl1v/faBafpfh4tZ4nIfNJJgHoX8A0zq069rQZeP9Gd3P0m0teYkuc40XVWs9ZMDW4UnERkvmkboMwsAIrA0cBRRAHnHncvzUDfRERkHmsboNw9NLOPuvvpRNtuiIiIzIhOkiS+b2avNV3Qk2q63pTqklWrtSslSojIfNPJGtR7gQVAxcyGqVbjcV/c1Z7NImk7a5jZhIkSjYdDh8CiIFW38y7jJZD0Z4KIzBedbLexaCY6MlvVAkZKNl8nQSq+a011n8OA5iBVvaEgJSLzQSfbbZiZ/b6Z/b/xz4eY2Ybud212MUuf7ptoZlQzdyIi6TpZg/oMcDrwxvjnvcCnu9YjEREROluDOtXdn2tmvwBw96fNrNjlfs170cygN43A4mZN84nInNdJgCqZWY74s9HMVgFhV3s1iyWWpPbr/GotPsMxaFqPio4pUInI3NXJFN8ngP8EDjSzvwVuAv6uq72apczSA4aZRV8tzm8XZJR+LiLzVSdZfJeb2W3AS+KmV7n73d3t1uxmcd55U+xos8tuWqp6/V3To5hGUCIyV3UyxQdR9fLqNF/He0GJiIhMVSdp5h8ELgOWAyuBL5vZpHbUFRERmaxORlAXAM9x9xEAM/sH4Hbgb7rZsflmoqWkxuoS4+2a5hORuamTJImHgP7Ez33AA13pzZwV7TrYLtGhXYwJWyRKVH9SooSIzEWdjKBGgV+Z2Q+IPhNfCtxkZp8AcPeLu9i/WataWaK6oWH6pobx9w4er74EElRDWrIEUvV5RUTmgk4C1H/GX1U3dKcrc1PLjL7Y1AY/yugTkbmvkzTzy2aiIyIiIkmdrEHJfprs6Mk9mtJrXFtyqu3N96xWlxARmSs6vQ5KpsBbTO1N1J6sIxVU94hi/Ct0yKkEkojMcRpBdUl9EBovddQuOFVoLnIYApWU0VTFNZoSkbmt7QjKzNYCbwDOBA4GhoHNwHeBa91dRWMnwwzDpxY8UkZEgaWXQNIISkTmgpYBysy+DKwBvgP8I/Ak0fVQzwLOBv7CzN7v7jfOREdFRGR+aTeC+qi7b05p3wx8M94T6tDudEtEROa7lmtQLYJT8viYu98//V2SydA2HCIyV02YxWdmZwB/BTwjPj9a63c/vLtdm7vSttYwizL20hb1nPE7JNeWKg45a67T54xXsdBalIjMVp2kmX8ReA9wG1GimXSgsdQRRD8HQVyiyJ1KWB+oAgPz9De59hhx9LH48UIfL3qkEkgiMpd0EqB2ufu1Xe/JHFUNVM3tRi5wymHz+bkWQQrGA1OS19pUAklE5o52WXzPjW9eb2YfBr5JVDgWAHe/vct9ExGReaxtFl/Dz+sTtx148fR3Z36ZSgHZ5HRhXXu8UNVqa3gRkdmmZYBy9xcBmNnh7v5g8piZKUFiP1Wz7OIlqfE9n+Lj1TCTDFaNbdVki+R6U7IEUu18JUuIyCzUSamjq1LavjHdHZkvaoVggWoJpOhWlMHn8Q9m4+tXta/EApQnz0+olkCy6jOYgpOIzE7t1qCOBo4DlpjZaxKHFlO/w67sJ2uzadRUgksuSC+BJCIym7RbgzoKOBdYCvxuon0P8Efd7JSIiEi7NahvA982s9Pd/Wcz2CcREZHWa1Bm9mozW+7uPzOzVWZ2mZn90syujKucS4+1KnM00TERkdmgXZLE37r7jvj2p4BNwDnAtcCXu92x+cbMUv9npO6u6457tG2Hx7eTymF1TygFKRGZvdoFqFzi9hHu/i/uvsXdLwVWdbdbc5dZXNIo5VguMPLJEkmJvAmvtTkh45mA1eCVDEYOlEIotwhgIiKzQbsAdYOZfcjMBuLbrwIwsxcBu2akd3NY2xJI1nB9UywZkBrb04TVlHNl9InILNQuQL2T6FKbXwOvI9oDqprB9wcz0DcREZnH2mXxlYi22fgrM1sC5N19+0x1TCanOo3XOFqqTv9pFCUis00nlSRw913J4BRfxCvTIK06OUAxgHwHMcXdKVecsQqMlZ0wrJ/wc6KkiVDrUCIyy3QUoFJ8f1p7MU+1KkNkFpVACiwKVNZ0LLpdCb2WsQfR97GKUyqHmHvd/SohlCtKmBCR2aNdqaNPtDpEVF1CpkmrSkfVaTlvPoJZ61FRPrDUKT3N8onIbNKu1NFbgD8lsQdUwgXd6Y50m9aiRGS2aBegbgU2u/tPGw+Y2V91rUciIiK0D1DnAyNpB9z9sO50R7pNGX0iMlu0TJJw9x3uPpRsS2wDLzOgWgUi13wkam/xf68UOmFKBYnapohKlBCRWWCyWXz/1pVeSEMJpPrgkgugUD0W1+ADCMzIB81p6qHDSNkphel1+ioqgSQis8BkA5TmhrqsdQkkKATpmX75XPr/lmowauQqgSQis8BkA9T/7UovREREGrRLkqgxszXAM4AdZvYCAHe/sZsdk+nRqtSRSiCJSNZNGKDM7B+B1wN3AZW42QEFqC4KqhXNG6boikE0PVeuJTxE5wVxnGmodETZoVKBYs4JEsGoWgIpF9S3i4hkRScjqFcBR7l72gW70gW1UY17XDGi/liAU7Roz6ew4T5BXHfCfbxyhAOjFciZUwjq154qIYREGYEaTYlIlnSyBvUgUOh2R6RZtSZfq/a0BIjq+dH3+vsGLeKP4pKIZFG7WnyfJPrjewjYZGY/JFH2yN0v7n73ZGpaRxyNkkRktmg3xbcx/n4bcHXDMV1AIyIiXdVuw8LLAMzsXe7+8eQxM3vXRA9sZl8CzgWedPdnpxx/IfBt4Ddx0zfd/UOdd11mQnItS0RkJnWyBvXmlLaLOrjfpcDZE5zzY3c/Kf5ScGqhVXwopByIqkN47XviCKVK5yWQojT08dsiIjOt3RrUBcAbgcPMLDnFtwiYcOt3d7/RzNbtbwcFgsBw96YU8nzOyLkzFu+YW60QURUFlvo7jZSdfEBdNp8RV5Ygeb7V3zP+QaMpEZkp7dagfgpsBVYCH0207wHunKbnP93M7gAeB/6Pu/8q7SQzeyvwVoBDDz10mp56dqmmlzcGKTOjLwf7xpq3NWylHEYZfYWGHX2T1101BiJHda5EZkryM2/lQWt63JvesYkKhprZgcAaos+ox939iY4fPBpBfafFGtRiIHT3vWb2cuDj7n7kRI+5fv1637hx40SnzVmVxggV2zcWTipzpZgzii1q+EF6tp+hEZTINOrot+nwY07wB++erjFBZnT02ttN8Z0EfA5YAjwWN681s53AH7v77fvTO3ffnbh9jZl9xsxWuvtT+/O4IiIyN7Sb4rsUeJu735JsNLPTgC8DJ+7PE5vZQcAT7u5mtoEoYWPCta35rrEEkntUTaKYN0KHUuLq3TCxJpVLTOe5w2jZGSs7AwUjl3IFb2OdPg2cRGSmtQtQCxqDE4C732xmCyZ6YDO7AnghsNLMtgB/SVyRwt0/R7Rj79vNrAwMA29wbVA0ITMbL4EUjhdHrK5R9eWNUsUZq9Tfr+Jg3rzmtK/kFAKnP99ctcLd42m9anWKrr0sEZEm7QLUtWb2XeArwKNx2yHAm4DvTfTA7n7BBMc/BXyqw35KQjVglBtWnWrtLdappvZcCkwi0hvtLtS92MzOAV5JlCRhwBbg0+5+zQz1T2aAApCIZFHbaubufi1w7Qz1RUREpKZlJQkzOyFxu2BmHzCzq83s78xscGa6JyIi81W7UkeXJm7/A3AE0QW7A0Tp59JjrS5jKrT4v5q2ASJEF+5Gx1qVOlLuiojMvHZTfMmPv5cAp7h7ycxuBO7obrekE/nACNxrAaaqkA/I5ZzRslMJ6++TrH6UTDsfKkWbGRZz1TOttutuYNFGiNF9tGAlIjOjXYBaYmavJhpl9bl7CSC+bkl/UmdEYEYhiOrxNbYPFIyh0ZCGQ7Ut4hsrQ5TCaLRUyFldexiPpnKdlBYWEZkm7QLUj4Dz4ts3m9mB7v5EfIGtqj1kSDSqSf+boc2h1Ow9/eUhIlnRLs38LS3af0s05SciItI17bL4nt/ujma22MyaisBKbxSC5qQJdyikFIV1d8Yq0RpV2JAAUfFoS47G9up6VNiQSFEtp5RMqhARmQ7tVhVea2Y/NbMPmtkrzGyDmb3AzP7QzL4KfIcoo08yIDAjZ1AMopJGocdrTUG0djRQMHLmlEOnHB8LgdFKVL+vLoOPKEiNlsOmDL5KCOXqxoeJ9ur2iApSItPv3295pNdd6Il2U3zvMbNlRDXzXgesJqqZdzfwr+5+08x0UTpVzbCr33hjPOGhGrQatYop7dq1FbyIdFu77TZOB2529y8AX5i5Lkm3tAo4yYKwje0iIr3SborvzcBtZvY1M7sozt4TERGZEe2m+P43gJkdDZwDXGpmS4DriaqZ/8TdK63uL7NL4/5PIiK9NuGll+5+j7v/i7ufDbwYuIloTappryjpPfcWJZDcW15oW4kzG9ISIhpLINU/pDIiRKR72lYzr4qTJQ4mSpL4nrbbyJ5a3TwgFxhWK4HkUQo4kA8CBgtxCSSPAox7lFo+Vok2LSwETnX1yYlKIBVzUXs0wLK6tSn3+nYNwkRkurRLklgCvAO4ACgC24B+4EAzuxn4jLtfPyO9lAk1jmWqJZBGK+PbvlfbBwrGvrFKUw2/kbJTievxJaf7xipOGEJ/oXUQUnASkenWbgR1FdFuume6+87kATM7GfgDMzvc3b/YzQ7K1LUtgdTiSLvUcgUnEZlJ7ZIkXtrm2G3AbV3pkYiICJ2vQa0BnpE8391v7FanZPJajYgKgRE6lBLzfO5OPogu4B0p17ePlqMpvQWFgFwQr0XF61R7x5zBAnXtTjSFmAucoINhVPU+1T4nt/yYTLuIzH0TBigz+0fg9cBdQDWt3AEFqIypfnYnA1Vg0Vc+MEbLYW3dKRfvJVUoGsOlkNEwqrVXfYCdoyEDeaM/b1T3hsJhz5hTzHmtvaoSRhUsckH6Rb/JQFNrI700Urv2xr2sRGTu6mQE9SrgKHcf7XZnZOqSH9jJD/fGYFFXBCk+ZmaUw+aIUKo4fQ17Q0Xt0JdrDhLtks6nKyF9PMdQROa6TragexAodLsj0n2TrbnXapTSy9GLpvlE5o92aeafJL4UBthkZj8EaqMod7+4+90TEZH5qt0U38b4+23A1Q3HVEJARES6ql3LltNzAAAgAElEQVSa+WUAZvYud/948piZvavbHZPplzOjnPK3RWppJKLMPWiu01dLeGhRvy/t/OlSTZTQNJ/I3NfJGtSbU9oumuZ+yDQKLD2RoJgPGCwEBLUU7mjjwSAwFhWNfEOiRTmEp0dCxlI2NNw75oxVquWSxo9Fu+6Ot6dl7+0PrUGJzB/t1qAuAN4IHGZmySm+RcD2bndM9k+1iERjcMgFxmAhYN9YWNtZt9q+sM8YGgsZSdSod2BvyekLo2ujkiOj0UpU72+wUB8xqtvA5+oz0ffv9aDAJDLftFuD+imwFVgJfDTRvge4s5udkulhlj69ZmYElj6yyQU2PreXkNIE1Nf5azRdKeEKTiLzU7s1qIeBh83sVcAaos+bx939iZnqnIiIzF/tpvhOAj4HLAEei5vXmtlO4I/d/fYZ6J/sp7QSSFGpI+jHGEkMjaqjrWLOGKtrdyoh7B2rMFjI1UodQbSX1J4xZ7AQkA/qEyNKDvmGEkjJ9apON0isVpZQCSSR+aXdFN+lwNvcvW5jQjM7DfgycGIX+yXTqPrZHSaCQz4XkA+cvryxrxQyVokSHMyi7TYKgTFacUpxoAodRiswWqkwmDf680HdGtZoJWQgD4OFAK+WRgJKIQREAbFRNVhNJlBVH9hbtCtQicwd7bL4FjQGJwB3vxlY0L0uyXQya7dFhhGYkTcbr8MXHcGsPqsvabjsjIXeNDIbLlMXtKpCJlirmkQeujO5bUJEZPZqN4K61sy+S7Qn1KNx2yHAm4DvdbtjMr3abA1FmN6Mt6iR3m6U0upQu4FNpyMoEZlf2iVJXGxm5wCvJEqSMGAL8Glt+S4iMrP+/ZZHeOOph/a6GzOqbTVzd78WuHaG+iKS0CpJvV3yevoxb1F5orHihYhkS8s1KDM7IXG7YGYfMLOrzezvzGxwZron06nVR3GuxZGgxR3CeEGpee3I44y7+nZ3x/GU9vHjTed7yuO3ak88Ruvn6Pz86peI9Fa7JIlLE7f/ATiC6ILdAaL0c5llgsDqgk71g7iYN5b2B001+XKBsaBY3+7uDJed3+4tM1KOSiURlzUqV5yn9pYZLoWJUkfVFHVntOKJ9mi8U453602eD+PVKJLtdZl7yfbEfeqP1SdV1N8nfvyG86k7X4FKpJfaTfElP65eApzi7iUzuxG4o7vdkm4xMwKioJH87C0ExrL+gD1jIaOJUkeBGQOFgJFyyL4xj+ryxcd2DFfoyxtL+4K6LL09oyFDpZCl/bm65xmrRJsgDhatKU28HJdGagwIoU/uOqfq+Y2Dwtq1VA3PUQ1irTZf1ASgSO+0C1BLzOzVRKOsPncvAbi7m5n+rpzFrEVKn5lRyEXXPzW2B2aMVcKme5UqTiVsXsuphNGxtB19wxZrQrXgsp/aVTxvNSJqtU4lIr3TLkD9CDgvvn2zmR3o7k+Y2UHAU93vmoiIJP37LY/Ubs+HjL52aeZvadH+W6IpP5nFkiWDqpwoMaI/b4yUm4cagwVjqFS/FlQOnaeGyiwbyNeVOiqHzpP7yizrz1FMlJEI3dk1EjJYDCjm6ttHSiGFwOhLnO/uDJWdXNyv6ojM3RkpR30ZaGgfrUTrYQv7cnUjuLFKNEU5WLC68kvl0Bkthww0lHEK3SmVQ4r5oK7d3SmFkAuiPbaS7dWRWOr+Wagsk8hktE0zb2Rmn3f3t3arMzJzap+3Fk3HhUQfooXAyAfRh/7esWgvqEoYrUUNFgIGCtEa03A5ZDgukQSwc6TEysGAJX0Be8acfaXoI/np4ZDlAwErBnOUwqgKBcC+UoWBfMiS/oByCHtHq9OHTrEUsqgvR8Vhz1hYW98aKjtL+gKcqA/VmcjhsrO4GGAGu0YqlMLq+WWW9OfIB7Bn1GtTl0MlZ2HR6M9Fj1ldcxuplXEySiG1eoRjYyHFHPTlAioerZcBVCoQmFOIq3WMJ2lE/wnixbC0skwKUiITm1SAAtZ3pRcy45IjjrrEAbPaX/n9+WjElLyPEY2kntgbNn3wPjUUMlRqKA4L7BgOwaCQq08aHS47pX2VutEJwFgI24crTetX5RC2D4dN6e+hw9MjlaZySqHDzuFKUzRwYM+YM5ySwzpUdsopa2pjFWqbOzY+RxnINT9UvNbWfG2WgpRIZzrZUTfpya70QnrGzFLSJSKtkhYqrRIQGr4n2xuD0Pjzt+hXiz61u0S3dXv60VYJE60u3p1sXzs5KiKtTSpAufvZ3eqIzE/6+BaRViac4jOzZwHvA56RPN/dX9zFfsk80a5wkYjMb52sQX2DqHLEF4DKBOfKLNQqQARB+vRfzmzSFRbCeOGlcfNCdyNMa48rqQfWuOYTpcQ1ThlWq040nl8OnVxgGPXrSpXQyeWi+zU+d3UtK62vntrX6sW+aXUAVe9PuiOZct7ObE5H7yRAld39s13vifRMLjAsTp1OKuaM5YM5do1UqMTHPM5OW7OkwLZ9ZcbK4/X3Qo8SIhYWA4q5KPW7mj6+adcYhywpsGpBvpbxNlwKeWDHGAcszHPI4gJB3D5Wce55apSFfQHPXFasJXGUQuf2x4fJGTz34IFaWnvFnZsfHWbvWMgZhw7Qnw8I4uzE2x4f5rE9ZV502AKW9OXIBdEa2gPbR3lwxxinrB1g5YI8ufi5tw2VefjpMY5c0ccBC/IEFgWf0bLz+O4xlg7mWNI/Xr2w4lFGYV/OWFAcD/bO+AaQhXgi3RIlLrz2n+QxEUmyVpvFmdny+ObFRMkR/wmMVo+7+46u9y7F+vXrfePGjb146jmven1PU5KDO7tHKuwrOckiE+7O7tGQJ/eWKTdsYFhNV981GtZdUzVYMA5dWuCpfRV2DI8PyPvzxmHLiuwaCdm6t5x4HDhkSYG9oyH3bh+tjW5yARy1og+A27eO1FLCDThmVZElfTlufmy4Lgvx6JVFjlhe5M4nRtkzOh6N1yzK8+wD+3h4V4ndifal/QHHH9DH3jGvew2FnLF6UZ5yIhUdotT9pX1RuntDMY4oYJM+ytK1UfNWR//XDz/mBP+bS78z5SfJ6Aiqo9febgR1G/VLBO9LHHPg8Kn1S7LKzMiZ03iNrpnRVwjYM1Zpah8oWFNwgmi0s324eSvEoZJz71NjTe0j5WjU1PjvthTCL3870lR+qRLCL7Y2tzvwyydGKYfNGzHet32MJ/Y2z1I/vqdMSPNU3M6RkCf3VZqmDUsVZ89omDrNOFpxCo1Vd+Nj+Vbl4UUkVbtKEocBmFm/u48kj5lZf7c7JiIi81snaeY/7bBN5oi0P/SrF+imWdrfvKNUGK9JNc4gV0Jn274yo+Ww6fx7nxpl+1C5qf2eJ/ay5enhpuf9zZatPPjIY017Oj3+6MPcfcdthGH9c+zYsZ1fbNpEuVz/HHuHx/jpr7cyPFbfPlZxbn1smN2j9aOu0J0n9pYZali0c/e4ukbDqM6j111uvJKY6tpd815ZItJmBBUXhV0DDJjZcxife1kMaMPCOSqw6K8Wt2iRv/qxWcwZhYEci/thx1ClVgcvHxgHLMyzckGe3+4tsXukwmjZScYZ86jgz+7RkEd3lXBgC3DgghwHLcrzxN4yP354KHpMhyNXFFm/ZoCnh0p8647fsn3fGO7wjOUDnPWslZRHR7jmxp+z5YmoZvHKZUt58Rkb6C8E/M93/5N7796MYdx64w/4nVddwIoDV3PrzT/hV7+8EwN+fNNPeMU557DusHX84oGt3HLf4xjw/Tsf4XdPPoznHraKh3eWuOO3IzjwPw/u5ax1Czh17SDD5Sg4hQ6P7Cpx4MI8hy4pANTekz1jFfrjPbYCMyoerUmVQidnTl8+Wo9KxqtoClAZfzL9Os32m6purnG1S5J4M3ARUXmjWxkPULuBy9z9m13rVRtKkpgZ1cKnpZR/HqVKyNY9laZ1p0oYcsdvRxktN69JPfz0GPtKYd2HcmCwZVeJ7cOVuqSCfAC79w2xa88QlcT6Vs6gtPdpRn77AHhIJX6wwIxwbIixx+7GPKRcGR/x5At95JccgJlRSoycCsU+Bo84BcsVKFXGR0LFfMAxRzyTgYE+koO8QgCnHzLImsXFutcWGBy4MMeqBYWm92kgbyzuC5qCTmBQDNITJvIt2mVOmpEkiW6bYoDavyQJd78MuMzMXuvu/zGVHsjsZRZdn5R2IZR7+r8tx1KDE8Desea9pEKPsvwaM97KIQyPlJqmxCoOpaHdVCrNU27loT3QEJyq9wnC5jp9Zbc4ANVP042VQwp99cEJomSNFYP51NewsJhWiY9aqn2jVvPqiksi9SZcg0oGJzP7n+52RySmD2uRea/dGtSdjU3As6rt7n5CNzsmIiLzW7sR1EPAncDvAb8bfz2ZuN2WmX3JzJ40s80tjpuZfcLM7jezO83suZPuvXRVy4rill4fvN2gJxekH89ZetYgWPrjWXPGIIAF1jT1B4BXqFSar8eysLqGVv9KDOLHaX6Fo5X0bLtyGGXiNaqWX2rqUkr/2x3wlGxIkfmg3RrUeWb2auDzwEfc/WozK7n7wx0+9qXAp4CvtDh+DnBk/HUq8Nn4u2REYEY+8Kb1mELOOGhhjqeGohJI1c/OfGCsP7ife56qT4jIGRx3QD+P7Bxj1+h4e+jOygU5dgxX2Dfm8VXhTrkSMlYq4R7G+1MZRkilUqH09GNU9m4nWLCMwAIcJ6yUGd1yN6NPPEDf6meRy+dxhzCsMLb1Xkq7n2LRSWeT6xvALQdhhfLOx9n+wK0se97ryS1aAbkCeIiXR9l084941rEnsmDxEoJcvhYQf/jgXk45eIDViwqJi3SdB58e4+BFeVYM5qO/+OLFpH2lEMdYENc6qq5HVRwsjLL2ku1OdAFyLkhGo+bSSFqrkiyZTJbgZBMqWmbx1U4wWwD8NXAE8Fx3X9vxg5utA77j7s9OOfavwA3ufkX886+BF7r71naPqSy+mZcsgWQk6s25s3Okwq5RJx+MF1d1d57cV+GebaMUclaXmbZzpMJdT44wXHJGEtkRYxXnsZ3D7Bse4ende+tGDDa2l8reHez7zSa8HFehsADLFwlH9rH3jv8mHNkTtQc5iqvWYYU+hu7/OeHQrlr74LOeR9/BRzHyyJ1U9lYrdRmDRz+fRSefS2XvDnx0qPa8y1YdxNHP2YBZ/UTDqsEcZ61b0FRMtz9vHLmij0JgBDYeSAKjVgfQsLoAkw+iaYy0TL+oDqBKI81hcyKLbzISAWq/Sx0B4O77gPea2YnA6VPvWpM1wKOJn7fEbU0ByszeCrwV4NBDM1lXak6zeCTVOFNmZizqyzFSrk85NzMOXJjnkV2lpsda2p9jqORNJYqKOWPXrt3sHk4pg7TlLka3NfyV5iEjv/kFo49ujkZaVWGFkUc3w8huwuQFuWGF4ft/TmXn1oaZNGf4gZ/Tf+izsaD+1+Hpbb/FQ8cakvS2DVUYKjsLivW/YyNlZ6wS0p+vv0O1AO6CXPOMelRNvam5NspUIJqfkp95Kw9a0+Pe9E7bLD4ze4GZHRX/uAhYaGavmKbnTvvVSx3Oufvn3X29u69ftWrVND29iEg2JT/zFi1dPvEd5qiWAcrMPgb8A/BVM/tr4J+AAeA9ZvbhaXjuLcAhiZ/XAo9Pw+NKl7T6Y74vn35k5WAudXRwyJICfSkFVQ9ZNsCi/uZB/dLlK1m4dEVTe3H5wfQffFRTe27JAfQdcnxTezC4hOLqo2h8JVYcIL9oZdNwJQiCllvVu5OaGGEt2qF1u/IfRNK1m+J7KfBsoqD0GLDG3YfM7B+AX1Bf3XwqrgbeaWZfI0qO2DXR+pP0jhFl4lUX8qsKAawYCKg47BiuMFaJKyXk4OhVRSoh3Lt9jO1DFYq5aIrvwAWDVBxu2TLEPU+NMZA3Dl9epG/1gYTubPzNDm75zQ7yuYCj1ixn4VEHg5/B4w/cxV23/piKG8uO2kBx5aF4WGH4wdt5/JpPURnZy5KTz2Xg6BeAO0ue+g3bvvtxyrueYNFzzmHJhtdiQUBlaBc7b76K8s6tLDjqDFa86C0Ehb7osbbeT2VoF8tXHsCxJ55MIZ6uK4VRbcElfQEnHtSPuzM0BsVcVL28Lx+wdnGUPDFahpw5hVy0LjeQt6gKRgUC89qaXHWdCYsuF7bkV4u/BjTjJ/NJu1JHm9392XHl8q3Awe4+bGY54JfufmzbBza7AnghsBJ4AvhLoADg7p+zaFX4U8DZwBDwFnefMPtBSRK9VS2BFG1SWP9BWglDdo1GKdeNu9c+vrvESLm+vRw6j+4q8fRI/ZYWlTBk6+4xtuyuEERpcdGBsMLT+0Z4cNs+glyAxxMA5hXGhofYt20LQS6Px2tJRkhYGiMc2kWuUMRzxfHXUSlTXLKK/OAiyPeNv4gw5JkHLGLpokGCXP3fb89cVmDN4kLdqMogqse3tBCPFseP9eWienzW0B4YDMSFdxuTI6qD0cb26k9ak5pTMpkkMUP7R+13ksR3zezHQD/wb8DXzexm4Czgxoke2N0vmOC4A+/opJOSHdEHpzcFp+oxT9lXKRcYY5WUD+PAGC43b9OeCwL2laIt5+sEOfaVwXL1JYfccrg7ueJA3TSaE0CuQK5/YdM0WlAcIL9oOTRk6BEELF+6OPW1H7So0DTl58CqBbmm1wBRVl9aqaNqrkTrIJT+u6vgJPNNu+ugLjGz06ObfrOZPRN4NVGwumqmOihzW7uLgSf9WFHsbG4PDE/Z6mLSjz+FIyIyde1KHZm7/6z6s7s/AHwk5Ryt8YqIyLRrl2Z+vZn9iZnVTUiaWdHMXmxmlwFv7m73ZDZpN45oNSIKLP1+tQSCBnmLEg4a5YL0f8o5s9RyQ1ZLS2jWKtuu0mIUFrWnbUaYXupIRDrTLkCdDVSAK8zscTO7y8x+A9wHXAD8i7tfOgN9lIypZqCltS/pC5qOGfDM5UUGClb7B2dEAei0tQMcsDBfCzoBUWmk/+eIhRx3QB/5+A6BRVUXXnHsMl76rKUUclYLboXAOOuY1bxmwzPpywfk4w4U8wEnHnYAv3fmsxko5inEiz/FfI7DDljM+SesZHF/jmL85H0546BFRU4+qMjivqD23Pkg2k14yUAQryuN9ylnUA7DplAXr9Slv4GKWSIdabcGNQJ8BviMmRWIsvGG3X3nTHVOssvMCPDxigdEo6S+fEAxZwyVQoZKUQp2dV+ko1YG7Biq8PjuEov6AlYuyBOYccDCAlv3lNj4+DBL+gKOOaCfYs541sp+NuwucfXduxkoBLzkmQtZ0p+DdYt40RFL+bdbfstY2bnw5ANYuyTKxDvr6NV8+Uf3sGXHPi488yiedfAyAE4/ei3/8ZO72PzwNl57xrGcePhBmBnPXbec/7l3Bz9/eCfnHLOK5x2+lMCM4w6KtqDf/OQox6zq4zmr+8kF0Whs75izfajM8oEc65YVawGx+l4UgmijwsZsP7MocaIp+YNq0VyVNRJJmrAWX9YozXz2GC6n/9saGmuuLg6wd6zStHkhwHAppJRyl7GyM5Zyh1LojJbTNkhsLtdU1ddic8H+fHpW3YqBXOpFvIv7AoopJY2KOaOYckFzQDzNqeA0H01LmvkMpYVPt45e+4QbFoqIiPSCApR0Tat/XPkWBwpB+n0KgaXep5AjtWRSdc2oUc5gUbH5gQxYmNIO0ZRcmrREjepjtTLLJitEem7CauYiU1XMRfkApcp4KZ9CAP19AQ7sHQspVaL2gYKxpC+6AHfHcIV9pWijiaX9AYPxfkrb9lXYPhxtSrhyMMfywagM0ba9ZbbsLhN61L5mcfTPevtQmft3jFEKo/bDl/URGOwerfDrp0YZKTtL+wOetaKPQs4YKUd7Ow2XnP68sWZxnmLOCB2e3FdmpByVKTpgQaEWuEbKIeVwPMjlU6b9CjmrlYlyH5+6q5Y1gvFsv+pUn2b2RBSgpIvMol1xizmv21ai2r6oGMQf7l7Xvnwgx5I+r/1c/dBetSDH8oEcbtU09Lh9YbRZYCms35dq5YI8ywfzDJecYs5qa0ZL+3OcsmaQobEKffnxZIaBPBy7qo/dIxXyufHnDgwOWpinHDqFXFDXp/58EL/G5tJFgY0X0k22O+O/eE3VOOJdt7T2JKIAJTOgmvHXVNrHjMCq2yCOC8ziub609uYaS9X2Ao0BIg48hfrHqvZnsFi/b5OZgTvFlPnEwMazERvvUwhITXLIB+mJF7VU+5R0/LR2kflKAUpmRKv6cm3u0eqBpvAc0/OJr7ghvTRLs/X2i5IkREQkkxSgRDo02SS8KCmi83spyU+kngKU9FSLDWtb/sMMSJ9qy5Ge+h2QntYexGtHjcyMwUKQ+hyFwFLboxjUHF7CuLREY5BqF4iqmX5KSRfRGpT0WDVhoboJIsSVFQIjR331h3wiGSF0r1WXKAQQxPs6he6MlqMqeIXAalUaqu0VjzYFjB4roBi3lz3ap6ma2NCXN4ZLISNlpxBEJZwsLjxbqjhjYbxzcNynakBxomWyYq6xdFGUDJJvUQS3+l6Mn9m8IaTIfKMAJT1n1jqtOjAjyKW397Vo789b8waFZvTnx+vlJZ+7L29NvwhmxmAxRyHnTe3FvJELG9uj78UgZaNFWo/kko+b1JzbKDL/aIpPRCTj5mMGHyhAibQ1naMYrSuJTI6m+GROqpUQamgPbDwRIdkWxOcmC7AHQF9crmmsMv5Yxvh0XSUcb88HRtq+idXEjsZSR419Teu/yHymACVzTvXDv3HEYokDcdGIuvNxKBhU6oJIdGJfDsbCZOJCdKdc4PF3a3oOJ0rIqGunfn2pqTJFw2sQmc8UoGTOsnjY0jiKqmbjNZUaioNW047AiQjWfB9L38/JjFxKeSdoHZzq+i0iWoMS6R5FGpH9oQAlIiKZpAAl0kXK3BOZOgUokYR26z+tLrQNPb3mXlzooilIea1d0UukHQUomdMsLi3UGHcCizYwTCs71LLdLKoUUfc40fmtEh4qRIGqljFIfV+iIOXjx7RsJVKjLD6ZF6oZfbXbtXYjIKrhlwwyZkbOoNJU0sgo5KAU1m/R3o4DWHoAqv7Yqj6fyHymACXzRrtYMtkNFSe/AWPr51dsEkmnKT6RGdBYvaLxmIg0U4ASoXltqCq6CLe5Pcckf3k8Xo9qkTARag8okSYKUDLvWSIINcYis2iTwubqEnGCRMp96s6rHo9PCqmv9xdNFY7vA6UgJTJOAUok1nKNyFpn6Zm1nqKrhp7Gu1Z/6aayjiUynyhAiXSTYpDIlClAiYhIJilAiYhIJilAiSS0vVaqRXurXyKv/afeeAkkZUSItKMAJdIgrTQSQJBSAqm6u27BGkogEWf4pZyfS/ysUkciramShEgKi4NU43VL1RJIYa3U0fj5eYu2gE+2Q2KDQtIDFqjUkUgaBSiRKWidkt6ifYJjItJMU3wiIpJJGkGJtFEd3UyUzjDVfAdP3KiOsNzrn09rUzJfaQQl0ka7EkhV1YCSFqNaVZNIHgfDida7Qm9+HJVAkvlKIyiRCSRHNvXt0YGwRfSoJUe0KaEkIq1pBCXSRYpBIlOnACUiIpmkACUiIpmkACWyn1rN4im5QWT/KECJdGgyJZBgguy9SSxOKc1c5itl8YlMghmQkgreqgRS8/0nF2lUAknmM42gRCapbcVzm77gpNgk850ClIiIZJIClMgUVCtENArajKAmLpjUfLaSLGQ+U4ASmYLGEkjjJY0MozmhYqrTddUSSApUMh8pSUJkimojJa//ubrW5DjmKAVPZIoUoET2U+v40yIvXUQ6oik+ERHJpK4GKDM728x+bWb3m9n7U45fZGbbzGxT/PW/utkfERGZPbo2xWdmOeDTwEuBLcCtZna1u9/VcOqV7v7ObvVDRERmp26OoDYA97v7g+4+BnwNeGUXn08kU1qVRmqlVeq6Sh3JfNXNALUGeDTx85a4rdFrzexOM7vKzA5JeyAze6uZbTSzjdu2betGX0W6wqyzckXV66eq59fS1lFwmo+Sn3mVoV297k7PdDNApf1aNV7N8V/AOnc/AbgOuCztgdz98+6+3t3Xr1q1apq7KdJbqaOmRMCS+UefeZFuBqgtQHJEtBZ4PHmCu29399H4xy8AJ3exPyIiMot0M0DdChxpZoeZWRF4A3B18gQzW5348Tzg7i72R6Rn2q0viUi6rmXxuXvZzN4J/DeQA77k7r8ysw8BG939auBiMzsPKAM7gIu61R+RLIh366jdFpHWulpJwt2vAa5paPtg4vafAX/WzT6IZEHdWpKntIlIE5U6EplhCkwinVGpIxERySQFKBERySQFKBERySQFKBERySQFKBERySQFKBERySRzbyyPl21mtg14uNf9aLASeKrXnWhBfZsa9W1q1LfOPeXuZ090kpl9r5Pz5qJZF6CyyMw2uvv6Xvcjjfo2Nerb1KhvMp00xSciIpmkACUiIpmkADU9Pt/rDrShvk2N+jY16ptMG61BiYhIJmkEJSIimaQAJSIimaQANQlmdraZ/drM7jez96cc7zOzK+Pjt5jZugz17SIz22Zmm+Kv/zWDffuSmT1pZptbHDcz+0Tc9zvN7LkZ6dcLzWxX4j37YNp5XerbIWZ2vZndbWa/MrN3pZzTq/etk7715L0zs34z+7mZ3RH37f+mnNOz31OZJHfXVwdfRLsCPwAcDhSBO4BjG875Y+Bz8e03AFdmqG8XAZ/q0Xv3AuC5wOYWx18OXEu0yexpwC0Z6dcLge/06D1bDTw3vr0IuDfl/2mv3rdO+taT9y5+LxbGtwvALcBpDef05PdUX5P/0giqcxuA+939QXcfA74GvLLhnFcCl8W3rwJeYjYj29N10reecfcbgR1tTnkl8BWP3AwsNbPVGehXz7j7Vne/Pb69B7gbWNNwWq/et0761hPxe7E3/rEQfzVmgvXq91QmSQGqc2uARxM/b6H5l7J2jruXgV3Aioz0DeC18dmCgrIAAASNSURBVFTQVWZ2yAz0q1Od9r8XTo+ni641s+N60YF4Cuo5RKOBpJ6/b236Bj1678wsZ2abgCeBH7h7y/dthn9PZZIUoDqX9hdW419mnZzTDZ08738B69z9BOA6xv+CzIJevW8TuR14hrufCHwS+NZMd8DMFgL/Abzb3Xc3Hk65y4y9bxP0rWfvnbtX3P0kYC2wwcye3XBKVv+9SQMFqM5tAZKjjrXA463OMbM8sISZmUKasG/uvt3dR+MfvwCcPAP96lQn7+2Mc/fd1ekid78GKJjZypl6fjMrEAWAy939mymn9Ox9m6hvvX7v4ufdCdwANBZa7dXvqUySAlTnbgWONLPDzKxItLh6dcM5VwNvjm+fD/yPu8/EX2YT9q1hbeI8onWDrLgaeFOclXYasMvdt/a6U2Z2UHVtwsw2EP2+bJ+h5zbgi8Dd7v7PLU7ryfvWSd969d6Z2SozWxrfHgB+B7in4bRe/Z7KJOV73YHZwt3LZvZO4L+Jsua+5O6/MrMPARvd/WqiX9qvmtn9RH+RvSFDfbvYzM4DynHfLpqJvgGY2RVEWV0rzWwL8JdEi9e4++eAa4gy0u4HhoC3ZKRf5wNvN7MyMAy8YQY/yM4A/gD4ZbyeAvDnwKGJ/vXkfeuwb71671YDl5lZjigoft3dv5OF31OZPJU6EhGRTNIUn4iIZJIClIiIZJIClIiIZJIClIiIZJIClIiIZJIClIiIZJIClGSWmVUS2zVsarUtgpkNmtnlZvZLM9tsZjfFZXiSj7HZzL5hZoNx+974+zozG47PucvMvhJXSUjbMmKTmf1Oiz603ObBzC41s98kHuOkuN2sB9tliMwWulBXsmw4rqk2kXcBT7j78QBmdhRQanwMM7sc+N9AY/WDB9z9pPjizh8AvwdcHh/7sbuf20EfRoEXu/veOMDdZGbXxlXGAd7n7lc13Occ4Mj461Tgs/F3EUEjKJkbVgOPVX9w918n6g4m/Rg4otWDuHsF+DlTqAje4TYPjXqyXYbIbKEAJVk2kJgW+882530JuMTMfmZmf2NmRzaeEBcFPQf4ZasHMbN+ohHM9xLNZzZM8T2zzf3bbfPwt/E03r+YWV/c1vPtMkSyTAFKsmzY3U+Kv17d6iR330S0m/CHgeXArWZ2THx4IA4aG4FHiOqwNXpmfM524BF3vzNx7MeJPpzk7g+06UerbR7+DDgaOCXu3yVxu7Z9EGlDa1AyJ8TTa98EvmlmIVER1bvpbB2ruga1GrjBzM6Li4pOtS87zewGom0eNicqjI+a2ZeB/xP/nMltRkSyQiMomfXM7AwzWxbfLgLHAg9P9nHiQPJ+ohHPZPvQcpuH6rpSvP3Eq4DN8d0yuc2ISFZoBCVzwTOBz8YBIAC+S7SZ3lR8C/grMzsz/vnMxJYSAH+Tko0HLbZ5iI9dbmariKb0NhFlEkLvtssQmRW03YaIiGSSpvhERCSTNMUns4aZvQz4x4bm37TL8OtCH1YAP0w59BJ3n5Ht4EXmC03xiYhIJmmKT0REMkkBSkREMkkBSkREMkkBSkREMun/B6/ou8d5ZI29AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 3 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "skew=(cat['FErr_SPIRE_350_u']-cat['F_SPIRE_350'])/(cat['F_SPIRE_350']-cat['FErr_SPIRE_350_l'])\n",
    "skew.name='(84th-50th)/(50th-16th) percentile'\n",
    "g=sns.jointplot(x=np.log10(cat['F_SPIRE_350']),y=skew, kind='hex')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For 350 $\\mathrm{\\mu m}$ depth is ~ 6mJy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/scipy/stats/stats.py:1713: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.\n",
      "  return np.add.reduce(sorted[indexer] * weights, axis=axis) / sumval\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n",
      "/Users/Steven/anaconda/envs/herschelhelp_internal/lib/python3.6/site-packages/matplotlib/axes/_axes.py:6462: UserWarning: The 'normed' kwarg is deprecated, and has been replaced by the 'density' kwarg.\n",
      "  warnings.warn(\"The 'normed' kwarg is deprecated, and has been \"\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa8AAAGoCAYAAADxbmq5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XmcZGV1+P/Pubeqept9YQZmYVhFkZlRB8HgAkYU1ESNJG5R9KXB+DUBY2LiLzES0W8SxcSEYDRoFPy55QtKYhRcgyJfZJnBAYZ9G5hmYJaerdequvee7x/33qpbVfdWV890dXd1n7evetH93KWe7hnrzPM8555HVBVjjDGmkzjT3QFjjDFmoix4GWOM6TgWvIwxxnQcC17GGGM6jgUvY4wxHceClzHGmI5jwcsYY0zHseBljDGm41jwMsYY03Fy092BSWblQowxnU6muwOdwEZexhhjOs5sG3kZMyHfvP2pCZ3/9jPWtqknxpiJsJGXMcaYjmPByxhjTMexaUMz6010atAYM/PZyMsYY0zHsZGXMROQNYqzRA5jppaNvIwxxnQcC17GGGM6jk0bmlnDEjOMmTts5GWMMabjWPAyxhjTcWza0HQcmx40xtjIyxhjTMex4GWMMabj2LShMZPAHl42ZmpZ8DIzlq1tGWOy2LShMcaYjmPByxhjTMex4GWMMabj2JqXmXa2tmWMmSgbeRljjOk4FryMMcZ0HJs2NFPGpgeNMZPFRl7GGGM6jgUvY4wxHcemDc2ks+lBY0y7WfAypo2s5qEx7WHByxw2G2EZY6aLrXkZY4zpOBa8jDHGdBybNjRN2dSgMWYmspGXMcaYjmPByxhjTMexaUMD2PSgMaaz2MjLGGNMx7HgZYwxpuPYtOEcY9ODxpjZwIKXMdPAykYZc2Rs2tAYY0zHsZHXLGXTg8aY2cxGXsYYYzqOBS9jjDEdx6YNO5xNDxpj5iIbeRljjOk4NvIyZgZpNpK2NHpjqix4dQibHjTGmCqbNjTGGNNxbOQ1w9gIyxhjxmfBy5gOYSWljKmy4DUNbHRljDFHxoKXMR3ORmRmLhJVne4+TBoR+SGwLOXQMmDvFHdnomZ6H2d6/8D6OFlmeh9nev/gyPq4V1XPm8zOzEazKnhlEZHNqrppuvvRzEzv40zvH1gfJ8tM7+NM7x90Rh87naXKG2OM6TgWvIwxxnScuRK8rpruDrRgpvdxpvcPrI+TZab3cab3Dzqjjx1tTqx5GWOMmV3mysjLGGPMLGLByxhjTMex4GWMMabjWPAyxhjTcSx4GWOM6TizKnidd955CtjLXvayVye/WjZLP/NaMquC1969M73cmTHGTJ65/Jk3q4KXMcaYucGClzHGmI5jwcsYY0zHsc0ojUlRLpfp7+9nbGxsurtiZqnu7m5Wr15NPp+f7q50JAtexqTo7+9n/vz5rFu3DhGZ7u6YWUZVGRgYoL+/n+OOO266u9ORbNrQmBRjY2MsXbrUApdpCxFh6dKlNrI/Aha8jMlggcu0k/39OjIWvIwxxnQcC17GmBlFVbn44os58cQTWb9+PXfddVfqed/61rc47bTTWL9+Peedd17lgd2tW7dy5plnsnHjRjZt2sQdd9wxld2vuOaaazjppJM46aSTuOaaa6alD7OZBS9j5hDP89pyX9/3J+1eN954I4888giPPPIIV111FR/4wAcazvE8j0suuYSbbrqJe+65h/Xr13PllVcC8Od//udceumlbN26lcsuu4w///M/n7S+tWrfvn184hOf4Pbbb+eOO+7gE5/4BPv375/yfsxmFryMmYG2b9/OKaecwoUXXsj69eu54IILGBkZAWDLli284hWv4EUvehGvec1reOaZZwD40pe+xOmnn86GDRt485vfXDn/3e9+Nx/+8Ic555xz+Iu/+At+8YtfsHHjRjZu3MgLXvACBgcHUVU+8pGP8PznP5/TTjuN//iP/wDg5z//OWeffTYXXHABp5xyCu94xzuId19ft24dl112GS996Uu59tprJ+1n/6//+i/e9a53ISKceeaZHDhwoPIzxlQVVWV4eBhV5dChQxxzzDFAuJZ06NAhAA4ePFhp37x5M+973/sm9Ls+XD/60Y8499xzWbJkCYsXL+bcc8/lhz/84RHd09SJ/xLMhteLXvQiNWYy3H///dP6/k888YQCesstt6iq6nve8x69/PLLtVQq6Ute8hLdvXu3qqp++9vf1ve85z2qqrp3797K9X/1V3+lV1xxhaqqXnjhhfq6171OPc9TVdXXv/71lfsODg5quVzW6667Tl/1qlep53n67LPP6po1a3Tnzp1600036YIFC3THjh3q+76eeeaZ+stf/lJVVY899lj99Kc/ndr/r3/967phw4aG15vf/OZxf/bXve51lfdQVX3lK1+pd955Z8N51157rc6fP19XrlypL3vZyyo/3/33369r1qzR1atX6zHHHKPbt29v+n5Zv+t6n/nMZ1J/pj/+4z9uOPfyyy/XT37yk5XvL7vsstR7Zvw9m+ufeS397DbyMmaGWrNmDWeddRYAv//7v88tt9zCQw89xLZt2zj33HPZuHEjn/rUp+jv7wdg27ZtvOxlL+O0007jG9/4Bvfdd1/lXr/7u7+L67oAnHXWWXz4wx/miiuu4MCBA+RyOW655Rbe9ra34bouK1as4BWveAV33nknAC9+8YtZvXo1juOwceNGtm/fXrnvW97yltS+v+Md72Dr1q0Nr+uuu27cn1u1sbB4fWZeuVzmC1/4Ar/+9a/ZuXMn69ev5+/+7u8A+MIXvsDnPvc5duzYwec+9zne+973jvueab/reh/5yEdSf6YrrrjisH4Gc2TaFrxEZI2I3CQiD4jIfSJySco5Z4vIQRHZGr0+njh2nog8JCKPishH29VPY2aq+g87EUFVOfXUUysfnPfeey8//vGPgXB68Morr+Tee+/l0ksvrXmGqK+vr/L1Rz/6Ub785S8zOjrKmWeeyYMPPpj6YRvr6uqqfO26bs26WfK+Sd/4xjcqU5PJ1wUXXNBw7uc///nK8Z07d7J69Wp27NhROd7f31+Z+ott3boVgBNOOAER4fd+7/e49dZbgTBR4nd+53eAMGi3krCR9ruud/nll6f+TBdffHHDua38DOYItTpEm+gLOBp4YfT1fOBh4Hl155wNfD/lWhd4DDgeKAB311+b9pqlQ2gzDWbKtOGtt96qqqrve9/79LOf/awWi0U94YQTKu2lUkm3bdumqqpLly7VXbt2aalU0le96lV64YUXqmo4bXjttddW7v3oo49Wvn7DG96g119/vX7nO9/RV7/61ep5nu7evVvXrl2rzzzzjN500036ute9rnL+Bz/4Qf3qV7+qquG04Z49eyb9Z//+97+v5513ngZBoL/61a/09NNPbzjn6aef1pUrV1amTz/2sY/phz/8YVVVPeWUU/Smm25SVdWf/vSn+sIXvlBVVW+//XZ95zvf2XCvrN/1kRgYGNB169bpvn37dN++fbpu3TodGBhoOM+mDVO19LO3rTyUqj4DPBN9PSgiDwCrgPtbuPzFwKOq+jiAiHwbeEOL1xozKzz3uc/lmmuu4f3vfz8nnXQSH/jABygUClx33XVcfPHFHDx4EM/z+NCHPsSpp57KJz/5Sc444wyOPfZYTjvtNAYHB1Pv+0//9E/cdNNNuK7L8573PM4//3wKhQK/+tWv2LBhAyLCZz7zGVauXMmDDz44xT81vPa1r+WGG27gxBNPpLe3l69+9auVYxs3bmTr1q0cc8wxXHrppbz85S8nn89z7LHHcvXVVwNh4soll1yC53l0d3dz1VVXAfDUU0/R09OT+p5pv+sjsWTJEv76r/+a008/HYCPf/zjLFmy5Ijuaeq0GuWO5AWsA54CFtS1nw0MEI6sbgROjdovAL6cOO+dwJUZ974I2AxsXrt27eT/G8DMSTNh5HXqqadOax9mmz/7sz/Tu+++u6F9On/XhzPymgOfeTMjYUNE5gHfAT6kqofqDt8FHKuqG4B/Af4zvizlVqmT8qp6lapuUtVNy5cvn6xuG2Nmmcsvv5z169dPdzeOmH3mhdoavEQkTxi4vqGq360/rqqHVHUo+voGIC8iy4B+YE3i1NXAznb21ZiZZN26dWzbtm26uzEn2O+6M7Uz21CAfwceUNV/zDhnZXQeIvLiqD8DwJ3ASSJynIgUgLcC32tXX41Jo00y8Iw5Uvb368i0cz+vswjXqu4Vka1R218CawFU9YuEa1sfEBEPGAXequGfqCcifwT8iDDz8Cuqel/9GxjTLt3d3QwMDNi2KKYtVMP9vLq7u6e7Kx2rndmGt5C+dpU850rgyoxjNwA3tKFrxoxr9erV9Pf3s2fPnunuipml4p2UzeGxnZSNSZHP522HW2NmMCsPZYwxpuNY8DLGGNNxLHgZY4zpOBa8jDHGdBwLXsYYYzqOBS9jjDEdx4KXMcaYjmPByxhjTMex4GWMMabjWPAyxhjTcSx4GWOM6TgWvIwxxnQcC17GGGM6jgUvY4wxHceClzHGmI5jwWsGsV3BjTGmNbYZ5QyQDFrx17bzvDHGZLPgNc3iYJUcdAkWxIwxphmbNpxmSm3gIvG9BS5jjElnwcsYY0zHseA1Q2lG9oaqJXYYY4yteU0T1cbpwrC92uoHiiMgIo3nq00rGmPmLgte06CVwBULFEQ1zOKgGq00uo9gQcwYM/fYtOEUywpcTa/JaLeYZYyZqyx4dYz0UGWjLmPMXNS24CUia0TkJhF5QETuE5FLUs55h4jcE71uFZENiWPbReReEdkqIpvb1U9jjDGdp51rXh7wp6p6l4jMB7aIyE9U9f7EOU8Ar1DV/SJyPnAVcEbi+DmqureNfWwrbXtShWKTh8aYuahtwUtVnwGeib4eFJEHgFXA/Ylzbk1cchuwul39mUr15Z7iANYsxT0rNb7+uEQ3U1UUcKhtN8aYuWBK1rxEZB3wAuD2Jqe9F7gx8b0CPxaRLSJyUZN7XyQim0Vk8549eyaju4ct+QyWJl917bXX6LiBq9n5gcbvMbH7GGM600z6zJtObQ9eIjIP+A7wIVU9lHHOOYTB6y8SzWep6guB84EPisjL065V1atUdZOqblq+fPkk937isso9pbVP2ntGAcxGX8bMfjPtM2+6tDV4iUieMHB9Q1W/m3HOeuDLwBtUdSBuV9Wd0X93A9cDL25nX40xxnSOdmYbCvDvwAOq+o8Z56wFvgu8U1UfTrT3RUkeiEgf8GpgW7v6aowxprO0M9vwLOCdwL0isjVq+0tgLYCqfhH4OLAU+NdoystT1U3ACuD6qC0HfFNVf9jGvk4aobXpwWa1C0nco9UKGqrga7WclDHGzGbtzDa8hXHyuFX1fcD7UtofBzY0XjGzVWJGRnJGrGnR3fo2JlYGKogusCBmjJnNrMJGG4hMPHDB5JWBspBljJntLHi1yWQHEBtEGWNMlQUvY4wxHceC1xxjzzEbY2YDC14zRGYFDqpJG1p3fvJVf014jjacn3wvY4zpVLYZZZtIImde0cyAEbcHpCdsJCtzOIkTkuc70YnJdbEg/r7yxtKwE3Oln8YY02EseLVRKwV504JW1ulBk3aHxiSR5GitPkhZPXpjTCezacMp0CxI2AyeMcZMnAUvY4wxHcemDdssrWpGs/bwmIZrVoTTfck9vILoIqemPdz5UxRykr6Opao1FTdsytAY08kseLVRNUBJWHUjWvxSBT/1/GpwgkRpqGjjyaQgag+jkFTOL2uY2OGmBLH4/R3L0jBmVtg3XJruLkwbmzZsk7SRlUR1o9ICF1ATuGrulfUeh9EvR+LRnGUaGmM6lwWvqTbpAaP1G1qsMsbMFi0FLxF5qYi8J/p6uYgc195uGWOMMdnGDV4icinwF8D/FzXlga+3s1PGGGNMM62MvN4E/DYwDKCqO4H57ezUbJCdSTix9qbvoXESRk3hKIKorX77lbSyUcYY04layTYsqaqKiAKISF+b+9SxsgNT9YAIuJqorKHVElAiMm5gSR6OsxGdyndVpSBMm3eoTZGvLxtlG1YaYzpRKyOv/yMi/wYsEpE/AH4KfKm93epc9eOgNCJhKrtLY3mo8YKJ0vgeQUpxXgBPwdfGkZZqdmajMcZ0gnFHXqr6WRE5FzgEPAf4uKr+pO09mwOSxXtr28cfgbXKYpQxZjZq6SHlKFhZwDLGGDMjZAYvERkk/R/u4YqJ6oK29aqDZQymUsUFMhqqyicqbcQlogACVcp+2OY6tWWjvKjkfN6lYY3LU8g5WlNOKl4vcx21dS9jTMfJDF6qahmFhykOBdXtSKpBJvxvda0rOXUY1JWHiu+hgeIrlQAF4PuQk7BslJ+4xvcg72pUHip834AwgcMVxakrJ+UF4U5fyWBojDEzXbOR1wJVPSQiS9KOq+q+9nWrM9VvBtl4XAgCbSgPVRlZZWzYVQ7S7+dlDPGCAFw3pV1rR3IxWxczxnSaZmte3wReD2yhce9CBY5vY786Xtb0YVaSRjMTzd2Y6ADKxlvGmE7TbNrw9dF/rRSUMcaYGaWV8lA/a6XNGGOMmSqZwUtEuqP1rmUislhElkSvdcAx491YRNaIyE0i8oCI3Ccil6ScIyJyhYg8KiL3iMgLE8cuFJFHoteFh/fjdZas6busdtXGElDQZGsVrV5X01533BhjZrpma17vBz5EGKi2UP0MPQR8voV7e8CfqupdIjIf2CIiP1HV+xPnnA+cFL3OAL4AnBEFzUuBTYSfrVtE5Huqur/1H216xWtbaXt65QizCuP8jDho1F8Tt7uuIKr4QdyuBMSZgpB3aks9BQrFKOswmV3oK3g+FJzaslGOhNc4UtcXY4yZoTJHXqr6z9F615+p6vGqelz02qCqV453Y1V9RlXvir4eBB4AVtWd9gbgaxq6jbAE1dHAa4CfqOq+KGD9BDjv8H7E6SNCFDzq2wXXEXKSeOYqanccqQSTZBkoR4ScA34QPtMVp80rYRp8Wop9yQ/PC1Rr0ulLieudKKVeRFoqbWWMMTNBK+Wh/kVEfgNYlzxfVb/W6ptEU40vAG6vO7QK2JH4vj9qy2pPu/dFwEUAa9eubbVLUyo781DQlCOSkZIoIgQoaRn1WUGnUog35fw4cBljOkfyM2/ZytSPxTmhlYSN/x/4LPBS4PTotanVNxCRecB3gA+p6qH6wymX1KflJ9sbG1WvUtVNqrpp+fLlrXbLGGM6UvIzb/6i1Mdw54RWahtuAp6nh1EpVkTyhIHrG6r63ZRT+oE1ie9XAzuj9rPr2n8+0fefSbLKQLmJKcJku5OYUqw/H2qnCVU1fJA5UPKu1FT0KPmKCHTlnJr2ACh6SiHn4CRGX6rhWNDBRmXGmJmrlS1RtgErJ3pjCT/5/h14QFX/MeO07wHvirIOzwQOquozwI+AV0dZjouBV0dtHSmOAUKlMGQl488h3HfLlbp2qb5UtbLW5Ui4VpZ3wnsFUV1DX8O1rJGy4gXh+WOeUgqg6MNgMaDsV+8faFi5Y7gUUPKChszFQMNgaBtXGmNmolZGXsuA+0XkDqAYN6rqb49z3VnAO4F7RWRr1PaXwNro+i8CNwCvBR4FRoD3RMf2icgngTuj6y7r9HJU1RJQmtrupBTpFQnbynX7b4lIFASDSgZiTIExTxsSRRQoegFuXhpGVCVfyaXUNsyavzXGmOnWSvD6m8O5sarewjiffdFU5Aczjn0F+MrhvPdMJol09Ib2jGuyBj8T3VAyaxbQZgeN6VzfvP0p3n7GzExWa6dWsg1/ISLHAiep6k9FpJdwE2BjjDFmWrSSbfgHwHXAv0VNq4D/bGenjDHGmGZaSdj4IOH61SEAVX0EOKqdnZqL4o0p02RN66U9AB3fLC3RImv6MatsVPW4JW0YY2aWVoJXUVVL8TciksMKMRw2EWlMpojLQEntH4hqtWxT8po4M1AEcjV/gmF7KVCKvhKkZBCOlMNMxOTGmL7CUCnMWgzvTc01Cg3ZiMYYM51aCV6/EJG/BHpE5FzgWuC/29ut2S0uD5V8lksT7bno2a/4VWl3qqWevER7wYVAw1T4MGiF1415ih8FnHgTSgVGvfD5r0Cr9T3i9tFyY5GoOIgaY8xM0Urw+iiwB7iXsFjvDcDH2tmpuSKuJ5janlHUN63+YJzi7qfcLB651U89+po+jZixmbMxxsworaTK9wBfUdUvAYiIG7WNtLNjxhhjTJZWRl4/IwxWsR7gp+3pztyTlnRRKZpb3x6VjWpcM6tW5ahtj7dHCRrWqwLV1HYIK2+kJnxge34ZY2aGVkZe3ao6FH+jqkPRs15mEsTBKIgSJ+LY4Drhw3SBUkmkAMg7EtYmjMo7+ap4PjjikHfCNaxytO6lhNd6UVt3PkwK8YPq+5T8gJ68Q96VSl98Bd+HnNRV3lBQqWZG2sPNxpjp0srIa7huh+MXAaPt69LcEu+lVT9qSrbXr4zFiRqoUvargSg8X8I6iXXvEwCjZQ3rGybaw0SNIAyWdSWivDCLpHLv+mhlozBjzHRpZeR1CXCtiOyMvj8aeEv7ujQ3NdvbKytFPS1BA0AzdhDL2r/LrTQ1HhOyq8vbyMsYM12aBi8RcYACcArwHMLPsgdVtTwFfTPGGGNSNZ02VNUA+AdVLavqNlW91wKXMcbMLN+8/anp7sKUa2XN68ci8maxnQnbqtn6UWY1+MzzM8o8kV4lI4i3Y2mxPFT4vVXcMMZMn1bWvD4M9AG+iIxS+ZzTBW3t2RwSx4Cc1GYcxgqu4EWbUSav6coJ4itFP9kelndyHWr2+lINq2+UROnNOTXrX4HCoTGfvi4XF61Z4xrzlIIrOFF7cjNLkWrn7d82xpip1MqWKPOnoiNzVbKShkS7JKtqmOlHNWGi4Ao5J9wdOS7VJCJ05YS8qwyVAvy4jJMIDiCOUg4UPwAvvkjhUCmgNxdeGwexgHC35S5X6MlXg5ECRT98vqyrbiOcuO+ZBYKNMaZNWtkSRUTk90Xkr6Pv14jIi9vftbkhqzyUQ/iHk4wLTsboJk6Pr68/GAcgL6UwYTlIzz4sZxQxDNSqyxtjZo5W1rz+FXgJ8Pbo+yHg823rkTHGGDOOVoLXGar6QWAMQFX3E6bPmzZKG2NptM6UVTbKzZi+S2sPVBnNKBs1Um5s17pKHzXHaEw4Ua2+jDFmsrWSsFGOivEqgIgsx4qPT5o4rjRWio8zY6rBwQdyTrjpV7yOpaoEQN51yLthUsdIuVoeKucIrkilxmG8XYoHlHxluOSzsNsl74RThvH1h4oBi7tdunJOtZRUEL4KTlrZKEVVUn8eKydljJlsrYy8rgCuB1aIyP8GbgH+tq29mkPiqkv1n+txeSiJ6x7WtbsOOKKVf0VItH6VjxIuqsEjPj9M+khmMyphduP+UZ+RclCTzRgo7Bv1U9fASgEE1JeNqiZ4pA22rKivMWYytZJt+A0R2QL8ZtT0RlV9oL3dmntE0j/cm5WNytog0g/Sy0NlDZczR0TjTEVOJD3eRl7GmMnUysgLoJewyLlD7fYopsNY/DDGzAatpMp/HLgGWAIsA74qIraTsjHGmGnTSsLG24AXqOoYgIj8PXAX8Kl2dsyMT9JnB5uUk0qfgswUP9esmjpFmNWeeTu1qUNjzORoZdpwO9Cd+L4LeKwtvZnjsj7X8076H1TBdcinHOhyoafunyVh2ShtWMNSDff4OjAW4Ae19QoDhT3DHl4QptAnFaNKH/bgsjFmOrQy8ioC94nITwj/LX4ucIuIXAGgqhenXSQiXwFeD+xW1eenHP8I8I5EP54LLFfVfSKyHRgkTLLzVHXThH6qDlWfHl9tF/JuGEDKUdaFAI4juI5LQZWRUlBTZqonL3TllENjYYmoUrQJpeMIomFKfBCEx+LrBkZ9+vJCb96J9gQLK3HsHvbozQuLe1yEaONMCVPrHQmDaziiytj3CxtxGWMmVyvB6/roFft5i/e+GrgS+FraQVW9HLgcQER+C/gTVd2XOOUcVd3b4nvNKlmZh2EZqMaJP0eEnEMlsNW0u8pwWRuCoUht4IqNekp3rjHYjJSVZb2N7dWMx6yNLFN/RGOMOSKtpMpfczg3VtWbRWRdi6e/DfjW4byPMcaYuafVVPm2EZFe4DzgO4lmJdxHbIuIXDTO9ReJyGYR2bxnz552dnVGE7Kru+dS/pSzkidUw8obaWtZRS+93dfG9niNzdbEjJlcyc+8wQP7xr9glpr24AX8FvB/66YMz1LVFwLnAx8UkZdnXayqV6nqJlXdtHz58nb3dcpUa1bUciXc96v+D66Qc+jNO3RFGRmq4PnhXmCLul0WdoV7eKkqfjTXl4+qbsTrbIGGFTcOFoOwuoYfnucIFFw4UAzYM+JT9GrnJ70g3DYlTviIA5YS3jMILIgZM1mSn3nzFy2Z7u5Mm5kQvN5K3ZShqu6M/rubcL1tzm3BMl7ZKKcugMXteVfozgnlIK6oIZX9wBZ2OfgBNfuBOSLkncaMe19h/5iPUFvHMFDYPxYwllLUN6yNmF241wKYMWayNF3zEpHVhMHlZcAxwCiwDfgBcKOqHlGBXhFZCLwC+P1EWx/gqOpg9PWrgcuO5H06WdOyUSkHRCS1HmFWOxBlFjYecwTyrqQ+y9WV0Z5VTiquvWiMMZMhM3iJyFeBVcD3gU8Duwmf9zqZcI3qr0Tko6p6c8b13wLOBpaJSD9wKZAHUNUvRqe9Cfixqg4nLl0BXB990OWAb6rqDw/3BzTGGDP7NBt5/YOqbktp3wZ8V0QKwNqsi1X1beO9uapeTZhSn2x7HNgw3rXGGGPmrsw1r4zAlTxeUtVHJ79L5kg1KRKf3p5xoNkSVbiOlbW+1Xq7McYcjlYK854lIj8RkYdF5HEReUJEHp+KzplQVtBxU3dVDtejutzG87tzwvwupzEJhDCDsV6gsGvIr2x6mbR/1MdP2VnZC8I1tNSdlameb4HMGHMkWqmw8e/AnwBbqO6JaKZQs7JROQnLRvmaOC5Cb8GlK1CGSkFNduH8LpfevMPAiEfJj+8n5FzBjeoc+lHavBJW3Og/6LGo22FRd1gHSgh3ct436kcB0a2O3kT3PchxAAAgAElEQVTwop2Tc9SO6uL+ORMpDmyMMSlaSZU/qKo3qupuVR2IX23vmWmQNb3niIA25gu6Tpg2X891woDTeP8wiPl191LC9Hgh/AuTvOOYp5Vrk9mEWTsqx8cs89AYcySaZRu+MPryJhG5HPguYZFeAFT1rjb3zRhjjEnVNNuw7vtkZXcFXjn53THt4EiygG5VWjHfrLUoVWXMC+jOOQ2jJj/Q1Oe44oeTs8pQ2eDLGHO4MoOXqp4DICLHR+nrFSJyfLs7ZtLFn/f1ISauX+hHJZ7iIJR3hHwBin44xReWbwofPl7W6zLmKQeL4bqYH4RbrsQPGvuJzSgBnh70Kbg+K/pydEVvmHOEobIiZaWv4JCPLo7X6HwBR+sqzCuoKKpileeNMYellTWv61Larp3sjpjWjFc2Kqx9qHXnh9mH8wrVpAmJrunOCct6HbwgqFTgiN/DoTHLsORD/yGP0XJA3hWcqBqwAkOlgGI5QAhHYnEnA+JSVYlOUb3OEg+NMRPVbM3rFOBUYKGI/E7i0AJqd1Y206B52SioH5uJCL7fuH+XiOD5ip9S6EvSq0ahwPxuNzXpopBLLxvlZJSHspGXMZPjm7c/xdvPyKwbMes0W/N6DuFOyIsIK7/HBoE/aGenTPtkxKPMYGiMMTNRszWv/wL+S0Reoqq/msI+GWOMMU1lrnmJyJtEZImq/kpElovINSJyr4j8R1Rt3nQcSR91MTWjrvSyUel7fanaSNAYk61Zwsb/TmwQeSWwlXBzyBuBr7a7Y2Z8zcpGpe2qXHChu26srdF+XWllo7LeQRWe3F+m7CtBXYQZLAWVTSmT4gefK+1azXysr3uYvNSCmDEmTbPglSzBcKKqfk5V+6NK8LNny+IOJlGQSss8LLgOhbrNtUSEeV05FvW4uNGuyl4AnoYVN1bMy9HlhgkfQRAGnJqqGYndlgdLAQ/uKbJ32K8JPH4AB8YCRsrpASxO46+p4BGXo4oDGtS84nOMMSbWLHj9XEQuE5Ge6Os3AojIOcDBKemdaUmzslFxdahk+aacI3Tno92WE0HBdYSF3S6B1qa2x9fFNQ9jCjw75FWeB0sa8zS1XUl/YJqoPe2YxS1jTL1mweuPCD/DHgJ+l3APrzjT8J1T0DczCeprDibbjTGmUzXLNiwDfwP8jYgsBHJWkHd2yUqbd4T0UZNqatAr+xpuz1JfHipa00q7ZqLloeIHq40xBlqrsIGqHkwGrugBZjODJApa1HAECk615BOEQSUnwvI+l/mFaqJGvEZ1VJ/Lou5qexBtlRJEX1f25Ir+9/BAkUcGiowmCiU6AiNlZbAY4NXNBVbWszKSMVTTMxADS94wxkRaCl4pfjypvTBHrNWyUUGgBFotG9WbD8tDxckb8TW9eYcV81xQpRQFLojWrKAhy3DUUx4ZKDIwXMYVrWQ7KjBcChgt+4hqTZBNBrE02TsyT+x3Y4yZfZqVh7oi6xBh1Q0zA41XNqq+QFSz8lCqMJY2f5hBgUU9GWWjnMZq9GBTgcaYw9OsPNR7gD8lsYdXwtva0x0zGbLWsqb7XpnvYRHMGDNBzYLXncA2Vb21/oCI/E3bemSMMcaMo1nwugAYSzugqse1pztmppmK5SXbmNIYM1GZCRuquk9VR5JtIvLC9nfJtIsrtVmHsYIr9BWkYf3JEalkHSaPqSplP0hNqNhxsIyXUjZq1GssGxUnf4T3nNjPYmWjjJnbJppt+OW29MJMqjjrMC3zsCvn0JVSNmpxT46j5uXIR38j4pqDi7pzrFqQpysnqCp+oIyUlTEfhkphCj3Re+Wi9PgH9hYZSJSNcqIzhstK0W8sAeUTBrH6YJT1ILWVjTLGTDR42eROh4hT59O4jlBI+ZMvuMKSHrcSuGI5R1jRl2OkrIx61XxFJUyRd6DmIWVVeGbIo+yHf8GSQajkQylonI5MBqSsqiD1LG4ZM3dNNHh9otUTReQrIrJbRLZlHD9bRA6KyNbo9fHEsfNE5CEReVREPjrBPpqErBCQFRyaBY2sYCFyePczxpjD1Sxho0JEVgHHAvtE5OUAqnrzOJddTbiVyteanPNLVX193Xu5wOeBc4F+4E4R+Z6q3t9KX83U04yMi7hSxpGWjYqnHy0QGmNi4wYvEfk08BbgfsLlCQj/Ed40eKnqzSKy7jD69GLgUVV9PHr/bwNviN7fHIa0Z7UcgZ6cUA6qlTVUw8oYx8x3GSwph4rhgUCVsXLA6vk5RsoBB8aCaI1KCQLYM6J0ucqCLgfXiSvXw45DHt05YeW8HF256iDfC8APlEKOyvmQXMfSxBRkcp+v2iBmocyYuauVkdcbgeeoatrDykfqJSJyN7AT+DNVvQ9YBexInNMPnJF1AxG5CLgIYO3atW3oYmerDFbq1rHiAJB3wBVlzNPKdiQiwvwC9OVd+g+VGSpVg0ZfwaUn77B72GOkXA0gRV/ZO+KzsMthUY+LE91/zFOePFBmaa/Dst5cNSgBRU9xRenOhQt01VqKgGrqhprxDxOfbYMxM9ckP/OWrVw1zb2ZPq2seT0O5Nvw3ncBx6rqBuBfgP+M2tM+jjLX5lX1KlXdpKqbli+3PTKzZH3Ih8FEGvbRioNMHLiSHBFGy3Fx3ioF5ndVA1eyvSfvkiYeeaX+oWf8qQvNE1KMmc2Sn3nzFy2Z7u5Mm2a1Df+F8HNnBNgqIj8jUSpKVS8+kjdW1UOJr28QkX8VkWWEI601iVNXE47MzBGaaKmneBuStGsykzeatFtShzFmsjSbNtwc/XcL8L26Y0ecpSwiK4Fdqqoi8mLCUeAAcAA4SUSOA54G3gq8/UjfzxhjzOzRbDPKawBE5BJV/efkMRG5ZLwbi8i3gLOBZSLSD1xKNP2oql8kLD/1ARHxgFHgrRquyHsi8kfAjwAX+Eq0FmamwUT/lWLPXhkzfb55+1O8/Yy5sfbfSsLGhcA/17W9O6Wthqo2rTyvqlcSptKnHbsBuKGFvpkjpKoI4U7I9bufOCiLuh0OjgUNQWlxt8OBunZV2DvisaIvhzi15ab2jHh05fI4Qs2amBcoOSd7WjGr7qHVQzRmbmu25vU2wum640QkOW04n3B6z3QYSSxgKVpJiBARevKCFyhFL6xLGJZtEo6al2NBd8Czgx4lv3qvxT055hWUPSMeY174PFfJD/f/OlgscfS8HAu6wnyggDDx4+GBIkf15Via2PMr0LCkVN6BgktN1mElMNYFKo3PqPR/Un9NxpgO0GzkdSvwDLAM+IdE+yBwTzs7Zdon/qBPy+TLOQKuMliube/OOaxdmOPRfV5Ne94Vjpmf56G9RbzE/XyF/kGP490chcTzXYHCs0MeBQcWdNf+1SsHAFoJYOl9b8xitLhlzNzUbM3rSeBJEXkj4bNXCuxU1V1T1TnTPllZhCLpRw6nbNRE240xplXNpg03Al8EFhJm/QGsFpEDwP9S1bumoH/GGGNMg2bThlcD71fV25ONInIm8FVgQxv7Zdqovmp8kiOwsMth1FNKUQZHfP5xi/McHPPZPxomaoTrXMrKedWyUXF7oPDw3hLzCg6rF+bJR9uwuAKDRWXMK7OkJ1dph3C6ccwLpw7dlPIaaWWjFBBtXhzYGDP7NKuw0VcfuABU9Tagr31dMu2iGq47xTtwpX3YiwiOCL05YUHBSUwvCq4jLO5xOXZxnpwDw+UwScMRYV7B5Zj5ObpcKPvhPl++wqFiwAN7igyMlOnNCX15ByXcGmXXkMeBUY9ECkZYNspXil5Qaa9J1tDaDS3jawKlod0YM3s1G3ndKCI/IKwKH9caXAO8C/hhuztmpoaIpH7oiwho0JA+LyLkhLA4b31qvQgjpdr2qEwhC7ty5JzagKmEI6z6AAWJDTVtMGWMSdEsYeNiETmfsKL7KsLPkn7g89FzWGaWy0regOy6g1ljH9fJGOmRPt1nMcsY00zTh5RV9UbgxinqizHGGNOSzDUvEVmf+DovIh8Tke+JyN+KSO/UdM/MGrYcZYyZRM0SNq5OfP33wImEDyv3EKbQm1lOCPf7qqeqLO11Uqb2wg0p06b8dg6W8YPGZIuRsqa2ewGVzMWJsJwNY+aGZtOGyc+g3wROV9WyiNwM3N3ebpl2iJewGndVDg/EyRXV84UFXS4lXxkqBZWSUn4AR8/Ps7A7oP+gR9EPg085UPoKLgXXYf+oT9GP8xRhYNRn885RTlxSYHFis8qirzx9yGNRt8P8Lqdms8pRLwyeead2B+UsVjbKmLmj2chroYi8SUTeDHSpahkgqvxu/77tUCLhs1zhFpTR11JNkU/7vC+4wqIuwfOrIyKA3rzDSUvzeH5A0a/uxJx3heV9LgU3vL8TJWuUA3hgb4mDY37N/RXYPxYwWAqiPlaTO8oBFGtPz/i5arMYjTGzW7OR1y+A346+vk1EVqjqrmgfrr3t75ppp+ydldOn3kQks5xUVrsjkvo+9en3sfrU+1bZw8nGzD3NUuXfk9H+LOE0ojHGGDMtmmUbvrTZhSKyQESeP/ldMtOtvjKTarityeIel65cYnpOw+1TTl5aYEWfWzPl6AVK3hW63Nqb+YGybdcYW3aOMuYFNfc6OObz9KFypSxV3O4HyqgX/jfZJ1/BUwhShorx+p0lcBgzOzWbNnyziHyGsJrGFmAP0E2YdXgOcCzwp23voZlS8RScE00G+gGES06C68C8gkNPHg6OeXh+GCS6cg4r5glLe3M8sb/EvlGfchAmgjgu5BwY9QKGS0FlT7DRsseeYY+TlxVYu7AQvptCOVCeOlBiUbfDkt5cpTyVKox6Ss5RCq4QJEKlHwUwVzLWvmzjSjOHzJXdlJtNG/6JiCwGLgB+FzgaGAUeAP5NVW+Zmi6a6SAiBIFSnyshIrgoZa82McIRoeDCYCmI9uaqvWa0rJQTN4u278Lzw8BTH3T8gNTtkv1oxFUfjGyAZczc0mxLlJcAt6nql4AvTV2XzEyRVR0qK0kDoqCTIe0a18koD5WVUJJ9e2PMHNIsVf5CYIuIfFtE3h1lGRpjjDHTrtm04R8CiMgpwPnA1SKyELiJcB3s/6pqC0/gGGOJE8aYydVs5AWAqj6oqp9T1fOAVwK3EK6BNez1ZWafrGm67lz6A83Let2GbMXw/PR77x3xU8tDjXmaukdXnHCYVTYquz212RjTocYNXgAislhETiVM2vihqv6xqm5qb9fMdBMR8k64+3G9xb05FvVUA5VqWH1jzaICz13eVQluqkrJC1CEvny4F1j1/nBgzOe2/lEGRv1KyrsQZh32D3oMlYKGVPiSH66taSIXXog3pGzcsDIuCaNq6fPGzBbNnvNaKCJ/KSL3ArcB/wb8H+BJEblWRM6Zqk6a6SMi5BwhnzaayjscNS+Hr1COngUDmN/lsmFlN4Eqw6WgUuPQdYSefBgQc051j6+Sr9y7q8iuIS8sWRW1Bwp7RwMOjAU1/YHoGa9o0joudQXVXZXTxM2WNm9M52v2nNd1hLsov0xVDyQPiMiLgHeKyPGq+u/t7KCZGSaaeShRncT6Y2EwDB8urlf0NLOclEKlmG/N/bAMRGPmomYJG+c2ObaF8MFlY4wxZsq1uua1SkR+Q0ReHr9auOYrIrJbRLZlHH+HiNwTvW4VkQ2JY9tF5F4R2Soim1v/cUw7aDTycaVxlKMaJmn01s0rqionLe1izcJ8zTWqSkDjfQJVHttX4lc7RhrKRg0X08tGeSllo6r3S0/esLJRxswOzaYNARCRTwNvAe6HSsEFBW4e59KrgSsJpx7TPAG8QlX3i8j5wFXAGYnj56iqVa+fZnHgAom2Ugn3/UpujdKVc8i7wvyCMjDiMRZV31jQ7TKv4LBqQZ4Hdo8xMOpXKso7ElfSUEbLAcOl8MC+UZ/+g2U2Ht3N8Uu6qg89+8rwvhJLe12W9Lg1/YvLRnW5UvPAcxBFKsfKRhkz64wbvIA3As9R1eJEbqyqN4vIuibHb018exuweiL3N+1XDVxVEu2ZkrahZd6FQOPwFrU7QrcTBr7kACn+Mg5c8feq4Xm7h3zWLNSadS4FhkoBi7udhqocfgC4TFhKBSpjTAdoZdrwcSDf5n68F7gx8b0CPxaRLSJyUbMLReQiEdksIpv37NnT1k6aUNb+WSKCnzEfVwoaAyFE03sp7fm0/HzCqcvJnPGzwGU6TfIzb/DAvtRzvnn7U1Pcq6nXrLbhvxB+TowAW0XkZ0Bl9KWqF09GB6KU+/cCyS1YzlLVnSJyFPATEXlQVVOnKVX1KsIpRzZt2mQrGcaYWS35mXf8c9fP2c+8ZtOGcaLEFuB7dccm5RcmIuuBLwPnq+pA5eaqO6P/7haR64EXM/4amzGTwqYSjZn5mqXKXwMgIpeo6j8nj4nIJUf6xiKyFvgu8E5VfTjR3gc4qjoYff1q4LIjfT8zuVIf+1KlKyf4pcY1scXdDgfG/IYHiHMZE9fxuULtdilFL4jW4mrbg7AFUrZXCbumqdOd9YEqnvW0AGbMzNbKmteFKW3vHu8iEfkW8CvgOSLSLyLvFZE/FJE/jE75OLAU+Ne6lPgVwC0icjdwB/ADVf1hC/00k0xSUuNjtWWjtJICv6w3x7K+ZNmosDzUUX15nrusiy5XcBLtg0VtCGgCHBzz+Z/Hh9gz7FVqH6oqw2Xlwb1FDkTlpOL2QMO9xIp+tS0uC+VpNXW+lbJRVk7KmJmv2ZrX24C3A8eJSHLacD4wkH5Vlaq+bZzj7wPel9L+OLCh8QozHeIA1hBgJKxT6KpWUuNjvXmH7gXCE/tKjHlaSY+f3+WycWU3tz89yr7RgGKizEagYfp8/BIRxjzltv5RTl5a4LjFhcp7+Ao7Bj0GSwGrFuRqhkhFP3z+qzdf+++yIHqPXEo0zopNSlTBw0Zgxsw4zda8bgWeAZYB/5BoHwTuaWenzMyTUR0qKg/VeMQRqex6XH++AzWBK/kebkpJ+kPFAD9QnLpjRT8MjPWJiXEGo8UcM5d98/anePsZa6e7G23TbM3rScIivG8EVhF+HuxU1V1T1TljjDEmTbNpw43AF4GFwNNR82oROQD8L1W9awr6Z2aQrNFXwRW8oHbtSlVZ1utyqBhwqBjUnL92YZ7unMNj+0o1CRXDnuJ4MK9QfQhZVdkz7HHH0wEbj+6mO8rwUFVGysojAyXWLc5X2gG8QDk45tNXcGueFwtUGfOh4NaO8FTDEZwj8TRh9b3jn8l1ku2JCvVUpxUn2m6MOXzNpg2vBt6vqjWbTorImcBXsXWpOSUqrFGz9UjMFXDdcAuTohfgRwkS8wsOfQWHpT0a1iaMYtjK+XmW9+U4YUmBzTtHeWbQY6QcVN5n1PNZ0OUghEkYQfS+T+wv8aJjejh+cZ7hsoZ7ehFuaLlmYY7VC/KU/epU5Zjn05MX5hUcAq2u2416iitKwZVKcgaE1wngRBkayZ/RC8B1FKG2ir4mTkxrrw/2cSKIBTFjjkyzbMO++sAFoKq3AX3t65KZqUTSP3BFwpqCcQmoyghMBEeErpzQU5dA4TpCX8FhSY/DaDmoBJH4+gNjAfvHArwg/N5XKAfw0N5ipT1+m0CJAqA2rLGNlpWS35jR6Ct4KWWulDC5I7saSOORZACsbzfGtEezkdeNIvIDwsK6O6K2NcC7AEtdn8PiUVhju6RuBCkilDN2iBwsKUHqkfRpypwr4XvXBVGR9HbSm5q2N2fDJWNmgmYJGxdH1d7fQJiwIUA/8HlVvWGK+meMMeYwJWsczrbMw6ZV5VX1RmoL5hpjJiQraT+93Sp7GNOazDWvqO5g/HVeRD4mIt8Tkb8Vkd6p6Z7pNG7G36junKR+KC/udhqe04I4069xqnGkFK2P1R0r+0G0ZlbbHmcMTqS9vhJHq+3pxxr7Gn/f2F77X2NMtmYJG1cnvv574ETCh5V7CFPozRyWNTjozTt015WxEOCY+TnWLswndmMOg9O6RQXOXtdHX14qQSxQZagUMBg9nBzzgzA1/upf72fHoTJlXytBbs+wz42PDNJ/qFwpJxWocmDM595dY+wZ9mrKSY2UA/oPlhvKTJV9Zf9YEK7FJdq9APaPBRyM+pQsSzVcChiOsiKTgazoK2OeNrT7Gm4RU9tem/hhZamMaa7ZtGHyE+g3gdNVtSwiNwN3t7dbZqaLy0bVb1gpUXZh3lWGiwHiEG0oKSzsDtPnH9xbZLhUzQw8al6O1z9nPv/z+BBPHfAYS6QMHioGCOG5A6N+5QP9uvsOcdqKLk4/podDUeAAuPPpUbb3ltiwspuDxTArEaD/kMeeEZ/jFuYZ8arth4ph4Fnal6Psa6W96CulUaU3H1YKKUV9KvnKwKhPXz7MpCwl+jpUCsLajQ6V+8T3ykmYYZlMTvEURBuLEyefCTPGpGsWvBaKyJsIR2ddqloGUFUVEfs3oQGyMw8dkTAzsL7dEXKONKS0OyIs6s7xsF9uuNeopwyXgob3efqQx3GLfXJ1ZaP2j/nsGfFrdmGGsCTVwaLfuAtzNHqqP1+BkbKmTneOeZpaEb8caLShc+N71Fe8j9/D1rmMmbhmwesXwG9HX98mIitUdZeIrAT2tr9rxhhjJst4uyt3WjZi5pqXqr6n7rUran9WVX9z6rpoZrqsQUPOgfo6u6rKwi6HxT1uQ/uSXpfnLivU3C9eD0pLBBkaGePOx3YzVvZr2g+M+vz40SEOjtW2D5UCbn5yhIERr6Z9zAvY/PQIu4Zq271AuX/3GM8M1o4GA1W27y+x42CpdosVVXYPlXn6ULkhEWSkHLA/Wl+rf4+hUmN7oErJD1ISSuK1ssbfhzFzSdNU+XoicpWqXtSuzpjOlVY2KucIOaLqGH6YABEoLO11WdzjsnKeyxP7y4yUwzWrNQvzHD0/x2kru/n5E8M8O+RRjuYXu1yHvBMmQJQ8n8ED+xkdK7JD4N4d+3nZc1ZwwoqF3Le7yPYDZRS4bccI5xw/j5es7uGhgRIP7CkSKGx9doz1K7r5jTU9PHmwzN3PjhEo3LuryHGLC5y5uof9YwF37RytJIws7c2xaVUPRS9g2+5ipV+P7y+zcWU3OUd4bH+Rsaha/tOHPE5e1sW8vMPAqF+pon9gLOCovhw9+XA6NI67Y164jtadq63G73tKztFoarR2k00rM2XmsgkFL2BTW3phOlqyAG1te3jARfGolmgSCTML+/IOy3tdHj9QTWPIOcLCbpcNK7t56qHB2qK2IrjqsfvZZxHCEZmvYRbiLx/axX37hZzjVD74A4WfPzHEUwfLlaAAYTLFPbvG2HGwTE9OSO7O8sT+ErsGy/TknZp1uT3DHv/z+BA9Oacm6eJQMeCOp0dY0OXWBO6irzwyUOSovhw1QUfDey3sdmrWvxQYitbX6tfFvCDe58yeDTMm1spOykm729ILMytkfYA2Kxs1krKvF4Q7KTvSWB6q5Hm4TuPmmD5OZdfkpLIfFg726mpQhYV2G8/3NcwKrE8oUcL7pJWycqQxMQXAjUtW1cn8PaU3J46ln2GBy8xFEwpeqnpeuzpiTOvs09qYuW7caUMRORn4CHBs8nxVfWUb+2WMaUlW+SljZrdW1ryuJayo8SXAH+dcY1JJxk6WBScc/tdPx/UV0icFcq6LH6RM3gXhA8yOKJr8MNewIkbOUZDEPVXxfCXnBDiOW9Ne9sPyF26u9v8eJT/cLCXn1mZKlv0AL5CG583CbEG3ofyVH0Qbq9Q99xVEZfG14XmwjPYoY6PxfGMmbrxU+lZNVcp9K8HLU9UvtL0nZlbI2m25O+dQ8rWmIoUqrJyXJ+86bD9QivbLCp24tItCzuGnjw4x5sX7eik+wvyFixkZOkTgB2EqfRBQOrCbp/vvZ8XJLyA/bxHiuqjvM/Ts4/z0l79mw+m/wVGr1+Hm8vhemX27nuFnd9zCCzZu4JTnnUrOdfE8jwP7B/iPH93Ic08+iZee9Rvkci6+HzA0Osb3b7qNo5cs4LUv30QhFyZojJUD/vueZynkhLdvOoaevIPjCH4Ad/SPUvKV1568gHl5B9cJf+Znhzwe3lvieSu66Ms7lccJRsvK/tEyy6NsREek8qjAkB+ElUvqq3HEyRqJxTULZGYukPrioJUDIkuiLy8mTNS4HijGx1V1X9t7N0GbNm3SzZs3T3c3DI1lo2KBhtUy/KD2uB8oDw2MMVQKWyUaPfmB8oOHD/Lg7iKDxXKieK0ycmAfIwf3Uh7YgZYrfzXpWbKC+ctXse+Ruyge3FNpX7LiaE5efzqPP3Avu3fuqLQvXrKEM844g8cefpCHH7y/0r5g/nzOPe+1PDtwkM33PVJ5pqunq8Brz3kpY5Ln5kcG8KLskbwrvPkFx7BiYQ/37ipWArUALzu2l5OWdrFz0KvZ22zVghzHLSowXA5/J7HevLBiXq7h9+hIWD8yTRwELXh1vJb/AI9/7nr91NXfb2dfJmwSRl4t/fzNRl5bqJ1Q/0jimALHH16/zFzQrGxU2ujMdYRFXTmGS15D+9qFBe5+erjmfiJCLufg792OerXXjOx9hqEntzW8975dz3Dbz35A/f839u/bxw//+/qGKu+HBgf5/s9uQXNdNe2jxRLf3/wozoLlNe1lX/nFYwc5YSU176HA/XuK5BzBqZta3Dvss6jLx61rHy0rftA4HRgX880KUBa4zFzRbDPK4wBEpFtVx5LHRKS73R0zxhhjsrSSKn9ri23GtCTnSGpR2+6cNJSNApjf5bDh6N6GuYRCPsfyVcciTu3NVANwGu8TeCWKzz5KUBqtPd8vM/rUNvzhA7XtQcDwY3dS2vtk3f2Vkce3MLp9a8N7DO/bxTNPPdEwihsqBTy8r5RaBmpgxCNIeRAurDzS2B6v/9VL2zvMmNmq2ZrXSmAV8HXg7VTnQRYAX1TVU6akhxNga14zS/Kvlta0h9+VfWWkHGbeuU74FyzeH+vRfSWKntKVk+ghY2XfqM/Xfz1A/1wjYn4AABpRSURBVMESS/ryzOvKgQaUSyUeuftO9u95FgnKBF4pmppU1Cujvkd573bGntqGqBIAPWueT27FCZR3P8HwvT8Gr4Qq9KxbT/fxm/AO7mb47h8RFMPpyr4TXsS8TW8kKI0wtPl7eIN7URF61zyfRed+AKd7HqUd2ygP7sN1hHkLFnHaS86md8Fi9o36DJYU1wnXq85e18fKeS5BEGZZSvTzH7uwwMJuh5wj4e8jmmKdX5AwEUQkqsABIOQkrB+ZNlXopFTqMB2jo9e8DkfdOllLP3+z4HUh8G7CklB3Jm54CLhGVb877s1FvgK8Htitqs9POS7APwOvBUaAd6vqXYn3/1h06qdU9Zrx3s+C18yUlbwRb/5Y9GvXcDTajPKJ/eWG8wdGPb55z8GGe5bGRtny0+sJynWjG1UGt/4QHRsk8KtrY46bo3xwN1oeQxPbsMTZiAQ+mjjfzeXRXDdOzzxEg0oAdlwXd8la5q1/TbTOFyVpCHQtWMqxL3kdjuvWVASZVxDedMoCCnVbxrgCG4/upqtuWCrAsl6XvCsNASnvhNfVt8c1Dy2AdSQLXi1otuZ1DXCNiLxZVb9zmH26GrgS+FrG8fOBk6LXGcAXgDOiTMdLCQOnAltE5Huquv8w+2GmUVbyRlg2qjH5QEQoetqQ2CEiDBbDIrWluvpNXrmIBo3V2VUD/OHGvzaB7xEUh0FrnxnzvTLiewSB39Ce75mHBn5NnwLfpzBvCUS1FqvvC4V5iyqjyaS+KFuw/lcS73VWTyE1cEF64AILXGb2G3fNKxm4ROR/JnJzVb0ZaJZS/wbgaxq6DVgkIkcDrwF+oqr7ooD1E8BKU81Gh/H5mnWJNDkyGW+QfX9jzFTLHHmJyD31TcDJcbuqrp+E918F7Eh83x+1ZbWn9fMi4CKAtWs7azM1Y4yZqORn3rKVqR+Lc0Kzkdd24B7g94Dfil67E19PhrR/ymYVa0tdnFPVq1R1k6puWr58edopZgbLGsvUP/cU68oJfsocZD6fb5jqA3AcBxynISNRREAcHDclKzHQ1HbP8xruE17gpc6LquelVqH3Am3YpBNAA82sEB/u35XyHhNsN50v+Zk3f9GS8S+YpZqtef22iLwJuAr4rKp+T0TKqvpk1jWHoR9Yk/h+NbAzaj+7rv3nk/i+ZopllY0quIKglOo+5Zf0ODiSo/+gV1M26vjFed522kL+84FBxryAcrTX1TFLFvDWC97Ej376PwwODlIqe4jjUpi3mGPechn7b76G4u4nCcpFHDeH07OARc97BcVHb2fs2cdQv4yby0OuQNexG/H3Pc3Ys48iGuC4OQJx6D7uBUjgUdz1GA4g4hCIQ27R0Uh3H1ocxXHCH1RFCPLdDA4NM3/evOjh5DBbsCfv8Oywx4q+HK4TTkYK0JUX9o54LO7O1WQROhJuETOv4JJzqmuE8XpavPVKMvCFCRs2zWmm1lTVNYQm2YaVE0T6gE8CJwIvVNXVE3oDkXXA9zOyDV8H/BFhtuEZwBWq+uIoYWML8MLo1LuAF41XksqyDWe+ms0lqX7gBhqlzSuVNHEIy0Nt319iqBQwr+BURmReoPzksUHuebbIcYvzzO8KR0pBELDl7m384o676Vn1HAqLVkTvqww/eicDP7+GwqpT6FpxQuU9Snu2M/jrG8gvWUVh5cmV0ZX//9q79yDJyrOO49/f6cvMzsze2AsQFthk2VzBAK4kkUqMuSDBFEm8FaRKwTKVihqv/4ha5QX9Ay+llqXmppTEiiQRc1kNSSQXCqIhYRM3gYDIconZbBJWVhJgh92Z6cc/zunhdPfpnp5l+3J6fp+qZrrfc+lnDtv9zHnPe573if9j/sF9VDdsZXrnBahSS9/j2FEWHnmAZN0GZs97DcnUbPoejQaiQVKdYnbn91GZnl3+vbdtmmP9uimet7W+HCvAKdMJW2aqbJ2ttJR9mq6KjVPpCMP8oIxqAuun0mK/7YM1mmdzHiZfeqUdbXiSktczLg8FQEQ8Cfy6pBcDL1tVBNKNpGdQWyUdJB1BWMv2+07gZtLEdYB0qPzPZsuOSPoD0iH6ANeOYy1FWz2p+F9mIjFV7Zw0spKkNf6qR1tH+VUT8dIds1QTtWyTJAkvetGLuGt+E0u5YX6SmNt9EQtPHOno4qtv28mGF19KY6l1aH5lbjMbL7ysY7RgMjXDhgsuQ3NbWn+3JKF+yhnMnfYcGm2/ZTUWueD0DVTauh2PLQXP2lDtGGV4bDGYnu0cYbjYSIfHJ763y9a4nslL0iuA70TEfcB6YE7Sj0bEx/rZeURcucLyAH6xy7Lrgev7eR8zM1tbug7YkPQXwHXAP2RnQX8MrAN+TdKfDCk+W2OKz8pgXbVzSTWBHRtqhe17dsxRSdrbxYufv4upeus2ktj9ghcyOzfXsa8zn7ObDQUXxU897TS2bN7U0b5lwyzb1k91tK+fqjBXMEdZ8ybjdvku1XbdykC5PJStJb3OvF4LnEuasL4JnBERRyVdB/wnrVXmzZ6xqqBagaWAhawrMAHW1xPm6glPLQaPPLnIYgPm6gmnr69yzpYpHj+2xK0PPcmj80ucOlvlnC0z/NDOOd58bIl3/cc3uec7Rzln2wyXn7edqcrZLCws8k8fv5X9997P6adu50cveTWzszNENPjsJz7Gl+74HJu2bOOSN13JKdtPA8RXvvA5vnjbp5memeG1r/8xduzcBRIPPPw//Pud/4kEr3zJhTxv106Q+Pb3jvG5B45wbLHBK8/ZxEt2bqKSiOOLwTcfX+T4UnDWxhrP3VJfTl5L2fXAuXrC9tnq8jWs5vQx1SRdthRBA5FELHcVNgfERIAifJOyTbxe5aHujohzswry3wKeFRHzkirAXRHxwmEG2g8P2JgMzQkYF7N/msq1NwKeWmz9co4IFhrBQ0eOZ3USn/7SPrawxIEjx0GiljsVO76wwLcfmycqtXRm5GxfiwvHefzoU0Rlimq1utzeWFxkcXGBerVKpVpB2azM0WjQaCwxVUuoVpKnZ2XOPldbZyrUq0nb0P9g+0yVekUt7RHBaetrTFXVck0rIlhXVWGVjaKBG81j5gRWWgMdsDHMEYEn6BkP2PiYpNuBaeBvgQ9KugP4IeC2Zx6fWbHlslEF7UR0zImVfkEHlUQdw/Hr1YRatfOerXqtRnWqc/h+tVantq7aOUijWmVdtdKZJJKEqWpCvb1MvsRMTUzXko5tpipJR4KCtFtzutqZoKTixAUuD2VrV6/7vH5D0svSp3GHpF3Am0gT2U3DCtCsX93uJTuR/RS2n0AyWO023e7QN7NWvcpDKSI+33wdEQ8Af1qwjq8Qm5nZUPUqD/VZSb8kqbVWvVSX9CpJNwBXDTY8W6tWe/YhVFiKqec2Pd5kNSV+ky73rnU7DYyCLtHm6r1/79WVgfLflTbJeiWvS4El4EZJhyTdI+kh4H7gSuDPI+LvhxCjrUEShfX/EkHBiHOmqmLnpqysUnMfwFQ14QfOmGa2puVRfYlgqiJecfYMW2cqy7M6V5SWq3r5zhnO2lhbbk+UjvS74FnTnHNKbXk/yrZ57pY6z9uajhpsvndFsH22yo711Zbklgjmpipsmq60tAuoJeq41tZcls4A03lAguKuSectm3S9rnk9BfwN8DeSasBWYD4iHuu2jdnJIqXJJolYHiqeiOXSTbUIji0GSwH1Spos1lWrbJiq8O0nFnn06BJbZypsnE4HWZw6W+Whx45z/6ML7NhQ5exNNRKJMzfWePixBb50aJ4zNtQ4/7R11Cri+VunOfi9BW57+Em2zFT4wTNnlss3PX/7Iv/+9aPUK+Lis2fYkJV72rW5zp2H5jm+FOx51gyb16Xt2+dqfP2x4zy50OCcU6aW29dPB4/NLzG/0GDLTIXZejq4o1lCS6S/W71gsEaCZ1G2tW3F2oZl4qHyk6nXTMzF6wdPdRaYB+DoQnHnYntZqqalolOhHrrFKqBadDcyzRGDne2zteIRhnUnrUl3QkPlSzAEvl99/f4rTkZpZmY2bpy8bCJ1mQ6M9tuxmmpJ8Z97lS7tSY9rct3KPXWLabUnS81KGv22m02iFavKm42D5SobubZ8Msj37lUSsS5J244tpVU5kmwwxrpqhUbA48eXWMzmAltfT6hmNzh/96kGx5bS0YCzdVGvJATwxLFGWtkDmKknTGe1Fo8uNJhfSN98XU3L18XmF4InjzcI0ulNZutpCaeFRrosSJNccz8RsNDItyeFSa2iLEFlN7W1rBMQCiLUszai2SRw8rKx1/wSjmg9C8qXh8onsuVJHAmmq+kIPuWWVQQbpyosNWK5PFOzPuDG6YTFpac3aLavn0qYqQFqnehxppYwXe0sx7SuBtO1Sjorc669lgTVejqNS359KahL6Zlbkr/WlZbKSkQ2ASY0g4vsP8vTzORLSuWWmU0iJy8rjZYkVpCsIqJlwEJaTgoSdd5BJall0sumREIFXYuSSJIu+1FnUM1SVpXCUlaQdLxHOstyZ0xp+9OJqz2uLuWhuqxvNimcvKx0un0pdx9pV9y++pF5PdZ3prARmaBRhqviARtmZlY6Tl5mJdf9frchB2I2RE5eNtGWBzOsQsGkzUB23anoPbq006N9tcPmu+Wh5ijL9gTmvGWTzte8bOJl0321fKE3BzS0V8QQkCTpLMVLkdYUFFkVjCzjNCKWZ3qu5u73auQm0Eyy/UCaWPITaybZoIyIYLGRblfJDdbIxyR6DNboaI+sXb4EZxPPycvWhG5nYN3blSamglOYRKJesCAR1Og862nuq6i91jlP5nJM3c7ouo0kzA/hN5t07jY0MyupU2brow5hZJy8zMysdNxtaNZDUVmqnutnXX7t19KafyUGtEyaWenSvtoBIE/fwB3uOrQ1wcnLrIdclabWgR25qh4UtEMsJ7H2/VR4er6ufHtC60CNjljaY8o0B4w4adla4uRl1ofmyMTOdhVkqKdHE7bnk5ayVu3tPU7vug3S8BxetlYN9JqXpEsl3SfpgKRrCpb/uaT92eO/JT2WW7aUW7Z3kHGaPSNpgcHVbTKgUMzWioGdeUmqAH8NvBY4CNwpaW9E3NNcJyJ+Lbf+LwEX5HYxHxHnDyo+MzMrr0GeeV0EHIiIByPiOPB+4A091r8SuHGA8ZiZ2YQYZPI6A/hG7vXBrK2DpLOBZwOfyTVPS9on6Q5Jb+z2JpLemq237/DhwycjbrNCg+7q69Xz2O1yWET32oY2mfydlxpk8ir6KHb7lF0B3BQRS7m2syJiD/Bm4C8k7SraMCLeHRF7ImLPtm3bnlnEZj0UVeMQ2USRBcsSiUqiwjqG3QZZFO2n+R5F6yYtE1faWuDvvNQgk9dB4Mzc6x3AoS7rXkFbl2FEHMp+PgjcSuv1MLORUC5RtY8AzJd1SnJjOJqzMS8/X549uds8Y0+/R34/y8tpn23ZbO0ZZPK6E9gt6dmS6qQJqmPUoKTnAZuBz+faNkuayp5vBS4G7mnf1mxUug0w7N6++mTTfdLNVe3GbCINbLRhRCxKejvwSdL7Mq+PiK9JuhbYFxHNRHYl8P5o7bh/AfAuSQ3SBHtdfpSi2VrWvMHZbC0b6E3KEXEzcHNb2++0vf69gu3+AzhvkLGZDdsK9yGvSiPCU5/YmuYKG2ZDsNoyUyvsrbkrIrpX3zCbZK4qbzZEveYPW3nbzutmTly2Vjl5mZlZ6Th5mZlZ6Th5mZlZ6XjAhlmJxfJ/fO3L1hafeZkNWb6CRmt79xuZew3oaI5TTOscnpwYzcadk5fZiHSvoKGer4ssz8Dssy9bI9xtaDZC3W5cdt1Cs9585mVmZqXj5GU2YoU3Lece3ZYV8TUvWyucvMxGqFmFPp+QlFtGW7tWaHdvo60VvuZlNgbytQ/b5whbqb3ltdka4eRlNkZWO4eXk5atVe42NDOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0nHyMjOz0hlo8pJ0qaT7JB2QdE3B8qslHZa0P3u8JbfsKkn3Z4+rBhmnmZmVy8AK80qqAH8NvBY4CNwpaW9E3NO26gci4u1t254C/C6wh7Ru9peybf9vUPGamVl5DPLM6yLgQEQ8GBHHgfcDb+hz2x8BbomII1nCugW4dEBxmplZyQwyeZ0BfCP3+mDW1u7HJX1V0k2SzlzltmZmtgYNMnkVzTTUPkn5vwA7I+L7gE8BN6xi23RF6a2S9knad/jw4RMO1sysDPydlxpk8joInJl7vQM4lF8hIh6NiGPZy/cA39/vtrl9vDsi9kTEnm3btp2UwM3MxpW/81KDTF53ArslPVtSHbgC2JtfQdLpuZeXA/dmzz8JXCJps6TNwCVZm5mZ2eBGG0bEoqS3kyadCnB9RHxN0rXAvojYC/yypMuBReAIcHW27RFJf0CaAAGujYgjg4rVzMzKRRGFl5JKac+ePbFv375Rh2Fm9kwUXfMvNKHfeX39/q6wYWZmpePkZWZmpePkZWZmpePkZWZmpePkZWZmpePkZWZmpePkZWZmpePkZWZmpTNRNylLOgx8vWDRVuB/hxzOao17jOMeHzjGk2XcYxz3+OCZxfi/EdHXFFCSPtHvupNmopJXN5L2RcSeUcfRy7jHOO7xgWM8WcY9xnGPD8oRY9m529DMzErHycvMzEpnrSSvd486gD6Me4zjHh84xpNl3GMc9/igHDGW2pq45mVmZpNlrZx5mZnZBHHyMjOz0pnI5CXpFEm3SLo/+7m5y3pLkvZnj71Diu1SSfdJOiDpmoLlU5I+kC3/gqSdw4hrFfFdLelw7ri9ZcjxXS/pEUl3d1kuSX+Zxf9VSRcOM74+Y3ylpO/mjuHvDDm+MyV9VtK9kr4m6VcK1hnpcewzxlEfx2lJX5T0lSzG3y9YZ6Sf54kWERP3AP4YuCZ7fg3wR13We2LIcVWAB4DnAHXgK8AL29b5BeCd2fMrgA+MWXxXA381wv+3rwAuBO7usvwy4OOks7G+FPjCGMb4SuBfR3gMTwcuzJ6vB/674P/zSI9jnzGO+jgKmMue14AvAC9tW2dkn+dJf0zkmRfwBuCG7PkNwBtHGEveRcCBiHgwIo4D7yeNNS8f+03AqyX1PS34EOIbqYi4DTjSY5U3AO+N1B3AJkmnDye6VB8xjlREfCsivpw9fxy4FzijbbWRHsc+Yxyp7Ng8kb2sZY/2EXCj/DxPtElNXqdGxLcg/RAA27usNy1pn6Q7JA0jwZ0BfCP3+iCdH8jldSJiEfgusGUIsbW8d6YoPoAfz7qSbpJ05nBC61u/v8OovSzrbvq4pBeNKoisG+sC0rOGvLE5jj1ihBEfR0kVSfuBR4BbIqLrcRzB53miVUcdwImS9CngtIJFv72K3ZwVEYckPQf4jKS7IuKBkxNhoaK/uNr/UutnnUHp573/BbgxIo5JehvpX5WvGnhk/Rvl8evXl4GzI+IJSZcBHwF2DzsISXPAPwO/GhHfa19csMnQj+MKMY78OEbEEnC+pE3AhyWdGxH5a51jcRwnUWnPvCLiNRFxbsHjo8B3ml0c2c9HuuzjUPbzQeBW0r/uBukgkD9T2QEc6raOpCqwkeF1Qa0YX0Q8GhHHspfvAb5/SLH1q59jPFIR8b1md1NE3AzUJG0dZgySaqRJ4X0R8aGCVUZ+HFeKcRyOYy6Wx0i/Q9qL5I7y8zzRSpu8VrAXuCp7fhXw0fYVJG2WNJU93wpcDNwz4LjuBHZLerakOukF3PZRjvnYfwL4TEQM6y+1FeNru+5xOem1iHGyF/iZbLTcS4HvNruQx4Wk05rXPSRdRPo5fHSI7y/g74B7I+LPuqw20uPYT4xjcBy3ZWdcSFoHvAb4r7bVRvl5nmil7TZcwXXAByX9HPA/wE8CSNoDvC0i3gK8AHiXpAbpP/rrImKgySsiFiW9Hfgk6ci+6yPia5KuBfZFxF7SD+w/SDpA+hfaFYOM6QTi+2VJlwOLWXxXDys+AEk3ko4y2yrpIPC7pBfKiYh3AjeTjpQ7ABwFfnaY8fUZ408APy9pEZgHrhjyF9rFwE8Dd2XXawB+CzgrF+Ooj2M/MY76OJ4O3CCpQvod8sGI+Ndx+TxPOpeHMjOz0pnUbkMzM5tgTl5mZlY6Tl5mZlY6Tl5mZlY6Tl5mZlY6Tl5mZlY6Tl5WOmqdymZ/t2kmJM1Iep+kuyTdLelzWbmh/D7ulvRPkmay9ieynzslzWfr3CPpvVnFh6KpOPZLek2PeB/OYtgvaV+uvXDqnuzG4JFO62I27py8rIzmI+L83OPhLuv9CvCdiDgvIs4Ffg5YaNvHucBx4G0F2z8QEecD55GWR/qp3LLb22L41Aox/3C23p5c2zXApyNiN/Dp7DXA60hr9O0G3gq8Y4V9m605Tl42yU4Hvtl8ERH35eoy5t0OnNNtJ1nx1S9y8quqd5u6Z+TTupiNOycvK6N1ue66D/dY73rgNyR9XtIfSuqoOJ4VS30dcFe3nUiaBl4CfCLX/PK2bsNdPeII4N8kfUnSW3Pt3abuGZvpSMzG1aTWNrTJNp915/UUEfuVTndzCWnR1DslvSwi7iVLgNmqt5PWoGu3K1tnN3BTRHw1t+z2iHh9n/FenE29sx24RdJ/ZRNWduNpNMxW4ORlEy2bMuNDwIeyIsyXkVbC7ycBPhAR52dddrdKujwrtrraGJpT7zySnSleBNxGNnVPRHxLrVP3jHw6ErNx525Dm1iSLs6N4KsDLwS+vtr9ZF161wC/eQIxzEpa33xOehbYnKyw29Q9Yz+ti9mo+czLJtku4B3ZnE8J8DHSyQ1PxEeA35P08uz1y3PdjgB/GBE3FWx3KukMu5B+3v4xIprXzgqn7mH005GYjT1PiWJmZqXjbkMzMysddxta6Un6EeCP2pofiog3DTGGLaQ3Grd7dUQMbWp6s7XC3YZmZlY67jY0M7PScfIyM7PScfIyM7PScfIyM7PS+X90nSOmEwn0vwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x432 with 3 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "skew=(cat['FErr_SPIRE_500_u']-cat['F_SPIRE_500'])/(cat['F_SPIRE_500']-cat['FErr_SPIRE_500_l'])\n",
    "skew.name='(84th-50th)/(50th-16th) percentile'\n",
    "g=sns.jointplot(x=np.log10(cat['F_SPIRE_500']),y=skew, kind='hex')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For 500 $\\mathrm{\\mu m}$ depth is ~ 6mJy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Add flag to catalogue"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "cat.add_column(Column(np.zeros(len(cat), dtype=bool),name='flag_spire_250'))\n",
    "cat.add_column(Column(np.zeros(len(cat), dtype=bool),name='flag_spire_350'))\n",
    "cat.add_column(Column(np.zeros(len(cat), dtype=bool),name='flag_spire_500'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "ind_250=(cat['Pval_res_250']>0.5) | (cat['F_SPIRE_250'] < 6)\n",
    "ind_350=(cat['Pval_res_350']>0.5) | (cat['F_SPIRE_350'] < 6)\n",
    "ind_500=(cat['Pval_res_500']>0.5) | (cat['F_SPIRE_500'] < 6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3492 12024 32257 112471\n"
     ]
    }
   ],
   "source": [
    "print(ind_250.sum(),ind_350.sum(),ind_500.sum(),len(cat))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "cat['flag_spire_250'][ind_250]=True\n",
    "cat['flag_spire_350'][ind_350]=True\n",
    "cat['flag_spire_500'][ind_500]=True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set XID+ cahtalogue back to orignal order of objects, as used in MF detection files\n",
    "use = cat['HELP_ID'].astype(int) -1\n",
    "use = np.argsort(use)\n",
    "cat = cat[use]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "# galaxies =  112471\n",
      "# galaxies =  112471\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING: MergeConflictWarning: Cannot merge meta key 'EXTNAME' types <class 'str'> and <class 'str'>, choosing EXTNAME='GAMA-12_SPIRE500_cat_MF0.fits' [astropy.utils.metadata]\n",
      "WARNING: MergeConflictWarning: Cannot merge meta key 'DATE-HDU' types <class 'str'> and <class 'str'>, choosing DATE-HDU='2018-06-09T12:09:15' [astropy.utils.metadata]\n",
      "WARNING: MergeConflictWarning: Cannot merge meta key 'STILVERS' types <class 'str'> and <class 'str'>, choosing STILVERS='3.1-' [astropy.utils.metadata]\n"
     ]
    }
   ],
   "source": [
    "# Reads MF table, removes duplicate RA and DEC\n",
    "cat2=Table.read('./data/GAMA-12_SPIRE_all.fits')\n",
    "print('# galaxies = ',np.size(cat2['RA']))\n",
    "print('# galaxies = ',np.size(cat['RA']))\n",
    "del cat2['RA']\n",
    "del cat2['Dec']\n",
    "cat_all = hstack([cat,cat2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Created HELP_ID, and changes HELP to HELP_BLIND to avoid confusion with HELP-Masterlist objects\n",
    "ID = gen_help_id(cat_all['RA'], cat_all['Dec'])\n",
    "ID_new = [IDs.replace('HELP','HELP_BLIND') for IDs in ID]\n",
    "ID_new = Column(ID_new,name=\"HELP_ID\")\n",
    "cat_all['HELP_ID'] = ID_new"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# all flux denisties are in mJy in the final BLIND catalogues\n",
    "cat_all['F_BLIND_MF_SPIRE_250'].unit = 'mJy'\n",
    "cat_all['F_BLIND_MF_SPIRE_250'] = 1000*cat_all['F_BLIND_MF_SPIRE_250']\n",
    "cat_all['FErr_BLIND_MF_SPIRE_250'].unit = 'mJy'\n",
    "cat_all['FErr_BLIND_MF_SPIRE_250'] = 1000*cat_all['FErr_BLIND_MF_SPIRE_250']\n",
    "\n",
    "cat_all['F_BLIND_MF_SPIRE_350'].unit = 'mJy'\n",
    "cat_all['F_BLIND_MF_SPIRE_350'] = 1000*cat_all['F_BLIND_MF_SPIRE_350']\n",
    "cat_all['FErr_BLIND_MF_SPIRE_350'].unit = 'mJy'\n",
    "cat_all['FErr_BLIND_MF_SPIRE_350'] = 1000*cat_all['FErr_BLIND_MF_SPIRE_350']\n",
    "\n",
    "cat_all['F_BLIND_MF_SPIRE_500'].unit = 'mJy'\n",
    "cat_all['F_BLIND_MF_SPIRE_500'] = 1000*cat_all['F_BLIND_MF_SPIRE_500']\n",
    "cat_all['FErr_BLIND_MF_SPIRE_500'].unit = 'mJy'\n",
    "cat_all['FErr_BLIND_MF_SPIRE_500'] = 1000*cat_all['FErr_BLIND_MF_SPIRE_500']\n",
    "\n",
    "cat_all['F_BLIND_pix_SPIRE'].unit = 'mJy'\n",
    "cat_all['F_BLIND_pix_SPIRE'] = 1000*cat_all['F_BLIND_pix_SPIRE']\n",
    "cat_all['FErr_BLIND_pix_SPIRE'].unit = 'mJy'\n",
    "cat_all['FErr_BLIND_pix_SPIRE'] = 1000*cat_all['FErr_BLIND_pix_SPIRE']\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAELCAYAAAAhuwopAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3XdgHMX1wPHv22uSLFu2XAA3bMB0g7FND4SWXwghQEIJKbRQQ+8t1NBD72CaIfQSWgATesfgQjEYg6nuBctF7dq+3x+7J51kldPpTqfyPiDpbm9vZ3Q+7buZnXkjqooxxhjTFk6hK2CMMabrseBhjDGmzSx4GGOMaTMLHsYYY9rMgocxxpg2s+BhjDGmzfIaPETkXhFZLCIz0raVi8grIvKN/7Ofv11E5CYRmS0in4nI2HzWzRhjTPby3fKYCOzeaNvZwGuqOgp4zb8P8BtglP91FHB7nutmjDEmS3kNHqr6NrCs0ea9gfv92/cD+6Rtf0A9HwJ9RWStfNbPGGNMdoIFKHMNVV0AoKoLRGSQv30IMCdtv7n+tgWNDyAiR+G1TujVq9e4DTbYML81NsaYPOjw/B7JJDJ3DvLzz0yFpao6MNtDFSJ4NEea2Nbka6uqE4AJAOPGj9f3PpySz3oZY0xedGR2KOeZpwmfeCwsX07yjLMIXX3Vj+06Xq4q1gaLUt1R/s/F/va5wLC0/YYC8zu4bsYY070sWkToTwcQOeAP6KA1iL03mcTlV7b7sIUIHs8Bh/i3DwGeTdt+sD/qahtgRap7yxhjTBupEvj3AxRtthGB558l/s/LiH7wMTpuXE4On9duKxF5BNgJGCAic4ELgSuBx0XkcOAnYH9/9xeBPYDZQDVwWD7rZowx3ZX89BOh444m8PIkkttuR/yOu9GNNsppGXkNHqr6p2Ye2rWJfRU4Lp/1McaYbs11Cdx5O6F/nA2qxK6/ieTfjwMn951MnemCeXZsORJjTJelND1WqO1k1ixCxxxB4L13Se72K+K3TUBHjMjJsZti6UmMMaZDKarel3+v7r+sPg3H4wT/dSWR8ZvjfDGD2N33EXvh5WYDR65GeHX9locxxnRq2uCE3dK5W9O+17dHmm+ZyPTphI8+HOeT6ST3+QOxm26FNddset/cNHDqWMvDGGNyrr51oUpdm6ItH/rrn9NEq6S2luD5/yCy3ZbIgvlEH32S2ONPrRY4ROq/cs1aHsYY007aqC8oH5diU60S5/33CB91BM7Xs0gcfCjxf10L5eVAfoJEc6zlYYwxOZBN66JNKisJnXIikZ13hGgt0RdeJnHPfUj/8ry1LlpiLQ9jjOnknP+9TPjYo5E5c0gcexyJSy5HevcubJ0KWroxxpjmLVtG+PDDKNrzN1BUTPSNt4lffxOUlha6ZhY8jDGmMwr85ymKN9+EwMMPEj/rHGqnTMfdbvtCV6uOdVsZY0xnsmAB4ZNOIPjMf3DHbEH0+ZfQMWMKXavVWMvDGGM6A1UCD0z0WhsvvUDssiuofe/DThk4wFoexhjTbu0dYSU//ED42KMJvPoKye1/QeyOu9ANNshJ3fLFWh7GGJOF1CRAtz35PpJJgrfcRNEWo3E+/IDYTbcQfe3NZgOHkD7fvLCJ/azlYYwxGajPRZUbMnMm4WOOJPDB+yR/vTuxW25H11579f3Sb0vjLYVjLQ9jjGlCesvC1azTFq4uHid4xWUUbbkFzqyviN57P9HnXqgLHJL25QiISN1XZwkcYC0PY4ypk+vWRWMyfRqRIw/H+exTEvsdQOz6G2GNNTpt66Il1vIwxpg0eQkcNTWEzj2bou22RhYvIvrEf4g9/CjiB47O2rpoSddveXSN19kY00M5775D+Ogjcb75msRhfyN25dXQr1/d49LRSalyxFoexhiTDytXEjrxOIp2+SXEY9S+9D9id97dIHB0ZRY8jDEmx5xJL1G0xWiCd95B/MSTqJ3+Oe6uuxW6WjnV9butjDGms/j5Z8Knn0rwoX/jbrgR0bfexd1m20LXKi+s5WGMMe2lSuCJxynebGMCjz1C/NzzqP14WrcNHGAtD2OMaReZP5/QiccTfO4ZkmPHEXvpFXSzzQpdrbzr+i2Pws7QN8b0VKoE7ruHos03IfC/ScSuuIroux9kETi65knMWh7GmE6npXRROR/Z6s8eB/803rjsJsqT774j/PejCLzxOskddvQSGY4a1fai/W9dcbRu1295GGO6PvUCRuqrxV3bsG9zhaVSj3jpR9KO1dTx0hcnTyQJ3nQDRWM3w5nyMbFbbiP6yuttChz1qUcEp25iYNdjLQ9jTEG0Jxltc8dpOrOHNtgn23Llyy8IH30EgY8mk9x9Dy+R4fBhrT8v9bOLBonmWMvDGNMh2t9iyKCMBuWkWhftLDcWI3j5JRRtNRbn29lE73+Q6DPPo8OGNWyVpElvXUgXbl20xFoexphuKRcBypnyMeGjj8CZ8TmJAw4kdt2NMHBgk/umclT1FNbyMMaYxqqrCZ19BpEdtoVlPxN98hli/3642cABPS/NnrU8jDEmjfPWm4T/fhTOt7NJ/O0IL5FhWVmhq9XpWMvDGGMAVqwgdNwxFP3fLuC61E56ldjtEyxwNMOChzGmx3NefIGiLTYleO/dxE86hdppn+HuvEuhq9WpWfAwxvRcS5YQPuSvFP3+d1DWl+hb7xH/17VQUlLomnV6FjyMMXmXr6G5zZfX+kzDwOOPUjxmEwJPPUHs/AupnTwVd6utO6aC3UDBgoeInCIiX4jIDBF5RESKRGSkiEwWkW9E5DERCReqfsaY9sn3nI6GZTWcNd4SmTuX8L77EDnoz7gj16F28lQS510I4bafbkTqv3racKuCBA8RGQKcCIxX1U2BAHAgcBVwvaqOAiqAwwtRP2NM23XEJMCG5WUWLOq4LoG7J1C0xaYEXn+V2FXXEH3rPXSTTdtUbnrA8L661trjuVLIbqsgUCwiQaAEWADsAjzpP34/sE+B6maMaUUhg0XGAcMns2cT2X03Iscdgzt2HLVTPyNx8qkQCLT+3LRA4Tj1waInBox0BQkeqjoPuAb4CS9orACmAstVNeHvNhcY0tTzReQoEZkiIlOWLF3SEVU2xvg6KljUl9f2YFEnmSR4/bUUjd8cZ/o0orfdSXTSq+i66zb7lAbBQixYNKdQ3Vb9gL2BkcBgoBfwmyZ2bfIdo6oTVHW8qo4fOKD5GZ/GmJ5LvphBZMftCJ99BslddqP2ky9IHn5kRvnP64JFD0o30laF6rbaDfheVZeoahz4D7Ad0NfvxgIYCsxv9Uj2b2uMSReLEbrkIoq2Hofz4w9EH3yE2FPPoEOa7MgwWSpU8PgJ2EZESsRrC+4KfAm8Aezn73MI8GyB6meM6YKcjyZTtPU4Qpf+k+R+B1DzyRck9/+jtSDyoFDXPCbjXRifBnzu12MCcBZwqojMBvoD9xSifsaYLqaqitAZpxLZcTtYuYLap58nNvHfMGBAoWvWbRUsMaKqXghc2Gjzd8BWBaiOMaaLct543Utk+P13xI88mvjlV0GfPoWuVrdnM8yNMV3T8uWE/34kRbvvBo5D7StvEL/ldgscHcSChzGmywk8/5w32W/ifcRPO4PaqZ/i7vjLQlerR+n663l0cM4cY3q2xn9w+b4Q3ai8xYsJn3oSwScew910NNEnn8EdNz7PdTBNyTh4iMgaeJP2FJivqovyVitjTOeh2uAU3nC+nnfHG8yUi0CiDcqoK0qVwKMPEzntZFi1ithF/yRx2plZ5aNqiQ3KylyrwUNExgB3AGXAPH/zUBFZDhyrqtPyWD9jTIfTBgEik8nd3j71O7YtmDQqr9GjMmcO4eP/TnDSiyS33obonXejG21c/3iGpTQnFTB60vrjuZBJy2MicLQ/vLaOiGwD3Adsnod6GWM6jKb+r9/Szu7g9GAiDb/RWrCo47oE77qT8D/OhmSS6DXXkzj2+NXyUbW1I82CRW5kEjx6NQ4cAKr6oYj0ykOdjDEdID1fVD5zVWnDbzS81TT55hsifz+SwDtvk9xlN6K33YmOHJl5edQHEQsW+ZFJ8HhJRF4AHgDm+NuGAQcDk/JVMWNM/nVogsNMdkokCN1wHaFLLoJIhOidd5M45LCsLkak0qWb/Gg1eKjqiSLyG7xEhkPwAvpc4FZVfTHP9TPG9BDOZ58SPvoIAtOmkthrH2I33oIOHlzoaplmZDTaSlVfAl7Kc12MMT1RNEroiksJXX0VlJdT+/BjJP+wnw196uQyGW21map+5t8O4eWf2gqYAVyqqtX5raIxprtyPvyAyNFH4Hw1k/hfDyb2r2uhf/9CV8tkIJMZ5hPTbl8JrAdcCxTjDeE1xpi2qawkfNrJFO30C6iqpPbZF4jdM9ECRxeSSbdVettxV2BLVY2LyNvAp/mpljGmu3JefYXIsUfj/PgD8WOOJXbpFdC7d6GrZdook+BRJiK/x2ulRPzFm1BVFRFLDmJMF+Qt7Vo/2irflxdUQSsqiJx9OqH778MdtT41r72F+4sd8luwyZtMgsdbwF7+7Q9FZA1VXSQiawJL81c1Y0yupOZ0uM2sP14XRPz77Q0mjdOLBJ59mshJxyFLlhA7/Sxi510IRUV1O+QieKXmITb8HZT2z0E3TclkqO5hzWxfiNeNZYzphNJbF5l2EaT2a2swWS0XVep5ixYRPvVEQk89QXLzMdQ+/V/cLcY2+/xMympw/PRg0dYnm3bJKCW7iPQRkXWb2L5Z7qtkjMmGquKqkkwqiaSSdP2WRnuO6X+lWiyNc16p1pfRMF+iEnzwAUrGbEzw+WeJXnwpNe9ObjJwrP57NPxKJ4AjDb9ExJsMaIGjQ7UaPETkAOAr4CkR+UJEtkx7eGK+KmaMyYxqWrBw87dKQXogaTJg+OTHHynaaw+KjjgUd4ONqP74E+JnnQuhUNZl1wULpz5YeLPHLWAUSiYtj3OBcao6BjgM+LeI/MF/zP7ljDEe1yV0x62UjBtN4P13iV5/EzWvv41usGG7DtugdWGnnE4jkwvmAVVdAKCqH4nIzsB/RWQonWEpJnsvGVNw8vUsio45ksD775LY7f+I3nIHOmJEoatl8iiTlseq9OsdfiDZCS/X1SZ5qpcxpiuIxwldfSUlW47BmfkFtXfdR+3zL1ng6AEyaXn8nUaf71V1lYjsDhyQl1oZYzo955PpRI45gsAn00n8fl+i19+MrrlmoatlOkirLQ9V/VRVZzexPa6qD6Xui8gHua6cMaYTqq0lfP65FG+/FbJgPjWPPEHtI09Y4OhhMl7DPANFOTyWMaYTct5/j6JjjsD5ehbxgw8leuU1UF5e6GqZAshonkeGCnPxvPCX7I3p/latInzKCRTvuiPU1lLz/CSiE+61wNGD5bLlYYzpaJq+HnjjBVhzI/DKy0SOPRqZO4f4sccTvfhSKC1FOij1h7eCrdqqgJ1MLoOH/csak2+qDRrbroK62qgB7t2TVLKnbC1bRuTMUwk9+ADJDTak9rW3SG63faNS/LLqtubuNCCAOP5PCxydTiYzzDdMux1p9Ng2aXcPymG9jDE+L0eV+ulHIOl6M8oTScVdLXCkP6/+uV7zJJM+Xm+/wNNPUrLFJgQfeYjoWedSPXlqg8DR9LNA/f+an3/eMkcg4EAwIAQCglM3OdB0Nplc83g47XbjEVW3pW6o6oyc1MiYni49WLjeV6IuX5Wulu8po0PSKJg0OLlr/Ul/wQIiB+5H8Z8OQNcaTPX7HxG7+BIvA25by6NxMFmd4AWLVMBwHAsWXUVbF4Nq/K9q/8rG5FCqJaGaXZDIVN2xU0vyqBJ8YCJFZ50ONTVEL72C2MmnQjA3PdupqyPWFdV9ZPLO0GZuN3XfGNMOSbfj/6Tkhx8oOu5ogq+9SmL7Hai9YwI6av2cl+MlNrSA0V1kEjyGishNeB8WUrfx7w/JW82MMfmVTBK68zYiF/wDRKi98RbiRx7tpa41phWZBI8z0m5PafRY4/vGmC7A+Wqml8hw8gck/u/X1N5yBzp8eKGrZbqQTFYSvL+5x0Rk7dxWxxiTV/E44euuJnz5JVBaSs09E0n8+a+2kJJps0xXEtxWRPYTkUH+/c1E5GHg3WwLFpG+IvKkiHwlIjP9MspF5BUR+cb/2S/b4xtjGnKmTaVk+62IXHQ+id/tTdW0GST+cpAFDpOVTOZ5XA3cC+wLvCAiFwKvAJOBUe0o+0ZgkqpuCGwOzATOBl5T1VHAa/59Y0x71NQQPu8cSnbcFlmyhJrHnqL2wUfRNdYodM1MF5bJNY/fAluoaq3fEpgPbKaq32RbqIj0AXYEDgVQ1RgQE5G98dYKAbgfeBM4K9tyjMlWc8Nk8/UhXVX9yX8uIlK3el57Bd59m6K/H4Uz+xtih/6N6BVXo337phVcfzOfDRCpO37HpDQx+ZdJt1WNqtYCqGoFMKs9gcO3DrAEuE9EpovI3SLSC1gjbdXCBcCgpp4sIkeJyBQRmbJk6ZJ2VsWY1AS6+q9M9mtfeV6wiMaTVMcSVMeSRBMuLpBUJe4qsaRLIuni1k3sa4OVK4mcfDwlv9oZEgmqXvwftbff1TBwtPC75WKOiTdbXAg6QsAR6x3rZjJpeawrIs+l3R+Rfl9V98qy3LHACao6WURupA1dVKo6AZgAMG7ceJtrYrLS/gDQ8H5LJ8fUpL+kqh8QMivDxZs46B/Fmysh0uIku8CkFyk64Vhk3lyix59E9KJLoFevzApsUOf625mc+EXAoWGQsImA3VcmwWPvRvevzUG5c4G5qjrZv/8kXvBYJCJrqeoCEVkLWJyDsowBcvNpum3H91oXSVdJ5Gjyn6vgaioZoeL4XVyIIEuXEjnzVMKPPERyo42peeNdkltv0/IBM9RUoEwFMMeCRY+UyVDdt3JdqKouFJE5IrKBqs4CdgW+9L8OAa70fz6b67JNz5TvwNFYLJHMWcBojuK1ZJKuUvzMUxSdeiJSUUH0nPOInnUuRCKtHiMbIhD0I4akNpgep9XgISKf00IaElXdLMuyTwAeEpEw8B1wGN41mMdF5HDgJ2D/LI9tTEG5HRStZMF8Sk85gfB/nyM5djzVL/wPd3S2f5IZlpn6aUGjR8uk22rPfBSsqp8A45t4aNc2Hcjev6YnUiXywH2UnHsmEo1Sc/lVxE84OWeJDI1pTSbdVj+KyD7AesDnqvpy/qtljGmO8/13lB5/DKG33iD+ix2pvPVOnFGjrCVgOlQmkwRvA04B+gOXiMj5ea+VMWZ1ySRFt9xA3623IDBtCpU33srKF1/BXXe9QtfM9ECZtHF3BDZX1aSIlADvAJfkt1rGmHSBL7+g13FHEfr4I2K770HVjbfiDhla6GqZHiyTSYIxVU0CqGo1dpXBmI4Ti1F8xaWUbb8lgW+/ZdW9D7DqiWcscJiCy6TlsaGIfObfFrxJg5/5t7Udo61yw6YImox1rdQYgakfU3rsUQS/mEF0/wOp+td16MCBha6WMUBmwWOjvNfCmLypX85V/e9pMxRyX1pqAp8I6rp1t9ukupqSyy6m6OYbcNdci5WPP018j5YHPapXuD/lIr8B0j6vGchwtFX6fRHpj3cd5CdVnZqvihmTvYbrfze9dnJqlnZKdifc9JxTyQb5qPyjir/8qqSV0kIwCb79FqUnHEPg29nUHnYE1ZdeiZaVNbt/+uxub+1zb+Z5Si6DSWoZWVt/3EBmkwT/C5ytqjP8lCHT8FYQXFdEJqjqDfmupDEtazlYtPzMhrcyaZWkAoarEPeDRXOzyVP5rJL+/VQgSZ17UydhWbGCkvPPoejeu0iusy4rXniFxC93Wu149dlpmz+Bp9ckPZi0NZCIQEC8XFUWLExjmXRbjVTVGf7tw4BXVPVgEekNvAdY8DAdxzsbZh0sWj184yM2Ong86RJPeskNsym3YW4qEFGKJr1I6SnH4yxcQM2Jp1B93kVQUlL3HKeVYNGaurSKrbRKUoEpkEo9YgHDtCCT4BFPu70rcBeAqq4SETcvtTKmGep/66h+93jSJekq8aTmPOWILF1Cn7NPo/jJx4hvvAnLH3wMttq6QdLBXEv/DcT/SqVLt2Bh2iKT4DFHRE7Ay4Q7FpgEICLFQCiPdTOmSR15wTbhusQSOS5RlaKnHqfPWafirFzBqnMuoPKUMyAcpqiDT+LBgFjQMFnJZJ7H4cAmeKv+/VFVl/vbtwHuy1O9jOmWnHlz6XfgH+h3xMEkR4xk6VuTqTzrHxAOF7pqxrRJJqOtFgPHNLH9DeCN1H0RuVlVT8ht9YzpJlyX4gfupc8F5yDxOCsv+xdVxxwPgUCha2ZMVnKZgnP7HB7LmG4j8N1syk48lsi7bxHdYSdW3HQbyZHrFrpaxrRLJt1WxphsJBL0uvk6Bm43jtBn01l+0+0se26SBQ7TLVjyf2PyIPjFDMpOOJrwtCnU/mZPVlx7E+7gIYWuljE5k8vgYUM2TN6oKgq4aYPD6yfL5b6spKtEEy7V0QSuesNZA47gtFZYNErptVdRet1VuH37UXHvg9T+fr9WK5makOeV33B7PjgiBKzfwbRDJjPMg6qayOBYN+agPsYADWdxNze9oi5nlbY/kLh+sIjGk0QT9elF6h5PenM9AAIOBByHQKNhtaEpH1F2wtGEZn5J9QF/ZuUVV6P9BzRZnuCl+kjN4K6ve8NfoPHvnu3vlwpOqfQi3jb7vGeyl0nL4yO8+R0tjqhS1Yk5rJfpgVS9NCPNZPpo5bkNf7YWTFS9YBCNJ6mNe+lFhMzmkCRdSKaSHgLBmmr6XflPSu+4GXfwEJY99jTRX++x2vO81CRSP4O7pQo2I9NWSSo4OSIN8l9ZwDC5kknwSH+32YgqkzOprih1cz/xr6lgknSVWCJJTTxJLLF6sMimDkXvvMnAU48j9OP3rDj0SJZfcAlOWRmOal03V6DBp33IVQ9v41ZJwGkYnLzyLFiY/MgkeFgGZpMXyQ5MblMbT7K8Ot4gYLTnje2sWE7/i8+lz4MTiY1cl3nPvEztdjt4DyaVXkGHcMDpkBTpAEUh7wJGNq0ZY7LRlsWg0heCgs6yGJT9nZgM+PkUc/JJqOSl5xl45kkEliym4vhTqTjjH2hxcYN9HOn4tB/WyjAdyRaDMiZDgSWL6f+P0+n9zJNENx7Nwn8/QXTMuEJXy5iCaPNiUMb0OKqUPvkoA847A6eqkp/PvoDlJ5wGIcsLanquTIbqrmL1TM6pxaBVVfvkqW7GFFxg3lwGnnECvV59mdpxW7H4htuJb2CNcWMy6bZ6DVgT+A/wqKr+lN8qGdMJuC597r+b/pecD26SpZdezYrDj7FEhsb4Mum22kdEyoA/AHeJSBHwGF4gWZbvChrT0ULffsPAU46l+MP3qN5xF5ZcewuJtUcUulrGdCoZpSdR1RXAfSJyP/BH4GagCLguj3XLjA0kNhlq9a2SSND39pvod/WlaKSIxTfewaoDD7Khr8Y0IaPgISLbAX8CdgDeBX6vqu/ks2KmAzWxrGu+zpeqDaflNZwRnvtC65aO1YZlNx7WGp7xGYNOPobIZ59QucdeLL3qepJrrJV1uYmkS9DxFnrtiNjjqtblxjKmI2RywfwHYDnwKHAUkPC3jwVQ1Wl5rJ/Jk9aW485VTqX0E7ar3izvpKsN02z4x3ekcVqRtheaSnGiQE0swaraBNWxJMm0nCdBRwgGhFDAQWJRyq+/in43X0uybzkL73mIqj33afMv7AiEAg6hgEPQaTjHIx+JDutSnQQkb+udG9OSTFoeP+D9Lf4a+D8a/kUrsEvuq2VyronWRZue3uhk39J5PRUwFEgmlaRqg2y4TVTNz2nlPU/8E6MjzbcU0stKJUaMxl1WRRNURRPEEs0XmHCVhKsEPvyQ4WcdT/HsWazY/8/8/M8rccv7Z3wiDjlCKOgFjFaz7dbV1/99/N8703O+l6vKS8jodPA658Y0JZML5jtlciAR+ZWqvtLuGpn20wY/8nN4XW0LkGpZeFlqsy1fFZKqJP37qRZJKmVT6rhJV6mqTVAZTVATT7bamkpxqioZfM2lDJx4B7HBQ5k98SlW7rSb91htwmtBBL28VOnpPgIihAJewAi0cwZ5XYqUFlp4jtSngrdgYTqbXK7ncRVgwaMT6OgxBPGE4qpmlQ03E67fNEkCyaTLipo41bEkiSwK7P3O6ww/5yQic39i8SFHMf+MC3BLezcoK5pwifqLEJSVhCiNBFfrisqXUMAPWta6MJ2cLQZl2i3ZjlZGW7kKq6KJjFsZKYEVFQy59DwGPPEgteuMYtYTk6jacttWnyd41zI6irUyTFeRy+Bhg2ZNp1Q26XmGn38awWVLWXjsqSw48Sy0qKjQ1TKmSyvoGuYiEgCmAPNUdU8RGYk3qqscmAYcpKqxQtbRdF3BxYsYdtEZ9HvxWao3Hs3s+x6nZtMxha6WMd1CLtvjP2TxnJOAmWn3rwKuV9VRQAVweA7qZXoaVcqfeoSNf7UVZa9OYt4ZF/DVs29Y4DAmhzIKHiLSX0ROEJFb/a/jRaR/+j6q+oe2FCwiQ4HfAnf79wVv2O+T/i73A/u05ZjGhOf+xLqH7MuI046hdr0NmPniuyw6zjLgGpNrrQYPEdkImAGMA74GvgG2BD4XkQ3bUfYNwJlAakB+f2C5qvrjXJgLDGmmTkeJyBQRmbJk6ZJ2VMF0G67LwPsnsNH/bUPplA+Zc/HVfP3EJKLrrV/omhnTLWVyzeMS4CRVfTx9o4jsC1wG7NvWQkVkT2Cxqk4VkZ1Sm5vYtcmL8Ko6AZgAMG7ceLtQT9rEvGZejVyP4PHWA3eJJpLEk4oIBEUIBnI/iU1VWRVNsKw6xpLKKDXxJJGgQ0kwQFEoQOn337L2WcdTOuVDVuy4K3Muv4HY0OFZlycCxaEAJeEAJeEgqprXEVCCP58jkCojteKBMZ1XJsFjtKru13ijqj4lIpdnWe72wF4isgdegsU+eC2RviIS9FsfQ4H5WR6/22stWDS3f0pbT4aqSizpEkt4X26jZV1VIa5K3PWm9gUEggGHgCMZz75OV5tIUuEHixW1cQTB1fohwbUJl1hNlPUeuJ0N77gOt6iYr6+6hZX7/wXHafulvHDQoTgUoFcCY76+AAAgAElEQVQkQCjg1J2+8xU0Ag44jvj5rzw2RNd0JZkEj6osH2uWqp4DnAPgtzxOV9W/iMgTwH54I64OAZ7N5vjdVVsDRibHgqZPWqpeDqqoHywSrq62BnhL1UgqJP0UIalP1sEWZksnXWV5TYyfq2Msq4oRc10c/zhNlVY283PGXnQqfb+awbxf7cmnZ19GfOAg3FW1hB2hOBSgOBQk3Ex5AX+fXpEAkVAA/BQnqX1zfRp3pD5YpFfHAobpqjIJHoNE5NQmtgswMMf1OQt4VEQuBaYD9+T4+F1KLoNFJuW46nVFxZJK3D/xZxosWjw+9fmkoD7tRizhUlETY0lVjKpoAkfSgwV16UnSOdFaNrzzOkZNvI1Y33ImX3s383f7bYMKxlwlHk2w0p8mXhR0KAkF6V8SplckSEk4gON/4q9rFeXhHO6lFvF+pliwMN1FJsHjLqB3M4/d3d4KqOqbwJv+7e+Ardp7zO4k34EjXVU0SbSFhIK54iq4SWX6vOUI9SMmkq38rv2nTWaLi0+j9w/f8uPeB/L56RcS79O3yX3TD1WTcOkVgTX6FHVY2o9I0MFxLFiY7iuTxIgXd0RFsmZ/mznT+LpI3ssjs9ZMsKqSjW+8nHUfu4+qwcN47/ZHWLzdTm0qSxAUxZEOSjViualMN5fJeh4XtPCwquolOayPMQ0Meu8NtrjkDIoXzmf2nw/nyxPOIVnSq9DVMqbHy/aCeS+82d/98YbyGpNToRUVbHb1hQx//glWjVyPtyc+y7IxWxa6WsYYXybdVtembotIb7yUIofhjYi6trnnGZOtwa/8l80vP4fwyuV8deRJzDryZNyIJTI0pjPJdA3zcuBU4C94aUPGqmpFPitmep7IkkVsfsW5DHntRSo2Gs37tz/Cig03LXS1jDFNyOSax9XAH/BmdI9W1cq818r0LKoMf/YxRl9zEYFYlBknn8fsg45GgwVN+myMaUEmf52nAVHgPOAfaSNIBO+CeZ881S0zlpykSyuZ+xNbXHIGgz58m6Vjt2b6hddSOWLdvJRlo5+MyZ1Mrnl03DJqnZU2HaOk7lsOi0obLqvqTdxLLyLXJ0DXzzOSTHoTBGsTLuGA4w01zXF5qt5EQUeEZZW1jHjwLkbffAXqOEw/9wp+2P9gyCK1SHOE1Lrn3qzuWMIlEnL8x/Kbq8p1XRxxLGCZbiuTbqvylh5X1WW5q07nkcmUB6375snqPJGWr0nxTubaQrAS6oOJZBG91F9rXICqaILl1XFW1SaIN5qhFwwIRaEAxWGnLteTU19oxhKuNwUwmnD5saKG7ytqiM2YwT7X/4O1v5rOV+N35JnjLkaHDaNfdZx+RSGCgUazv9sgIN4kxF7hAANKI5SXhOkVDqx2Ek9fODcXCUkchxbTrxjT3WTSbTWV5tN8KrBOTmtUKM2csNt0iEYHaO4ckt66cP3WRabBKn0/URBpOZjUBQuBWNytCxbVsaaSf9RLJJXKZILKWu9+JORQFHQoCgfq0m00dZJ0/daFAAtXRZm9rJo5y2uojCVxEnF++eTd7PrIrUSLe/Hoaf/ik5338ioXTbIimuQHooQcoawoyICSIKVh7y2anuIjXWpz0HHoXxKivFeEvsWhZvdviqZ9h8yCiUh9vq6OmrVuTGeSSbfVyOYeE5Em19voCjpiMnV9Gat3Rbm5SG6YVoZQ3zJJL2dlTYKVNV7AaE+Z0bhLNO6yosbLQeW1SgKEQ47ftSasqI3z3bJqflxew+LKWINgPOSbGex347ms9f0sPt3hNzx3zPlU9e3fZFlxV1laHWdpdRyAXiGHsqIg/UtCFAWduky0ZUUhBpRG6FccoigUyP6XayQ9mKR3bwUcL0hZ2hFj2r+G+QdA9gsn9BCqXtbYfMar9EBSURXj51WxvOWpchWqY8m61svUhStZXBUl1kRyqmC0lt0eupkdnr6Pqr7lPHDerXy57W5tKq8q7lIVjzF/VYzy4hD7jl6DsqJQh5zAAw6Egk5WXWjGdGftDR72F5WhjhwU5qp3cbijLK+JNxk4Rs74mD/ceB4D5//AR7/enxf/dia1pe0bnKdASSjYYZ/8RbJbj8SY7q69wcMGyprVRKor2f2+a9j2xUf4eY2h3HXZRL4ds22hq2WMyaFMRlvdTPODf5rOh216rA0+fovf33ohfZYu5J19DuV/B51EvKik0NUyxuRYJi2PKVk+ZnqQkhXL2POuKxj7xnMsGr4et1/zKHM2HFPoahlj8iST0Vb3d0RFTBelykZvvcCvb7mY4sqVvPqn43jjj8eQDIULXTNjTB5l0m01ADgOqADuBa4GdgC+BU5T1dl5rWE3oKp1czs64kJvPOkSTbiEg/lNDhBatIB1LjqdbV99kTmjNuXuy+5j4cgN81pmPmb1G2PaLpNuq4fxuqdGAR8B9wE34gWQu4Gd8lW5rkxVSfrrdruaPgFNcRzJaeoPVWX+yigzF63ik/krmb+iFgXKIkFG9C1mWN9iBvUK52zUUCyexLnvbsbcehlOLMr1vz6aR7fZl7J4CQOXVDKgTxHFkdwlNSwNBxhZXsKGA3sxtKy4bp5HvqTWWA84jj8Jsbk5ssb0XNLa0qMi8qmqbi7eme5HVR2e9tgnqlrQju1x48bre5Pbf+mlvZMGvZncXsBoy5wOgayCycraBLOWVPLZ/JV8vaSqbqZ6otFMQEcgIIKrsFbvMCP7lTC0TxG923ByV1XmLa/hywUrWTDtCw7/92Vs+d10pozYnEv2Oo05/evnigYDgusqoaDDwL7FDCgrorx3hGAg81ZQ0BGG9S1mvf4lrNe/F8UhB5RmjqF1p/Zs/wmFVLCQBjPTbSKg6c6KQzJVVcdn+/xMziBJ8NLnisjSRo913GSCPGt8nmgtmKifkyqZdEm2Y8a44k0gBL9VIk2nu0gkXb5bVs2XCyv5bMFKltcmCApEm5hfkS4VVADmroyyqDKGq0ok6DC8rJi1+xYzuHeEUKMT86raODMXrOLTucv5auFKJJlkv/ef4sxX7yXpBLj0d6fw9Ng90EaJDBN+faJxl3lLqli0rJqEq/QuDjHIDya9S1af4DewlxfYNhxUysBeYZKqhDLKE1U/B1z8VzSTf4qA+MEikPZ8CxbGZCyT4LGOiDyH97eZuo1/v9nUJV1dg/NIWt6rpOuSSCrJPOQ3UfCO62dcTCRdJv+0nE/nr+Kn5TUEHSGacOvrkkUZcT9QVcddZi2t4rsK7+ReXhyibzjA0opavlywgorqOMGAUBt3WWfxD1z47DWMnjuTd9bfmsv2PIXFZQMz+n1SCRdXVsepqo3zw8JVKDCgT4RdNlmTHdbtz9r9ihGRuhM6QCDrbqL0hCINA0ko4B0/vfvOAoYx2ckkeOyddvuaRo81vt89SX2Pd75SfjRlzvJa/vvl4roTfuMuqfZSqJsZvrQ6zrTvfqZiZbR+h1iMI999hCPefojKSC/O2e8fvLzpzlmmD4akmyoVFi2v5ZAth7WpO6vt6gNJwIFQwFKkG5MrmQzVfSuTA4nIU6q6b/ur1Ebd+FygeH3x8RwHjeakv5Qbz/uKC5+5hlGLv+fF0btyzW+OZXmvHM8J7dB/u278RjGmAHK5zmf3SM3ewxXFajnmjYn85YOnWNq7nJP/fAlvb7BdoatljOlkchk8LM9VF7f57Gmc/MiVDFs2nyfH78lNvzqSyqLSQlfLGNMJ5TJ4mC6quHoVBzx8PTu99iRzygdz1KHXMGXkFoWuljGmE8tl8LBO5S5o82lvcfA9l9C3YimP7/wnbtjur9SGiwpdLWNMJ5dJepLhqvpTBsc6Kwf1MR2k98pl/Pn+q9jm/ZeYO2w9bjnlet4uXZva9NFWxhjTjExaHs8AY6HlEVWq+r9cVixjHXilpbXZ+Dkvz3WJJ3M8NFiVbd5/kT9PvIrimkqe3u9YXtj7cJLBEO6ClbktqwVCaiJmR6X+sBQjxuRSJsEj/S+ux42oSgUMBRL+HA/X1bqpDrmcN6CqVNfGCAUD/DB3KS+8+SnffLOYSP/+DFp7COIIwUAAnOzmRvT7eSEH33spY6a9zbfrjebeoy7ip7W8f9Jli1eyZEEFK6ti9O5TQigS8n7HHP5+RSEH14XexUG2Wac/3y2pZIM127eyYCZSM/YtfBiTO5kED23mdreU3rpIuko86c30Xm2qhX9fREk/x7Y1mERjCVzXpaY2xpuTZ/LSW5/z3rRvWFlZs9q+ZQPLWWvkUNbeaD36DCgHdQkEg62e4MV1+eXrT3HAw9fhJJM8/NfTeGb7fflxwUrmv/01FRVVDdKxLF20AhGhuFeEPn1L6NW7GCcQIOCAtuH0G3SEcNAh6Sqjh5ax4/oDGb92P9YoW/2aSq4bdQEHgo6D49gscmPyIZPEiEmgCu9DWzFQnXoIL+VV/j86tiAXiRFTr4GrXjrzRFLr8k21lQikcutJakMa13Wpro0RDAT47Ks5PPfaNN76eBbfz1nSpnICwSADh63F0FFrM2S9tQkXF3lJFoMNPw8MnP8jh959MRvPnMK0dbfgjO2PZGqsF/F425KbBEMBSkqL6NuvF5GSiJ/Q0Vnt00RJOEAs4TKkXzE7jBrAVuuUs+GafRokHGxNNoHEES8po+PYmuPGZKK9iRFbDR6dXbbBI5XYMOG3LBLJzDPhZiq9RaIocxcs46W3P+O1979kyuffE09kk52qacW9e7HmiKGM2Ggd+g8dTEBddnx6In/+7wRiEuDs0X/k/hE75Kwbqqg4TK/eRZT1LSVcHKJPUZCtRpaz3XoDGDu8H6VFuRvI19RbVKQ+E67TRCJJY0zLOiKrbs6JyDDgAWBNvMy8E1T1RhEpBx4DRgA/AAeoakU+6lAVTWbdusiU4p34VlXVss2+F1JVE8tbWTWrqvj+81l8//ksNg243PHdJMZVfM9/1xrDyVscxILifjktr7YmRm1NjJ8Xr+T1K/Zi1OCyvJ3A0w8rQDgoOFle9zHG5EahJgkm8FYhnCYivYGpIvIKcCjwmqpeKSJnA2fTTYYAx9rYTZSNsJvkzEVTOHPxNJaHe3HQ1sfwnyFb5vSid1M6OuGgtTKMKbyCBA9VXQAs8G+vEpGZwBC8DL47+bvdD7xJNwke+bZV1UJun/M6G9dW8PCgzThn6yNYGrbUIsaY/Ch4ehIRGQFsAUwG1vADC6q6QEQGNfOco4CjAIYNH97ULj1GSTLOhQsnc/yST5kfKmWfdX7Ly0O3JBgp7QFj44wxhVLQ4CEipcBTwMmqujLT7ghVnQBMAO+Cef5q2LnttGoOt815k5GxldzZf1POH7wtqwLhQlfLGNMDFCx4iEgIL3A8pKr/8TcvEpG1/FbHWsDiQtWvMytLRLl8/vv8bdmXfBMpY7f1fs97pYMLXS1jTA9SkCEr4jUx7gFmqup1aQ89Bxzi3z4EeDYf5aeG6XaUZauiaHFvcALtPtaeK75j+lcPc/CymVwzaCxbbXDgaoGjT78+9O1bku/r5AAUhwMUhQve+2mM6WCF+qvfHjgI+FxEPvG3nQtcCTwuIocDPwH756IwVcXV1Izx+tni6SfXXE53icWTTPtuKW9+Pp/XPpvPgmXVhIatj6PgxKMkViwjWVmB1lRmfMxB8WqunfcO+y2fzadF/dl3nd8yvcS7JBQuLmLQ8CEMWX8E/YcORrxp1YgIyyuqmDtvOUuWrKSmJp6T32/jYf3YefPB7LX1SDYe1g/Em00O3kTLXAdnx5/T4Thi6UWM6SQKNdrqXZpPM7Rrjsog4ao3CbCZ+Rzp11hE0vJYtfHMp6p8v2gV7365kJenz+WzH34mHHSoiSZI1h3L8coIFxMaOITggDVRBamtIrZ8KW7Vcog3MQ9ElQMrvuaaee9Q6sa5cM2tuX6tcZQNGcym6wxnrXXXpqi0BFFFgqv/c/Yf0Jvy8l64OoREPMnixSuZv2AFP/9cSTLDpIsD+hTxy00H89uthrPDJoMJ+GlHQsHVW1IBUVS911WpD9xtIQIBSc0WT99uocOYzqLb9DekTlLxpEvC1axaEqmTk/gpX1OHaOpYK6pifDBrEa9/Oo+3v1xAbSyJKtT68zliieZPzC6ABEBAS/oQKemN6yq4CbRyOYmVy3CrVjIsuoKb5rzJ7qt+4qM+Q/jnzn8jvs12/HrggIzzWgGI4xAAAgGHYcP7M2SI11qoqowyd14FS5asYsWK+lxaoYDDVhsM4ldjhrLH+LVZs7yEpKuURDLLo1mX58v/5mj96+jS9L9NaqZ4ehoTCxbGdF5dPngoUN3gE36OSH0XSXqr5IE3vuaRt77lpyWVREIOlbWJdhfl4ifEcsJI30EU9SnnkNlv8M+vHyMgcPf+J/PKbn+kNBRKy6ib/fUTJ+Ado3efYjbqXcT6o9ZAgSHFQX6/2WBGjygnnnApjgTblJOqaWnBREDSUtuK0KAryoKFMV1Hlw8eQO4DRxNSJ7ZbX/iSikpvwaScr7UBrLtyATdPvoftlszi3aGb8cjxl/DzmsPy19cvQsDvftpr3HDGj/DSmERC7b+433Rx9b9JwDLeGtNldYvg0R0E3CTHffUS53z+H6JOiOO2PpIXx/6aLQcO6bAhcXYeN8ZkyoJHJ7BpxY/c8uHdjKn4geeHjuf0LQ9hUXFf+tnZ3BjTSVnwKKBIMsbpM57l5C9fYFmklIN/cSLPDd+y0NUyxphWWfAokK2WfM3Nk+9hg5XzeXjkLzh37F9YHrFEhsaYrsGCRwfrFa/l/E+f4KivX2FeSTl/2OkMXh+8WaGrZYwxbWLBowPttOBzbvzoXtauWsqE9X/FPzffn8pQcaGrZYwxbdblg0dHrqLrqlIcDlATDlAby3xxp7JYFZdNe5i/fvc2X/dZi913O48PB23Q6vNCAS/FSNCRZmfJ50pAIJr0Vlds/9yOtkib+GGM6TK6fPAQgaKgQ8J1Sbq5XcLCm7WuiAhV0QTLqmJccug2TP1mCVO+XsTnPyzDdb194s3MKP/dnI+5+uP7GRBdxbUb/45/jd6HaDNp0yOhAAIUhQNsu+Ea7LjpYEasVcb8yihfLqpk4cpaggGHWMLNye/pvW7KoN4Rxg0tY1hZSYcFjvpiLHAY0xWJduRH9zwYO268vvvBx3X3VZWk6+WzymbyoOsqiJdeZFl1jBU1cSprE02erFWVOUsq+fS7pXz81SJ+XLyKcDBATSzBwOrlXD3lAfae8zGf9Vub47c+gs/KRzR4viNQHA4SS7psOrycXTYbwjYbDmLYgNImJ8/Fki4/Lqth1pJVzFpSRXUsiSMQy/AXDQUEQQgHhc3W6sMWQ/uy8ZqllDSTFTeV4DAXbxGhfnK8TQw0pvCKQzJVVcdn+/xuFzzSpWfTTTaToC+1j6qysjZORXWcFTXxrLqJovEkX/24jN6PPshBz9xMUTzGtWP+wPXr707C8U7QJZEg8YTLGv2K+eWmg9l+ozUZPaKccBNJBluzoibOtz9XM3PRKr6vqEYQkq5bFzQdIBJyiCeVdfqXsOXwfmy6Vm/W6B1p8wk8PWlkJq+MAOL4Py1YGNPptDd4dPluq5aICIG0ZHuqXmsk6brUxF2iiSTLqrxgURPP/BpGc/oumc8R15zBoPfeYNkWWzHpjCtY7vZj9KxFzFlSyWYjytlp9BC2HDWQ8t5F7S6vrDjE2KFljB1ahqvKgpW1fL2killLKkm6yhZD+zJmSB/WHdCLUKB989QbJI2kPpikx1iR+u4oCxjGdG/duuXRkulzKlhe0/6khgC4Lms/eh8b3nApiPDVyefx4x8PreunEYF1ynt12Am1tChIv5JQh57ALVYY07VYy6PAen33DZtddCrl0z9i8fY7M+OCq6kZPKzQ1TLGmLyy4JElicdZZ+KtjLr9WpLFJXxy2c3M+93+9hHcGNMjWPDIQp+Zn7PZBSdT9tUMFvzqd8w493JiAwYVulrGGNNhLHi0gROtZdTt17DOxNuI9evP1OvvZeFuvy10tYwxpsNZ8MhQv2mT2ezCUyj94Vvm7HMgX55+MYmyvoWuljHGFESPCx6peR2lkSCqwqpoAreFEWeBqko2vOEyRjx6L9VDhjH5zsdZut0vMy5PgD5FQUoiARJJJZ7nZQ8dgXCwo5aPMsb0VD0ieLiqJJJKLOEST3qzpvuXROhXTN1s8uU1MVbVJqhOm+8x8N3X2fSfZ1C8cB7f//VIZp1wDsmSXq2WVxxy6FscZkBp2AtSKI5IXRanWEKJxpNEE26TExfbKhJyKA4FKIkECDrS4dmibIyAMT1Ptwweql56knjCJZZsemY5InUT2opCAdYIFjGw1Dvp1i5cRN/zz2LQ04+yap31ef+B51k+pvlFmgKO0Lc4RHmvEP1KwjgifjqO1FlV0r5DJOSlCOmNN8nOCyRefTOJJcGA1AWLSNBB1TuB103ka/JFyV3eLwsWxphuEzySrtclFEu4WaUWERECKMXPPc2QM0/CqVjGitPOZtkpZ1PqOiRr41TFEnWBqHfEm4jXvzRMUSiAqtdllOnEvNR+AYHicICisHfST7hKbdwlFq//PRzxAlxx2PsSaZj2I6MipWFQacvcUKn7Zowxni4fPFShojre7uR9zsIF9DvzJIr/+yyxMWNZ+tQLxDfdjBBQrkrfkjAiEAxAJBgAvyuqLli04+QqfksFvOSFQUfQiJdhNxRwCKR1ReVq1njjwzR+/ax1YYxpSZcPHqDtCxyqlDz8AH3POxOJ1rL8osuoPPZkCNa/NKkcWQBFIQdHhHx+FE8FE/HzcqUHl/yVmecCjDHdSjcIHtkL/Pg9/U4+lqK3Xie63Q5U3HAbifXWb/E5do41xpieGjySSUrvuo0+l14AToCKa2+m6pAj6hecMMYY06IeFzyCX82k30nHEPn4Q2p2+zXLr7uV5FBLZGiMMW3Rc4JHPE7vG6+hzzWX45b25uc7J1Kz34HW2W+MMVnoEcEj9Mk0+p1wFOEvPqf69/uz/MrrcAdaIkNjjMlW9w4eNTWUXXUJpbdcjztoDZY++AS1e+xV6FoZY0yX1+WDR3OjdMPvv0O/k44h9O1sKg/+GysuvgJtZyJDwZsRHuionq6uvcijMaYb6/LBwxGhrCRIPOml90hUrKDsn/+g9N4JJEaMZMkzk4juuHM7y4BwUPwJe4AIruvNL8nH+V0EAiJp6U2MMaZz6fLBA7wAEgkKJa9MInzcMciC+dQcfzIrzrmARHHriQwbE7z8UeGgQyjQ9Enc8Zsfql4QcbWZHFoZCgg4jtOmFCfGGFMonS54iMjuwI1AALhbVa9s9UlLlxI5/RSCjzyEu9HG1D76BO5WW9MbP0liUoknvSSJzc1GDzj4wcJrXbQlR5UIONQHE9cPJi3NfHf82eMNUpwYY0wX0amCh4gEgFuBXwFzgY9F5DlV/bLZ51Qso2TMJlBRQey8C4ifeQ5EIunHJBQUQkGHEsBNJVBMugQEQn7rInc5o7xUJgHEa5X4Zap6WXYDjuQ0R5UxxhRCpwoewFbAbFX9DkBEHgX2BpoPHt9/jztuPNFJr6Kbjm61AMcRIo4QCTl5n+KRyknlBCxgGGO6l84WPIYAc9LuzwW2bryTiBwFHOXfjQanTpnBuM07oHqd3gBgaaEr0UnYa1HPXot69lrU26A9T+5swaOpj+arXTlQ1QnABAARmaKq4/Ndsa7AXot69lrUs9einr0W9URkSnue39kyAc4F0hNNDQXmF6guxhhjmtHZgsfHwCgRGSkiYeBA4LkC18kYY0wjnarbSlUTInI88DLeUN17VfWLVp42If816zLstahnr0U9ey3q2WtRr12vhWh71281xhjT43S2bitjjDFdgAUPY4wxbdZlg4eI7C4is0RktoicXej6dCQRGSYib4jITBH5QkRO8reXi8grIvKN/7NfoevaUUQkICLTReS//v2RIjLZfy0e8wdgdHsi0ldEnhSRr/z3x7Y99X0hIqf4fx8zROQRESnqKe8LEblXRBaLyIy0bU2+D8Rzk38u/UxExmZSRpcMHmlpTH4DbAz8SUQ2LmytOlQCOE1VNwK2AY7zf/+zgddUdRTwmn+/pzgJmJl2/yrgev+1qAAOL0itOt6NwCRV3RDYHO816XHvCxEZApwIjFfVTfEG4BxIz3lfTAR2b7StuffBb4BR/tdRwO2ZFNAlgwdpaUxUNQak0pj0CKq6QFWn+bdX4Z0ghuC9Bvf7u90P7FOYGnYsERkK/Ba4278vwC7Ak/4uPeK1EJE+wI7APQCqGlPV5fTQ9wXeaNJiEQkCJcACesj7QlXfBpY12tzc+2Bv4AH1fAj0FZG1WiujqwaPptKYDClQXQpKREYAWwCTgTVUdQF4AQboKWvt3gCcCbj+/f7AclVN+Pd7yvtjHWAJcJ/fhXe3iPSiB74vVHUecA3wE17QWAFMpWe+L1Kaex9kdT7tqsEjozQm3Z2IlAJPASer6spC16cQRGRPYLGqTk3f3MSuPeH9EQTGArer6hZAFT2gi6opfn/+3sBIYDDQC697prGe8L5oTVZ/L101ePT4NCYiEsILHA+p6n/8zYtSzU3/5+JC1a8DbQ/sJSI/4HVf7oLXEunrd1dAz3l/zAXmqupk//6TeMGkJ74vdgO+V9UlqhoH/gNsR898X6Q09z7I6nzaVYNHj05j4vfp3wPMVNXr0h56DjjEv30I8GxH162jqeo5qjpUVUfgvQ9eV9W/AG8A+/m79ZTXYiEwR0RS2VJ3xVvOoMe9L/C6q7YRkRL/7yX1WvS490Wa5t4HzwEH+6OutgFWpLq3WtJlZ5iLyB54nzBTaUwuK3CVOoyI/AJ4B/ic+n7+c/GuezwODMf749lfVRtfNOu2RGQn4HRV3VNE1sFriZQD04G/qmq0kPXrCCIyBm/gQBj4DjgM70Nij3tfiMjFwB/xRidOB47A68vv9u8LEXkE2AkvBf0i4ELgGZp4H/jB9Ra7evMAAALhSURBVBa80VnVwGGq2mrG3S4bPIwxxhROV+22MsYYU0AWPIwxxrSZBQ9jjDFtZsHDGGNMm1nwMMYY02YWPIwxxrSZBQ/T7YlIUkQ+Sfsa0cx+JSLykIh87qfxftdPAZN+jBki8oSIlPjbK/2fI0Skxt/nSxF5wM8CgIjsJCIrGtVht2bq0GS6ff+xi0RkXtox9kh77Bw/pfYsEfl1rl47Y5rTqdYwNyZPalR1TAb7nQQsUtXRAP5M7XjjY4jIQ8AxwHWNnv+tqo7xlwx4BTgAeMh/7B1V3TODOqTS7U8Tkd7AVBF5RVW/9B+/XlWvSX+Cn47/QGATvDxOr4rI+qqazKA8Y7JiLQ9j6q0FzEvdUdVZzcw+fgdYr7mD+Cftj8giY2sL6fZbsjfwqKpGVfV7YDbesgXG5I0FD9MTFKd19Tzdwn73AmeJyAcicqmIjGq8g59U7zd4qWGaJCJFwNbApLTNOzTqtlq3tUo3Srefcry/2tu9Ur8ioC1RYDqcBQ/TE9So6hj/6/fN7aSqn+CtiXE1Xu6jj0VkI//hYhH5BJiClxfoniYOsa6/z8/AT6r6Wdpj76TVYYyqfttShZtJt387sC4wBm+NimtTuzf167R0fGPay655GJNGVSvx0nf/R0RcYA+8rqNMrpukrnmsBbwpInupapuzPTeTbh9VXZS2z13Af/27PX6JAtPxrOVhjE9Etk91Bfmp/jcGfmzrcfx01mcD52RRh+bS7afWYEj5PTDDv/0ccKCIRERkJN5a1B+1tWxj2sKChzH11gXeEpHP8dJ1T8FrAWTjGaBERHbw7ze+5rFfM8/bHjgI2KWJIbn/8ocRfwbsDJwCoKpf4KXa/hLvOstxNtLK5JulZDfGGNNm1vIwxhjTZnbB3PQ4/gzsqxpt/r6lkVh5qEN/4LUmHtpVVX/uqHoYky3rtjLGGNNm1m1ljDGmzSx4GGOMaTMLHsYYY9rMgocxxpg2+38SVcJCxgsKFwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# XID+ flux density vs. MF flux densities\n",
    "plt.hexbin(cat_all['F_SPIRE_250'],cat_all['F_BLIND_MF_SPIRE_250'], cmap=plt.cm.Blues,gridsize=(600,600))\n",
    "plt.plot([0,100],[0,100], color = 'red')\n",
    "plt.xlim(0,100)\n",
    "plt.ylim(0,100)\n",
    "plt.xlabel('F_SPIRE_250')\n",
    "plt.ylabel('F_BLIND_MF_SPIRE_250')\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Add field name\n",
    "cat_all.add_column(Column(['GAMA-12']*len(cat_all),name='field'))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING: UnitsWarning: 'mJy/Beam' did not parse as fits unit: At col 4, Unit 'Beam' not supported by the FITS standard. Did you mean beam? [astropy.units.core]\n"
     ]
    }
   ],
   "source": [
    "cat_all.write('./data/dmu22_XID+SPIRE_GAMA-12_BLIND_Matched_MF.fits', format='fits',overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "*This is a default HELP jupyter notebook *\n",
    "\n",
    " ![HELP LOGO](https://avatars1.githubusercontent.com/u/7880370?s=75&v=4)\n",
    "\n",
    "**Authors**: S. Duivenvoorden\n",
    "\n",
    " \n",
    "For a full description of the database and how it is organised in to `dmu_products` please the top level [readme](../readme.md).\n",
    " \n",
    "The Herschel Extragalactic Legacy Project, ([HELP](http://herschel.sussex.ac.uk/)), is a [European Commission Research Executive Agency](https://ec.europa.eu/info/departments/research-executive-agency_en)\n",
    "funded project under the SP1-Cooperation, Collaborative project, Small or medium-scale focused research project, FP7-SPACE-2013-1 scheme, Grant Agreement\n",
    "Number 607254.\n",
    "\n",
    "[Acknowledgements](http://herschel.sussex.ac.uk/acknowledgements)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}