{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# A quickstart of the GraphEM solver"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Expected time to run through: 3 mins**\n",
    "\n",
    "This tutorial demonstrates a quickstart of the GraphEM solver."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as  plt\n",
    "import GraphEM"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.set_printoptions(precision = 4) # set print forma"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Original Temp:\n",
      " [[0.3745 0.0206 0.6119 0.6075 0.122 ]\n",
      " [0.9507 0.9699 0.1395 0.1705 0.4952]\n",
      " [0.732  0.8324 0.2921 0.0651 0.0344]\n",
      " [0.5987 0.2123 0.3664 0.9489 0.9093]\n",
      " [0.156  0.1818 0.4561 0.9656 0.2588]\n",
      " [0.156  0.1834 0.7852 0.8084 0.6625]\n",
      " [0.0581 0.3042 0.1997 0.3046 0.3117]\n",
      " [0.8662 0.5248 0.5142 0.0977 0.5201]\n",
      " [0.6011 0.4319 0.5924 0.6842 0.5467]\n",
      " [0.7081 0.2912 0.0465 0.4402 0.1849]]\n"
     ]
    }
   ],
   "source": [
    "#Generate data\n",
    "np.random.seed(42)\n",
    "# Temperature\n",
    "temp = np.random.rand(5,10).T   # 10 years, 5 grid points\n",
    "target = np.copy(temp)\n",
    "print(\"Original Temp:\\n\",temp)\n",
    "\n",
    "temp[4:10,:] = np.nan  # Insert missing values after instrumental period"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAN8AAAGdCAYAAACSHqb/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAxA0lEQVR4nO3deXwb5Z0/8M/MSBrJliX5jC/Zju0E5+BKKAUKgSWvDZAssITSZpNybdkWdtOFbhfYtLT5NT1ZyrZdaDgXKEcphXAVaAOFLWkSoCUh5QpJiO3E8hVfki3rGGlmfn94bTKRLI2wrWee5Pt+vfLyK/NIz3z1lb965HlmnhF0XddBCMk7kXUAhByrqPgIYYSKjxBGqPgIYYSKjxBGqPgIYYSKjxBGqPgIYSSvxReJRLBz505EIpF87pYQS7Llc2cfffQRFi9ejKbV/wbXrNop9eVzObC0pRavfhRAMKpMqa/Kn22f0vNnQsdTC6elH5/kxDLfbLwcbENQjU2pL10XpiWm6TT7m6PT0k+Rz4XTl83BGy/vw0gwOqW+fr/nVlOP4/Zrp64DUSUJOjkuMx06ImoSOihRmei6jlg0gXyebZnXkW86hWIKXvrgIOswLC+kxvFCcB/rMCwvHIphy29353Wf3I58hQ4bPttQgUIHt58feVEo2nG6uwaFop11KJbmKnTghNPr4Cp05G2f3BZfeZELN5x7IsqLXKxDsbQKeyG+UX0aKuyFrEOxtJIKN67497NRUuHO2z65LT5CeEfFRwgjVHyEMMJt8SlJDW39w1CSGutQLE3RVbTGhqDoKutQLC2hJBFoHUBCSeZtn9weKuwKjeKbz7/FOgzL61RGcNPB11iHYXmHOofxsxtfyus+uR35COEdt8XXUFKEh69cioaSItahWNps2YvHm/8es2Uv61AsrXp2MX7869Wonl2ct31yW3wQALskAtY73dBiBNhFCZSozAQIsNklCHnMU87Fp2ka1q9fj5qaGhQWFmL58uVoa2ubidgIOarlXHzf+973sHHjRtx7773Yvn07VFXFeeedB0WZ2pUFhBxrcio+RVFw++23Y8OGDVixYgVOPPFEPPHEEwgEAti0adNMxUjIUSmn4tu1axdGRkawdOnSiW0+nw+LFi3Cli1bpj24TDqDo7jx6e3oDE7P9VxHq05lGDe0v4xOZZh1KJbW2xnCbTc8j97OUN72mdM8XyAQAAD4/X7D9urqanR0dKQ8vrGx0fD/eDw+9nhfAbylnxylHI0n0BeOwS6JqPGlngDcPjACAKjyFEC2S4Y2hyQioWooctpRWug0tMUSKnqGIxAEoD7NUdGOwTBUXUdpbTFcRbKhLdg7gvDgKFxFTpTW+gxtSiyJQ239AICalkoIR/yN3tvaj0Q8ieIqLwp9xhO/RwZGETo0ArnAgfL6EkObmtTQve8QAKDO4YEkGD8bu5UwYnoSJTYnvJLxtY6oCvqTETgEETUOzxGvVIeij52MUOMogkMw5vBQYhSjWgJeSUaJzRhvREugNzEKCQL8sjEPAHAgHoQGoMpeCOcRV04MJCMYVhW4RTvKjzixO64l0ZUIAwBmp+k3oAwjoWuosBegUDReaTCUjCGoxuASbaiZbcxhMqmit2OsgKrqfRBFYw77uoehxJLwlhTA7TXmcHQkhqSiwu6QUFFjPDqs6zq62ocAALNqvbAd8Xs40DuCWCSBIp+xz0xyKr7x5R9k2fiL6nQ6MTg4aLqftWcfj/nz50/8f+vH3fjFlvdRUiDjRxeflvL4f3jgFQDAtUsWYG6Fz9D2yzc/wu8/7MBps2fhH0+fZ2j7a6AfP375HThtUtp+v/KrP2IklsDf33QeFp5znKHt2ds24/VH3sTc0xpx1e2XGdoCu7tx+xfvBQDc8OiXYTvisqZbL9mInv19WPbVJTht5SJD2x/u34oX//tV1M6vwtoHrjK0BXuH8d2//SkA4Fs1n0OpvcDQvr7jdXwQ7cf5viasLGkx9htqw929O1FhL8Rt9UsNbZqu45/bfo/+ZATXV34GjU7j4fTbu97EG+FOnFVUh6sqTjC0/SXchVu73kChZMdPGs7Fkb6073lEtSSumXUSTiqcZWi7r3cXfh9sxSJ3Ja6v+oyhbU90AN88+DoApO33X1o3oycxilVl83G2p87Q9kT/bvxmYDeOc5bg6z8509DW3z2MH699DgDw1fV/m1Jgd6z7HQ7s7ceSi+bh7AvnG9p6O4K4/wevocAt4+s/WWFoi0UU3HL5EwCAy/99CSr9PkP7Az/6X3z4dgCf+ZvmlNcyGSGXG6Vs2rQJn//85xGJROByffIJ+YUvfAHxeBzPPfdcxufv3LkTixcvxln/fAu81fUT2z/NyFftLcDXzjkBG176C3b3BKc08i14arflRj7h+TOmZeSrcRThhqpTceOBV9EWD05t5HMY8wCwH/kW32N8LZ925Kuo8WDN18/CT//9RRzqDE1p5Hvyre+mvJ50chr5xr9udnV1oampaWJ7V1cXTjjhhMmelqIrGMGQPJKyPaFqE4WWTvdw6sJL0cTYOYsjsQRGYom0z9N1ZOx3IDA0aVt0JIbA7p5J2zs/mrxtqDuEoe70f0PEI0rGfg9m+BttMBnDYDL9eiyKrqEtHpz0uZ3K5HkIqXGE1HjaNhV6xn67E5P/7R3WEghneG6mfg8lIgDSL7gV1ZLobJv8b7TuA5P3GxqMIDSYvt+EoqKzbfJvcr2Byfc5EjS/Tk5OB1xOPPFEeDwe/PGPf5zYFgwGsXPnTixZsiSXrgg55uU08smyjLVr1+Lmm29GeXk5GhoacOONN8Lv9+PSSy+dqRgJOSrlfFXDhg0bkEwmcc011yAajWLJkiXYvHkz7Pb8rhESiip47q9tCE1x2cCjXUiN4enBjxCa4rKBR7uRUBSvPv0+RkJTWzYwFzkdcJmq8QMu07Fu53Q6mtftnE5H87qd0+moX7fTaZMwr7IYTpuU/cHHMKdgwwJXGZwCt5du5oXstKFpwSzIzvzlidviq/QW4DvLT0GltyD7g49hVQ43vus/G1WO/K3KxaOyKg+u27AMZVVHnqAwc7gtPkJ4R8VHCCNUfIQwwm3xqZqOgdEYVI1uAJKJqmsYSESg6rTKWyaqqiE4MApVzV+euD0E1jEUxton/sQ6DMs7qAzjq22/Yx2G5fUcDOL7X3k6r/vkduQjhHfcFp+/2I07v3gW/MV0CD2TOocH98y+AHUp1/iRw1XW+XDLvStRWefL2z65LT5JFFBa6IQkWu+sCyuRBBGl9oKUy5OIkSSJ8JUWQpLylyd6RwhhhMkBlweueBQnHT+1mxAqySb0DJ+GBy9/FA7b/in15b7J/KX/+bL8b4/L/iATqhuSwI+A6ruT0NvTX+9ollaQvxtHmrXnuopp6SdZNHax9YHLirF/JD9lQSMfIYxwW3w2qQsVRTfBJnWxDsXS+ntGcNeG19DfM/kV7AToikSwbsef0RVJf3X7TOB2nk8UonDa32MdhuXFY0ns332IdRiWF1VVvD80+XIiM4HbkS+plSIYuQpJrZR1KJbmKXbhglUnwFNM967PpESWcUXTHJQcsTLfTOK2+DTNh+HYF6FpPtahWFqR14mlF89Hkdd6B5WspNjhwGWzG1HsyN9BJW6LjxDeUfERwggVHyGMcFt8ojCCQvn3EAU6hJ5JJBzHW6/tRyScfjFcMmY4kcDLnQEMJ6Z2IkIuuJ1qsEmHUFr4c9ZhWN5QfwRP3vcX1mFYXl8shjt2f5DXfXI78mm6A0qyDppuvVOerMRmlzCr1pNybwFi5BBF1BUWwiHSidVZJVU/eobvQVL1Z3/wMWxWjQc33rYcs2rokqJM/IWF+MXpZ8JfmHqjnpnCbfERwjsqPkIYoeIjhBHOiy9/h4V5pUNHMqFCB63ylokOIKFpec0St1MNDtt+1JVcxDoMy+tqD+I/rniSdRiW1zoygpWvvZLXfXI+8hHCL26LL6H60R26AwmaasiootqDG364DBXVNNWQSW1BIX526umoLaCphqx03YGE2gydJtkzsjsk1M4ugd1Bk+yZyJKIJo8HMq1eRsjRj4qPEEao+AhhhNvis4m9KHP/ADaxl3UoljbYF8bDP9uGwb4w61AsrScaxY/f3YWeaDRv++R2nk8UwyhwbGUdhuVFRxN4960O1mFY3mgyiW2H8vtBzu3Ip2o+DEcvgUoLKGXk9spYsvw4uL35W5WLRz6HAxfX1cNHCyhlp2qlCEa/ApWWDszIW1yAiy4/Gd7iAtahWFqpLOOauS0opaUDCTn6UfERwggVHyGMcFt8ghCBy/4mBCF/N7bgUSySwAc7OhGL0OVXmYwmk3ir7xBGk8m87ZPbqQa71I3you+yDsPyBg6F8eBP/sQ6DMvriUbx/b++k9d9cjvy6boEVfNC1+mE4UxESUBhkQxRottnZyIJAjx2OyQhf3nitvgSagM6g79GQm1gHYqlVfl9+O69l6DK72MdiqU1uN147Oxz0eB2522f3BYfIbyj4iOEESYHXH7ccxpKfGVT6mOWWIGrXcCPek9Hr9Y0pb7OLf5oSs+fCUMnFk9LP0WzigAAw3OLMOSdWl/FL+2ehoim15ybRqeln5qWSuBXZ8D/890QP+qZWmc3mHsYjXyEMMLtVMMhrQ//NfoLJGj5wIw6Do3gX/7rD4gnVNahWFrX3l6sO+s2KFG6S1FWOnQoUFiHYXm6DsQUKrxsdE1HfDS/v0/cfu0sFnz4gnwJigUf61AsraK4AF//wmJU0FUNGZX5i/GVX6xCmX96/tY2g9vicwgONNoa4BBo9bJMnA4JCxvL4KTVyzKSC2W0nN4EuZAuKSLkqEfFRwgjVHyEMMJt8Q1rI3g5/hqGtRHWoVja4HAMj778IQaHY6xDsbRgzzA2/fj3CPYM522f3E41RBHFzuRfWYdheeFoAv+7k1Yvy2Y0GMG23+zI6z65HfmckLFAaoETtCpXJoVOO05bUIVCp511KJZW4HFi8fKFKPA487ZPbovPK3pxofMCeMUpnrB4lCv1OvFPF56AUm/+fql4VFztw5rvX4zial/e9slt8RHCu5yLb3BwENdeey1qa2vh8Xhw5plnYutWWjmakFzlXHyrVq3C9u3b8fjjj+Ptt9/GSSedhGXLlmHPnj0zER8hR62ciu/jjz/GK6+8grvuugtnnXUW5s6dizvuuAPV1dV47LHHZirGtBJ6Ap1qFxI6XdWQSTyhYn9nkK5qyEKJKmh/NwAlmr+Tq3OaaigrK8OLL76IU045ZWKbIAgQBAFDQ0PTHlwmg/oQHok9kdd98qh3MIIfPvIW6zAsr+/AIP77ql/mdZ85FZ/P58Py5csN2zZt2oSPP/4Y559/fsrjGxsbDf+Px+MAgBKxBLPEiontMT2GkD4MCRLKxNR7L/Rqh8aeJxTDLhgPmYe0EGKIwwUXPGKRoU3RFQzpQQgQUCGWp/Tbp/VDgwZZL4Oku4zPFYJICiOQ9ALIujEmTUggJoxd7Vygpd4TPir0QBcScOglsOnGe3wnhBEkhCBEXYZTrzC06VARFbsAADXlbkiicSWtQ0MRxBQVPrcMT6HxhPLRaAIDwzHYbSKqSlPvK36wd+xkhKrSQthtxi88/cEoIvEkPAUO+IqMUzfReBJ9wSgkUUDN3NQcdn3cD13TUVbrg1xgfG+CvWGMhqJwFckoqTLeE16JJtHXMfaBna7f3vZBJBUVJVUeuI6IaXhgFCMDEcgFDpTUGt/zpJJEb2v/2GudUwHxiNs89x0YgBJNwFtRBHeJMU+jwQiCPcOwO22oaDCutKBrOrr2jt3FaNbsMthkY+kMBIYQC8dRlCb3k5nSJPv27dtx9dVXY+XKlVixYoXp510kX4D5rvkT/38/uRsvxH+PIsGNq11rUh7/49GfAgBWyMtQI1Ub2l5XtuGNxJ8xzzYXy+RzDW2tyXb8Jv4M7LCn7ffno3cjiijqkl9AsX6ioe2g9CR6pD/Ao7dgTvKrhrZR4SA+sP8AADA/eTNEGH/p3rP9P0SFblSrK1ChnWlo6xJ/h4DtWRTq9ZiX/IahTcEQdjn+AwBww2WLUXLEnNN//urP2HNwCOcursOK040fbFv+GsAvf/cByn0urL/6jJTX+t0Ht+Ng7wj+6cITUF9pLIS7ntmFt/f04rQFVfji0hZD2659h3DHpnfgkm34+v2pOfzW+RsRjyi45Ia/wXGn1hvanv7pa9j+zLuYd9psrP628cO5/f1u3PnPY99c0vX7o394EAOdIZz35dOxeNk8Q9vLD76Jlx98E/ULqvCV2y8xtPV3DOKHF98FALju7tVwFxuL4edXPYQD73bi7C99Fud86bMp+7199f8AAL7xqy8btsfCcXxzyU8AAFf+50pUNhk/MP7nht/ggy37cOpFxt+jTARd13XTjz7Mc889h9WrV+Nzn/scnn/+eTid2eeRdu7cicWLF2PNw1dhVkvlxPZPM/KVCsW4yLkcj0WfRIcWmNLId4Gv33Ij3wt3XTEtI19laSG+etGJE8U3lZFv3t6+lNfKfOQrMbZ92pGvoqEUl//w73H76v/Bofb+KY18m7rvT3k96Xyqke/OO+/E9ddfj8suuwwPP/wwHDne02xQGwS01GM9KtSJQkv7PH0IGP+o+L+nj1/NHkUUUS39XUV16Bn7jQv9wCRrpapCBJEMS9JHxMlP3VKEQSjCYNo2TYgjIkz+3M4Md5INhuMIhuNp2xJJbeIrZjrdA5MvODQcUTAcSX/AQdV0dKYpvnH9geCkbdGRODpHJn9upn4Hu4eB7vRt8YiCzkOTH2vo3jf5ex46NILQofR5SsSS6MywiFJvW/+kbSMZ8nuknKca7rrrLnzta1/D2rVr8fjjj+dceISQMTmNfHv37sX111+PSy65BOvWrUNv7ye30XW5XPB66VQvQszKaeR76qmnkEgk8Mwzz6Cqqsrw7/rrr5+pGNPq1wZwd+QB9GsDed0vb7r6R7Hu7j+hq3961rc8WvW29uGHF29Eb+vkX4GnW07F981vfhO6rqf999BDD81QiOmpUBHUQ1BBk8eZJFUNh4IRJFWNdSiWllRU9HcMIZnHld64PbHaK3jwd/L58Aqe7A8+hpV5Xbjm745HmdeV/cHHsJJqL9Z8/yKUVOfvTydui88pOLHQNg9OgS6VyaTAacPpC6tR4OT2uum8cHlcWLz8eLg8+fuQ4rb4COEdFR8hjFDxEcIIt8UX1kexVXkDYZ0OoWcSCsfx3NaPEZrkjBgyZrg/jM33bMFw/+RnFk03bv8KH9VHsTXxJuswLC80quD5rftZh2F5I/1hbL7nT3ndJ7cjnwMOzJbq4QCd3paJ0yFhwexSuldDFnKhA8ed3gi5MH+/T9wWX7HowxedK1Es+liHYmkVxQX4ty+eQncpyqLMX4Kv/uIfUOYvyds+uS0+QnhHxUcII0wOuLzXXQ3ZVZ39gRk0Or1AI/BedxVaY1P7StX5rTlTev5MiC6cns/FuE+c+BlVp9bnK+9vno6QplVHcnrOWdXVZmD0y/jGczshSB9PS5/ZcDvyJXQN3UoYCZ1OGM4koWroCUWQoBOrs0gAQtfYzzzhdqqhIz6C6z5+hXUYltcZHMW/PbmNdRiWJ0gHgKKr87pPbkc+QnjHbfHVyx78cu5y1Mt0SVEm/hI37l5zNvwlbtahWJquzoY+8gR0dXbe9slt8UmCAK9NhiRMsvIRATCWJ4/LQXnKSgJ039jPPOG2+AjhHRUfIYxQ8RHCCLfF1xUP4+a219EVz98lIDzqDkWw/vk/ozs0+cK/BIAYAApuGPuZJ9zO88V0FXui6VeDJp+IJ1XsOxRiHYblCUIMsO3O6z65HflKbU5cPet4lNpoAaVMSgpkfOmzc1FSIGd/8DFM18qgx74CXSvL/uBpwm3xeW0yLi5thtdGv1SZeFwOLD++Hh4XXfeYke4DlEv/b7ohP7gtPkJ4R8VHCCNUfIQwwm3xDScVvDTYiuFk/m5gz6ORWAIvf9iBkVj+LpXhkhAC7M+P/cwTbqca+pNR3NvzV9ZhWN7AaAwPbf+IdRiWJ4h9gOsXed0ntyOfQ5DQ6PTCIdCqXJk4JBENpUVwSNy+1Xmh6zJ0tRm6nr+j59y+I7WyG//VeC5qZbpUJpNqXyF+eMlpqPYVZn/wsUzzA6O/GPuZJ9wWHyG8o+IjhBEqPkIY4bb4dB2IqAnoOutIrE3XdUSUJHRKVBYagNH/+5kf3E41tMVDWL3nBdZhWN6BwTCuefh/WYdheYLUCnhW5nWf3I58hPCO2+KrdRThvxuXotZRxDoUS6vxFeI/Lz0dNTTVkJGu1kEP3wtdrcvbPrktPocoos7pgUPk9iXkhV0SUVvshp0m2bNwAFr92M88oXeEEEao+AhhhIqPEEa4Lb5eJYIfHHwDvQqtypXJoZEofvLyLhwaibIOxdrEbsC1fuxnnnA7zzeqJfCXcA/rMCwvoiSx82Af6zAsTxBGAfubed0ntyOfT5Jxaelc+CRaQCkTr8uBi05sgJcWUMpI14qhx78IXSvO2z65Lb4SuxOXz1qAEjstHZhJcYGMVZ+Zg2JaOjAzvRSI/+PYzzzhtvgI4R0VHyGMCHoeT3ffuXMnFi9ejB07dmDRokVT6isajWL//v1oamqCy+WapgiPPpQnc1jkiduRT5IkeDweSBKt4ZIJ5ckcFnnidqrB4XCgri5/J8HyivJkDos8cTvyaZqGRCIBTcvfxY88ojyZwyJP3BZfPB7Hnj17EI/HWYdiaZQnc1jkidviI4R3VHyEMELFRwgjVHyEMMLtVIPT6cT8+fMhCALrUCyN8mQOizxxW3yCINAvlAmUJ3NY5Inbr53xeBytra10CD0LypM5LPLEbfFpmoZIJEKTx1lQnsxhkSdui48Q3n3q4tu7dy/cbjceeuihaQyHkGPHpyq+RCKBNWvWYHR0dLrjIeSY8amKb/369fB4PNMdS07sdjuqq6tht9uZxmF1lCdzWOQp56mGLVu24J577sGuXbuYXqpis9lQUlLCbP+8oDyZwyJPORVfMBjE5ZdfjjvuuAN+f/Z7Vzc2Nhr+P34YNxaLIRr9ZB1JSZLgcDigaVraQ73jVxbH4/GJo1GqqiIcDqO4uBiyLCOZTCKRSBieJ4oiZFmGruuIxWIp/TqdTgiCYOh3nN1uh81mg6qqUBTF0CYIApzOsYWbDn8d42RZhiiKUBQFqqoa2mw2G+x2e9Z+Y7FYyj31HA4HJElCIpFAMpk0tE2WQ1VVMTo6itLSUthstoz9zlQOM/ULWCOHqqoiEomgpKQEoihm/D3M1q/Z0TOn4rvuuutwxhlnYPXq1bk8LUUgEDB8bfV6vfD7/UgkEti/f3/K4xcuXDjxvCPfKFEUUVFRgVAohO5u44KnbrcbDQ0N0DQtbb8tLS2w2Wzo6enByMiIoa2yshJlZWUIh8Po6OgwtDmdTjQ3NwMAWltbU96I5uZmOJ1O9PX1YWhoyNBWVlaGyspKRKNRtLe3G9psNhtaWloAAO3t7SkF1tDQALfbjYGBAfT39xvaiouLUVNTA0VR0r5Wj8cDm82GQCCQUkR+vx9erxfBYBA9Pca1UIuKilBfXw9VVdP2O2/ePEiShO7uboTDYUNbVVUVSktLEQ6HEQgEDG0ulwtNTU0AkLbfOXPmQJZl9Pb2IhQKGdrKy8sxa9YsRCIRHDhwwNDmcDgwd+5cAEBbW1tK4TY2NqKgoAD9/f0YGBhI2W9RUVHamERRxPz58wEAHR0dKcVZV1cHj8eDYDCI8vLylH7TMb2GyyOPPIJ169bhvffeQ3Hx2NqGgiDgwQcfxFVXXWVqZ+NruGzbtg0nn3zyxPZPM/LF43EEAoGJX0ga+dLncDxP42uT0MiXPoeH50mWZWuNfA888AB6e3tTvm5ee+21eOKJJ/C73/3ObFdwOp1pF6kRRTHj4jXjb9bhxtfcsNlssNnSvxxBEHLu9/D+Mz03U5vDMflCtdn6Hf8FSsdut0/6BmfLYaZ+ZyqHmfoFjq4c5nLAxnTxPfrooymfUHPmzMGGDRuwZs0a0zskhIwxXXw1NTVpt1dUVEzaNpNEUURBQQFEujlmRpQnc1jkidurGmRZTjmaSlJRnsxhkacpFV8e19tNu29d1+mSmSwoT+awyBO330VisRg+/PDDtEfgyCcoT+awyBO3xUcI76j4CGGEio8QRqj4CGGE66mG4447ju6+kwXlyRwWeeK2+ERRpIljEyhP5rDIE7fviqIoOHjwYMqJtcSI8mQOizxxW3yqqmJ4eDjljHdiRHkyh0WeuC0+QnhHxUcII1R8hDDCbfHZbDbMmjUr40WahPJkFos8cfuO2O1202tlHMsoT+awyBO3Ix8dxTOH8mQOHe3MAc1fmUN5Mofm+Qg5hlDxEcIIFR8hjHBbfIIgQJZlWpckC8qTOSzyxO1Ug9PpxJw5c1iHYXmUJ3NY5InbkY8Q3nFbfNFoFB9++GHadf7JJyhP5rDIE7fFByCvN6/nGeXJnHznieviI4RnVHyEMELFRwgj3BafLMsTNzIkk6M8mcMiT9zO82W7gSEZQ3kyh0WeuB35FEVBV1cXna2fBeXJHBZ54rb4VFXF4OAgXaeWBeXJHBZ54rb4COEdFR8hjFDxEcIIt8Vns9lQWlpKq3JlQXkyh0WeuH1H7HY7qqqqWIdheZQnc1jkiduRT1VVRCIROoqXBeXJHBZ54rb4FEVBa2srzV9lQXkyh0WeuC0+QnhHxUcII1R8hDDCbfEJggBJkmhVriwoT+awyBO3Uw1OpxPz5s1jHYblUZ7MYZEnbkc+QnjHbfHFYjHs3bsXsViMdSiWRnkyh0WeuC0+XdehKAp0XWcdiqVRnsxhkSdui48Q3lHxEcIIFR8hjHBbfA6HA/X19XA4HKxDsTTKkzks8sTtPJ8kSSgqKmIdhuVRnsxhkSduR75EIoHe3l4kEgnWoVga5ckcFnnitviSyST6+vqQTCZZh2JplCdzWOSJ2+IjhHdUfIQwQsVHCCPcFp8kSfB6vZAkiXUolkZ5ModFnridanA4HPD7/azDsDzKkzks8sTtyKdpGuLxON3yOAvKkzks8sRt8cXjcezbtw/xeJx1KJZGeTKHRZ64LT5CeEfFRwgjn6r4Hn74YcyfPx9OpxMLFizAk08+Od1xEXLUy7n4Hn30UXz5y1/G2rVr8cEHH2D16tVYtWoV3njjjZmIj5CjlqDncN28rutobGzE5z//edx2220T28877zycc845WLduXcbn79y5E4sXL8aOHTuwaNGiTx81IUeBnOb59u7di/b2dqxevdqwffPmzdMaFCHHgpyKb8+ePQCA0dFRnHfeeXjnnXcwe/Zs3HLLLbjwwgtTHt/Y2Gj4//hh3Fgshmg0OrFdkiQ4HI6JuZYjuVyuieePz8MoioJDhw6hpqYGBQUFSCaTKZeDiKIIWZah63raVamcTicEQUg7v2O322Gz2aCqasrNMwRBgNPpBADD6xgnyzJEUYSiKCl3vbHZbLDb7Vn7jcViKYv5OBwOSJKERCKRcvb9ZDkcz1NdXR1kWc7Y70zlMFO/gDVyeHie7HZ7xt/DbP3a7faU56aTU/ENDw8DAK644gqsX78et956KzZt2oSLL74Yr7zyCpYuXWqqn0AgAI/HM/F/r9cLv9+PRCKB/fv3pzx+4cKFE8878o0Kh8MoKChAKBRCd3e3oc3tdqOhoQGapqXtt6WlBTabDT09PRgZGTG0VVZWoqysDOFwGB0dHYY2p9OJ5uZmAEBra2vKG9Hc3Ayn04m+vj4MDQ0Z2srKylBZWYloNIr29nZDm81mQ0tLCwCgvb09pcAaGhrgdrsxMDCA/v5+Q1txcTFqamqgKEra1zpeGIFAIKWI/H4/vF4vgsEgenp6DG1FRUWor6+Hqqpp+503bx4kSUJ3dzfC4bChraqqCqWlpQiHwwgEAoY2l8uFpqYmAEjb75w5cyDLMnp7exEKhQxt5eXlmDVrFiKRCA4cOGBoczgcmDt3LgCgra0tpXAbGxtRUFCA/v5+DAwMpOx3/MPryJhEUcT8+fMBAB0dHSnFWVdXB4/Hg2AwiPLy8pR+08npb74nnngCq1atwsaNG3HddddNbF+xYgUA4MUXX8z4/PG/+bZt24aTTz55YvunGfni8TgCgcDELySNfOlzOJ6npqYmuFwuGvkmyeHheZJl2XojX21tLQDg+OOPN2xfsGABXnjhBdP9OJ3OiRdyOFEU024fN/5mHW78RFibzTbpLX0FQci538P7z/TcTG2Z1gPJ1u/4L1A6drt90jc4Ww4z9TtTOczUL3B05dBs4QE5TjUsWrQIRUVFePPNNw3b33vvvYmvYYQQc3Ia+VwuF2666SZs2LABNTU1OPXUU/HrX/8aL7/8Ml599dWZijEtu92O2tranD5pjkWUJ3NY5CnnS4puueUWFBQU4Fvf+hY6Ozsxb948PP300zjnnHNmILzJ2Ww2+Hy+vO6TR5Qnc1jkKacDLlM1nZPsyWQSoVAIXq83498TxzrKkzks8sTtidWJRALd3d20JF4WlCdzWOSJ2+IjhHdUfIQwQsVHCCPcFp8oinC73RBFbl9CXlCezGGRJ24Pf8myjIaGBtZhWB7lyRwWeeL241DXdaiqSrc7zoLyZA6LPHFbfLFYDLt3787rDex5RHkyh0WeuC0+QnhHxUcII1R8hDBCxUcII9xONTidTrS0tNDdd7KgPJnDIk/cFp8gCHSWvgmUJ3NY5Inbr53xeBwHDhygG4BkQXkyh0WeuC0+TdMwMjJCt77KgvJkDos8cVt8hPCOio8QRqj4CGGE2+Kz2+2orKykVbmyoDyZwyJP3B6DttlsKCsrYx2G5VGezGGRJ25HPlVVEQqFUpYSJ0aUJ3NY5Inb4lMUBR0dHSlr9RMjypM5LPLEbfERwjsqPkIYoeIjhBFui2/8PmyCILAOxdIoT+awyBO3Uw2H3x2WTI7yZA6LPHE78hHCO26LLxqN4oMPPkh7S2HyCcqTOSzyxG3xAaC1KE2iPJmT7zxxXXyE8IyKjxBGqPgIYYTbqQZZltHc3AyHw8E6FEujPJnDIk/cFp8oinA6nazDsDzKkzks8sTt105FUdDZ2Uln62dBeTKHRZ64LT5VVTE0NETXqWVBeTKHRZ64LT5CeEfFRwgjVHyEMMJt8Y0veEP3IciM8mQOizxx+46ML/VGMqM8mcMiT9yOfKqqIhwO01G8LChP5rDIE7fFpygK2tvbaf4qC8qTOSzyxG3xEcI7Kj5CGKHiI4QRbotv/Da+tCpXZpQnc1jkiduphvEb2JPMKE/msMgTtyMfIbzjtvhisRg++ugjxGIx1qFYGuXJHBZ54rb4dF1HMpmklbmyoDyZwyJP3BYfIbyj4iOEESo+QhjhtvgcDgcaGhpoVa4sKE/msMgTt/N8kiTB7XazDsPyKE/msMgTtyNfIpFAT08PEokE61AsjfJkDos8cVt8yWQS/f39SCaTrEOxNMqTOSzyxG3xEcI7Kj5CGMm5+JLJJL7zne+gvr4eRUVFWLJkCd58882ZiI2Qo1rOxff9738f9913H+677z688847OO6443D++eeju7t7JuKblCRJKC4uhiRJed0vbyhP5rDIU87F9+yzz2L16tVYtmwZmpubcfvttyMUCuGNN96Yifgm5XA4UFNTQ/NXWVCezGGRp5yLr6KiAi+88ALa29uhqiruvfdeyLKME088cSbim5SmaYjFYtA0La/75Q3lyRwWecp5kv3nP/85LrvsMsyePRuSJEGSJGzatAlNTU0pj21sbDT8Px6PAxi7fOPwG89LkgSHwwFN0yYecziXyzXx/PHkxONxBAIBNDQ0wO12I5lMpszRiKIIWZah63raS0WcTicEQTD0O85ut8Nms0FV1ZQVrQRBmLid1OGvY5wsyxBFEYqipCxFZ7PZYLfbs/Ybi8VSzrB3OByQJAmJRCLlkPhkORzPU1NTE1wuV8Z+ZyqHmfoFrJHDw/Mky3LG38Ns/drt9pTnppNz8X344Yfw+Xx49tlnUVNTg/vuuw9r1qzB66+/jpNOOslUH4FAAB6PZ+L/Xq8Xfr8fiUQC+/fvT3n8woULJ5535BsViUTgdrsRCoVS/u50u91oaGiApmlp+21paYHNZkNPTw9GRkYMbZWVlSgrK0M4HEZHR4ehzel0orm5GQDQ2tqa8kY0NzfD6XSir68PQ0NDhraysjJUVlYiGo2ivb3d0Gaz2Saupm5vb08psPEPmoGBAfT39xvaiouLUVNTA0VR0r7WcYFAIKWI/H4/vF4vgsEgenp6DG1FRUWor6+Hqqpp+503bx4kSUJ3dzfC4bChraqqCqWlpQiHwwgEAoY2l8s18YGdrt85c+ZAlmX09vYiFAoZ2srLyzFr1ixEIhEcOHDA0OZwODB37lwAQFtbW0rhNjY2oqCgAP39/RgYGEjZLzBWiEfGJIoi5s+fDwDo6OhIKc66ujp4PB4Eg0GUl5en7fdIgp7DBUwdHR1obm7GH/7wB5x11lkT28866yyUlpbi2Wefzfj8nTt3YvHixdi2bRtOPvnkie008qX2SyMfjXwGb731FhRFwWc+8xnD9tNOOw0vvfSS6X6cTufECzmcKIppt48bf7MON350ymazTbrOviAIOfd7eP+ZnpupLdMf79n6zXSXVLvdPukbnC2HmfqdqRxm6hc4unJotvCAHA+41NbWAgDeffddw/Z33313YqjPJ1qRyxzKkzn5zlNOI9+pp56KM888E1deeSU2btyI2tpaPPzww3j11Vexbdu2mYoxLZfLhQULFuR1nzyiPJnDIk85jXyiKOL555/Hueeei6uuugqLFy/Ga6+9hldffRWf/exnZypGQo5KOR1wmarxAy47duzAokWLptRXLBZDIBBAbW1txu/gxzrKkzks8sTtidXjR99oVa7MKE/msMgTt8VHCO+o+AhhhIqPEEa4LT6HwwG/309n62dBeTKHRZ64Xr3M6/WyDsPyKE/msMgTtyMfLQxkDuXJHFpAKQe0JJ45lCdzaOlAQo4hVHyEMELFRwgj3BafKIooKiqCKHL7EvKC8mQOizxxO9UgyzLq6+tZh2F5lCdzWOSJ249Dut2xOZQnc+i20DlgcQN7HlGezGGRJ26LjxDeUfERwggVHyGMUPERwgi3Uw1OpxPz5s2j+assKE/msMgTt8UnCALd9soEypM5LPLE7cdhPB5He3t72mW9yScoT+awyBO3xadpGsLhMN36KgvKkzks8sRt8RHCOyo+Qhih4iOEEW6Lz263o6qqKqdbMh2LKE/msMgTt1MNNpsNpaWlrMOwPMqTOSzyxO3Il0wmEQwGaVWuLChP5rDIE7fFl0gkEAgEaFWuLChP5rDIE7fFRwjvqPgIYYSKjxBGuC0+URThcrnobP0sKE/msMgTt1MNsiyjqamJdRiWR3kyh0We6OOQEEa4Lb5oNIr3338f0WiUdSiWRnkyh0WeuC0+QnhHxUcII1R8hDBCxUcII1xPNcyZM4culcmC8mQOizxxW3yiKEKWZdZhWB7lyRwWeeL2a6eiKOjo6ICiKKxDsTTKkzks8sRt8amqilAoBFVVWYdiaZQnc1jkidviI4R3VHyEMJLXAy7jp+7s3r17yn3FYjEEAgEMDw/D6XROub+jFeXJnOnOU0tLCwoKCjI+Jq/F197eDgD40pe+lM/dEpJ3O3bswKJFizI+RtDzeBPq/v5+bN68GQ0NDXC5XFPq68ILLwQA/Pa3v52O0I5alCdzpjtPZka+vBbfdGpsbAQAtLa2Mo7E2ihP5rDIEx1wIYQRKj5CGOH2aychvKORjxBGqPgIYYSKjxBGqPgIYYS74tM0DevXr0dNTQ0KCwuxfPlytLW1sQ7LcgYHB3HttdeitrYWHo8HZ555JrZu3co6LEvbu3cv3G43Hnroobzsj7vi+973voeNGzfi3nvvxfbt26GqKs477zy6Xu0Iq1atwvbt2/H444/j7bffxkknnYRly5Zhz549rEOzpEQigTVr1mB0dDR/O9U5Eo/H9aKiIn3jxo0T24aGhnSXy6X/6le/YhiZtezbt08HoG/dunVim6ZpelNTk/7tb3+bYWTWtW7dOv3cc8/VAegPPvhgXvbJ1ci3a9cujIyMYOnSpRPbfD4fFi1ahC1btjCMzFrKysrw4osv4pRTTpnYJggCBEHA0NAQw8isacuWLbjnnnvy9nVzHFdruAQCAQCA3+83bK+urkZHRweLkCzJ5/Nh+fLlhm2bNm3Cxx9/jPPPP59RVNYUDAZx+eWX44477kj5vZppXI18kUgEAFIWunE6nYjFYixC4sL27dtx9dVXY+XKlVixYgXrcCzluuuuwxlnnIHVq1fnfd9cjXzjlyHF43HDJUmxWAyFhYWswrK05557DqtXr8bnPvc5PPbYY6zDsZRHHnkEf/rTn/Dee+8x2T9XI9/414Kuri7D9q6uLtTU1LAIydLuvPNOrFy5EhdeeCFeeOEFupL9CA888AB6e3vh9/vhdrvhdrsBANdeey0uuOCCmQ8gL4d1pkksFtM9Ho9+//33T2wbP9r5+OOPM4zMejZu3KgD0P/1X/9V1zSNdTiWFAgE9H379hn+AdBvvfVWPRAIzPj+ufraKcsy1q5di5tvvhnl5eVoaGjAjTfeCL/fj0svvZR1eJaxd+9eXH/99bjkkkuwbt069Pb2TrS5XC54vV6G0VnHZN+WKioq8vJNiqviA4ANGzYgmUzimmuuQTQaxZIlS7B582ZaDv0wTz31FBKJBJ555hk888wzhrYrr7wy74fUSXp0PR8hjHB1wIWQowkVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGM/H+9TjBQsUOT/AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(temp)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "calib\n",
      " range(0, 4)\n"
     ]
    }
   ],
   "source": [
    "# Proxies\n",
    "proxy = np.random.rand(3,10).T # 10 years, 3 proxies\n",
    "\n",
    "# calibration period\n",
    "calib = range(0,4)\n",
    "print(\"calib\\n\",calib)\n",
    "G = GraphEM.solver.GraphEM()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJUAAAGdCAYAAAAbhjGWAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAs8UlEQVR4nO29eXwV1f3//5y5e+6SPYEkbAE3FEVwF1CxgkLdsP6KoFVbrdqqtFqsW6XVqh9/autWrEvVUi1SxaLiQhV3cakGFEVBZDELCWS5ufs2M98/UuPjkpvlmsncO3qejwePS+Y9c97vO3nlPTNnznkfSdM0DYFAR+RcByD47iFEJdAdISqB7ghRCXRHiEqgO0JUAt0RohLojhCVQHcMFVUkEqGuro5IJGKkW4HBWI109vnnnzN58mTefHEYEyfYB9VWSi0mGDsBr/MFrHLHoNo69I5fDer4oWDGGe/q0o4LD3vLk/hcrSNKaFBt/fnAZQPaz7SXPwkNWYogId4y9YWGRpI4moHnydBMpScW2U+ha0Wuw8h7YoT5RH3PUJ+mzVSK6iGcmIKienIdSl5jx8loaR/sOA3zaVpRpdRKdoWuIqVW5jqUvMZLEcdaf4SXIsN8mlZUgvxFiEqgO0JUAt0xragkKYHdshlJSuQ6lLwmRYpWbQcpUob5NG2Xgt1ST1XhglyHkfd00srTqQcN9WnaTCXIX0wrqniqlm3tK4inanMdSl5TyjDOsV5FKcMM82laUYEE2P73KegLi2TsXU7WolJVlUWLFlFdXY3b7WbWrFls3bp1KGITmJSsRXXDDTewePFi7r//ftasWYOiKMycOZNEQjyFCbrISlSJRILbb7+d66+/ntmzZ3PAAQewbNkyGhoaWL58+VDFKDAZWYlq3bp1BINBjj322O5tRUVFTJo0iTfeeEP34PrCZqmnyncRNku9oX7Nhp9Wlif/ip9Ww3xmdQfX0NAAwIgRI9K2V1VVUV/f85dbW5v+ZBaPxwFIpGqIp7zd22UphM3SgqrZSCoje7TjsH4JQFKpRtXS37Zrmh2kBIrqI6WWp9lkKYrN0oSmySSUMT3atVu2IUkKlR4XBfb0U9EeidMZS+C2W6nwuNJsCUWlsTMMwOgSb49HhYbOMElFpcztxOuwpdn80QQd0ThOq4XhvoI0W0rVqPd3DaQrpgJ5t7/5AO0kSVCAFxfuNFucKCE6sWCliLIe31X5X+dnIWVYd/u1B/GTIIYTN268abYEcYJ0IGWRf7IS1dfDgB0OR9p2p9NJe3v7gNtpDS9kR2B8989u+6uUe25DUcvYEbirx/6jS2Z3HRf6NXFlnzRbses+Cl3PEE5MpT3yi/S4rB8yzHcdGo6M7Y4oOgOLFOAnB+/FQSPSBbnkvxt57rOvmDC8hF8fdUCabWtbgCuf6xqj9McTDsFmST/hlz+9hobOMKftX8v0ParTbCvWb2Xp2s3UlvpYNPOgNFtbOMYvlr8JwEzrGbglX5r9udQSmrXtjJcP5gDLkWm2jepa3lJW4qWYU2znp9lUTeWJ1D2E6ORo6ymUScPT7KtTT7JN+4xx8n4capmRZtuubuJlZRmOLIbOZCUql6vrLzYej3f/HyAWi+F2u3vsv2XLlrSf6+rqmDx5MmXuWxnuS89UABa5leG+S3v1X+b5c3emSio1tIavwG7tevJ029/EYf08bX9ZigIgEc/Yrix1ZZsl/93Ikx99mWZrj3Rl1fU72rlyZfrQ3oSidv//2hfe75GpWkJdfpd/vIX/bEzP4P5o1wPNlrZAj3ZT6jejM1ellmbMVAAb1P+yVd2QZovT5TNIByuSD3RvL6SUY2xzcOAiRCevpVZkzFQAm9VP2KFuT/+uxP/XfoyBkpWovr7sNTU1MXbs2O7tTU1N7L///gNux25twGHtOUZdlpLdl7pM2CyNGY7pyp4WOYBFDmQ8TpLUPtv9WgSZCCdSbG0P9mrf1oetNRyjNZz5lxFLKX2228HOXm0RgkTIfKxCijaaez22s497qxhhYoQz2jTUjNszkdWN+gEHHIDP5+O1117r3ub3+6mrq2PatGnZNCX4DpNVpnI4HFx88cX89re/pby8nNGjR7Nw4UJGjBjBaaedNlQxCkxG1v33119/PalUivPOO49oNMq0adNYtWoVNput/4N1xCL7KXT+C4vsN9Sv2YgS5iPlbaK9XNaGgqxFZbFYuOWWW7jllluGIp4BY5XbKC74e05jMAMRgnygvmKoT9O+UFY1F9HkBFTN1f/O32Ns2BkmjcLG4CbvZoNpRZVUqmgJ/h9JpSrXoeQ1PkqYbf0JPkoM82laUQnyFyEqge4IUQl0x7SikkhhkVqRDJwlYkZUVMJaADWLHvHBYt7ZNNbtjCg+O9dh5D0d7OTx1J2G+jRtphLkL6YVVSI1ivqOv5NIjcp1KHlNMRXMtS6gmArDfJpWVBpWFK0MzbxXcEOQkXFLvh7DaIbWp0CgMzn5M7+3bQIVLcWDaqOUYZxi62qrjfL+D+gD31fGPRkNlI9+M1GXdmpqCmEhbLpjLxoaBjmhdPXAdhOZSqA7phVVgHaeSy3pHmIryMyuXSHuunsNu3YNrjJxNpj2LjdJgmZte/87fs+JxxU2b24z1KdpM1UBXg6Sp1Ow25QiQTqFhU5O/OHeFBaKQrL94sLNAZYje8x/E6Tj9To47rg98Hod/e+sE6YVlSB/EaIS6I4QlUB3TCuqOFE2qmu7Z+YKMhMOJ3jnne2Ew8aVejJtl0KITt5SVuY6jLynoyPK0sc/NtSnaTNVV3WTcizm/bswBJtNZtgwDzabeKHcL0WUcZrtwoxlcwTfUFnp5eqrjqGy0rj+PNOKSpC/CFEJdEeISqA7phaVoomZNP2iaaRSCmhiudt+aaOZR1I35zqMvKehMcBllz9vqE9TZypBfmJaURVSxsnW8ygUXQp9UlnpYeFvplFZadxa06YVlRUrZdLwHkVRBenYbBZGjCjEZrMY5tO0ohLkL0JUAt0RohLojmlFFcTP6tST3YXlBZlpa4vw0MMf0NYWMcynae9yE8TYpn2W6zDynmg0ybp1Owz1adpM5cTNfvKhOMXEhz7xeu0cc3QtXq8oJNsvbrwcapnRY9UnQTqFhS5OPXVfCguNq+JsWlEJ8hchKoHuCFEJdMe0okoQZ7u6qXs9OkFmYtEk69c3E4smDfNp2i6FIB28rCzLdRh5T2tbhAce/K+hPk2bqSRknBRktbbv9xFZlvC47cjy7uunDqFPwzzpTAkVzLddTomBBVLNSFWVj5tumklVla//nXXCtKIS5C9CVALdycmN+hnFdexXPriVTJXUOKKhrrYs1s2Daqv297sGdfxQ8Mz4Ul3asewzDJiG5YPPsHzW+4LdeiIylUB3TNulIFu24C48GYjlOpS8pnFjC1cefjMJ0U/VP5KkAsaNETIrmqoRN7CMEJj48qcq1URDN6Mq1bkOJa8pG1nCBffOp2ykWO62XzTNhZI6CE0szN0nTredvY8ch9MtxlMJTIwQlUB3hKgEumNaUcnyLuyuu5Hl/Ou4zCc6mgM8eePzdDQHDPNp3i4FuRO745lch5H3hDsivL1MDH0ZEJrqJZk4Fk0VEx/6osDnZPLsCRT4xNo0/aKqlcQjV6KqlbkOJa8pqS7izJvnUFJdZJhP04pKkL9kLar29nYuvPBCampq8Pl8TJkyhbfeemsoYhOYlKxFNXfuXNasWcPSpUv54IMPmDhxIjNmzGDjxo1DEZ/AhGQlqs2bN/PSSy9x7733MnXqVPbcc0/uvvtuqqqqeOyxx4YqxoxIUgzZsgFJEqMU+iIeTbLto3ri+TpKoaysjOeee46DDjqoe5skSUiSREdHh+7B9YVsaaDAu8BQn2Zk17Y27jzrIUN9ZiWqoqIiZs2albZt+fLlbN68meOPP77H/rW1tWk/x+Ndc/QUZQRK6psalJIUQrY0o2k2VGVUj3a+HtmpKjVoWvqjsSy3IMlBNLUQVS1Ps0lSFNnSiKbJqEp6LACyZSuSpODUyrCS3m6cTpJSEKvmwkn6KEyVJBGpBQCPVtOj3QgtqFISh1aMbbcCIgmCJKROLJoDF+nxaiiEpa4KLcP3qMBiTb+Q7NreTjySoLDCi7c0vd2wP0rHjk5sDiuVtel1UDUNGj/vGvVZMaYMuzP9197W4CcajOEpcVO023Ij0WCctoYOZOvAL2qD6vxcs2YN5557LnPmzGH27NkDPi4euYpoaPw3Qdhexum+BU0tJxq6t8f+nqLjAIhFFqIq49NsNsfDOFz/JJk8ikT0kjSbxfoBLs9VgDNju27fj0DqZA/1VMq0CWm2zfK/qZdepVjbi/3Un6bZgtTzgfVWACYrlyHvdhrfs9xEhGZGq8dTpR2eZtsuvcQWy7N4tREcqF6aZovh5x3rdQD8fPF8iirTZ8Dc89NH+PKD7UyZewg/OG9Kmu3dp+pY9vtnKa0p5vJlF/T4rrf/+D4aPmvmrP+bQ80+w9Nsj1z+BB+9tIHJsydwysKZabZPXtvI3y59HJd34P1ckqZ9u6rtTz/9NPPmzePII4/kmWeewens32ldXR2TJ09mxbMT2W+/wWUqVRlBPHo1TvflWG0fDypTvRScmneZ6r8n76NLpqoYU8ZZ/3dat6gGk6lWJQY2efdbZap77rmHBQsWcPrpp7NkyRLs9uzG6lgs9VisPSc+SFKyz0kMsqUhwzFdoz8luROL3JnxOElS+2w3JrX2aktJUUL09Ps1Ial3W1zqIE7me01FivfZ7o4vdvZq69wZpHNnMKMtGU/R0McEh51be/+uofYwofZwRpuaUns9bney7lK49957ueSSS7j44otZunRp1oISfPfJKlNt2rSJBQsWcOqpp3LVVVfR0tLSbXO5XBQWFuoeoMB8ZJWpnnzySZLJJP/+978ZPnx42r8FC4x9vJct2ynwno1s2W6oX7PR/OUubpx9F81fGjdEKCtRXX311WialvHfI488MkQhZkaSksiWJiTJuE49M5JKKLTWd5BKKIb5NO0LZVUZRiz8W1RlWK5DyWtKqouYf9OpYpTCQNA0D6nkD9A04xbyMSMFPicH/XB/MZ5KYG6EqAS6I0Ql0B3TikqS27E5liDJ7bkOJa8J7Arx4r2vEdgVMsynaWfTyHI7Dtc/ch1G3hNoDbHq3tcN9WnaTKVpBaSSB6FpBbkOJa9xuO3sdcRYHKKWQv+oShWx8M2oSlWuQ8lrykeWcOFfz6RcVH0RmBkhKoHu5ORGvUlx4UkNrofXohRQDDQpBSgMrlf9wcUnDur4oeCsT17UpR2f5gUNDv9XmIBkTD0FE2eqJIq0AxAvlPtCJUWYVlRShvk0bZeCYqmnw/2LXIeR94SkFt6UbjbUp4kzlSBfMa2oLMooSkKPYMkwUULwDR5tONPVP+DRhve/s06YVlRgQaYQsOQ6kLxGRsaOB9nAX7WJRSXIV4SoBLojRCXQHdOKSpGb8LuuRJGbch1KXhNmF+9KdxHGuNk0pu2nQoqRsoiaWP2hSAn8GDuNzbSZSlZLccfPRVb1WRfvu4pDK2Rv9SQcmnETfU0rKkkrxJU8CcnAk2VGHHgYzVE4Bvl+NBtMKypB/iJEJdAdISqB7phWVJoUIGp7Ac2gMUJmJUGY7bxNgsx1p4YC03YpqHIrYcf9uQ4j74lJfj6TnjLUp2kzFZodi1ILmii61heyZsOnVSNrPSsXDplPwzzpjEWtoTh6Oxa1Z81NwTd4qOAI7TI8VBjm07SiEuQvQlQC3RGiEuiOiUWlohIBBl6K+fuIhkaKGBrfqlz+t8K0XQqKZRvtnvm5DiPvCUpNvCxdY6hPE2cqQb5iWlFZ1BqKIneKLoV+cGuVHKkuxK0ZtyywaUWFZseqjhSdn/1gwYqXYVgMvNMxr6gEeYsQlUB3hKgEumNaUalyMwHnTahy78uQCSBCGx9KDxGhzTCfpu2n0qQICet/cx1G3pOSYuziU0N9mjZTSWoRrsQcJLUo16HkNXbNS602Hbvm7X9nnTCtqGStBHfiLGTNuAKpZsSJjz212Tjx9b+zTphWVIL8RYhKoDs5uVG/4ZMTcYQHV/+81lnIn8bCDZ/8kC2xzAtyD5Sy1vwb6fBK6166tFNtKeWIQvivfxSNyuAmlP5mgPuZNlOFlSRvdzYSVkQh2b6Iagk+TmwhqiUM82naLoWWZIRbG0SXQn+0q0EeDb1sqE/TZiqrJFFqdWKVpFyHktdYkCmU3FhEecb+Genw8be9jmekw7hHZTMyzFLMNcXzGWYpNsynaUUlyF+EqAS6I0Ql0B0hKoHumLZLYWuskx9teAZFy7+Oy3yiSWnjqvYHUQ2cymZaUWlASgiqXzRAMXhupGkvf1V2N38cPYUquzvXoeQ1ZXIhF3h/SJksCsn2i1O2sp+7DKds2mRrCA7JylhbFQ5JzKYRmJhvLapNmzbh8Xh45JFHdAxH8F3gW4kqmUwyf/58wmHj6kgKzMO3EtWiRYvw+XL7zm1XMso9jWvZlYzmNI58p0MN8WT4dTrUkGE+s757e+ONN7jvvvtYt24dI0eOHIqYBkRQSfCy39g1V8xIRIvzftzYNXyyEpXf7+ess87i7rvvZsSIEf3uX1tbm/ZzPB4HoNruwef85hE3pCTZmYxgk2RGOHrO+vh6ZGeV3YNT7lqJtEC2McFdxqv+epqTYXwWO2U2V9pxUTXFjkQYGRjt7PlIvT0WQEGjvMiFy5F+KvyhOIFwggKnlbLC9HaTKZUdbV2X/hGVXnYffLOjLUwypVLic+JxpRdwDYQT+ENxHHYLlcUFaTZF1Wjc1ZVRhlmKewxXaVUCxEnikwrwyukxhbU4fjWEFQuVlqLu7U7JzlhbFW/FPiGixamQi7BJ6au5tqtBoloCj+SiUE6PKaolaFeDyD2+Ze9kJaqLLrqII444gnnz5mVzWA8uqzmI8WPHd//8mr+eOxo/pNTm4k9jj+mx/ymfrgDg0upJ7F2QPnsmpak80bqJIwuruWD4AWm2taEW/rD9HRyyNWO7P/n8eQJKgh//YG8m7pleaHXZy5/z0vvb2Wd0KRfNmZhm294c4IaH3gHg6rMPw2ZN/+Vfd/9bNLWGOXHKWKZOTK9K8/yaLTz12heMHuZj4ZmHpNnaAzGuuOd1AH7qPYEiOX34718Dz7IltYMjnPsy3XVgmu392Oc8GXmDUouPBYWn9fiunya2EVHinOGZTrW1LM32j+BLrE9u5UDHOE4sODzNtiGxnUdCq3BJjh5t9oakadqASqz94x//4KqrrmL9+vUUF3eNzZEkiYcffphzzjlnQM7q6uqYPHkyh96+EN/YbzLdt8lU1XYPl484mGu3vsknkbZBZap9Xi/Ju0w1/pcBXTJVuVzEPO+x3Nm5nEalbVCZauW0OxkIA85UDz30EC0tLT0uexdeeCHLli3jhRdeGGhTNCZCtGaYrJDU1D4nMTQlet5sRtQUAAElQUDJPA5bhT7b3eXv/WY/EkvxVSzYq72+pXdbeyBGeyCW0RZPKHzVx7HNSkevtoAWIaBEMtpSKDQqvU9x36n6e7WFtCghJfO5ULMo7zhgUT366KNEo+kO99hjD66//nrmzxdlEgXfMGBRVVdXZ9xeUVHRq20oiakpPgm3EvtfphJkJq6l+DLZRFwz7jyZ9sVZUyLMtdveynUYeU+r2sl9wZWG+hyUqAZ4jz8kSIBFklG0bK723z8kQEZGzequaHCY9oXyGGchT44/iTEZnuoE31BlKeXmkvOoshi31rRpRSXIX4SoBLojRCXQHSEqge6Ytkvhq3iAn218kU4lnutQ8ppmpYMbOx4jpBk3RMi0okppGm2pzK9ABN+goNKpGTuY0rSXv0pbAQtrDqbSVtD/zt9jSmQvZ3p+QIksCsn2i9ti48jCatwW4xacNiMuyc7+9lpcknFr+JhWVIL8RYhKoDs5uVH/04QnGD9hcOlYU8ZB+Bj+tP+TSJbNg2prVnDBoI4fCkba9HmqdcldNVFd1iQe2ZgnZfNmKqkNHA91fQp6JaiFWR1fQ9DAJ0DTdilIcgc4luU6jLwnrEV4O/GhoT5Nm6k0zY2WPAxNEwU6+sKBnT0tY3Agnv76Rx0O0T90fQp6pVguZG7BDykWVV8EZkaISqA7QlQC3TGxqBIgb+/6FPRKihQ7lTZSiNk0/SJZvgLPz3MdRt7Tqnbw18g/DfVp4kwlyFdMKypNqUULPIWm1Pa/8/eYSrmM33ouoFIu639nnTCtqLpCd2Pqr2AAEhIOyY6URSmgwSJ+IwLdEaIS6I4QlUB3zCsquR7cv+z6FPRKq9rB/eHHaVV7r3elN+btp5LiMMjBed8HUqRoVncZ6tO0mUpTy9Giv0RTy3MdSl7jkzyc4DgKn+Tpf2edMK2o0AoheVLXp6BXCiQXB9v3p0By9b+zTphXVIK8RYhKoDtCVALdMa+oJD/Yl3d9CnolrEV4N7GWsJa5RPZQYN4uBbkVnPfnOoy8J6iF+U/c2IK7ps1UmuZES+2DpjlzHUpeY8NGjTwMG8bVnDCtqFBrIHJH16egV0rlIn7qPp1Sucgwn+YVlSBvEaIS6I4QlUB3TCwq5X/dCUquA8lrVFTCahQV1TCf5u1SsGwF749zHUbes1Nt4/bwg4b6NHGmEuQrphWVpoxCCz6MpozKdSh5TblcwsXusyiXS/rfWSdMKyqwgVbV9SnoFQsWSuQiLFj631knTCwqQb4iRCXQnZw8/f35vBkU2isG1Ub12FJ+dRfc8bPjaPxy0qDa2ruk54LfuWbYQwFd2imla/GCUnsY0KfN/jBtpmptCvDA716ktcmYE2VWArTzfGopAdoN82nafqp4NMmmusZch5H3JEnQoG0x1KdpM5W32MVx8w7EW2zcgH4z4sLDZHkqLsRsmn7xlRQwY/4kfCViwaO+KMDDZMs0CoSoBGZGiEqgO0JUAt0xragioTh1r24mEhLL3fZFgihfqOtJIJa77ZeOlhBLb3s912HkPUE6eVV5xlCfps1UVpuF0uFerDbjXpSaEQsWfBSLF8oDoXJkEVc++P9RObIo16HkNUWUM9f2C4owrjqOaUUlyF+EqAS6861EtWTJEsaPH4/T6WTffffliSee0DsugYnJWlSPPvooP/vZz7j44ov59NNPmTdvHnPnzuWdd94ZivgEJkTSNE0b6M6aplFbW8uPfvQjbr311u7tM2fO5Oijj+aqq67q8/i6ujomT57M4cPmDXo8lZ6oJd5ch9CDvR76Itch9OCuA5cOaL+s+qk2bdrEtm3bmDdvXtr2VatWZdOM4DtOVqLauHEjAOFwmJkzZ7J27VrGjBnDtddey4knnthj/9ra9HVj4vGu3u+KmkLKvKXd2yOhOB0tIaw2S8YugsYvu1Z0L68uxO7sCrm4wsOM+ZN44s43qf+iFbfPSVF5+nrK8WiS1qYAkixRNabnbJId29pRFY3SCjdOV/oEis6OKKFAHFeBjZLd2k0mFXY2BQGoHtUz3pamIKmkQnFpAQWe9LWLg4E4gY4oDqeVssr0kQOKotHc0Nn1/ShH3q1vKUA7SRIU4MXFbt+VCCECWLBSxDfr0HjxMdkyjZeVp+iknSJKsew2WSSInwQxXLgpID1rJ4gRxI+UxZ1SVqIKBLpGWf7kJz9h0aJF3HLLLSxfvpyTTz6Zl156iWOPPXZA7cxbeDTjx4/v/rnu1c0sve11CssK+NVdp/TYf+HsvwHw419PY9Q+6ZfNPSdVU/9FKwdMHcOpvzgizbbxwwYevG4Vdoc1Y7u/P+MxwoEYJ8+fyL6TqtNsTz+2jjde3MQe+1Vy9iXp7TZs6+DPv3sJgEsXHdujA/b/v/JFWhoDHHfKeA49Ov0Pa/Wzn/H8v9ZTM7qYX1xzTJrN3x7hhgUrATjBOheP5EuzP5v6Bzu0r9hXnsyBliPTbJ+r63hDeQ4fRZxm+1mP72pVusQ93XoKZdKwNNtLqeVs1T5nnLwvh1uOS7NtVzexSnkCBwMv2ZTVPdWyZcuYO3cuixcv5qKLLurePnv2bACee+65Po//+p7q5IMuosxb1b3922Sq8ppC5l9xDPde+Rxb1jcPKlMV7z0s7zLVYQ/5dclURZRyrO0Ulif/RhvNg8pUdx74WI/vmomsMlVNTVctqAkTJqRt33fffVm5cuWA29nZ0Enc7uixPZVUugWUiV2NnT22xSNJAMKBGOFALONxmqr12W7bznCvtmgkSeN2f6/2vmwdbRE62jKXRYzHUn0fS+8F9SMEiRDMaFNI0UZzr8f66f08RAkTJfO50LKoxZBVl8KkSZPwer28++67advXr1/PuHHjsmlK8B0mq0zlcrm44ooruP7666muruaQQw7h8ccf5z//+Q+rV68eqhgz0t4SZOltr9HekvkvVtBFCD+vpJ4mhN8wn1kPfbn22mspKCjgmmuuobGxkX322YennnqKo48+egjC651oKEHdq18a6tOMxImxWfvEUJ/f6jXNZZddxpYtW4jH46xbt46TTz5Z77j6xe1zcsTsfXD7RCHZvnBSwHh5Mk6MmyBi2hfKReVuTv3FET2e+ATpuPExxXI8bnz976wTphWVIH8RohLojhCVQHdMK6p4NMnGDxuIR5O5DiWvSRKnXv2SJMbNOjLtbJrWpgAPXidGR/RHgA5eUB431KdpM5UkSzhcNiRZynUoeY2EhA07EsadJ9OKqmpMCX988icZXxQLvqGESs61LaSESsN8mlZUgvxFiEqgO0JUAt3JydNfsqaYuHdwEx8SNYX/+ywhPsgp3UfftWZQxw8F78zQZ9GByLhCWAz1Cwpp3DzI7pemge1m2i6FpqYAV1+9ikhE9FP1xY6t7fz+9EeIhhKG+TStqFRVI2TgiTIrqqIS7sw8InaoMO09VVlZAeeffzBlZaLmZ1+UDvdxzh+Op3S4GKXQL06njQkThuF0irVp+sLptrPv4aNxuu3976wTphWVIH8RohLojhCVQHdMK6rOzij//vendHYaVyDVjHS2hnn2vjV0tvY+t1FvTNulEAwmePVVY9dcMSMhf5Q3ln9sqE/TZiqXy8bEicNxucTTX1+4PHb2n1qLyyOe/vqltLSAn/70IEpLRT9VX5QM83HW72ZQMkz0UwlMjBCVQHeEqAS6Y1pRJZMK9fWdJJNKrkPJa5LxFA1f7CIZTxnm07RdCi0tIW699Y1ch5H37Kz3c+cvlxvq07SZSpC/mFZUNTU+/vSnWdTUGPeobEaqxpZy88rzqRpb2v/OOmFaUYGE1WoBA+ezmRFJkrDaLUiSmPcnMDFCVALdEaIS6I6JuxSC3HTTq7T1UlJa0EXLVx3cdv4y2nYEDPNpWlElkyrNzaFch5H3pBIKLds7DPVp2stfcbGLM87Yn+JiV65DyWuKKjz86NdHUVTh6X9nnTCtqNxuO4cfPgq3gbNEzIjb5+TQE4yt4mxaUQnyFyEqge4IUQl0x7SiCgbjvPTSFwSDxhVINSPBjiivPF5HsMO4WUem7VLo7Izx7LOf5zqMvCfQFuaFh9431KdpM5XDYWHcuFIcjsHVpvqu43DZqN2/CoeBs45MK6rycg+XXnoE5eXG9b+YkbLqQi667STKqgsN82laUQnyFyEqge4IUQl0JydPf7bP67FbBveIK/tL8e88EPmTbdi39L7Y9EB4Z+boQR0/FDSeXqtLO3Kxh7ZQjB3Tq2k80Jj7KtN2KTRvaeOPcx7MdRh5T31HiEuWvmmoT3H5E+iOaUU1rLaUa586j2G1xs0SMSMjij3cfcZURhSLoS/9YrFaKKrwYrGKzs++sMgSpR4nFgNXGzOtqAT5ixCVQHeEqAS6Y1pRtdZ3cO8lT9Bab+ygfrPRHIhww8oPaA4YN+vItP1U8WiSL9c25DqMvCeWVPhsh5hNMyB8ZW5OuOBIfGXuXIeS1xQXOPjxweMoLnAY5tO0ovKWuDn2rEPwlghR9UWhy87JE8dQ6BLViQUmRohKoDtZiyqVSnHdddcxatQovF4v06ZN49133x2K2AQmJWtR/fGPf+SBBx7ggQceYO3atey1114cf/zx7NixYyji65VIZ5T3nl1PRKxN0yeheJJXP28kFDduWeCsRbVixQrmzZvHjBkzGDduHLfffjudnZ288847QxFfr3S0BHnilpfpaAka6tdstIZiPPDmBlpDxi15m7WoKioqWLlyJdu2bUNRFO6//34cDgcHHHDAUMTXK1a7hcoxpVjt4oVyX9gsMtXFbmwW426fs+78vPPOOzn99NMZM2YMFosFi8XC8uXLGTt2bI99a2vTRy/G410TP8tHllDmrujeHg3GaN8R6BLK6J5DWRo37ew6bkQx9v9NNaoYWcz8389i8SVPsGVtA+4iF0UV3nR/kQStDX4kWaJqXHmPdnd82YqqqJQO9/VYDrazNUzIH8XlsfdY1yUZT7Gz3g90FWrdvZ5my1cdpBIKRRWeHoUxgh0RAm0RHC5bjxkuSkqleVs70DVkZfeRBc2BCLGkQnGBo0cXQSiepDUU6xJR0TfdLFVFBVw8fX+ufupdtrUFqSpyY99NYLuCUcKJFD6XnZLd+rMiiRQ7g1EsWdQMzVpUGzZsoKioiBUrVlBdXc0DDzzA/Pnzef3115k4ceKA2pi/6ATGjx/f/fOHqz5j6Q0vUlju4dcPze+x/2+m/BmAH18zg9H7VaXZxkyoYsvaBg6YvidzLpueZtv43jYeuPzf2J22jO0u+uFfCfujnHjB4ex7+Og027P3reGNp9azx4E1nHXtcWm2hi92cefFTwFwyR2n9siWt/38X7Rs7+AH8yZx6An7pNleeXwtLzz8PtV7lHHRrSel2fy7Qtx45mMAXHH8gZR60gV5w8oP+GxHBzP2HcHJE8ek2V79vJEH3txAhc/FTXMO6/Fdv+aXx+zHmLL0P5I7X/6I97buZMq4YZx52F5ptg+37+L2/6yjwDFwqUiapmkD3bm+vp5x48bx8ssvM3Xq1O7tU6dOpbS0lBUrVvR5fF1dHZMnT+ak8edQ5h7WvT3Xmaps3LC8y1Ty+fvnXaZae8evdj+FGckqU7333nskEgkOPvjgtO2HHXYYzz///IDb2fVVO/EM1/hUQukWUMbjMrw8jocTAIT9UcL+zE+Cmqr12W5fpQujoQSNm1t7tTd92fukC//OEP6dmav9xaPJPtut7+i9SmBHJE5HJHMNiaSisq2t94eXJn/vK5QGogkC0URGmzLw3JPdjXpNTQ0AH3+cvtLlxx9/zJ577plNU4NG0zRSiRRZJNrvJRpdQjPyLGWVqQ455BCmTJnC2WefzeLFi6mpqWHJkiWsXr2at99+e6hizEjTF7u4cvrdhvo0I9vbgpz90GpDfWaVqWRZ5plnnmH69Omcc845TJ48mVdeeYXVq1dz6KGHDlWMApORdedFcXExf/nLX9i+fTuBQIC3336bo446aihi65OKUSX86m/zqBhVYrhvM1FV5ObGUw+lqsi40RymfaFsc1ip2asSWxaPut9H7BaZMWW+Hk98Q4lpRSXIX4SoBLojRCXQHdOKqn1HJ0t+t5L2HZ25DiWv2RWMcufLH7ErKArJ9ks0GOfjV7/IdRh5TziR4r2tvb9NGApMm6k8xQVM+/EkPMUFuQ4lr/G57MyaMBKfmPjQP4XlHk665CgKRSHZPikpcHDmYXv1eFE8lJhWVIL8RYhKoDtCVALdycnTX/txY4mV1QyukUIXdVt20njUSHZO7DkALxs+uP7ewcUyBMw6+jRd2vFWhvl0RAPe1c1UtRizkqtpuxR2dka5Y+W6XIeR97S1hHj4ltcM9Wnay59FlvC6bIaWHTQjskXC7XMgW0R5xn6pKfXwl/OPoaZUdCn0xfCRxfzhodMZPrLYMJ+mFZUgfxGiEuiOEJVAd4SoBLpj2i6Fr1qD/Pyvq4knlVyHktc0be/gmrOWkYinDPNpWlFpGsQSQlD9oaka8ahxZYTAxJe/ysICFp48icpCMfSlL8qGeTn/2umUDfP2v7NOmFZUTruFCaPKcIpSQn3icNnYa6JYmFtgcoSoBLojRCXQHdOKqj0U4++vfUa7gbUszYi/LcxTD76Pv633EkJ6Y9ouhWA0yeqP63MdRt4TDsRZ8+ImQ32aNlO5HVaO2Gs4blFLoU9cHjuTpo7B5RGzafqlzOfiwpkTKPO5ch1KXlNS7mHegiMpMXDWkWlFJchfhKgEuiNEJdAd04oqnlTYvMMvRin0QyKeYtvGXWKUwkBo9ke4/on3cx1G3rOrKcA916wy1KdpM5UgfzGtqEaVe1ly6QxGlRs3pMOMVI8p4bYnz6R6jHEFd00rKkH+IkQl0B0hKoHuCFEJdMe0XQpN7WF+8/c36QhlXk1K0EVLg5+bL15BZ1vEMJ+mFVVSUdkpFuXul1RSpa3ZmBJCX2Pay1+Zz8UFM/YToxT6oaTCzRmXHklJhVibpl/cDitH7l0lxlP1g8vtYPK0MbjcopCswMQIUQl0x9BrRzTadWMd87cMui2/UsCGDYX4m7cT6Rjck03dx/k3eSIQa9alHUdnlA0bNtDa2UggNrglV+rq6th7770pKOhnVrhmII8++qhG17K+4p9J/3344Yf9/p4lTTNuZevW1lZWrVrF6NGjcbkG99R24oknAvDss8/qEdp3Fr3P00AylaGi0pPa2loAtmzZkuNI8ptcnCdxoy7QHSEqge6Y9vInyF9EphLojhCVQHeEqAS6I0Ql0B3TiUpVVRYtWkR1dTVut5tZs2axdevWXIeVd7S3t3PhhRdSU1ODz+djypQpvPXWW4b4Np2obrjhBhYvXsz999/PmjVrUBSFmTNnkkgkch1aXjF37lzWrFnD0qVL+eCDD5g4cSIzZsxg48aNQ+98SF/26Uw8Hte8Xq+2ePHi7m0dHR2ay+XS/vnPf+Ywsvziiy++0ADtrbfe6t6mqqo2duxY7Xe/+92Q+zdVplq3bh3BYJBjjz22e1tRURGTJk3ijTfeyGFk+UVZWRnPPfccBx10UPc2SZKQJImOjo4h92+qYZMNDQ0AjBgxIm17VVUV9fWiVOPXFBUVMWvWrLRty5cvZ/PmzRx//PFD7t9UmSoS6Ro35XCkD411Op3EYvk3JipfWLNmDeeeey5z5sxh9uzZQ+7PVKL6erhMPJ4+LSsWi+F2Gzew30w8/fTTHHfccRx22GE89thjhvg0lai+vuw1NTWlbW9qaqK6ujoXIeU199xzD3PmzOHEE09k5cqVOJ1OYxwP+aOAjsRiMc3n82kPPvhg97avn/6WLl2aw8jyj8WLF2uAdumll2qqqhrq21Si0jRNu/rqq7XS0lLt6aef1j766CNtxowZ2p577qklEolch5Y3bNy4UbPZbNqpp56q7dixI+2f3+8fcv+mE1UqldKuuOIKrby8XPN4PNqsWbO0rVu35jqsvOLGG2/sdYz52WefPeT+xXgqge6Y6kZdYA6EqAS6I0Ql0B0hKoHuCFEJdEeISqA7QlQC3RGiEuiOEJVAd4SoBLojRCXQHSEqge78PyJj3zR5DEZqAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(proxy)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Estimating graph using greedy search\n",
      "Iter    TT      TP      PP\n",
      "\n",
      "001  10.000   0.000   0.000\n",
      "002  10.000  13.333   0.000\n",
      "003  10.000  13.333   0.000\n",
      "004  10.000  13.333   0.000\n",
      "005  10.000  13.333   0.000\n",
      "006  10.000  13.333   0.000\n",
      "007  10.000  13.333   0.000\n",
      "008  10.000  13.333   0.000\n",
      "009  10.000  13.333   0.000\n",
      "010  10.000  13.333   0.000\n",
      "011  10.000  13.333   0.000\n",
      "012  10.000  13.333   0.000\n",
      "013  10.000  13.333   0.000\n",
      "014  10.000  13.333  33.333\n",
      "Running GraphEM:\n",
      "\n",
      "Iter     dXmis     rdXmis\n",
      "\n",
      "001     0.0544     0.1120\n",
      "002     0.0742     0.1516\n",
      "003     0.0347     0.0694\n",
      "004     0.0175     0.0345\n",
      "005     0.0089     0.0173\n",
      "006     0.0045     0.0087\n",
      "007     0.0023     0.0045\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/fzhu/Github/GraphEM/GraphEM/GraphEstimation.py:270: RuntimeWarning: invalid value encountered in reciprocal\n",
      "  dinv = np.sqrt(d)**(-1)\n"
     ]
    }
   ],
   "source": [
    "#Fit test\n",
    "G.fit(temp,proxy,calib,sp_TT=3,sp_TP=3,sp_PP=3,graph_method = 'glasso')\n",
    "temp_r = G.temp_r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Computing verification statistics\n",
      "\n",
      "Mean MSE = 0.0832, Mean RE = -0.0908, Mean CE = -1.0142, Mean R2 = nan\n",
      "temp:\n",
      " [[0.3745 0.0206 0.6119 0.6075 0.122 ]\n",
      " [0.9507 0.9699 0.1395 0.1705 0.4952]\n",
      " [0.732  0.8324 0.2921 0.0651 0.0344]\n",
      " [0.5987 0.2123 0.3664 0.9489 0.9093]\n",
      " [   nan    nan    nan    nan    nan]\n",
      " [   nan    nan    nan    nan    nan]\n",
      " [   nan    nan    nan    nan    nan]\n",
      " [   nan    nan    nan    nan    nan]\n",
      " [   nan    nan    nan    nan    nan]\n",
      " [   nan    nan    nan    nan    nan]]\n",
      "temp_r [[0.3745 0.0206 0.6119 0.6075 0.122 ]\n",
      " [0.9507 0.9699 0.1395 0.1705 0.4952]\n",
      " [0.732  0.8324 0.2921 0.0651 0.0344]\n",
      " [0.5987 0.2123 0.3664 0.9489 0.9093]\n",
      " [0.6632 0.2401 0.353  0.6908 0.3905]\n",
      " [0.6631 0.217  0.353  0.7117 0.3905]\n",
      " [0.663  0.173  0.3531 0.7514 0.3905]\n",
      " [0.6652 0.8992 0.3516 0.0956 0.3905]\n",
      " [0.6643 0.603  0.3522 0.3631 0.3905]\n",
      " [0.6651 0.8556 0.3517 0.1349 0.3905]]\n",
      "proxy:\n",
      " [[0.9696 0.3887 0.7722]\n",
      " [0.7751 0.2713 0.1987]\n",
      " [0.9395 0.8287 0.0055]\n",
      " [0.8948 0.3568 0.8155]\n",
      " [0.5979 0.2809 0.7069]\n",
      " [0.9219 0.5427 0.729 ]\n",
      " [0.0885 0.1409 0.7713]\n",
      " [0.196  0.8022 0.074 ]\n",
      " [0.0452 0.0746 0.3585]\n",
      " [0.3253 0.9869 0.1159]]\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/fzhu/Apps/miniconda3/envs/LMRt/lib/python3.8/site-packages/numpy/lib/function_base.py:2642: RuntimeWarning: invalid value encountered in true_divide\n",
      "  c /= stddev[:, None]\n",
      "/Users/fzhu/Apps/miniconda3/envs/LMRt/lib/python3.8/site-packages/numpy/lib/function_base.py:2643: RuntimeWarning: invalid value encountered in true_divide\n",
      "  c /= stddev[None, :]\n"
     ]
    }
   ],
   "source": [
    "print(\"Computing verification statistics\\n\")\n",
    "V = GraphEM.solver.verif_stats(temp_r,target,calib)\n",
    "print(V)\n",
    "print(\"temp:\\n\",temp)\n",
    "print(\"temp_r\",temp_r)\n",
    "print(\"proxy:\\n\",proxy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAN8AAAGdCAYAAACSHqb/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAxA0lEQVR4nO3deXwb5Z0/8M/MSBrJliX5jC/Zju0E5+BKKAUKgSWvDZAssITSZpNybdkWdtOFbhfYtLT5NT1ZyrZdaDgXKEcphXAVaAOFLWkSoCUh5QpJiO3E8hVfki3rGGlmfn94bTKRLI2wrWee5Pt+vfLyK/NIz3z1lb965HlmnhF0XddBCMk7kXUAhByrqPgIYYSKjxBGqPgIYYSKjxBGqPgIYYSKjxBGqPgIYSSvxReJRLBz505EIpF87pYQS7Llc2cfffQRFi9ejKbV/wbXrNop9eVzObC0pRavfhRAMKpMqa/Kn22f0vNnQsdTC6elH5/kxDLfbLwcbENQjU2pL10XpiWm6TT7m6PT0k+Rz4XTl83BGy/vw0gwOqW+fr/nVlOP4/Zrp64DUSUJOjkuMx06ImoSOihRmei6jlg0gXyebZnXkW86hWIKXvrgIOswLC+kxvFCcB/rMCwvHIphy29353Wf3I58hQ4bPttQgUIHt58feVEo2nG6uwaFop11KJbmKnTghNPr4Cp05G2f3BZfeZELN5x7IsqLXKxDsbQKeyG+UX0aKuyFrEOxtJIKN67497NRUuHO2z65LT5CeEfFRwgjVHyEMMJt8SlJDW39w1CSGutQLE3RVbTGhqDoKutQLC2hJBFoHUBCSeZtn9weKuwKjeKbz7/FOgzL61RGcNPB11iHYXmHOofxsxtfyus+uR35COEdt8XXUFKEh69cioaSItahWNps2YvHm/8es2Uv61AsrXp2MX7869Wonl2ct31yW3wQALskAtY73dBiBNhFCZSozAQIsNklCHnMU87Fp2ka1q9fj5qaGhQWFmL58uVoa2ubidgIOarlXHzf+973sHHjRtx7773Yvn07VFXFeeedB0WZ2pUFhBxrcio+RVFw++23Y8OGDVixYgVOPPFEPPHEEwgEAti0adNMxUjIUSmn4tu1axdGRkawdOnSiW0+nw+LFi3Cli1bpj24TDqDo7jx6e3oDE7P9VxHq05lGDe0v4xOZZh1KJbW2xnCbTc8j97OUN72mdM8XyAQAAD4/X7D9urqanR0dKQ8vrGx0fD/eDw+9nhfAbylnxylHI0n0BeOwS6JqPGlngDcPjACAKjyFEC2S4Y2hyQioWooctpRWug0tMUSKnqGIxAEoD7NUdGOwTBUXUdpbTFcRbKhLdg7gvDgKFxFTpTW+gxtSiyJQ239AICalkoIR/yN3tvaj0Q8ieIqLwp9xhO/RwZGETo0ArnAgfL6EkObmtTQve8QAKDO4YEkGD8bu5UwYnoSJTYnvJLxtY6oCvqTETgEETUOzxGvVIeij52MUOMogkMw5vBQYhSjWgJeSUaJzRhvREugNzEKCQL8sjEPAHAgHoQGoMpeCOcRV04MJCMYVhW4RTvKjzixO64l0ZUIAwBmp+k3oAwjoWuosBegUDReaTCUjCGoxuASbaiZbcxhMqmit2OsgKrqfRBFYw77uoehxJLwlhTA7TXmcHQkhqSiwu6QUFFjPDqs6zq62ocAALNqvbAd8Xs40DuCWCSBIp+xz0xyKr7x5R9k2fiL6nQ6MTg4aLqftWcfj/nz50/8f+vH3fjFlvdRUiDjRxeflvL4f3jgFQDAtUsWYG6Fz9D2yzc/wu8/7MBps2fhH0+fZ2j7a6AfP375HThtUtp+v/KrP2IklsDf33QeFp5znKHt2ds24/VH3sTc0xpx1e2XGdoCu7tx+xfvBQDc8OiXYTvisqZbL9mInv19WPbVJTht5SJD2x/u34oX//tV1M6vwtoHrjK0BXuH8d2//SkA4Fs1n0OpvcDQvr7jdXwQ7cf5viasLGkx9htqw929O1FhL8Rt9UsNbZqu45/bfo/+ZATXV34GjU7j4fTbu97EG+FOnFVUh6sqTjC0/SXchVu73kChZMdPGs7Fkb6073lEtSSumXUSTiqcZWi7r3cXfh9sxSJ3Ja6v+oyhbU90AN88+DoApO33X1o3oycxilVl83G2p87Q9kT/bvxmYDeOc5bg6z8509DW3z2MH699DgDw1fV/m1Jgd6z7HQ7s7ceSi+bh7AvnG9p6O4K4/wevocAt4+s/WWFoi0UU3HL5EwCAy/99CSr9PkP7Az/6X3z4dgCf+ZvmlNcyGSGXG6Vs2rQJn//85xGJROByffIJ+YUvfAHxeBzPPfdcxufv3LkTixcvxln/fAu81fUT2z/NyFftLcDXzjkBG176C3b3BKc08i14arflRj7h+TOmZeSrcRThhqpTceOBV9EWD05t5HMY8wCwH/kW32N8LZ925Kuo8WDN18/CT//9RRzqDE1p5Hvyre+mvJ50chr5xr9udnV1oampaWJ7V1cXTjjhhMmelqIrGMGQPJKyPaFqE4WWTvdw6sJL0cTYOYsjsQRGYom0z9N1ZOx3IDA0aVt0JIbA7p5J2zs/mrxtqDuEoe70f0PEI0rGfg9m+BttMBnDYDL9eiyKrqEtHpz0uZ3K5HkIqXGE1HjaNhV6xn67E5P/7R3WEghneG6mfg8lIgDSL7gV1ZLobJv8b7TuA5P3GxqMIDSYvt+EoqKzbfJvcr2Byfc5EjS/Tk5OB1xOPPFEeDwe/PGPf5zYFgwGsXPnTixZsiSXrgg55uU08smyjLVr1+Lmm29GeXk5GhoacOONN8Lv9+PSSy+dqRgJOSrlfFXDhg0bkEwmcc011yAajWLJkiXYvHkz7Pb8rhESiip47q9tCE1x2cCjXUiN4enBjxCa4rKBR7uRUBSvPv0+RkJTWzYwFzkdcJmq8QMu07Fu53Q6mtftnE5H87qd0+moX7fTaZMwr7IYTpuU/cHHMKdgwwJXGZwCt5du5oXstKFpwSzIzvzlidviq/QW4DvLT0GltyD7g49hVQ43vus/G1WO/K3KxaOyKg+u27AMZVVHnqAwc7gtPkJ4R8VHCCNUfIQwwm3xqZqOgdEYVI1uAJKJqmsYSESg6rTKWyaqqiE4MApVzV+euD0E1jEUxton/sQ6DMs7qAzjq22/Yx2G5fUcDOL7X3k6r/vkduQjhHfcFp+/2I07v3gW/MV0CD2TOocH98y+AHUp1/iRw1XW+XDLvStRWefL2z65LT5JFFBa6IQkWu+sCyuRBBGl9oKUy5OIkSSJ8JUWQpLylyd6RwhhhMkBlweueBQnHT+1mxAqySb0DJ+GBy9/FA7b/in15b7J/KX/+bL8b4/L/iATqhuSwI+A6ruT0NvTX+9ollaQvxtHmrXnuopp6SdZNHax9YHLirF/JD9lQSMfIYxwW3w2qQsVRTfBJnWxDsXS+ntGcNeG19DfM/kV7AToikSwbsef0RVJf3X7TOB2nk8UonDa32MdhuXFY0ns332IdRiWF1VVvD80+XIiM4HbkS+plSIYuQpJrZR1KJbmKXbhglUnwFNM967PpESWcUXTHJQcsTLfTOK2+DTNh+HYF6FpPtahWFqR14mlF89Hkdd6B5WspNjhwGWzG1HsyN9BJW6LjxDeUfERwggVHyGMcFt8ojCCQvn3EAU6hJ5JJBzHW6/tRyScfjFcMmY4kcDLnQEMJ6Z2IkIuuJ1qsEmHUFr4c9ZhWN5QfwRP3vcX1mFYXl8shjt2f5DXfXI78mm6A0qyDppuvVOerMRmlzCr1pNybwFi5BBF1BUWwiHSidVZJVU/eobvQVL1Z3/wMWxWjQc33rYcs2rokqJM/IWF+MXpZ8JfmHqjnpnCbfERwjsqPkIYoeIjhBHOiy9/h4V5pUNHMqFCB63ylokOIKFpec0St1MNDtt+1JVcxDoMy+tqD+I/rniSdRiW1zoygpWvvZLXfXI+8hHCL26LL6H60R26AwmaasiootqDG364DBXVNNWQSW1BIX526umoLaCphqx03YGE2gydJtkzsjsk1M4ugd1Bk+yZyJKIJo8HMq1eRsjRj4qPEEao+AhhhNvis4m9KHP/ADaxl3UoljbYF8bDP9uGwb4w61AsrScaxY/f3YWeaDRv++R2nk8UwyhwbGUdhuVFRxN4960O1mFY3mgyiW2H8vtBzu3Ip2o+DEcvgUoLKGXk9spYsvw4uL35W5WLRz6HAxfX1cNHCyhlp2qlCEa/ApWWDszIW1yAiy4/Gd7iAtahWFqpLOOauS0opaUDCTn6UfERwggVHyGMcFt8ghCBy/4mBCF/N7bgUSySwAc7OhGL0OVXmYwmk3ir7xBGk8m87ZPbqQa71I3you+yDsPyBg6F8eBP/sQ6DMvriUbx/b++k9d9cjvy6boEVfNC1+mE4UxESUBhkQxRottnZyIJAjx2OyQhf3nitvgSagM6g79GQm1gHYqlVfl9+O69l6DK72MdiqU1uN147Oxz0eB2522f3BYfIbyj4iOEESYHXH7ccxpKfGVT6mOWWIGrXcCPek9Hr9Y0pb7OLf5oSs+fCUMnFk9LP0WzigAAw3OLMOSdWl/FL+2ehoim15ybRqeln5qWSuBXZ8D/890QP+qZWmc3mHsYjXyEMMLtVMMhrQ//NfoLJGj5wIw6Do3gX/7rD4gnVNahWFrX3l6sO+s2KFG6S1FWOnQoUFiHYXm6DsQUKrxsdE1HfDS/v0/cfu0sFnz4gnwJigUf61AsraK4AF//wmJU0FUNGZX5i/GVX6xCmX96/tY2g9vicwgONNoa4BBo9bJMnA4JCxvL4KTVyzKSC2W0nN4EuZAuKSLkqEfFRwgjVHyEMMJt8Q1rI3g5/hqGtRHWoVja4HAMj778IQaHY6xDsbRgzzA2/fj3CPYM522f3E41RBHFzuRfWYdheeFoAv+7k1Yvy2Y0GMG23+zI6z65HfmckLFAaoETtCpXJoVOO05bUIVCp511KJZW4HFi8fKFKPA487ZPbovPK3pxofMCeMUpnrB4lCv1OvFPF56AUm/+fql4VFztw5rvX4zial/e9slt8RHCu5yLb3BwENdeey1qa2vh8Xhw5plnYutWWjmakFzlXHyrVq3C9u3b8fjjj+Ptt9/GSSedhGXLlmHPnj0zER8hR62ciu/jjz/GK6+8grvuugtnnXUW5s6dizvuuAPV1dV47LHHZirGtBJ6Ap1qFxI6XdWQSTyhYn9nkK5qyEKJKmh/NwAlmr+Tq3OaaigrK8OLL76IU045ZWKbIAgQBAFDQ0PTHlwmg/oQHok9kdd98qh3MIIfPvIW6zAsr+/AIP77ql/mdZ85FZ/P58Py5csN2zZt2oSPP/4Y559/fsrjGxsbDf+Px+MAgBKxBLPEiontMT2GkD4MCRLKxNR7L/Rqh8aeJxTDLhgPmYe0EGKIwwUXPGKRoU3RFQzpQQgQUCGWp/Tbp/VDgwZZL4Oku4zPFYJICiOQ9ALIujEmTUggJoxd7Vygpd4TPir0QBcScOglsOnGe3wnhBEkhCBEXYZTrzC06VARFbsAADXlbkiicSWtQ0MRxBQVPrcMT6HxhPLRaAIDwzHYbSKqSlPvK36wd+xkhKrSQthtxi88/cEoIvEkPAUO+IqMUzfReBJ9wSgkUUDN3NQcdn3cD13TUVbrg1xgfG+CvWGMhqJwFckoqTLeE16JJtHXMfaBna7f3vZBJBUVJVUeuI6IaXhgFCMDEcgFDpTUGt/zpJJEb2v/2GudUwHxiNs89x0YgBJNwFtRBHeJMU+jwQiCPcOwO22oaDCutKBrOrr2jt3FaNbsMthkY+kMBIYQC8dRlCb3k5nSJPv27dtx9dVXY+XKlVixYoXp510kX4D5rvkT/38/uRsvxH+PIsGNq11rUh7/49GfAgBWyMtQI1Ub2l5XtuGNxJ8xzzYXy+RzDW2tyXb8Jv4M7LCn7ffno3cjiijqkl9AsX6ioe2g9CR6pD/Ao7dgTvKrhrZR4SA+sP8AADA/eTNEGH/p3rP9P0SFblSrK1ChnWlo6xJ/h4DtWRTq9ZiX/IahTcEQdjn+AwBww2WLUXLEnNN//urP2HNwCOcursOK040fbFv+GsAvf/cByn0urL/6jJTX+t0Ht+Ng7wj+6cITUF9pLIS7ntmFt/f04rQFVfji0hZD2659h3DHpnfgkm34+v2pOfzW+RsRjyi45Ia/wXGn1hvanv7pa9j+zLuYd9psrP628cO5/f1u3PnPY99c0vX7o394EAOdIZz35dOxeNk8Q9vLD76Jlx98E/ULqvCV2y8xtPV3DOKHF98FALju7tVwFxuL4edXPYQD73bi7C99Fud86bMp+7199f8AAL7xqy8btsfCcXxzyU8AAFf+50pUNhk/MP7nht/ggy37cOpFxt+jTARd13XTjz7Mc889h9WrV+Nzn/scnn/+eTid2eeRdu7cicWLF2PNw1dhVkvlxPZPM/KVCsW4yLkcj0WfRIcWmNLId4Gv33Ij3wt3XTEtI19laSG+etGJE8U3lZFv3t6+lNfKfOQrMbZ92pGvoqEUl//w73H76v/Bofb+KY18m7rvT3k96Xyqke/OO+/E9ddfj8suuwwPP/wwHDne02xQGwS01GM9KtSJQkv7PH0IGP+o+L+nj1/NHkUUUS39XUV16Bn7jQv9wCRrpapCBJEMS9JHxMlP3VKEQSjCYNo2TYgjIkz+3M4Md5INhuMIhuNp2xJJbeIrZjrdA5MvODQcUTAcSX/AQdV0dKYpvnH9geCkbdGRODpHJn9upn4Hu4eB7vRt8YiCzkOTH2vo3jf5ex46NILQofR5SsSS6MywiFJvW/+kbSMZ8nuknKca7rrrLnzta1/D2rVr8fjjj+dceISQMTmNfHv37sX111+PSy65BOvWrUNv7ye30XW5XPB66VQvQszKaeR76qmnkEgk8Mwzz6Cqqsrw7/rrr5+pGNPq1wZwd+QB9GsDed0vb7r6R7Hu7j+hq3961rc8WvW29uGHF29Eb+vkX4GnW07F981vfhO6rqf999BDD81QiOmpUBHUQ1BBk8eZJFUNh4IRJFWNdSiWllRU9HcMIZnHld64PbHaK3jwd/L58Aqe7A8+hpV5Xbjm745HmdeV/cHHsJJqL9Z8/yKUVOfvTydui88pOLHQNg9OgS6VyaTAacPpC6tR4OT2uum8cHlcWLz8eLg8+fuQ4rb4COEdFR8hjFDxEcIIt8UX1kexVXkDYZ0OoWcSCsfx3NaPEZrkjBgyZrg/jM33bMFw/+RnFk03bv8KH9VHsTXxJuswLC80quD5rftZh2F5I/1hbL7nT3ndJ7cjnwMOzJbq4QCd3paJ0yFhwexSuldDFnKhA8ed3gi5MH+/T9wWX7HowxedK1Es+liHYmkVxQX4ty+eQncpyqLMX4Kv/uIfUOYvyds+uS0+QnhHxUcII0wOuLzXXQ3ZVZ39gRk0Or1AI/BedxVaY1P7StX5rTlTev5MiC6cns/FuE+c+BlVp9bnK+9vno6QplVHcnrOWdXVZmD0y/jGczshSB9PS5/ZcDvyJXQN3UoYCZ1OGM4koWroCUWQoBOrs0gAQtfYzzzhdqqhIz6C6z5+hXUYltcZHMW/PbmNdRiWJ0gHgKKr87pPbkc+QnjHbfHVyx78cu5y1Mt0SVEm/hI37l5zNvwlbtahWJquzoY+8gR0dXbe9slt8UmCAK9NhiRMsvIRATCWJ4/LQXnKSgJ039jPPOG2+AjhHRUfIYxQ8RHCCLfF1xUP4+a219EVz98lIDzqDkWw/vk/ozs0+cK/BIAYAApuGPuZJ9zO88V0FXui6VeDJp+IJ1XsOxRiHYblCUIMsO3O6z65HflKbU5cPet4lNpoAaVMSgpkfOmzc1FSIGd/8DFM18qgx74CXSvL/uBpwm3xeW0yLi5thtdGv1SZeFwOLD++Hh4XXfeYke4DlEv/b7ohP7gtPkJ4R8VHCCNUfIQwwm3xDScVvDTYiuFk/m5gz6ORWAIvf9iBkVj+LpXhkhAC7M+P/cwTbqca+pNR3NvzV9ZhWN7AaAwPbf+IdRiWJ4h9gOsXed0ntyOfQ5DQ6PTCIdCqXJk4JBENpUVwSNy+1Xmh6zJ0tRm6nr+j59y+I7WyG//VeC5qZbpUJpNqXyF+eMlpqPYVZn/wsUzzA6O/GPuZJ9wWHyG8o+IjhBEqPkIY4bb4dB2IqAnoOutIrE3XdUSUJHRKVBYagNH/+5kf3E41tMVDWL3nBdZhWN6BwTCuefh/WYdheYLUCnhW5nWf3I58hPCO2+KrdRThvxuXotZRxDoUS6vxFeI/Lz0dNTTVkJGu1kEP3wtdrcvbPrktPocoos7pgUPk9iXkhV0SUVvshp0m2bNwAFr92M88oXeEEEao+AhhhIqPEEa4Lb5eJYIfHHwDvQqtypXJoZEofvLyLhwaibIOxdrEbsC1fuxnnnA7zzeqJfCXcA/rMCwvoiSx82Af6zAsTxBGAfubed0ntyOfT5Jxaelc+CRaQCkTr8uBi05sgJcWUMpI14qhx78IXSvO2z65Lb4SuxOXz1qAEjstHZhJcYGMVZ+Zg2JaOjAzvRSI/+PYzzzhtvgI4R0VHyGMCHoeT3ffuXMnFi9ejB07dmDRokVT6isajWL//v1oamqCy+WapgiPPpQnc1jkiduRT5IkeDweSBKt4ZIJ5ckcFnnidqrB4XCgri5/J8HyivJkDos8cTvyaZqGRCIBTcvfxY88ojyZwyJP3BZfPB7Hnj17EI/HWYdiaZQnc1jkidviI4R3VHyEMELFRwgjVHyEMMLtVIPT6cT8+fMhCALrUCyN8mQOizxxW3yCINAvlAmUJ3NY5Inbr53xeBytra10CD0LypM5LPLEbfFpmoZIJEKTx1lQnsxhkSdui48Q3n3q4tu7dy/cbjceeuihaQyHkGPHpyq+RCKBNWvWYHR0dLrjIeSY8amKb/369fB4PNMdS07sdjuqq6tht9uZxmF1lCdzWOQp56mGLVu24J577sGuXbuYXqpis9lQUlLCbP+8oDyZwyJPORVfMBjE5ZdfjjvuuAN+f/Z7Vzc2Nhr+P34YNxaLIRr9ZB1JSZLgcDigaVraQ73jVxbH4/GJo1GqqiIcDqO4uBiyLCOZTCKRSBieJ4oiZFmGruuIxWIp/TqdTgiCYOh3nN1uh81mg6qqUBTF0CYIApzOsYWbDn8d42RZhiiKUBQFqqoa2mw2G+x2e9Z+Y7FYyj31HA4HJElCIpFAMpk0tE2WQ1VVMTo6itLSUthstoz9zlQOM/ULWCOHqqoiEomgpKQEoihm/D3M1q/Z0TOn4rvuuutwxhlnYPXq1bk8LUUgEDB8bfV6vfD7/UgkEti/f3/K4xcuXDjxvCPfKFEUUVFRgVAohO5u44KnbrcbDQ0N0DQtbb8tLS2w2Wzo6enByMiIoa2yshJlZWUIh8Po6OgwtDmdTjQ3NwMAWltbU96I5uZmOJ1O9PX1YWhoyNBWVlaGyspKRKNRtLe3G9psNhtaWloAAO3t7SkF1tDQALfbjYGBAfT39xvaiouLUVNTA0VR0r5Wj8cDm82GQCCQUkR+vx9erxfBYBA9Pca1UIuKilBfXw9VVdP2O2/ePEiShO7uboTDYUNbVVUVSktLEQ6HEQgEDG0ulwtNTU0AkLbfOXPmQJZl9Pb2IhQKGdrKy8sxa9YsRCIRHDhwwNDmcDgwd+5cAEBbW1tK4TY2NqKgoAD9/f0YGBhI2W9RUVHamERRxPz58wEAHR0dKcVZV1cHj8eDYDCI8vLylH7TMb2GyyOPPIJ169bhvffeQ3Hx2NqGgiDgwQcfxFVXXWVqZ+NruGzbtg0nn3zyxPZPM/LF43EEAoGJX0ga+dLncDxP42uT0MiXPoeH50mWZWuNfA888AB6e3tTvm5ee+21eOKJJ/C73/3ObFdwOp1pF6kRRTHj4jXjb9bhxtfcsNlssNnSvxxBEHLu9/D+Mz03U5vDMflCtdn6Hf8FSsdut0/6BmfLYaZ+ZyqHmfoFjq4c5nLAxnTxPfrooymfUHPmzMGGDRuwZs0a0zskhIwxXXw1NTVpt1dUVEzaNpNEUURBQQFEujlmRpQnc1jkidurGmRZTjmaSlJRnsxhkacpFV8e19tNu29d1+mSmSwoT+awyBO330VisRg+/PDDtEfgyCcoT+awyBO3xUcI76j4CGGEio8QRqj4CGGE66mG4447ju6+kwXlyRwWeeK2+ERRpIljEyhP5rDIE7fviqIoOHjwYMqJtcSI8mQOizxxW3yqqmJ4eDjljHdiRHkyh0WeuC0+QnhHxUcII1R8hDDCbfHZbDbMmjUr40WahPJkFos8cfuO2O1202tlHMsoT+awyBO3Ix8dxTOH8mQOHe3MAc1fmUN5Mofm+Qg5hlDxEcIIFR8hjHBbfIIgQJZlWpckC8qTOSzyxO1Ug9PpxJw5c1iHYXmUJ3NY5InbkY8Q3nFbfNFoFB9++GHadf7JJyhP5rDIE7fFByCvN6/nGeXJnHznieviI4RnVHyEMELFRwgj3BafLMsTNzIkk6M8mcMiT9zO82W7gSEZQ3kyh0WeuB35FEVBV1cXna2fBeXJHBZ54rb4VFXF4OAgXaeWBeXJHBZ54rb4COEdFR8hjFDxEcIIt8Vns9lQWlpKq3JlQXkyh0WeuH1H7HY7qqqqWIdheZQnc1jkiduRT1VVRCIROoqXBeXJHBZ54rb4FEVBa2srzV9lQXkyh0WeuC0+QnhHxUcII1R8hDDCbfEJggBJkmhVriwoT+awyBO3Uw1OpxPz5s1jHYblUZ7MYZEnbkc+QnjHbfHFYjHs3bsXsViMdSiWRnkyh0WeuC0+XdehKAp0XWcdiqVRnsxhkSdui48Q3lHxEcIIFR8hjHBbfA6HA/X19XA4HKxDsTTKkzks8sTtPJ8kSSgqKmIdhuVRnsxhkSduR75EIoHe3l4kEgnWoVga5ckcFnnitviSyST6+vqQTCZZh2JplCdzWOSJ2+IjhHdUfIQwQsVHCCPcFp8kSfB6vZAkiXUolkZ5ModFnridanA4HPD7/azDsDzKkzks8sTtyKdpGuLxON3yOAvKkzks8sRt8cXjcezbtw/xeJx1KJZGeTKHRZ64LT5CeEfFRwgjn6r4Hn74YcyfPx9OpxMLFizAk08+Od1xEXLUy7n4Hn30UXz5y1/G2rVr8cEHH2D16tVYtWoV3njjjZmIj5CjlqDncN28rutobGzE5z//edx2220T28877zycc845WLduXcbn79y5E4sXL8aOHTuwaNGiTx81IUeBnOb59u7di/b2dqxevdqwffPmzdMaFCHHgpyKb8+ePQCA0dFRnHfeeXjnnXcwe/Zs3HLLLbjwwgtTHt/Y2Gj4//hh3Fgshmg0OrFdkiQ4HI6JuZYjuVyuieePz8MoioJDhw6hpqYGBQUFSCaTKZeDiKIIWZah63raVamcTicEQUg7v2O322Gz2aCqasrNMwRBgNPpBADD6xgnyzJEUYSiKCl3vbHZbLDb7Vn7jcViKYv5OBwOSJKERCKRcvb9ZDkcz1NdXR1kWc7Y70zlMFO/gDVyeHie7HZ7xt/DbP3a7faU56aTU/ENDw8DAK644gqsX78et956KzZt2oSLL74Yr7zyCpYuXWqqn0AgAI/HM/F/r9cLv9+PRCKB/fv3pzx+4cKFE8878o0Kh8MoKChAKBRCd3e3oc3tdqOhoQGapqXtt6WlBTabDT09PRgZGTG0VVZWoqysDOFwGB0dHYY2p9OJ5uZmAEBra2vKG9Hc3Ayn04m+vj4MDQ0Z2srKylBZWYloNIr29nZDm81mQ0tLCwCgvb09pcAaGhrgdrsxMDCA/v5+Q1txcTFqamqgKEra1zpeGIFAIKWI/H4/vF4vgsEgenp6DG1FRUWor6+Hqqpp+503bx4kSUJ3dzfC4bChraqqCqWlpQiHwwgEAoY2l8uFpqYmAEjb75w5cyDLMnp7exEKhQxt5eXlmDVrFiKRCA4cOGBoczgcmDt3LgCgra0tpXAbGxtRUFCA/v5+DAwMpOx3/MPryJhEUcT8+fMBAB0dHSnFWVdXB4/Hg2AwiPLy8pR+08npb74nnngCq1atwsaNG3HddddNbF+xYgUA4MUXX8z4/PG/+bZt24aTTz55YvunGfni8TgCgcDELySNfOlzOJ6npqYmuFwuGvkmyeHheZJl2XojX21tLQDg+OOPN2xfsGABXnjhBdP9OJ3OiRdyOFEU024fN/5mHW78RFibzTbpLX0FQci538P7z/TcTG2Z1gPJ1u/4L1A6drt90jc4Ww4z9TtTOczUL3B05dBs4QE5TjUsWrQIRUVFePPNNw3b33vvvYmvYYQQc3Ia+VwuF2666SZs2LABNTU1OPXUU/HrX/8aL7/8Ml599dWZijEtu92O2tranD5pjkWUJ3NY5CnnS4puueUWFBQU4Fvf+hY6Ozsxb948PP300zjnnHNmILzJ2Ww2+Hy+vO6TR5Qnc1jkKacDLlM1nZPsyWQSoVAIXq83498TxzrKkzks8sTtidWJRALd3d20JF4WlCdzWOSJ2+IjhHdUfIQwQsVHCCPcFp8oinC73RBFbl9CXlCezGGRJ24Pf8myjIaGBtZhWB7lyRwWeeL241DXdaiqSrc7zoLyZA6LPHFbfLFYDLt3787rDex5RHkyh0WeuC0+QnhHxUcII1R8hDBCxUcII9xONTidTrS0tNDdd7KgPJnDIk/cFp8gCHSWvgmUJ3NY5Inbr53xeBwHDhygG4BkQXkyh0WeuC0+TdMwMjJCt77KgvJkDos8cVt8hPCOio8QRqj4CGGE2+Kz2+2orKykVbmyoDyZwyJP3B6DttlsKCsrYx2G5VGezGGRJ25HPlVVEQqFUpYSJ0aUJ3NY5Inb4lMUBR0dHSlr9RMjypM5LPLEbfERwjsqPkIYoeIjhBFui2/8PmyCILAOxdIoT+awyBO3Uw2H3x2WTI7yZA6LPHE78hHCO26LLxqN4oMPPkh7S2HyCcqTOSzyxG3xAaC1KE2iPJmT7zxxXXyE8IyKjxBGqPgIYYTbqQZZltHc3AyHw8E6FEujPJnDIk/cFp8oinA6nazDsDzKkzks8sTt105FUdDZ2Uln62dBeTKHRZ64LT5VVTE0NETXqWVBeTKHRZ64LT5CeEfFRwgjVHyEMMJt8Y0veEP3IciM8mQOizxx+46ML/VGMqM8mcMiT9yOfKqqIhwO01G8LChP5rDIE7fFpygK2tvbaf4qC8qTOSzyxG3xEcI7Kj5CGKHiI4QRbotv/Da+tCpXZpQnc1jkiduphvEb2JPMKE/msMgTtyMfIbzjtvhisRg++ugjxGIx1qFYGuXJHBZ54rb4dF1HMpmklbmyoDyZwyJP3BYfIbyj4iOEESo+QhjhtvgcDgcaGhpoVa4sKE/msMgTt/N8kiTB7XazDsPyKE/msMgTtyNfIpFAT08PEokE61AsjfJkDos8cVt8yWQS/f39SCaTrEOxNMqTOSzyxG3xEcI7Kj5CGMm5+JLJJL7zne+gvr4eRUVFWLJkCd58882ZiI2Qo1rOxff9738f9913H+677z688847OO6443D++eeju7t7JuKblCRJKC4uhiRJed0vbyhP5rDIU87F9+yzz2L16tVYtmwZmpubcfvttyMUCuGNN96Yifgm5XA4UFNTQ/NXWVCezGGRp5yLr6KiAi+88ALa29uhqiruvfdeyLKME088cSbim5SmaYjFYtA0La/75Q3lyRwWecp5kv3nP/85LrvsMsyePRuSJEGSJGzatAlNTU0pj21sbDT8Px6PAxi7fOPwG89LkgSHwwFN0yYecziXyzXx/PHkxONxBAIBNDQ0wO12I5lMpszRiKIIWZah63raS0WcTicEQTD0O85ut8Nms0FV1ZQVrQRBmLid1OGvY5wsyxBFEYqipCxFZ7PZYLfbs/Ybi8VSzrB3OByQJAmJRCLlkPhkORzPU1NTE1wuV8Z+ZyqHmfoFrJHDw/Mky3LG38Ns/drt9pTnppNz8X344Yfw+Xx49tlnUVNTg/vuuw9r1qzB66+/jpNOOslUH4FAAB6PZ+L/Xq8Xfr8fiUQC+/fvT3n8woULJ5535BsViUTgdrsRCoVS/u50u91oaGiApmlp+21paYHNZkNPTw9GRkYMbZWVlSgrK0M4HEZHR4ehzel0orm5GQDQ2tqa8kY0NzfD6XSir68PQ0NDhraysjJUVlYiGo2ivb3d0Gaz2Saupm5vb08psPEPmoGBAfT39xvaiouLUVNTA0VR0r7WcYFAIKWI/H4/vF4vgsEgenp6DG1FRUWor6+Hqqpp+503bx4kSUJ3dzfC4bChraqqCqWlpQiHwwgEAoY2l8s18YGdrt85c+ZAlmX09vYiFAoZ2srLyzFr1ixEIhEcOHDA0OZwODB37lwAQFtbW0rhNjY2oqCgAP39/RgYGEjZLzBWiEfGJIoi5s+fDwDo6OhIKc66ujp4PB4Eg0GUl5en7fdIgp7DBUwdHR1obm7GH/7wB5x11lkT28866yyUlpbi2Wefzfj8nTt3YvHixdi2bRtOPvnkie008qX2SyMfjXwGb731FhRFwWc+8xnD9tNOOw0vvfSS6X6cTufECzmcKIppt48bf7MON350ymazTbrOviAIOfd7eP+ZnpupLdMf79n6zXSXVLvdPukbnC2HmfqdqRxm6hc4unJotvCAHA+41NbWAgDeffddw/Z33313YqjPJ1qRyxzKkzn5zlNOI9+pp56KM888E1deeSU2btyI2tpaPPzww3j11Vexbdu2mYoxLZfLhQULFuR1nzyiPJnDIk85jXyiKOL555/Hueeei6uuugqLFy/Ga6+9hldffRWf/exnZypGQo5KOR1wmarxAy47duzAokWLptRXLBZDIBBAbW1txu/gxzrKkzks8sTtidXjR99oVa7MKE/msMgTt8VHCO+o+AhhhIqPEEa4LT6HwwG/309n62dBeTKHRZ64Xr3M6/WyDsPyKE/msMgTtyMfLQxkDuXJHFpAKQe0JJ45lCdzaOlAQo4hVHyEMELFRwgj3BafKIooKiqCKHL7EvKC8mQOizxxO9UgyzLq6+tZh2F5lCdzWOSJ249Dut2xOZQnc+i20DlgcQN7HlGezGGRJ26LjxDeUfERwggVHyGMUPERwgi3Uw1OpxPz5s2j+assKE/msMgTt8UnCALd9soEypM5LPLE7cdhPB5He3t72mW9yScoT+awyBO3xadpGsLhMN36KgvKkzks8sRt8RHCOyo+Qhih4iOEEW6Lz263o6qqKqdbMh2LKE/msMgTt1MNNpsNpaWlrMOwPMqTOSzyxO3Il0wmEQwGaVWuLChP5rDIE7fFl0gkEAgEaFWuLChP5rDIE7fFRwjvqPgIYYSKjxBGuC0+URThcrnobP0sKE/msMgTt1MNsiyjqamJdRiWR3kyh0We6OOQEEa4Lb5oNIr3338f0WiUdSiWRnkyh0WeuC0+QnhHxUcII1R8hDBCxUcII1xPNcyZM4culcmC8mQOizxxW3yiKEKWZdZhWB7lyRwWeeL2a6eiKOjo6ICiKKxDsTTKkzks8sRt8amqilAoBFVVWYdiaZQnc1jkidviI4R3VHyEMJLXAy7jp+7s3r17yn3FYjEEAgEMDw/D6XROub+jFeXJnOnOU0tLCwoKCjI+Jq/F197eDgD40pe+lM/dEpJ3O3bswKJFizI+RtDzeBPq/v5+bN68GQ0NDXC5XFPq68ILLwQA/Pa3v52O0I5alCdzpjtPZka+vBbfdGpsbAQAtLa2Mo7E2ihP5rDIEx1wIYQRKj5CGOH2aychvKORjxBGqPgIYYSKjxBGqPgIYYS74tM0DevXr0dNTQ0KCwuxfPlytLW1sQ7LcgYHB3HttdeitrYWHo8HZ555JrZu3co6LEvbu3cv3G43Hnroobzsj7vi+973voeNGzfi3nvvxfbt26GqKs477zy6Xu0Iq1atwvbt2/H444/j7bffxkknnYRly5Zhz549rEOzpEQigTVr1mB0dDR/O9U5Eo/H9aKiIn3jxo0T24aGhnSXy6X/6le/YhiZtezbt08HoG/dunVim6ZpelNTk/7tb3+bYWTWtW7dOv3cc8/VAegPPvhgXvbJ1ci3a9cujIyMYOnSpRPbfD4fFi1ahC1btjCMzFrKysrw4osv4pRTTpnYJggCBEHA0NAQw8isacuWLbjnnnvy9nVzHFdruAQCAQCA3+83bK+urkZHRweLkCzJ5/Nh+fLlhm2bNm3Cxx9/jPPPP59RVNYUDAZx+eWX44477kj5vZppXI18kUgEAFIWunE6nYjFYixC4sL27dtx9dVXY+XKlVixYgXrcCzluuuuwxlnnIHVq1fnfd9cjXzjlyHF43HDJUmxWAyFhYWswrK05557DqtXr8bnPvc5PPbYY6zDsZRHHnkEf/rTn/Dee+8x2T9XI9/414Kuri7D9q6uLtTU1LAIydLuvPNOrFy5EhdeeCFeeOEFupL9CA888AB6e3vh9/vhdrvhdrsBANdeey0uuOCCmQ8gL4d1pkksFtM9Ho9+//33T2wbP9r5+OOPM4zMejZu3KgD0P/1X/9V1zSNdTiWFAgE9H379hn+AdBvvfVWPRAIzPj+ufraKcsy1q5di5tvvhnl5eVoaGjAjTfeCL/fj0svvZR1eJaxd+9eXH/99bjkkkuwbt069Pb2TrS5XC54vV6G0VnHZN+WKioq8vJNiqviA4ANGzYgmUzimmuuQTQaxZIlS7B582ZaDv0wTz31FBKJBJ555hk888wzhrYrr7wy74fUSXp0PR8hjHB1wIWQowkVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGMUPERwggVHyGM/H+9TjBQsUOT/AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAN8AAAGdCAYAAACSHqb/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA+G0lEQVR4nO3dd5hU1fkH8O+d3tv2xu7SFFBQFxUVUCGigiWgRgJqMBK7QU3UGDX8QkyiiZoYFWvUWIOKIrEBioqA2BBFQEBpW2d3dnZ6u+33x0bMZcruusst8H6eh2cf7rlz7jtn5p0zc+895zCiKIoghMhOp3QAhBysKPkIUQglHyEKoeQjRCGUfIQohJKPEIVQ8hGiEEo+QhQia/IlEgmsX78eiURCzsMSokoGOQ/29ddfo6GhAUNmXQ9rWXW/6vJYTZh8aDXe+boJoWSmX3WV/31tvx6/PzS+dNiA1OPRWzDFU4/loZ0I8al+1SWKzIDENJDqfxsfkHqcHiuOmzIMHy7fjmgo2a+63tp6Z6/20+zXTlEEkhkOdHNcYSJEJHgOIqihChFFEakkCznvtpS15xtI4VQGb2zao3QYqhfm03gttF3pMFQvFk5h1X+2yHpMzfZ8dpMBx9aVwm7S7OeHLOw6I45zVMGuMyodiqpZ7SaMPm4QrHaTbMfUbPKVOK24dtIYlDitSoeiaqVGO35VOQ6lRrvSoaiar9SBi359InylDtmOqdnkI0TrKPkIUQglHyEK0WzyZTgBOwMRZDhB6VBULSPy2JHqQkbklQ5F1dgMh6YdnWAznGzH1OypwpZwHL9d+pHSYahecyaKG/esVDoM1WtvjuDvN7wh6zE12/MRonWaTb46nxNP/Wwy6nxOpUNRtXqzG88P/THqzW6lQ1G1ynov7vj3LFTWe2U7pmaTDwxg1OsA9d1uqDIMjDo9qKEKY8DAYNSDkbGd+px8giBg/vz5qKqqgt1ux9SpU7Fz5879ERshB7Q+J98f/vAHLFy4EI888gjWrl0Lnudx6qmnIpPp38gCQg42fUq+TCaDu+++GwsWLMC0adMwZswYLFq0CE1NTVi8ePH+ipGQA1Kfkm/Dhg2IRqOYPHny3m0ejwdHHXUUVq1aNeDBFdIciuOGl9eiOTQw47kOVM2ZCK7dtRzNmYjSoaiavzmMv167FP7msGzH7NN1vqamJgBATU2NZHtlZSUaGxuz9h88eLDk/+l0unt/jw3uou/PUsbTLDpiKRj1OlR5sm8A3tUZBQBUuGwwG/WSMpNeB5YX4LQYUWS3SMpSLI+2SAIMA9TmOCvaGIyBF0UUVXthdZolZSF/FLFgHFanBUXVHklZJsWhfWcAAFB1aDmYfX6j+3cEwKY5eCvcsHukN35HO+MIt0dhtplQUuuTlPGcgNbt7QCAQSYX9Iz0s7E1E0NK5OAzWODWS59rlM8gwCVgYnSoMrn2eaYiMmL3zQhVJidMjLQN29k44gILt94Mn0Eab0Jg4Wfj0INBjVnaDgCwOx2CAKDCaIdln5ETnVwCET4Dh86Ikn1u7E4LHFrYGACgPke9TZkIWFFAqdEGu0460qCLSyHEp2DVGVBVL21DjuPhb+xOoIpaD3Q6aRt2tEaQSXFw+2xwuKVtGI+mwGV4GE16lFZJzw6LooiWXV0AgLJqNwz7vA87/VGkEiycHmmdhfQp+b6b/sFslr5RLRYLgsFgr+u5+sTDMXLkyL3/X/1NKx5Y9RV8NjP+fPa4rP1/+vgKAMDlE0dheKlHUvavdV/jrc2NGFdfhp8fN0JS9kVTAHcs/xwWgz5nvZc+9x6iKRY/vvFUHHbSIZKyJX9dhvefXofh4wZjzt3nScqatrTi7vMfAQBc+8wlMOwzrOnO6QvR9m0Hplw2EeNmHCUpe/ux1Xj9H++gemQFrn58jqQs5I/g96f8DQBwS9UJKDLaJOXzG9/HpmQAp3mGYIbvUGm94Z14yL8epUY7/lo7WVImiCKu3PkWAlwC88qPxmCL9HT63S3r8GGsGROcgzCndLSk7JNYC+5s+RB2vRF31U3Cvi7YvhRJgcPcsiNwhL1MUvaofwPeCu3AUY5yzKs4WlK2NdmJ3+55HwBy1nvVjmVoY+OYWTwSJ7oGScoWBbbghc4tOMTiw3V3jZeUBVojuOPqVwEAl80/JSvB7rv5TezeFsDEs0bgxDNHSsr8jSE89seVsDnMuO6uaZKyVCKDWy9cBAC48NcTUV7jkZQ//ud3sfnTJhx98tCs55IP05eFUhYvXoxzzz0XiUQCVuv3n5A/+clPkE6n8eqrrxZ8/Pr169HQ0IAJV94Kd2Xt3u0/pOerdNtwzUmjseCNT7ClLdSvnm/US1tU1/MxS48fkJ6vyuTEtRXH4Ibd72BnOtS/ns8kbQdA+Z6v4WHpc/mhPV9plQuzr5uAv/36dbQ3h/vV87340e+znk8ufer5vvu62dLSgiFDhuzd3tLSgtGjR+d7WJaWUAJd5mjWdpYX9iZaLq2R7ImXkmz3PYvRFItois35OFFEwXo7m7ryliWjKTRtactb3vx1/rKu1jC6WnP/hkgnMgXr3VPgN1qQSyHI5Z6PJSMK2JkO5X1scyZ/O4T5NMJ8OmcZD7Fgva1s/t/eMYFFrMBjC9XbziYA5J5wKylwaN6Z/zda6+789YaDCYSDuetlMzyad+b/Judvyn/MaKj38+T06YTLmDFj4HK58N577+3dFgqFsH79ekycOLEvVRFy0OtTz2c2m3H11VfjpptuQklJCerq6nDDDTegpqYG55xzzv6KkZADUp9HNSxYsAAcx2Hu3LlIJpOYOHEili1bBqNR3jlCwskMXv1iJ8L9nDbwQBfmU3g5+DXC/Zw28EAXDSfxzstfIRru37SBfdGnEy799d0Jl4GYt3MgHcjzdg6kA3nezoF0wM/baTHoMaLcC4tB3/POBzELY8AoazEsjGaHbsrCbDFgyKgymC3ytZNmk6/cbcPvpo5FudvW884HsQqTA7+vOREVJvlm5dKi4goXrlgwBcUV+96gsP9oNvkI0TpKPkIUQslHiEI0m3y8IKIzngIv0AIghfCigE42AV6kWd4K4XkBoc44eF6+dtLsKbDGrhiuXvSB0mGo3p5MBJftfFPpMFSvbU8It1/6sqzH1GzPR4jWaTb5arwO3H/+BNR46RR6IYNMLjxcfzoGZY3xI/+rfJAHtz4yA+WDPLIdU7PJp9cxKLJboNep764LNdEzOhQZbVnDk4iUXq+Dp8gOvV6+dqJXhBCFKHLC5fGLnsERh/dvEcIMNwRtkXF44sJnYDJ826+6HDf2fui/XKaeckjPO/VCZR0H/BmofIiDuCv3eMfeEmzyLRzZW1uvKB2Qejhn92Dr3ed58W1UnrSgno8QhWg2+Qz6FpQ6b4RB36J0KKoWaIviwQUrEWjLP4KdAC2JBG7+7GO0JHKPbt8fNHudT8ckYTFuVDoM1UunOHy7pV3pMFQvyfP4qiv/dCL7g2Z7Pk4oQigxB5xQpHQoqubyWnH6zNFweWnt+kJ8ZjMuGjIMvn1m5tufNJt8guBBJHU+BMGjdCiq5nRbMPnskXC61XdSSU28JhPOqx8Mr0m+k0qaTT5CtI6SjxCFUPIRohDNJp+OicJufgs6hk6hF5KIpfHRym+RiOWeDJd0i7Asljc3IcL270aEvtDspQaDvh1F9nuVDkP1ugIJvPjoJ0qHoXodqRTu27JJ1mNqtucTRBMy3CAIovpueVITg1GPsmpX1toCRMqk02GQ3Q6Tjm6s7hHH16At8jA4vqbnnQ9iZVUu3PDXqSiroiFFhdTY7XjguPGosWcv1LO/aDb5CNE6Sj5CFELJR4hCNJ588p0W1ioRIjiWhwia5a0QEQArCLK2kmYvNZgM32KQ7yylw1C9ll0h/OaiF5UOQ/V2RKOYsXKFrMfUeM9HiHZpNvlYvgat4fvA0qWGgkorXbj2T1NQWkmXGgqpttnx92OOQ7WNLjX0SBRNYPmhEOkie0FGkx7V9T4YTXSRvRCzXochLhfMNHsZIQc+Sj5CFELJR4hCNJt8Bp0fxY4/wqDzKx2KqgU7Ynjq72sQ7IgpHYqqtSWTuOPLDWhLJmU7pmav8+l0MdhMq5UOQ/WScRZfftSodBiqF+c4rGmX94Ncsz0fL3gQSU4HTxMoFeRwmzFx6iFwuOWblUuLPCYTzh5UCw9NoNQzXihCKHkpeJo6sCC314azLjwSbq9N6VBUrchsxtzhh6KIpg4k5MBHyUeIQij5CFGIZpOPYRKwGteBYeRb2EKLUgkWmz5rRipBw68KiXMcPupoR5zjZDumZi81GPWtKHH+XukwVK+zPYYn7vpA6TBUry2ZxO1ffC7rMTXb84miHrzghijSDcOF6PQM7E4zdHpaPrsQPcPAZTRCz8jXTppNPpavQ3Po32D5OqVDUbWKGg9+/8h0VNR4lA5F1eocDjx74iTUORyyHVOzyUeI1lHyEaIQRU643NE2Dj5Pcb/qKNOV4mIr8Gf/cfALQ/pV1yTv1/16/P7QNcY7IPU4y5wAgMhwJ7rc/avL+8aWAYhoYA27MT4g9VQdWg48dzxq7t0C3ddt/avs2t7tRj0fIQrR7KWGdqED98QfAEvTBxbU2B7FVfe8jTTLKx2KqrVs8+PmCX9FJkmrFPVIhIgMMkqHoXqiCKQylHg9EQUR6bi87yfNfu30Mh78xDwdXsajdCiqVuq14bqfNKCURjUUVFzjxaUPzERxzcD81u4NzSafiTFhsKEOJoZmLyvEYtLjsMHFsNDsZQWZ7WYcetwQmO00pIiQAx4lHyEKoeQjRCGaTb6IEMXy9EpEhKjSoahaMJLCM8s3IxhJKR2KqoXaIlh8x1sItUVkO6ZmLzUkkcR67gulw1C9WJLFu+tp9rKexEMJrHnhM1mPqdmezwIzRukPhQU0K1chdosR40ZVwG4xKh2KqtlcFjRMPQw2l0W2Y2o2+dw6N860nA63rp83LB7gitwW/OLM0Shyy/em0iJvpQezbz8b3kqPbMfUbPIRonV9Tr5gMIjLL78c1dXVcLlcGD9+PFavppmjCemrPiffzJkzsXbtWjz//PP49NNPccQRR2DKlCnYunXr/oiPkANWn5Lvm2++wYoVK/Dggw9iwoQJGD58OO677z5UVlbi2Wef3V8x5sSKLJr5FrAijWooJM3y+LY5RKMaepBJZrDryyZkkvLdXN2nSw3FxcV4/fXXMXbs2L3bGIYBwzDo6uoa8OAKCYpdeDq1SNZjapE/mMCfnv5I6TBUr2N3EP+Y8y9Zj9mn5PN4PJg6dapk2+LFi/HNN9/gtNNOy9p/8ODBkv+n02kAgE/nQ5mudO/2lJhCWIxADz2KddlrL/iF9u7HMV4YGekp87AQRgppWGGFS+eUlGXEDLrEEBgwKNWVZNXbIQQgQIBZLIZetEofy4TAMVHoRRvMojQmgWGRYrpHO9uE7DXhk0wbRIaFSfTBIErX+GaZKFgmBJ1ohkUslZSJ4JHUtQAAqkoc0OukM2m1dyWQyvDwOMxw2aU3lMeTLDojKRgNOlQUZa8rvsfffTNCRZEdRoP0C08glEQizcFlM8HjlF66SaY5dISS0OsYVA3PbsOWbwIQBRHF1R6YbdLXJuSPIR5Owuo0w1chXRM+k+TQ0dj9gZ2rXv+uILgMD1+FC9Z9Yop0xhHtTMBsM8FXLX3NuQwH/45A93MdVgrdPss8d+zuRCbJwl3qhMMnbad4KIFQWwRGiwGlddKZFkRBRMu27lWMyuqLYTBLU6ezqQupWBrOHG2fT78usq9duxYXX3wxZsyYgWnTpvX6cWeZT8dI68i9//+K24LX0m/ByThwsXV21v53xP8GAJhmnoIqfaWk7P3MGnzIfowRhuGYYp4kKdvB7cIL6VdghDFnvffGH0ISSQzifgKvOEZStkf/Itr0b8MlHoph3GWSsjizB5uMfwQAjORugg7SN91Gw/8hybSikp+GUmG8pKxF9yaaDEtgF2sxgvuVpCyDLmww/QYAcO15DfDtc83pL899jK17ujCpYRCmHSf9YFv1RRP+9eYmlHismH/x8VnP9fdPrMUefxS/OHM0asulifDgKxvw6VY/xo2qwPmTD5WUbdjejvsWfw6r2YDrHstuw1tOW4h0IoPp156MQ46plZS9/LeVWPvKlxgxrh6zbpN+OO/6qhX3X9n9zSVXvX/+6RPobA7j1EuOQ8OUEZKy5U+sw/In1qF2VAUuvXu6pCzQGMSfzn4QAHDFQ7Pg8EqT4d45T2L3l8048YJjcdIFx2Yd9+5Z/wQA/Oq5SyTbU7E0fjvxLgDAz/4yA+VDpB8Y/7z2BWxatR3HnCV9HxXCiKIo9nrv//Hqq69i1qxZOOGEE7B06VJYLD1fR1q/fj0aGhow+6k5KDu0fO/2H9LzFTFenGWZimeTL6JRaOpXz3e6J6C6nu+1By8akJ6vvMiOy84aszf5+tPzjdjWkfVcFe/5fNKyH9rzldYV4cI//Rh3z/on2ncF+tXzLW59LOv55PKDer77778f8+bNw3nnnYennnoKpj6uaRYUgoCQfa6HB7830XI+TuwCvvuo+O/DvxvNnkQSSSH3qqIixIL1ppkAkGeuVJ5JIFFgSvqELv+tWxkmiAwTzFkmMGkkmPyPbS6wkmwolkYols5ZxnLC3q+YubR25p9wKJLIIJLIfcKBF0Q050i+7wSaQnnLktE0mqP5H1uo3mBrBGjNXZZOZNDcnv9cQ+v2/K95uD2KcHvudmJTHJoLTKLk3xnIWxYt0L776vOlhgcffBDXXHMNrr76ajz//PN9TjxCSLc+9Xzbtm3DvHnzMH36dNx8883w+79fRtdqtcLtplu9COmtPvV8L730EliWxSuvvIKKigrJv3nz5u2vGHMKCJ14KPE4AkKnrMfVmpZAHDc/9AFaAgMzv+WByr+jA386eyH8O/J/BR5ofUq+3/72txBFMee/J598cj+FmBsPHiExDB508bgQjhfQHkqA4wWlQ1E1LsMj0NgFTsaZ3jR7Y7WbceEM82lwM66edz6IFbutmHvG4Sh2W3ve+SDmq3Rj9u1nwVcp308nzSafhbHgMMMIWBgaKlOIzWLAcYdVwmbR7LhpWVhdVjRMPRxWl3wfUppNPkK0jpKPEIVQ8hGiEM0mX0yMY3XmQ8REOoVeSDiWxqurv0E4zx0xpFskEMOyh1chEsh/Z9FA0+yv8LgYx2p2ndJhqF44nsHS1d8qHYbqRQMxLHv4A1mPqdmezwQT6vW1MIFubyvEYtJjVH0RrdXQA7PdhEOOGwyzXb73k2aTz6vz4HzLDHh1HqVDUbVSrw3Xnz+WVinqQXGND5c98FMU1/hkO6Zmk48QraPkI0Qhipxw2dhaCbO1sucdCxhscQODgY2tFdiR6t9XquZbhvXr8ftD8rCB+VxMe3R7/yb5/tW54qtlAxHSgGrkBuaeVZEfCsQvwa9eXQ9G/82A1NkTzfZ8rCigNRMDK9INw4WwvIC2cAIs3VjdAxZgWrr/ykSzlxoa01Fc8c0KpcNQveZQHNe/uEbpMFSP0e8GnBfLekzN9nyEaJ1mk6/W7MK/hk9FrZmGFBVS43PgodknosbnUDoUVRP5eojRRRD5etmOqdnk0zMM3AYz9EyemY8IgO52cllN1E490gOip/uvTDSbfIRoHSUfIQqh5CNEIZpNvpZ0DDftfB8tafmGgGhRaziB+Us/Rms4/8S/BICuCbBd2/1XJpq9zpcSeWxN5p4NmnwvzfHY3h5WOgzVY5gUYNgi6zE12/MVGSy4uOxwFBloAqVCfDYzLjh2OHw2c887H8REoRhi6lKIQnHPOw8QzSaf22DG2UVD4TbQm6oQl9WEqYfXwmWlcY8FiR4gc85/LzfIQ7PJR4jWUfIRohBKPkIUotnki3AZvBHcgQgn3wL2WhRNsVi+uRHRlHxDZTSJCQPGpd1/ZaLZSw0BLolH2r5QOgzV64yn8OTar5UOQ/UYXQdgfUDWY2q25zMxegy2uGFiaFauQkx6HeqKnDDpNftSy0IUzRD5oRBF+c6ea/YVqTY7cM/gSag201CZQio9dvxp+jhUeuw973wwE2qA+APdf2Wi2eQjROso+QhRCCUfIQrRbPKJIpDgWYii0pGomyiKSGQ4iNRQPRAAxP/7Vx6avdSwMx3GrK2vKR2G6u0OxjD3qXeVDkP1GP0OwDVD1mNqtucjROs0m3zVJif+MXgyqk1OpUNRtSqPHX855zhU0aWGgkR+EMTYIxD5QbIdU7PJZ9LpMMjigkmn2acgC6Neh2qvA0a6yN4DEyDUdv+VCb0ihCiEko8QhVDyEaIQzSafP5PAH/d8CH+GZuUqpD2axF3LN6A9mlQ6FHXTtQLW+d1/ZaLZ63xxgcUnsTalw1C9RIbD+j0dSoehegwTB4zrZD2mZns+j96Mc4qGw6OnCZQKcVtNOGtMHdw0gVJBouCFmD4fouCV7ZiaTT6f0YILy0bBZ6SpAwvx2syYefQweGnqwMLEIiD98+6/MtFs8hGidZR8hChEkRMuXmcCDk//1ljwGbt/w/icScQs/auLt8r3VaO3GG6A6uG//9vfOs/ddk7/AxpgzW/UDkg99T4n7jgD+M2iC7EzGO1XXZv+0rv9NNvzJYQMPkvuRkKg2csKiac5fPStH/H0AGXzASqe4bBulx/xjHztpNlLDQE+hoe7Vikdhup1RJP4x4ovlQ5D9dpjSfxtlbztpNmeTw8dPDob9Np9CrLQ6xj47GbodbQsdCF6HQOfTd520uw7t8rowV/Kz0GV0aN0KKpW43PgvgsnosZHs7wVMsjjwIPnTsQgj3ztpNnkI0TrKPkIUQglHyEKoeQjRCGavdTQyAZxZcuz4GWc6k2Ldgei+Nkjb4MXaOrAQnYFo5j9jLztpNnkEwFwlHg9EgFwlHg9UqKdNPu1s1TvxK+KTkGpnmYvK6TcbcMtZzWg3G1TOhRVq3Da8LspDahwytdOmk0+i86IQ8zlsOiMSoeiahajHiMrfbAYaSm1QixGPUaVy9tOmk0+QrTuByfftm3b4HA48OSTTw5gOIQcPH5Q8rEsi9mzZyMejw90PIQcNH5Q8s2fPx8ul2ugY+mTIB/HU6EPEeTpA6CQQCyFR9/bjEAspXQoqhaIp/Dw2s0IxOVrpz5fali1ahUefvhhbNiwAYMGyTev/b5iQhqrE98odnytiKVYvPd1s9JhqF40zWLlN/K2U5+SLxQK4cILL8R9992Hmpqe164ePHiw5P/pdBoAUG5ww2f07d0eFzLo5GMwQIfKHKMU9rBBAECZ3gWzrjtkK9N9tnNdYgfa+SgcOjN8euliICmBRTsfBQMGNcbsWama2S7wEFHss8FqkTZFKJJCNJaB1WJAsU96+pllebR1dPe41RUuMPuMQmlrj4HlBPg8Ftht0lnDIrE0wpE0zCY9Soul8fK8iBZ/9yjqGp8ja3hLWziBFMvDazdnzUYWS7EIxFIw6nWo8n5fr9VkwMhKH5Z9tQexFItKjx0mg/QLT3skiUSGg8tqgs8unWgpkeHQHklCr2NQpc8e8d/CByFCRLHOBTMjPfMcEuKIiynYGDO8OulogQw4dPBhAMhZr58PgQMPn84BKyONKSIkEBWTMDNG1Pukl5pYXkBTuPu1GeTJbsOWSAJpjofPZobb8n0b2kwGDPG58O63LUhzPKrc0tdGEEXs7uqeMaHKbYdpn7Uv/NEkEiwnqbMnfUq+K664AscffzxmzZrVl4dlmesdj5ElI/f+f11iBx4PrYFXb8etJdOy9r+05WkAwBzv8RhiKpGUcaKAN2IbMdZSi1meYyVlm1ItuDf4DsyMIWe917e9gJiQxrnTRmD0yDJJ2Uuvb8HK1Ttx6NBi/GL2UZKyPc1h3HH/GgDADVceB6NBenr6D39bhdb2GE6fNAwnHC39kFr23rd4ddlWDKpy47pLx0nKusJJ3HJH91p6N049Ej6HdGa225d+ii0tXZgyqgZnHVUvKXt3SzMee38zSl1W/PFcab0A8NmudsRSLK6cfBjqS6Q/Ge5d/gU+3tGOE4aV44LjD9nncR24560NsJkM+JV3ela9N3f+C2mRxQzH8TjUVC0pWxxbgzWpLRhhqsFs50mSsl2sH/8I/wcActb7p+ALCAgRnG4biwbLUEnZssR6LEusR52hFJedIX2ubZEE5i3pfm1um9IA1z7JcOsbH2N7IIxpI2txxsjsKSg2tnV/0N+xT72JDIeL/9392lx34mjU7DP06C8rN+Czpg6cPLQyq858GLGXS5Y+/fTTuPnmm7Fx40Z4vd29CMMweOKJJzBnzpxeHWz9+vVoaGjAKQ9dA9/wqr3bf0jPV25wYa53Au4KLMe2jL9fPV/5szWq6/lcR3oHpOer9Nhw1Y9G45aX1mFXINqvnm/Medlzmyjd84kfSj8sfmjPV+my4ZcTR+M3r61Dczjer55v7YKrsp5PLr3u+R5//HH4/f6sr5uXX345Fi1ahDfffLO3VaGNCyPGWrO2cxD2Jloufj4C8NJtKZEF0P0bMCakcz5OhFiw3kAw/5TzyRSHxpZI3vKm1vxlwVAKwVDuH/DpDF+w3sZg/kmhuuJpdMVzP1eWF7ArkH8CoJZQ/hNUkWQGkWTuOXF4QUQz35n3sQEh/3NJiGkk+NzxAihYb1CIAcjdFmmRRXOByY72hPK3YTCRRjCRO6YMLxScRKk5nL8Nw6nezynU6+R75plnkExK5/sfNmwYFixYgNmzZ/f6gISQbr1OvqqqqpzbS0tL85btTymBxdZ0G1ICK/uxtSTF8tjcEkSK5Xve+SCWYnlsapO3nTQ7qqGdj+LuzhVKh6F6beEE/rj0M6XDUL3WaAILlsvbTv1Kvl6eq9kvGHTPYMZDAA2YyY9B98xcvCBSOxWgRDtp9sbqGqMPCytno+Z/rheSbLXFTvzr0h+htpiGXhVS53Pi2Qt+hDqffO2k2eQjROso+QhRCCUfIQqh5CNEIZq91NDMhnBj22JEBRoqU0hjMIZrnl6FcJ47V0i3PaEYrnhpVZ/uUOkvzSYfDwEhIf9tYaQbL4gI5rkVjXyPF8S8t5vtL5r92lmsd+Ay70QU62kBkEJKnFb88pTRKHFm30tLvlfqsOK6iaNR6pCvnTSbfDadCQ3WWth0vR8/dTCymw04dkgZ7GbNfsmRhd1kwLi6MthN8rWTZpOPEK2j5CNEIfRdBIA+RdPO94bDSCduBpJme74Qn8DLkc8R4umMZyFd8TQWfbQ97+Bb0i2YSOP59dtlPeOp2Z4vIqTwVuwrpcNQvXAyg6Wf71I6DNULpzJY8tUuWY+p2Z7PyhgxxlwNK0NrNRRiMxlwVG0JbDKexdMim9GAhuoS2Ix0trNHJQYnrio6GSUGGipTSKnLil+dfgRKXXSdr5AypxU3TjoCZTJeD9Vs8hGidZR8hCiEko8QhWg2+ViRRwsbAivSrFyFsLyApmAMLE/XMgvJ8AIaQzFkZGwnzZ4Ca+XC+L+O/ygdhuo1d8Vx0wsfKh2G6jWH4/j1UnnbSbM9HyFap9nkqzZ4cW/5+ag2ZK/BQL5XW+TAYz8/GbVFNPSqkFqvA0/MPBm1XvnaSbPJp2MYWHUm6PZdpYRIMAwDq8kAhtqpIB3DwGYyyPp+0mzyEaJ1lHyEKISSjxCFaDb52rgwbu94HW1cWOlQVK0lFMctL60ruC4f6b7U8N3CmHLR7HW+jMgXXPCSdMtwhRfLJN16WhBzf9Bsz+fT2/BT9zHw6W0973wQK3JYMGf8oSjaZ313IlVkt+DnxxyKIrt87aTZ5HPoLDjZfggcOnpTFeK0GHHKYTVwWmjcYyEusxGnHloDl1m+dtJs8hGidZR8hCiEko8QhWg2+SJ8CitimxHhaaGUQiLJDN74YjcitFBKQeFUBq9t3k0LpfRGSEjgxYi8C9hrUTCexrMfblM6DNULJtJ4+lN520mzPZ+ZMWCwsRhmRrOfH7IwG/QYWuaG2aBXOhRVMxv0GFYsbztpNvnKDC78puR0lBlcSoeiahUeG34//RhUeOh6aCGVLhtun3oMKl3ytZNmk48QraPkI0QhlHyEKESzyceLAqJ8CrxIs3IVwgsiIskMeEFUOhRV4wURkZS87aTZU4XNXAi/8r+odBiq1xiM4Yp/va90GKq3JxTDL16Qt5002/MRonWaTb4Kgxu3l56NCoNb6VBUrcprx90/PQFVXrvSoahatduOe398Aqrd8rWTZpPPyOhRanDByNDF40KMeh3K3TYY9Zp9qWVh1OtQ7pK3negVIUQhlHyEKESzZzsHkiGqxjv+zUoHkKXcor65YLYqHUA/aLbna+ei+Hvn22jn1PeGUBN/OIE7XlsPfzihdCiq1hZN4E9vr0dbVL520mzPlxJZbE63Kh2G6iVZHhubOpUOQ/WSLI8vWuRtJ832fG6dFWc6R8Oto7XGC/HYTJgxdjA8NpPSoaiax2rCuWMGw2OVr520m3x6K850joFbT8lXiMdmxjljh8BjU99vSDXxWs04b8wQeK3ytZNmk48QraPkI0QhlHyEKESzyRcXMliX2IG4oMZrdOoRT7NYva0V8TSrdCiqFsuw+GBHK2IZ+dpJs5caOvkYHg+tUToM1euIpvDgyq+UDkP1OmIp3L9a3nbSbM9ngA4leicM2n0KsjDqdShzWenG6h4YdTqUOa0w6ujG6h5VGj34Y9mPUWn0KB2KqlV57bhn1ngaUtSDao8d/5g+HtUeGlJEyAGPko8Qhfyg5HvqqacwcuRIWCwWjBo1Ci++SHOpENJXfU6+Z555BpdccgmuvvpqbNq0CbNmzcLMmTPx4Ycf7o/4CDlg9elSgyiKuO2223DttdfiyiuvBADccsstWLVqFd577z0cd9xx+yXIXPawQVza8rRsx9OqXYEoZj+0QukwVG9nMIrzn5K3nfqUfNu2bcOuXbswa9YsyfZly5YNaFCEHAz6lHxbt3aPG47H4zj11FPx+eefo76+HrfeeivOPPPMrP0HDx4s+X86nQYAlBvc8Bl9e7fHhQw6+RgM0OW8dLCHDQIAyvQumHXdIft0dpzpGoOnuz7ELq4TDp0ZPr30NHFKYNHOR8GAQY3Rm1VvM9sFHiKKSuywWqVrcYdCScQiaVhtRhQVS+vNsDzaW7sH8VbVuMEwjKTc3xoBywrw+qywO6R3yUcjKYRDKZjNBpSUOSRlPC+gtTkCAKjxOaDXSettCyeQYnl47Wa49xn6EkuxCMRSMOp1kssKxQ4Lzhk7BP9Y8SVawwlUeuwwGaS/NtojSSQyHFxWE3x2abyJDIf2SBJ6HYNipiyrDTvFdogQ4Wa8MGKfmMQIUkjCDAucjHSWOQ4sQmL365qr3i4xAB48nHDDzFikMYkxJBCHESbU+5ySMpYX0BSOAwAGebLbsCWSQJrj4bOZ4bZ8H2+x3YIZo+vxjw++Qmc8hap9ZjETRBG7u2IAgCq3HaZ9rpv6o0kkWE5SZ0/6lHyRSPcb46KLLsL8+fNx5513YvHixTj77LOxYsUKTJ48uVf1zPWOx8iSkXv/vy6xA4+H1sCrt+PWkmlZ+3/39XKO93gMMZVIykZaKrEr1omxllrM8hwrKduUasG9wXdgZgw5672+7QXEhDR+fP5oHHZEhaRsyaIv8f6KbzB8RCnmXCGtt2l3CHf/YSUA4NrfngSDUTqD2p2/W4G2liimnDkC4ybUScrefmMrXn95E6rrPLj6homSslAwid/f+CYA4MapR8LnkL7pbl/6Kba0dGHKqBqcdVS9pOzdLc147P3NKHVZ8cdzx2U9V/N/Y7xy8mGoL5Gu7HTv8i/w8Y52nDCsHBccf4ik7LNdHbjnrQ2wmQz4ienirHofTd8DFhlMMJyCQTrph+0qdjm+EtajVjcEPzJKP5zbhGa8zHa/rrnqfSb9ECII4VjDBAzXHyYp+4RbjU/41ShnqvCLM6TPtS2SwLwl3Xc+3TalAa59kuHWNz7G9kAY00bW4oyRtVnHtRj0qHLbccc+9SYyHC7+97sAgOtOHI0aj/SD8y8rN+Czpg6cPLQyq858GFEUez0/9qJFizBz5kwsXLgQV1xxxd7t06Z1v7Fff/31go9fv349GhoacMpD18A3vGrv9h/S85UbXJjrnYC7AsuxLePvV8836L5S1fV85ikVA9LzVXpsuOpHo3HLS+uwKxDtV893xpzGrDZUuuf7dsXxkrIf2vNVumz45cTR+M1r69Acjver51u74Kqs55NLn3q+6upqAMDhhx8u2T5q1Ci89tprva6njQsjxmYPguUg7E20XPx8BOCl21Ji942wMSGNmJDO+TgRYsF6OzviecuSCRZNe0J5y5sbw3nLuoJJdAWTOcvSaa5gvY3BWP5642l0xXM/V5YXsCuQf16bllD+5xpJZvIuH80LIgKiP+9jw2JX3rI0UkiL+ZfvLlRvFGFExdxtzCKDncH8z3VPKH8bBhNpBBO52zDDCwXrbQ7nb8O+LCvdp0sNRx11FJxOJ9atWyfZvnHjRgwdOrQvVRFy0OtTz2e1WnHjjTdiwYIFqKqqwjHHHIN///vfWL58Od555539FWNOAS6Gf3atRoDL/+lGgI5oEgvf2YiOaO4emHRrjyVx3wcb0R6Tr536PKTo1ltvhc1mwy233ILm5maMGDECL7/8Mk466aT9EF5+CTGDj5I7ZT2mFsXTHNZsb1M6DNWLZzis3ilvO/2g28uuv/567NixA+l0Ghs2bMDZZ5890HH1yKEz4yTbcDh0NDFQIU6LEaeMqobTYux554OY02zElEOq4TTL106avbHap7djlufYrDOcRKrIYcGcCSNQtM9lCyJVbLfgkmNHoNguXztpNvkI0TpKPkIUQslHiEI0m3wpgcWmVAtSAs3KVUiS5fFlYwBJlu9554NY91oN8raTZmcva+ejuDco77VFLfKHE7jz9c+VDkP1ulcpkredNNvzMWBgYYxgwPS880GMYQCrUQ+GmqkgJdpJs8lXY/TiHxUzc94wTb5XW+TEY5dMQm2Rs+edD2J1Xiee/Okk1HnlayfNJh8hWkfJR4hCKPkIUYhmz3YOpDsXPaZ0CFku/ut1SoeQZef0YqVDyDZT6QB+OM0mXzPbhevbXkCSVikqqDEYw+VPvodEhlM6FFXb0xXD3EXytpNmk4+HmHfkOvkeL4iIpuhGhJ7wooiozMuoafY3X4negat8J6FE7+h554NYqcuK6087AqUuWru+kDKHFTecfATKHPK1k2aTz6ozYYylBlZd76dqOxjZTAY01JXAZtLslxxZ2EwGjK2Rt500m3yEaB0lHyEKoeQjRCGaTb4uPoEXwp+ii08oHYqqBeNpPLN2K4J55vkk3YKJNJ76ZGveuTz3B83+Co8KKbwd36J0GKoXSWbw5pd7lA5D9cKpDF7fIm87abbnszEmNFgGwcbQ2c5CbCYDjhlcSmc7e2A3GTCuthR2OtvZs2KDA5f5TkSxga7zFVLqsmLelDF0na8HpQ4rrjtxDErpOh8hBz5KPkIUQslHiEI0m3wZkcOeTCcyIt2tX0iGE7CzI4IMJygdiqpleAE7OyPI8PK1k2ZPgbVxEdweeEPpMFSvJRTHrYs/UjoM1WsOx/Gb1+VtJ832fIRonWaTr8bgxQMVs1BjoNnLCqktcuLJX0ym2ct6UOdz4pnZk1Hno9nLesQwDIyMPms9dCLFMIBRr6N5O3vA4L/tJOMxNZt8hGgdJR8hCqHkI0Qhmr3U0MqGMb99KQJcTOlQVK25K44bF61FeySpdCiq1hSO41evroU/Jl87aTb5WPBo5cJKh6F6LC+guSuudBiqx/ICmsLytpNmv3b69HZc6B5Ha7L3oNhhwdwTR6KY1mQvqNhuwWXHjaQ12XvDoTNjgn0YHDqz0qGomsNixMkjquCwGJUORdWcZiMmDauC0yxfO2k2+QjROko+QhRCyUeIQjSbfBE+iTejXyHC0yn0QsLJDJau34lwkhaUKSSUzGDJxp0IydhOmr3UEBKSeCUq7wL2WtQVT2PRx98oHYbqdSXTeP5zedtJsz2fmTFguKkMZkaznx+ysBj1GFHphcWoVzoUVbMY9BhZ5oXFIF87aTb5ygwu/Lp4CsoMLqVDUbVytw23njUW5W6b0qGoWoXLhvmnjkWFS7520mzyEaJ1lHyEKISSjxCFaPZsBS8K6OLj4MX+zzY185O5AxDRwBqoeZN5QUQwlgIviP2uq3l67QBEpE6cIKIzngI3AO3UW5pNvmYuhJv8Lysdhuo1BmO45pkPlA5D9RpDMVy5WN52oq+dhChEs8lXZfDgzrIZqDJ4lA5F1Wp8Dtx3wQTU+GhBmUJqPA4sPGcCajzytZNmk0/P6ODV26FnNPsUZKHXMfA5LNDraPqyQgw6BkV2CwwythO9cwlRCCUfIQqh5CNEIZpNPj8XwV2B5fBzEaVDUbW2cAK3L/0UbeGE0qGoWmskgd8v+xStEfnaSbPX+dIih20Zv9JhqF6K5bGlpUvpMFQvxfHY7Je3nTTb83l0Vkx3HgmPjtYaL8RrN+P8Y4bCa6eJpgrxWs346ZFD4bXK106aTT6X3orTnYfBpafkK8RtNeGso+rhtpqUDkXVPFYTfnx4PTwytpNmk48QraPkI0QhfU4+juPwu9/9DrW1tXA6nZg4cSLWrVu3P2Ij5IDW5+S7/fbb8eijj+LRRx/F559/jkMOOQSnnXYaWltb90d8ecWEND6Ib0dMSMt6XK2JpVi8u6UZsRSrdCiqFk2zWLm9GdG0fO3U5+RbsmQJZs2ahSlTpmDo0KG4++67EQ6H8eGHH+6P+PIK8nE8HV6HIE+LgBQSiKXw2PubEYillA5F1QLxFB7+cDMCcfnaqc/JV1paitdeew27du0Cz/N45JFHYDabMWbMmP0RX15G6FFhcMMImpWrEKNehyqvHUY9/bwvxKjXodotbzv1+SL7vffei/POOw/19fXQ6/XQ6/VYvHgxhgwZkrXv4MGDJf9Pp7u/IpYb3PAZfXu3x4UMOvkYDNCh0ujJqmcPGwQAlOldMOsM/63DhbneCbgrsBzbMn44dOasFYtSAot2PgoGDGqM3qx6m9ku8BBRZrTDrpMukNHJJRHm07DrjCgzSutNizyaM1EAQL3ZDWaflbybMhFkRAElBiuceul1oxCfQpBLwaIzoNIoHb7CQcCedPcdOzU+R9ZIhLZwAimWh9duzrp0EEuxCMRSe5PtO5UeG6760Wjc8tI67ApEUemxw2SQvsHaI0kkMhxcVhN8+1wPTGQ4tEeS0OsYVBdlD7fZHYxCFIFyly1resLOeArRFAu72YASh/SSUJrj0frfu27qipxZ9TaH4mB5ASVOC+wm6WsTSqQRSmZgNepR4pPONva/S30N8mS3YUskgTTHw2czw235vg0rXTb8cuJo/Oa1dWgOx1Hllr7mgihid1f3WpBVbjtM+ySpP5pEguUkdfakz8m3efNmeDweLFmyBFVVVXj00Ucxe/ZsvP/++zjiiCN6Vcdc73iMLBm59//rEjvweGgNvHo7bi2ZlrX/pS1PAwDmeI/HEFOJpGyoqRTbMn6MtdRiludYSdmmVAvuDb4DM2PIWe/1bS8gJqTx89LRONpZISl7wv8l/tP1DUbbS3FDlbTeHakQfr1rJQDgjtqTYNRJ33TzdqxAYyaK84pH4EeeOknZ4s6teLZjE4ZYPPjDoImSsk42iV98+yYA4MapR8K3z7Jety/9FFtaujBlVA3OOqpeUvbulmY89v5mlLqs+OO547Ke63eunHwY6kuk0y3eu/wLfLyjHScMK8cFxx8iKftsVwfueWsDbCYD/jQju95LnlyJJMtjzgmHYHR1saTsiTVbsGJzE46oKcZVJx8uKdvuD2H+0k8AIGe91y1aDX8kifMahmL8MOlrs/izb7F4/Q4MK/PgN6cfJSlriyQwb8kaAMBtUxrg2icZbn3jY2wPhDFtZC3OGJl7Wowqtx13nCGNKZHhcPG/3+2O7cTRWeP+/rJyAz5r6sDJQytz1pkLI4piryetaGxsxNChQ/H2229jwoQJe7dPmDABRUVFWLJkScHHr1+/Hg0NDTjloWvgG161d7vSPZ8+Waq6nm/4lxXq6/mKVdjzFauv51u74Kqs55NLn3q+jz76CJlMBkcffbRk+7hx4/DGG2/0up42LowYm31nCgdhb6Ll4ucjAC/dlhK7z07FhHTeM58ixML1svlP2sQFFjvSobzlO9P5V8ft4JLo4HKvJZESuIL1NgbzL3fdFU+jK577ubK8gF2BaN7HtoTyP9dIMoNInrUKeEHErs789bYVuCE5nuYQT+d/bKF6O6IpdCD3SZAky2NnMP9j94Tyt2EwkUYwkbsNM7xQsN7mAivYhlO9X+uhT78uq6urAQBffvmlZPuXX36J4cOH96WqfhNFEazIow8d90FJFLsTkpqpMBH/bScZj9mnnu+YY47B+PHj8bOf/QwLFy5EdXU1nnrqKbzzzjtYs2bN/ooxp0auC1e1PifrMbVod2cUcx59R+kwVG9XMIoLnpW3nfrU8+l0OixduhSTJk3CnDlz0NDQgJUrV+Kdd97Bscce23MFhJC9+nxRw+v14oEHHsDu3bsRiUSwZs0anHjiifsjtoLKDS7cWjwV5bRQSkGVHjtuP+dYVHrsPe98EKty23HHtGOzTrTsT5q98mpiDBhkKoKJlggryGTQob7ElXWGk0iZ9DrUF7myzmLuT/SKEKIQSj5CFELJR4hCNJt8AS6Gh4PvI8Dlv5BKuu9euXf5F2iP5L7YT7q1x5L42/tfoD0mXztp9mxFQszgs9QepcNQvUSGw8c72pUOQ/XiGQ7rdsvbTprt+Zw6C35kHwGnztLzzgcxl9WE00cPgosmUCrIbTFh2ohBfRqV0F+aTT6v3oafuMfCq5dvAXst8tnNuOD4Q7JumCZSPpsZFx19CHw2mjqQkAMeJR8hCqHkI0Qhmj3bmRQy+CLViKTQ+/FT+fzlyMUDENHAmv/enAGpJ5Hh8NmuDiQyXL/rqny7YwAiGliN00p63qkXEhkOnzYOTDv1lmaTr4OP4YHge0qHoXrtkSTueWuD0mGonj+WxF/f3SDrMTX7tVMPBg6dGXrQcseF6HUMnBYjLQvdAz3DwGk2Qs/QstA9qjJ6cU/5T1CVY24W8r0anwMPzTkJNb7s+VfI9wZ5HXjs/JMwyCtfO2k2+QjROko+QhRCyUeIQij5CFGIZi81NLJd+GXrv5EW5bsuo0W7O6OY+8+VSHF8zzsfxHZ1RTHneXnbSbPJJ0LcO2EuyU8UuyeXJYUp0U6a/dpZqndinm8ySvXZU42T75W5bbhp2pEoc9Poj0LKnTb89kdHotwpXztpNvksOiNGWSph2WeNBSJlNeoxuqYYViMtpVaI1ajHmEp520mzyUeI1lHyEaIQSj5CFKLZ5AvycTwX+ojWZO9BZyyFJz/Ygk5ak72gQDyFf360RdY12TV7qSEmpPFeYpvSYaheNMVixaYmpcNQvWiaxfKt8raTZns+G2PCsdZ62BialasQu9mAE4aVw27W7OesLOwmA8bXl8Nukq+dNJt8xQYHLvGOR7GBhsoUUuK04srJh6PEmb0SMPleqcOKayYcjlKHfO2k2eQjROso+QhRCCUfIQrRbPKlBQ7fZjqQFmhUQyFplsf2thDSdHN1QSmOx7aOEI1q6A0/H8GdgbeUDkP1WsMJ/N+ST5QOQ/VaIwnc9qa87aTZno8QrdNs8g0y+vBI5YUYZPQpHYqq1RU78ezlp6CumIZeFVLvc2LRRaeg3idfO2k2+QjROko+QhRCyUeIQij5CFGIZi81tLAh3OJfgi4aUlRQc1cc1z+3GsF4WulQVK0pFMcvX5G3nTSbfBwEdPBRpcNQPZYX4I8klQ5D9VhBgD8qbztp9mtnkd6Bn3tOQJGeRjUUUuK04IpJh6HEaVE6FFUrcVhw9fjDUOKQr500m3x2nQnjbINh19F4vkLsZiPGD6+A3UyzvBXiMBkxYXAFHCb52kmzyUeI1lHyEaIQWU+4JJPdP2iTjZ39ritoSGJzx2YEu5oR48L9quvbZKLf8Qy0RMfAzCcSYm3YvNmFUOtuJEL9e56RZNeAxDSQkv6BOTsZTna3U7h5N5KR/rXT+vXrceihh8Jm62H2a1FGzzzzjAiA/tG/A/7fZ5991mM+MKIoipBJIBDAsmXLUFdXB6u1f3NlnHnmmQCA//znPwMR2gGL2ql3BrqdetPzyZp8A2nw4MEAgB07digcibpRO/WOEu1EJ1wIUQglHyEK0ezXTkK0jno+QhRCyUeIQij5CFEIJR8hCtFc8gmCgPnz56Oqqgp2ux1Tp07Fzp07lQ5LdYLBIC6//HJUV1fD5XJh/PjxWL16tdJhqdq2bdvgcDjw5JNPynI8zSXfH/7wByxcuBCPPPII1q5dC57nceqppyKTySgdmqrMnDkTa9euxfPPP49PP/0URxxxBKZMmYKtW7cqHZoqsSyL2bNnIx6XcWaE/Xoz5wBLp9Oi0+kUFy5cuHdbV1eXaLVaxeeee07ByNRl+/btIgBx9erVe7cJgiAOGTJEvO222xSMTL1uvvlmcdKkSSIA8YknnpDlmJrq+TZs2IBoNIrJkyfv3ebxeHDUUUdh1apVCkamLsXFxXj99dcxduzYvdsYhgHDMOjq6lIwMnVatWoVHn74Ydm+bn5HU3O4NDV1D7OpqamRbK+srERjY6MSIamSx+PB1KlTJdsWL16Mb775BqeddppCUalTKBTChRdeiPvuuy/rfbW/aarnSyS6x1mZzWbJdovFglRKvoXstWbt2rW4+OKLMWPGDEybNk3pcFTliiuuwPHHH49Zs2bJfmxN9XzfDUNKp9OSIUmpVAp2u12psFTt1VdfxaxZs3DCCSfg2WefVTocVXn66afxwQcfYOPGjYocX1M933dfC1paWiTbW1paUFVVpURIqnb//fdjxowZOPPMM/Haa6/BYqEZzP7X448/Dr/fj5qaGjgcDjgc3TPhXX755Tj99NP3fwCynNYZIKlUSnS5XOJjjz22d9t3Zzuff/55BSNTn4ULF4oAxF/+8peiIAhKh6NKTU1N4vbt2yX/AIh33nmn2NTUtN+Pr6mvnWazGVdffTVuuukmlJSUoK6uDjfccANqampwzjnnKB2eamzbtg3z5s3D9OnTcfPNN8Pv9+8ts1qtcLvdCkanHvm+LZWWlsryTUpTyQcACxYsAMdxmDt3LpLJJCZOnIhly5bBaKR5Kb/z0ksvgWVZvPLKK3jllVckZT/72c9kP6VOcqPxfIQoRFMnXAg5kFDyEaIQSj5CFELJR4hCKPkIUQglHyEKoeQjRCGUfIQohJKPEIVQ8hGiEEo+QhRCyUeIQv4f5GQaKX2YYcEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(temp)\n",
    "plt.show()\n",
    "plt.imshow(temp_r)\n",
    "plt.show()"
   ]
  }
 ],
 "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}