{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "# Ctypes Example\n", "Author(s): Paul Miles, Joel Kulesza | Date Created: June 28, 2018\n", "\n", "**Note, this tutorial is useful for illustrative purposes. However, compiling the C++ code in Binder is non-trivial, so this is recommended as a Read-Only example.**\n", "\n", "In this example, we demonstrate how to use a model written in C++ within the context of using [pymcmcstat](https://github.com/prmiles/pymcmcstat/wiki). In this case we consider a linear model.\n", "\n", "$$y(x;q) = m x + b + \\epsilon, \\quad \\epsilon\\sim N(0,\\sigma^2), \\quad q = [m, b]$$\n", "\n", "but note that the model can be arbitrarily complex. Similar work can be performed with Fortran using the `iso_c_binding` module." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.9.0\n" ] }, { "data": { "text/plain": [ "{'divide': 'warn', 'over': 'warn', 'under': 'ignore', 'invalid': 'warn'}" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy as np\n", "from pymcmcstat.MCMC import MCMC\n", "import matplotlib.pyplot as plt\n", "import ctypes\n", "from numpy.ctypeslib import ndpointer\n", "import pymcmcstat\n", "print(pymcmcstat.__version__)\n", "np.seterr(over='ignore')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Define data set\n", "To define our data set, we generate a set of points along the line, $y = 2x - 3$, and add some random noise to the response. We have plotted the data with and without the random noise added to it." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWkAAAFZCAYAAABaAd7jAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3deXhTVfrA8e9J2MsqIlvZlIIbKiIKiLTgUkAUHa24DIpToTKOyyioOC7gLCg64IgwAu1PcButgqBssgVFdkQBEbEghUYqiEDLDm3O74+bxDRN0iTN3vfzPPfR3NybnAvl5XDOe96jtNYIIYSITaZoN0AIIYR3EqSFECKGSZAWQogYJkFaCCFimARpIYSIYRKkhRAihkmQFkKIGFYt2g3wh1JKAS2AI9FuixBChFA9YK/2sWAlLoI0RoC2RrsRQggRBsnAz97ejJcgfQSgoKCA+vXrR7stQghRacXFxbRq1QoqGCGIlyANQP369SVICyGqFJk4FEKIGCZBWgghYpgEaSGEiGFxNSbti81m4/Tp09FuhohDNWrUwGSS/oqITQkRpE+fPs2uXbuw2WzRboqIQyaTiXbt2lGjRo1oN0WIcuI+SGutKSwsxGw206pVK+kRiYDYbDb27t1LYWEhrVu3xlg3JUTsiPsgXVJSwvHjx2nRogV16tSJdnNEHGrSpAl79+6lpKSE6tWrR7s5QpQR993O0tJSAPmnqgia42fH8bMkRCyJ+yDtIP9MFcGSnx1RGVarFYvFgtUansoVCROkhRAi0nJycmjTpg19+vShTZs25OTkhPw7JEgLIUQQrFYrw4YNc2aV2Ww2srKyQt6jliAdZ5RSzJ492+c1Q4YM4ZZbbolQi363fPlylFIcPnw44t8tRCRZrVZyc3PLpf2WlpayY8eOkH6XBOkoePPNN6lXrx4lJSXOc0ePHqV69eqkpaWVudYR+Hbu3AlAYWEh/fr1AyA/Px+lFN9++23E2u5Ljx49KCwspEGDBtFuihBh4xjieOKJJ8q9Zzabad++fUi/T4K0i3BPADj07t2bo0ePsmHDBue5FStW0KxZM9auXcvJkyed5y0WC61bt+a8884DoFmzZtSsWTOs7QtWjRo1aNasmUzEiYTlPsThymw2M2XKFJKTk0P6nRKk7SIxAeDQsWNHmjdvzvLly53nli9fzsCBA2nXrh1r1qwpc753797O167DHe3atQOgc+fOKKXK9cJfffVVmjdvTuPGjXnooYc4c+aM1zaNHj2ayy67jHfeeYe2bdvSoEED7rzzTo4c+b3U7alTp3jkkUc455xzqFWrFj179mT9+vVl2uo63LF7925uuukmGjVqRFJSEhdddBHz5893Xv/dd9/Rr18/6tatS9OmTRk8eDAHDhwI4FdSiMjKy8vzGKAnTJhAfn4+mZmZIf9OCdJEbgLAVe/evbFYLM7XFouFtLQ0UlNTnedPnDjB2rVrywRpV+vWrQNgyZIlFBYWMmvWrDKft3PnTiwWCzNmzGD69OlMnz7dZ5t27tzJ7NmzmTt3LnPnzuWLL77gpZdecr7/5JNPMnPmTGbMmMHGjRtp37496enpHDx40OPnPfTQQ5w6dYovv/ySLVu28PLLL1O3bl0ADh8+TJ8+fejcuTMbNmxg4cKF7Nu3jzvuuKPiXzwhoiQlJaXcqmaz2cztt98e8h60k9Y65g+gPqCLioq0uxMnTujvv/9enzhxotx7/lq2bJkGyh0WiyXoz6zItGnTdFJSkj5z5owuLi7W1apV0/v379fvv/++7tWrl9Za66VLl2pA796923kfoD/55BOttda7du3SgP7mm2/KfPZ9992n27Rpo0tKSpznMjIy9KBBg7y254UXXtB16tTRxcXFznMjR47UV111ldZa66NHj+rq1avr9957z/n+6dOndYsWLfS4ceO01lpbLBYN6EOHDmmtte7UqZMePXq0x+/7+9//rm+44YYy5woKCjSgt2/f7rWd4RCKnyFRdWRnZ2uz2awBbTabdXZ2dlCfU1RU5Ig19bWP+Bf3y8JDwfG3o+s/Y8IxAeAqLS2NY8eOsX79eg4dOkSHDh1o0qQJqamp3H///Zw8eZLly5dz7rnn0rp164A//6KLLsJsNjtfN2/enC1btvi8p23bttSrV6/MPfv37weMXvaZM2e4+uqrne9Xr16dK6+8km3btnn8vEceeYThw4ezaNEirrvuOm677TYuueQSADZt2oTFYnH2rF3t3LmTDh06+P+wQkRQZmYm6enp7Nixg/bt24evB20nwx1AcnIyU6dOdQa1cE0AuHL85losFiwWC6mpqQC0aNGCVq1asWrVKiwWC3369Anq891rUCilKqwSGMw9vjzwwAP89NNPDB48mC1btnDFFVcwceJEwMhmuemmm/j222/LHHl5efTq1Svo7xQiEpKTk0lLSwt7gAYJ0k6ZmZnk5+djsVjCNgHgrnfv3ixfvpzly5eXmfTr1asXCxYsYN26dV7HoyGyNSfOO+88atSowcqVK53nzpw5w/r167nwwgu93teqVSsefPBBZs2axRNPPMG0adMAuPzyy9m6dStt27alffv2ZY6kpKSwP48Q8UKCtItI/u0IRpD+6quv+Pbbb509aYDU1FSmTJnC6dOnfQbpc845h9q1azsn3YqKisLW1qSkJIYPH87IkSNZuHAh33//PUOHDuX48eNe/0J77LHH+Pzzz9m1axcbN27EYrFwwQUXAMak4sGDB7nrrrtYv349O3fu5PPPP+f++++XQkci4iKVfhsMCdJR1Lt3b06cOEH79u1p2rSp83xqaipHjhxxpup5U61aNV5//XWmTJlCixYtGDhwYFjb+9JLL3HbbbcxePBgLr/8cnbs2MHnn39Oo0aNPF5fWlrKQw89xAUXXEDfvn3p0KEDkydPBoxhnZUrV1JaWsoNN9xAp06deOyxx2jYsKHUBBcRFcn022AobWRPxDSlVH2gqKioiPr165d57+TJk+zatYt27dpRq1at6DRQxDX5Gaq6rFYrbdq0KZc0kJ+fH/Z/URcXFztW5zbQWhd7u066LEKIKsvT4pRw1N+oDAnSQogqy9vilHCm3wZKgrQQosqKRvptoAIK0kqp4UqpzUqpYvuxWinVr4J7MpRSPyilTiqltiil+leuyUIIETrRSL8NRKArDq3A00AeoID7gDlKqc5a663uFyulegD/A0YBc4G7gdlKqcu11t9VquVCCBEiycnJwfeeS0vBZXVvqAXUk9Zaf6a1nq+1ztNa/6i1/htwFOjm5ZZHgYVa61e01tu01s8BG4G/+PoepVRNpVR9xwHU83W9EEKEm8dc6vx86NwZFi0K2/cGPSatlDIrpe4EkoDVXi7rDixxO/e5/bwvo4AilyP2MsyFEAnPEZhfffXV8rnU334L3bvDli3w2GNGjzoMAi6wpJTqhBGUa2H0om/VWn/v5fJmwD63c/vs530ZC4x3eV0PCdRCiAjKycnxWODfZrPx0bBh3F+7NqZjx6BTJ5g/P2xDHsH0pLcDlwFXAf8FZiilvBdvCILW+pTWuthxAEcqvKmKiOU9DsPJsSmBEJHgaweWwcBnNpsRoHv3hhUrIIzZIAEHaa31aa31Dq3111rrUcAmjLFnT34Bmrqda2o/X2Ulyh6HkfzLYMSIESxdujQi3yWEtx1YRgFvA9WB4wMHwoIFEOY9PUORJ20CvG26txq41u3c9Xgfw64SEnWPw3CqW7cujRs3jnYzRBXhvsjFBEwG/mV/vTk9nTqzZkEE/iwGmic9VinVSynVVinVSSk1FkgD3rO//7b9nMN/gL5KqSeUUucrpUYDVwBvhKb58SkW9jg8dOgQ9957L40aNaJOnTr069ePvLw85/uehhdee+012rZt63x/xowZzJkzB6UUSqkyz+MqLS2NRx55hCeffJKzzjqLZs2aMXr06DLX7Nmzh4EDB1K3bl3q16/PHXfcwb59v09nuLdn+fLlXHnllSQlJdGwYUOuvvpqdu/e7Xx/zpw5XH755dSqVYtzzz2XMWPGlPmXixC+uC5yqQ3MAoYDWikOjxnDJQsXQoQKgQX6Ledg9Pa3A0uBrkC61nqx/f3WgLNsm9Z6FUZu9DCMYZHbgVvCmiOtNRw7Fp0jgGJV0d7jcMiQIWzYsIFPP/2U1atXo7Wmf//+PjerdTVixAjuuOMO+vbtS2FhIYWFhfTo0cPr9TNmzCApKYm1a9cybtw4XnzxRRYvNn5sbDYbAwcO5ODBg3zxxRcsXryYn376iUGDBnn8rJKSEm655RZSU1PZvHkzq1evZtiwYc5dylesWMG9997Lo48+yvfff8+UKVOYPn06//znP/16NiHAWOSyZ+NGCi+8kIEANWuiPvqIhs8/H9mG+NpbK1YOAtnj8OhRrY1wGfnj6FH/NjfT0d3j8Mcff9SAXrlypfP9AwcO6Nq1a+vc3FyttbHn4aWXXlrmcydMmKDbtGlT5nsGDhxY4bOmpqbqnj17ljnXtWtX/dRTT2mttV60aJE2m816z549zve3bt2qAb1u3bpy7fntt980oJcvX+7x+6699lr9r3/9q8y5d955Rzdv3tzj9bLHofBo506tU1KMP9uNGmn91Vch/Xh/9ziU2h1R4rrH4YoVK8rscegYlw71HoeO/Qq3bdtGtWrVuOqqq5zvN27cmI4dO3rdr7CyHHsbemtPq1ataNWqlfP9Cy+8kIYNG3psz1lnncWQIUNIT0/npptu4j//+Q+FhYXO9zdt2sSLL75I3bp1ncfQoUMpLCzk+PHjYXk+kWA2bDByoPPyoE0bWLkSXPb3jKTE24i2Th04ejR63+0n1z0ODx06FBN7HLoymUyOf8U4+TsUEo72uHvrrbd45JFHWLhwIR9++CHPPvssixcvplu3bhw9epQxY8bwhz/8odx9Ui9aVGjBAsjIMIYwL7vMyIH2sflGuCVekFYK4mSPPMceh4cOHWLkyJHO8657HA4fPtzr/cHucXjBBRdQUlLC2rVrnePIv/32G9u3b3fuV9ikSRN++eUXtNbOsV73VL8aNWqEZKurCy64gIKCAgoKCpy96e+//57Dhw/73D+xc+fOdO7cmVGjRtG9e3fef/99unXrxuWXX8727dtjqtykiBM5OZCVZawevP56+PhjcNtoJNJkuCOKorXHYUpKCgMHDmTo0KF89dVXbNq0iT/+8Y+0bNnSuQVXWloav/76K+PGjWPnzp1MmjSJBQsWlPmctm3bsnnzZrZv386BAweC7mlfd911dOrUiXvuuYeNGzeybt067r33XlJTU7niiivKXb9r1y5GjRrF6tWr2b17N4sWLSIvL8+5f+Lzzz/P22+/zZgxY9i6dSvbtm3jgw8+4Nlnnw2qfaIK0BrGjIEHHjAC9L33wty5UQ/QIEE6qqK5x+Fbb71Fly5dGDBgAN27d0drzfz5853DEhdccAGTJ09m0qRJXHrppaxbt44RI0aU+YyhQ4fSsWNHrrjiCpo0aVJmJ/FAKKWYM2cOjRo1olevXlx33XWce+65fPjhhx6vr1OnDj/88AO33XYbHTp0YNiwYTz00ENkZWUBkJ6ezty5c1m0aBFdu3alW7duTJgwgTZt2gTVPpHgSkpg2DBwpIU+8wxMnw72f6lGm+xxKKo8+Rmqwo4ehUGDjHFnkwkmTzaGOyLA3z0OE29MWggh/LFvHwwYYGRy1K4NH3wAN98c7VaVI0FaCJGQrFYreXl5pKSklC/on5cHffvCTz/B2WfDZ59BN29l8aNLxqSFEAknJyenfP1nhzVrjBzon36Cc8+FVavKBGiPxf2jSIK0ECKhuJcZtdlsZGVlGUH300+hTx/47Te44gojQKekOO/1GdyjRIK0ECJhWK1WcnNzyy2UKi0t5fj48XDrrXDiBPTvDxYLuGRV+QzuUZQwY9LxkKUiYpP87CQGbzupAPxLKTpMmGC8yMyEN9+EamXDn6ca0qWlpezYsSP4TWpDIO6DdPXq1VFK8euvv9KkSRPn6jgh/KG15tdff0UpVW7puogf3nZSqQZkK8V9jr+IR4+G5583Via7cdSQdv0Ms9kc9ZWrcR+kzWYzycnJWK1W8vPzo90cEYeUUiQnJ5cpSCXii6decF3guw4daPPjj8b+g1OmGL1oLxw1pLOysigtLcVsNjNlypSo9qIhARazOJSWllaqAJCouqpXry4BOgJ8psSF4LPbtGnjDNTNgPlAZzAKn330kTEO7edn7dixw1kELVz8XcySMEFaCBG7XMeLTSYTU6dOJdNHrzbY78jKyiKltJQFQFuAc86BefOMTI4YI0FaCBET3Hu5YAxT5ufnh7ynun/WLBoNGUL1I0eM1LqFC41c6Bjkb5CWFDwhRFj5ypoIqZkzOefuu40A3a2bkQNtD9CxtkAlEBKkhRBh5b7zNoQha+L1141C/adOGfU3li41lnsTmwtUAiFBWggRVq47bwOhzZqw2WDkSHj0UaMm9IMPwsyZzl2SPC1QGTZsGLm5uXHTq5YxaSFERIQia6JMhkiTJnD//fC//xlv/vOfMGpUmRxoX1vQhWsC018ycSiESCiuGSINgA2tWtG+oMBYOZiTY+ym4sbTpKWrcE1g+kMmDoUQCcN12KIlsAJoX1DAEWDhww97DNBQfqjFXVgmMENMgrQQIuY5MkQuAlYDnYBCoBcw4PXXfY4vZ2Zmkp+fT25ubvgnMMNAgrQQIualpKSQphRfAa2AbUA34Fv86w0nJyeTkZERvgnMMIr72h1CiMSXvHIlS8xmzCUlrAAGAofs7wXSG87MzCQ9PT0iy75DRYK0ECK2jR8PTzyBGTjevz/zO3SgeOJEsBdBGjt2LHl5eQB+Bd3k5OS4CM4Okt0hhIhNNhs88QS89prx+uGHYcIEMJud6XwbNmzgqaeeCromSDiLPlVEUvCEEPHr5EkYPBg+/th4/corRsB2yYGubE2QSBR98kVS8IQQ8engQbjhBiNAV6/ObxMnYunSBevPP5e5rDI1QWJ1qyxPJEgLIWLH7t3QsyesWAH16zP/0Uc559FHPdbdqExNkIgVfQoBCdJCiNiwaRN07w7btkHLlvzy8cfcNH68195uZWqCRKToU4hIkBZCRJTHsqFLlsA110BhIVx0EaxezbZq1Srs7ToWqlgsFvLz8/0eUw5r0acQk4lDIUTYuGdPeJysq1nTKJRUUgKpqTB7NjRsGJHNAiK1VZYn/k4corWO+QOoD+iioiIthIgP2dnZ2mQyaUCbTCY9btw452vHMUoprY0io1rfeafWJ0+W+wyz2awBbTabdXZ2dpSeJvSKioocvw71tY/4Jz1pIYRHlckh9tQLNplMztcmYCLwZ8ebI0bAyy+DqfwIbDR7u+EkKXhCiKBVdjcTT9kTjiGO2sBMjABtAw6PHm3kQXsI0GCMH6elpSVUgA6EBGkhRBmhyCH2lj3x+vPPswy4BTgJWB58kIYvvBC6xicgCdJCiDJCkUPsKXvi/X/8g4fef59uwJl69Sj++GM6/u1vHjeIjeeNY0NNgrQQooxQ5RC7psf9/Omn3DFhAvz4I7RuTfW1a/ns8GGPQyrxvnFsyPmaVYyVA8nuECKiQppVsWCB1klJRgbHpZdq/fPPuqCgoFymh8lk0pMnTy533mw264KCgtA9XIyQ7A4hRKWEJKvirbdg6FAoLYXrrjN28q5f3+cGsZ5YLBbS0tKCa0OMkuwOIUSlOLIqgMDHh7WGF1+EP/3JCNCDB8O8eWDvZHkaUvEmVpdrR4oEaSGqIH8n5oIaHy4pgawscGRtjBoFM2ZAjRrOSyraINYhlpdrR4yvsZBYOZAxaSFCxn0loLfxZm/jxh9++KH3MeKjR7W+8UZj/Nlk0nryZJ9tKSgo0Lm5uR6/Jzc3NyHHoh38HZOOegD255AgLURoeAq83ibmli1bVuY69yBaLrjv26d1165GWKldW+vZs/1uVyIv//bG3yAtwx1CVCGB5ED7Gjcut8Blxw7o0QPWr4fGjWHpUhg40O92BVvNriqQIC1EFRJIDnRF48bO4L52rVEHeudOaNcOVq0yXgeoqi//9kaCtBBVSKB1lB093NzcXI/B/eJdu6B3bzhwALp0gdWroUOHsD9HVSJ50kJUQcHkQOfk5JCVlUVpaSlms5nld99Nz/feA5uNE717s37ECM695BKAqO3AHU9kt3AhRMhZrVZ25OVx+SefUH/iRAC29+zJJStXclprlH03b611VHbgjidhWcyilBqllFqvlDqilNqvlJqtlOpYwT1DlFLa7TgZyPcKIWJDctOmpM2Y4QzQxX/9KxfaAzTgmpEV0ztwx5NAx6RTgUlAN+B6oDqwSCmVVMF9xUBzl6NNgN8rhIi2I0dgwABjYYrZzMFXXiE7ORmbj3+Nx+oO3PGkWiAXa637ur5WSg0B9gNdgC9936p/Cbh1QojYUFgIN94I33wDdeqw6IEH6PfUU+XS+dxV9SXdoRBQkPaggf2/Byu4rq5SajdGz30j8IzWequ3i5VSNYGaLqfqVaqVQojg/fAD9O0Lu3dDkybs+7//o9/AgR4DtFIKpRQ2m02WdIdI0EFaKWUCXgNWaq2/83HpduBPwGaMoD4CWKWUukhr7W2wahQg2zUIEaTK7E9YxsqVcPPNcPAgtG8PCxfy/Z49HgP0hAkTuP322wESck/CqPG1HNHXAfwXyAeSA7yvOrAD+LuPa2piLAV3HC2RZeFC+MXf2hwOBQUFetmyZeWXhs+cqXXNmsYy76uu0nr/fuf1VaXmcziFtXYH8AZQALQL8v6PgP8FcL3U7hDCD4EGUK8BfeJErZUyQsRNN2l97Fi5+6parY1QC0vRf2UkQU4EbgXStNZ5gfbclVJmYCswX2v9uJ/3SJ60EH7wVkzfU9F8q9VKmzZtygxdVDOZOJiVRb3//tc4kZUFb7wB1cqPjIZkU4AqzN886UDHpCcBdwMDgSNKqWb280Va6xMASqm3gZ+11qPsr58H1mAMcTQERmKk4GUH+N1CiAo4anO4Bl5vGRbuxZaqA2/ZbL8H6H/+06gFbV+g4i45OVmCcwQEmic9HGPybzlQ6HIMcrmmNUYutEMjYBqwDZiPMXTRQ2v9fXBNFkJ4E0htDtdiS/WBBcA9gK5WDaZPh2ee8RqgReTIsnAhEpC/QxE5OTmMGTaMz2w2LgVO16xJjTlzID09co2toqR2hxCiYlu3UnLDDVTbu5fSc87BvHAhdO4c7VZVCbIRrRDCty++gJ49qbZ3L5x/Pua1ayVAxyAJ0kJURbm5cMMNcPgwXH21sWilbVu/N6gVkSNBWoiqZsIEGDQITp+GP/wBFi+Gs84KbmdwEXYyJi1EDAjZMm5f37FnD/rxx2k1c6Zx4uGHjYBtNnvMmTabzeTn50uaXZjImLQQcSISPdjpb77JqjZtnAF6bUYG/Oc/YE/VC2SDWhFZ0pMWIkqsViurVq3irrvuqrAHW5me9s/ffceOTp1IBU4DQ4Bct++QnnTkSU9aiBjm6D0PGjSowh5spXrae/bQcMAAUoEioC/wPw/fEegGtSJypCctRIR56rW6cu3BVqqHu3kz9OsHe/diBfoDWyr4DKnHETnSkxYiRnka/3Vw78EGPVa8bBlccw3s3QsXXcRX48bxvR+95OTkZNLS0iRAx5DK7swihAiQpyJIJpOJDz74gO7du5cJkBUVTPI4Vv3++zBkCJw5A6mpMHs2dzZsSM+77pJecjzyVcc0Vg6knrRIMIHUY3a/dty4cXrZsmX6lVdeKVsLeto0rV96yagBDVrfcYfWJ09G8KlEIMJSTzpaZExaJKJAxn8d127YsIGnPGwAawL+oxR/cfx5/utf4dVXwSQjmrFKCiwJkWC8TTjWAt4D/gBopVDjx8Njj0WjiSIA4Sr6L4SIEk+TiGcBnwJXA6eAo5Mm0Xj48Ci0ToSL/FtIiDjhWqQfoC2wEiNAHwKWPPmkBOgEJEFaiDjhuuCkM7AaOB842bQpp5Yu5caXX45yC0U4yHCHEHEkMzOTm2vVotHQoVQ7cQIuuYRaCxbQrEWLaDdNhIn0pIWIJ9On02TIECNAX3stfPklSIBOaBKkhYgHWsM//gH33w8lJXDPPTB/PhjZAV5JEf/4J0FaiFhXUgLDh8Nzzxmvn34a3n4batTweZsU8U8MkictRBAiUaQfgGPH4M47Ye5cUAomToSHHvKrfVJ6NLZJgSUhwiScPdQywxO//gp9+hgBulYtmDnTrwANUsQ/kUhPWogAhLOHmpOTw7Bhw7DZbKQoxYazz6b+r7/CWWfBZ59Bjx4x0U4RGtKTFiIMwtVDtVqtzgDdFfhKa+r/+islrVrBqlUBBWiQIv6JRPKkhQhARaVDg+UI/v2BXCAJ+Bo49dpr9OjYMajPzMzMJD09XcqTxjnpSQsRgHD1UFNSUhiqFHMwAvRC4FqTidZXXlnp9koR//gmPWkhAuRPDzWg7A+tSZ42jan2+aG3gOEmE5OmTpXgKiRICxGM5ORkrwHUdQLQZDIxdepUMjMzPX/QmTOQlQVvvQVA8aOP0m7gQHYEkdoXsbRAEVGS3SFECAWUVXH0KGRkwMKFRnH+//4Xhg0L6nsD+otBxATJ7hAiCvzO/vjlF2P/wYULoU4dmDMn6ADtmhkCYLPZyMrKkqXgCUKCtBAh5F7zGTxkf2zfDt27w8aNcPbZYLHAgAFBf6csXElsEqSFCKEKsz8cOc/5+XDeebB6NVQyg8OvvxhE3JIgLUSIZWZmkp+fj8ViIT8///ex4dmz0X36wMGDnL70UiNg+wik/lawk4UrCc7XVuKxcgD1AV1UVBSq3dSFiKxJk3SpUlqD/gx0XaV0dna218uzs7O1yWTSgDaZTD6vdSgoKNAWi0UXFBSEsuUiTIqKijSggfraR/yT7A4hKsln6pvW8Mwz8NJLAEwBHgJK8Z71IXU3qgbJ7hAiSIEUyvdZEe/0abj3XmeAfhZ4ECNAg/fJPZkIFK4kSAvhIpAypD5T34qL4cYb4d13wWzm4L//zVg/J/dkIlC4kiAthF2g+cbeerx71qyBXr1gyRJISoJ58zjr8cf9ntyTiUDhSsakhbCzWCz06dPH41SQEmcAAB16SURBVPm0tLRy5z2NHV9sMvFNs2ZU27sXmjaFefOgS5cy9/hblS6Qa0X88XdMWoK0EHbBTNjl5OSQlZVFaWkpqSYTn9eqRc3jx6FDB2M1Ybt2kWq+iDMycShEgIIZZnDkRH/3wgssM5upefw4p7p0gZUrJUCLkJCetBBuAh5meO019OOPo7TmE+CPSvH6tGlS4Ej4JMMdQoSbzQYjR8L48QC8ATwK2JC8ZlExGe4QIpxOnYK773YG6KeAhzECNEheswgdCdJCBOrwYUhPhw8/hOrV+e3113lV8ppFmEiQFlVCIKsIfSoogJ494YsvoH59WLCAxg8/LHnNImxkTFokvJDtWrJlC/TrBz//DC1awPz5cOmlzrclr1kEQiYOhSCExYqWLYNbbzWWe19wgZED3bp1GFosqgqZOBSCyhcrslqtbH32WXTfvkaAvuYa+OorCdAiYmS3cJHQHMWK3HvS/kzq5WRn8+OwYbxs/9fmT1dcwbmLFkGtWmFrrxDupCctElqwxYqsu3dzfOhQZ4CeAHTcuBHrgQPhbrIQZQQUpJVSo5RS65VSR5RS+5VSs5VSHf24L0Mp9YNS6qRSaotSqn/wTRYiMF63s/LmxAlqDB7Mw/aXfwUeB0psNsl9FhEX6HBHKjAJWG+/91/AIqXUhVrrY55uUEr1AP4HjALmAncDs5VSl2utvwu65UIEIDk52Wfv2bG7Ssezz6bF8OGcs3Ilp4DBwEf2ayT3WURDQEFaa93X9bVSagiwH+gCfOnltkeBhVrrV+yvn1NKXQ/8BWOjinKUUjWBmi6n6gXSTiEC4UjRa2WzsRBoAdCgAUuGDWPW+PFQWiq5zyJqKjtx2MD+34M+rukOjHc79zlwi497RgEvVKJdQvjFUej/EpuN+UBzYA9QY+ZMbrz2WvIfeaTSuc8+90AUogJBTxwqpUzAa8DKCoYtmgH73M7ts5/3ZizGXwCOQ36yRVjk5eVxrc3GlxgBehNGr+IH+0RjcnIyaWlpQQfXQLbjEsKTymR3TAIuBu4MUVuctNantNbFjgM4EurvEALg0m+/ZR7GeNpSoBewL0Rjz4FuxyWEJ0ENdyil3gAGAL201hX9xP0CNHU719R+XoiIs1qt5P34IykffUTym28C8L5SDNEaWwjHnn0tpJFhD+GvgIK0UkoBE4FbgTSt9S4/blsNXIsxNOJwvf28EBGVk5PD8KFDeV1retvPjQPU2LEsuuqqkNbdqMxCGiGctNZ+H8Bk4DBGKl4zl6O2yzVvA2NdXvcAzgBPAOcDo4HTwMUBfG99QBcVFWkh/FVQUKCXLVumCwoKnK+TlNJzQGvQpaAfAg1os9nsvC6UsrOztdlsdn5HdnZ2yL9DxKeioiKN8fNXX/uIf4EOdwy3/3e52/n7gen2/2/N77XP0VqvUkrdDfwDI686D7hFS460CCNPle86nnUWS7SmG3ACuAf4xH59uIYhMjMzSU9Pl+p4ImhSBU8kHE+V71JMJja3bEmtggJ+A24GVrncI9tdiUiTKniiynKfsLsC+Mpmo1ZBAUfOPpteJlO5AC0LVUSskip4Iu65LxZxnbDrD+QCScDpiy+m3uLFfF5Swo4dO0hKSuLYsWN+D0PIohQRFb4GrGPlQCYOhRfZ2dnaZDJpQJtMJufEXHZ2tn5AKX3GPklYcNFFWhcXh/x7hAiWvxOHMiYt4oprbxbwvOvKrl0k5+TAmDEAHMvIIOm996B69aC/MyS7uwjhwt8xaRnuEHHDNWNDKUVGRka5xSKqtBRzVhYsWGCcePZZkl58EZQK+ntlUYqIJulJi7jgqTfrLgmjrGg/AJMJJk+GrKywfLf0pEVlSXaHSCieerOuzsFI3u8HULs2zJ4dkgANwe/uIkQoSE9axAVfPekUYP1ZZ9Hg4EE4+2yYOxeuuiosbZBFKSJUpCctEoqjN2sylf2R7YaxKKXBwYNw3nmwalVYArSjDZUpWypEMCRIi7iRmZnJ7t27GTFiBGazmZuBZcDZAF27GgHanvUhRKKQ4Q4Rlw796180fO45lM0GN94IH34ISUnRbpYQfpPhDhF3rFYrFovFd1F8reGZZ2j0t78ZAXroUGOSUAK0SFASpEVM8GubqdOn4b77YOxY4/WLL8KUKVBN0v1F4pLhDhFVVquVVatWcdddd/nOQy4uhttugyVLwGyGadPg/vuj1GohKk9WHIqY57qC0F2ZFX1790L//rBpkzGs8fHH0LdvFFosRORJkBYR5+g9ewvQ4LLN1LZtRkDesweaNoV586BLlwi3WIjokSAtIspX79nBuaJv1y64+WY4fBg6dICFC6Fduwi2Vojok4lDETFWq9VngDaZTOTm5pKfn09mw4Zw/fVGgO7WDVau9CtA+5UhIkQckSAtIsZX/Q2z2czUqVPJyMggeeZMyMiAU6fglltg6VJjuXcF/MoQESLOSHaHiBhP9TdMJhMffPAB3bt3J7lFC3jySfj3v403//xneP11I5sjiM+WSnUilsliFhE13oYcPFWTc/aemzSBe+75PUC/9BK88QaYzX4NYfiq+SxEXPO1bUusHMj2WXHDn22mCgoKtMVi0QUFBcaJQ4e0TkvTGrSuXl3rd94J6PMcn+m4znGYzebfv0OIGOPv9llRD8D+HBKk40NQgXLPHq0vvtj4UaxXT+vFi4P+vOzsbG02m53XyT6EIpb5G6QlBU+ETMDbTG3ZAv36wc8/Q/PmMH8+XHZZ0J+XmZlJenq61HwWCUWCtAiZlJQUTCZTucm79u3bl7/YYjEyN4qL4fzzjRzoNm2C/zy75ORkCc4iocjEoQgZv7eZ+uADSE83AnTPnkYOtFuADujzhEhgkoInQs7rNlNaG9kbI0car2+7Dd59F2rVCu7zhIhj/qbgSZAWkVFaCo8/buQ9AzzyCIwf71cOtBCJSKrgidhx4gQMHgwzZxqvX33VCNhKRbddQsQBCdIivA4ehIED4auvoEYNmDED7rwz2q0SIm5IkBbhs3u3kWK3bRs0aGBsc5WWFu1WCRFXJEiL8Pj2W6NQf2EhtGxppNhdfHG0WyVE3JEUPBF6ixdDr15GgL74YlizRgK0EEGSIC1C6513jB70kSPG0MaKFSBpc0IETYK0CA2tjV28770XSkqMycGFC7EePSpF+IWoBAnSovJKS+Evf4FnnjFejxgB771HzrvvShF+ISpJFrOIyjl+HO6+G+bMMfKeX3sNHnnEa4H/NWvW0LVr1yg2WIjYIEX/RfgdOADXXmsE6Jo1ITfXWEmI5wp2NpuNbt26SY9aiABIkBbB+eknzlx1FaxZg61BA1iyBG6/3fm2o4KdO5vNRlZWVoVj1LKhrBAGCdLCK6+B8uuvOXHZZVT/6Sd2A52Ki8nZvr3MJY4Kdp4CdUXbWsmGskK48LUjQKwcyM4sEed126oFC3RpnTpag94IulkFO6asW7cuoN1VZBssUVX4uzOL9KRFOVarlWHDhjnHlB1DFAdffRUGDMB0/DiLgFTgF/s93nrHXbt2DagmtGwoK0RZsixclOMpUP6ttJSz7HWgj912GzfPmsUpl8wgXzumBLKtVTC7sQiRyKQnLcpxnfQzA1OBMY43n3mGpI8+YtK0aQHtmJKcnExaWlqFRftlNxYhypI8aVGG1WolLy+Pr7/+mhefeor3bTYGADalME2aBMOHl7k2XDumyG4sItFJ0f8E4wieKSkpYQtaOTk5zrHopkqxvWVLmlutlNasyffPPkujm27C9ZvDuemrbCgrhEGGO+JAJFLSXCcL2wMrtaa51crxOnW45vRpLnnuOUmHEyIKZLgjxnlaXm02m8nPzw9pT9NisdCnTx+uBOYCTYCfgP5Ksd1tgjDU3y1EVSTLwhNEpFLSUlJSuFkpLBgBej1wtVuAdnz3Rx99JCsBhYgQCdIxztPy6nCkpCXPm8cnQB1gPnCdycTjL7/sccXg448/LkMfQkSIBOkYF/aUNK3huefgwQcxac2xO+8kackStu7ezciRI8t8tyt/a3AIISrJ13JETwfQC/gM2IuxpPGWCq5Pw2WJr8vRLIDvrPLLwgsKCrTFYgnt8ujTp7W+7z6tjVCt9ejRWttsHr97/Pjxnn4PtcViCV17hKhCwrksPAnYBDwU4H0dgeYux/4gvrvK8ncxiN+OHIEBA2DGDDCbITsbXnjBqAnt4bszMjJ8DrtI1TohwiPgIK21XqC1flZr/UmAt+7XWv/ictgqvkWERWEhpKbCokVQpw58+ilkZvq8xdewi1StEyJ8KpWCp5TSwK1a69k+rkkDLMBuoCbwHTBaa73Sxz017dc61AOsVTEFL+R++AH69oXdu6FJE5g3DwLYKcV9JWCkUgSFSDSxlIJXCDwI3GY/CoDlSqnLfdwzCihyOarMv6HDOmywciVcfbURoFNSYPXqgAI0lB92kap1QoRX2IO01nq71nqK1vprrfUqrfWfgFXAX33cNhZo4HJUiS5ZqIYNPAb6WbPguuvg4EG46iojYJ93XqXbHKkUQSGqqmil4K0DvP4p1lqf0loXOw7gSOSaFh3eajgH2qP2GOgnTjS2tjp5Em6+GZYtM4Y6QkCq1gkRXtEqsHQZxjCIsPM1bOBvwHMP9Npm4+ADD/x+wYMPGgG7Wmh/2wOpFy2ECEzAf1qVUnUp2wtup5S6DDiotd6jlBoLtNRa32u//jFgF7AVqAU8APQBbqhs4+ORt2p2oSh27xroawBvAXfb3/ubUpzbpQuZHgJ0KCrsSdU6IcLEVxK1pwPvi1Om29+fDix3uf5JYAdwAvgNI9Ojd4DfmRCLWbzuG+jyvtlsdu7r5/5+RRz7A9YHvdS+QOU06ME+9gqsqE1CiPDwdzGLVMGLEH9T1QItdu/eC/7fK69w0ZNPcgnGQP5twGKX6y0WC2lpaQG1SQgRelL0P8b4O+YcyLCBa5F+k8nERy+8wF3TpgFwomFD0oqK2OhjH8JQjIMLIcJLCixFSKhT1dwnCXvabPR54QWwWuH886n9zTf8uYJ9CCV9TojYJ0E6QkKZqma1WsnNzXUG6DuARUBDoOjii40c6LZtyczMJD8/H4vFQn5+PpluS78lfU6I2Cdj0hFW2Q1WXYc4wFgRNN7+3izgyh9/JDklJaJtEkIEzt8xaQnSccR1ok8B/+b3ZZuvK8XpsWPpcuWVYd2sVggRGjJxmIAcE301gXeADPv5FTfdxOlrruGpp592TiJOnTq13PCGECL+SE86jlitVi5t3ZpPtKYXcBq4XykeW7uWbt26SSqdEHEklqrgCS8CrXiXbLOxo1kzemGUBuxvMtFn2jSOHj0qleiESFASpKMk4Ip3mzZBt240KiykpFkzfszJYfru3WRmZkoqnRAJTIJ0FARc8W7pUrjmGmNHlYsuotq6dXT905+cQxmSSidE4pKJwygIaKXfe+/B/ffDmTOQlgaffAING5b7TKlEJ0RikiAdQv5Wk/Or4p3WMG4cPP208XrQIGPT2Jo18UYq0QmReGS4I0QCGWOucHiitBT+8pffA/QTT8D77/sM0IGS3b2FiBO+SuTFykGMlyp1lAjFpXSrp7Kgnu6zWCxlrzt+XOtbbtEatFZK6wkTQt5eKU8qRPT5W6pUhjtCINhqcuWGJ377DW66ydggtmZNePddY9urELFaraxatcrjpGV6eroMlQgRg2S4IwRCkgK3axf06GEE6IYNYfHikAZox3DMoEGDJKdaiDgiQToEKp0Ct3EjdO8OP/4IrVsbVeyuuSZk7XNP+XMnOdVCxC4Z7giRoFPgFi40eszHjsEll8CCBdCiBRCavQfB83CMg+RUCxHbJEiHUMApcG+9BUOHGtkc110HM2eCvTaJ+64rlSmY5Cnlz2Qy8cEHH9C9e3cJ0ELEMBnuCEDI0ta0hhdfhD/9yQjQgwfDvHnOAB3wisQKeBqOmTp1KhkZGRKghYhxEqT9FHCtDW9KSiArC154wXg9apSxSKVGDeclvrJFglXRLi1CiNgkpUr9ELJdtY8dM1YOzpsHJhO88QYMHx6+7xNCxCwpVRpC3nq2H330kf9DEPv3Q+/eRoCuVQtmzfIYoEEKJgkhfic9aT946tk6+DWpt2MH9O0LO3dC48bw2WdGyp0f3ysFk4RITLLHYYjl5OSQlZVFaWlpufd8DkWsWwc33ggHDkC7dkaKXceOEWixECKWyXBHiDkm3saPH1/uPa+TenPnGuVFDxyALl1g1SoJ0EKIgEiQDkBycjIZGRn+LQGfOhUGDoQTJ4yhDosFmjWLYGuFEImgSgbpyuQ7VzippzU895yRZmezGQX7P/0U6tUL5SMIIaqIKjcmHaqVfB4n9c6cMVYQzphhvH7+eRg9GpSqVJuFEIlHJg49CEf+saO+RofmzWn52GPw+edgNsObb8IDDwTdViFEYvM3SFep2h3B1n32xtErb2KzMQ9oCVCnDuTmGhkdQghRSVVqTDokdZ/tHPU12ttsrAa6APuBfR9+KAFaCBEyVSpIh3IlX15eHlfZbKwC2gF5QHdgW9265a6V/QSFEMGqUkEaQldoqNOOHSwFGgNrgB7Abg+98pAVZhJCVElVauIwZCZNgocfBq35DBgEnLb3yh1B37Gf4F133SWFkoQQ5cjEYTjYbPDMM/Dyy8brYcPoPGoU8/Pzy6Tiuab5uavMRKUQouqRnrS/Tp82ivS/957x+h//MAK2Ww60r2JMID1pIYRBetKhVFQEf/gDLFsG1apBdjbcd5/HS2U/QSFEKEmQrsjPP0P//rB5M9StCx9/DOnpXi+X/QSFEKFU5bI7ArJ1q1H3efNmozjSl1/6DNAg+wkKIUJLxqS9+fJLo4rd4cNGedGFC6FtW79vl4L9QghfZEy6En77739p9OijmM6cgR49jCp2jRsH9BnJyckBBWdHDZCUlBQJ6kIIJxnucLNm0CAa//nPmM6c4RNgwo03Ytm8OayrBWXBixDCGxnucLDZOJKVRb3sbABeB/4KOKb/KlPW1BfZGVyIqkm2zyKAmhknT8JddzkD9EjgUX4P0AA2m42srKyQ96h9VeYTQoiEDdJ+DyEcOmRkbOTmoqtX549K8aqXzwxH8AxlZT4hROJJyCDtKCPq6KF67QUXFEDPnvDll9jq1ePbsWO5bNw4Z/qcu3AEz1BW5hNCJJ6EzO7wq7j/5s3Qrx/s3cuxhg25uqiITSNGYDKZeOmll+jatSsbNmzg6aefprS0NKzBMzMzk/T0dEnZE0KUk5AThxVOxi1bBrfeCsXFnOnQgZS8PHa7/Dq4Xiv5zkKIcKjSE4c+hxDefx/69oXiYkhNZc0rr5QJ0FB27Dk5OZm0tDQJ0EKIqEjInrRDmV5wy5Ywbhw8/bTxZkYGvP021gMHJAVOCBFxYetJK6V6KaU+U0rtVUpppdQtftyTppTaqJQ6pZTaoZQaEuj3BsPZC27e3CjS7wjQf/0rfPAB1KolE3dCiJgWzMRhErAJ+D9gVkUXK6XaAfOAN4F7gGuBbKVUodb68yC+PzAnTsA998Annxivx483grQLmbgTQsSqgIO01noBsABAuRW89+JBYJfW+gn7621KqZ4YC/rCG6R/+w1uvhlWrULXqMHWp5+mYUYGnkJwoLU2hBAiEiIxcdgdWOJ27nP7eY+UUjWVUvUdB1Av4G/dtQuuvhpWreJUnTr0OXOGTi++KLUxhBBxJRJBuhmwz+3cPqC+Uqq2l3tGAUUuR+BrsTdtgh9/pKRFC644cYLl9gnScC3vFkKIcIjVFLyxQAOXI/BxiFtugXffZe2ECXznI8VOCCFiWSRWHP4CNHU71xQo1lqf8HSD1voUcMrx2s+x7/Luvps2Vmu57aykNoYQIl5Eoie9GiOjw9X19vNhJyl2Qoh4FvBiFqVUXcDRDf0GeBywAAe11nuUUmOBllrre+3XtwO+AyZhpO31wSjXfKO/KXihqCftaXm37IYihIiWcC4LvwIjOH9jfz3e/v8v2l83B1o7LtZa7wJuxOg9bwKeAB6ISI60C/fl3bIbihAiHiT0snBvZDcUIUS0VekCSxWR3VCEEPGiSgZp2Q1FCBEvqmSQlowPIUS8qJJj0g5S0F8IES3+jkkn5PZZ/pKiSkKIWFclhzuEECJeSJAWQogYJkFaCCFimARpIYSIYRKkhRAihkmQFkKIGCZBWgghYlhc5UkXF3vN9xZCiLjibzyLlxWHLQlmn0MhhIh9yVrrn729GS9BWgEtgCMB3loPI7gnB3FvPEj054PEf0Z5vvhXmWesB+zVPgJxXAx32B/A69803rjsjXjE19r4eJXozweJ/4zyfPGvks9Y4fUycSiEEDFMgrQQQsSwRA/Sp4Ax9v8mokR/Pkj8Z5Tni39hfca4mDgUQoiqKtF70kIIEdckSAshRAyTIC2EEDFMgrQQQsQwCdJCCBHD4jZIK6VGKaXWK6WOKKX2K6VmK6U6+nFfhlLqB6XUSaXUFqVU/0i0N1DBPJ9SaqhSaoVS6pD9WKKUujJSbQ5UsL+HLvffqZTSSqnZ4WxnsCrxM9pQKTVJKVWolDqllPoxFn9OK/F8jymltiulTiilCpRSE5RStSLR5kAppYYrpTYrpYrtx2qlVL8K7glpjInbIA2kApOAbsD1QHVgkVIqydsNSqkewP+AHKAzMBuYrZS6OPzNDVjAzwekYTxfb6A7UGC/p2V4mxq0YJ4RAKVUW+BVYEUY21dZwfyM1gAWA22B24GOwFCCKIsQAcE8393ASxh5xRcAmcAg4F9hb21wrMDTQBfgCmAZMEcpdZGni8MSY7TWCXEATQAN9PJxzYfAXLdza4A3o93+UDyfh3vMGLUB7o12+0P5jPbnWonxB3w6MDvabQ/V8wEPAjuB6tFub5ie7w1gqdu5fwNfRbv9ATznQSDTy3shjzHx3JN218D+34M+rukOLHE797n9fKzz5/nc1cHo3QRyTzT5+4zPA/u11jlhbk+o+fN8NwOrgUlKqX1Kqe+UUs8opczhb16l+fN8q4AujmE4pdS5QH9gfpjbVmlKKbNS6k4gCeP3yJOQx5i4qIJXEaWUCXgNWKm1/s7Hpc2AfW7n9tnPx6wAns/dy8Beyv/QxBx/n1Ep1ROjB31ZpNoWCgH8Hp4L9AHewwhe7YHJGH/Zjgl3O4Pl7/Nprd9XSp0NfGUvQVwNo5cZq8MdKKU6YQTlWsBR4Fat9fdeLg95jEmIII0xLnYx0DPaDQmTgJ9PKfU0cCeQprU+Ga6GhVCFz6iUqge8AwzVWh+IVMNCxN/fQxOwHximtS4FvrbPKYwkhoM0fj6fUioNeAb4M7AW4y+h/yilntNa/z3cjQzSdoxOQQOMeYIZSqlUH4E6pOI+SCul3gAGYIyDVbR7yy9AU7dzTe3nY1KAz+e4ZwTGZMd1WuvN4WxfKATwjOdhTKh95lLD12T/jBKgo9Z6ZxibGpQAfw8LgTP2AO2wDWimlKqhtT4drnYGK8Dn+zvwjtY62/56i32icapS6p9aa1s42xoM+6/5DvvLr5VSXYFHgSwPl4c8xsTtmLQyvAHcCvTRWu/y47bVwLVu567H+/hS1AT5fCilngSeA/pqrTeEs42VFcQz/gB0wujVOI5PAYv9/wvC2NyABfl7uBJobx8+cOgAFMZagA7y+eoA7oHY8ReSIj6YgJpe3gt9jIn2TGklZlgnA4cx0oCauRy1Xa55Gxjr8roHcAZ4AjgfGA2cBi6O9vOE6PmewiiXeJvbPXWj/TyhekYPnzGdGM3uCPL3sBVGRs5EjOB8I8aY5t+i/Twher7R9ue7E2iHEcB2AB9G+3m8PONYoBfGv+A62V/bgOu9PF/IY0zUfxEq8YunvRxDXK5ZDkx3uy8DY4zpFPAd0D/azxKq5wPyvdwzOtrPE8rfQ7fPiOUgHezPaHeMtK2TGOl4zwDmaD9PiH5GqwEv2APzCWAPxnh2w2g/j5dnzLH/uTqFMVewxBGgffz+hTTGSD1pIYSIYXE7Ji2EEFWBBGkhhIhhEqSFECKGSZAWQogYJkFaCCFimARpIYSIYRKkhRAihkmQFkKIGCZBWgghYpgEaSGEiGESpIUQIob9P2NI6txY/SE3AAAAAElFTkSuQmCC\n", "text/plain": [ "