diff --git a/simulations/03_emitter_receiver_simulation.ipynb b/simulations/03_emitter_receiver_simulation.ipynb index bda4f08..a5e5b39 100644 --- a/simulations/03_emitter_receiver_simulation.ipynb +++ b/simulations/03_emitter_receiver_simulation.ipynb @@ -9,21 +9,23 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ - "%matplotlib inline\n", - "\n", "import numpy as np\n", - "import scipy.fft as ft\n", - "import matplotlib.pyplot as plt\n", - "import matplotlib.gridspec as gridspec\n", - "import matplotlib.ticker as tck\n", + "from functools import partial\n", "\n", - "rng = np.random.default_rng()" + "import matplotlib.pyplot as plt" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, { "cell_type": "markdown", "metadata": {}, @@ -33,156 +35,42 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ - "# copied from 02_discrete_signal_translation.ipynb #ae7aba6\n", - "class TravelSignal:\n", - " \"\"\"\n", - " Model an arbitrary digitised signal that can be translated to another position and time.\n", - " \"\"\"\n", + "from lib.TravelSignal import TravelSignal\n", "\n", - " def __init__(self, signal, sample_rate, t_0 = 0, x_0 = 0, periodic=True, interp1d_kw = None):\n", - " \"\"\"\n", - " Initialise by saving the raw signal\n", - " \n", - " Parameters\n", - " ----------\n", - " signal : arraylike\n", - " The raw signal to wrap.\n", - " sample_rate : float\n", - " Sample rate of the raw signal.\n", - " t_0 : float, optional\n", - " Time that this signal is sent out.\n", - " x_0 : float, optional\n", - " Location that this signal is sent out from.\n", - " periodic : bool, optional\n", - " Translated signal is 0 if it is not periodic\n", - " and the time/distance is outside the samples.\n", - " interp1d_kw : bool or dict, optional\n", - " Use scipy.interpolate's interp1d_kw for interpolation.\n", - " Set to True, or a dictionary to enable.\n", - " Dictionary will be entered in as **kwargs.\n", - " \"\"\"\n", + "####\n", + "from scipy.stats import norm\n", "\n", - " self.raw = signal\n", - " self.periodic = periodic\n", + "sample_rate = 3e2 # Hz\n", + "interp_sample_rate = sample_rate * 1/10 # Hz\n", "\n", - " self.sample_rate = sample_rate # Hz\n", - " self.sample_length = len(self.raw)\n", - " self.time_length = self.sample_length*sample_rate # s\n", - " \n", - " self.x_0 = x_0\n", - " self.t_0 = t_0\n", + "t_offset = 8\n", + "periodic = False\n", "\n", - " # choose interpolation method\n", - " if not interp1d_kw:\n", - " self.interp_f = None\n", + "time = t_offset + np.arange(0, 1, 1/sample_rate) #s\n", + "time2 = t_offset + np.arange(-1.5, 1, 1/sample_rate) #s\n", "\n", - " else:\n", - " # offload interpolation to scipy.interpolate\n", - " import scipy.interpolate as interp\n", + "signal = norm.pdf(time, time[len(time)//2], (time[-1] - time[0])/10)\n", "\n", - " interp1d_kw_defaults = {\n", - " \"copy\": False,\n", - " \"kind\": 'linear',\n", - " \"assume_sorted\": True,\n", - " \"bounds_error\": True\n", - " }\n", + "if False:\n", + " mysignal = TravelSignal(signal, sample_rate, t_0 = t_offset, periodic=True)\n", + " mysignal2 = TravelSignal(signal, sample_rate, t_0 = t_offset, periodic=False)\n", "\n", - " if self.periodic:\n", - " interp1d_kw_defaults['bounds_error'] = False\n", - " interp1d_kw_defaults['fill_value'] = (self.raw[-1], self.raw[0])\n", - " \n", - " # merge kwargs\n", - " if interp1d_kw is not True:\n", - " interp1d_kw = { **interp1d_kw_defaults, **interp1d_kw }\n", + " fig, ax = plt.subplots(1, 1, figsize=(16,4))\n", + " ax.set_title(\"Raw and TravelSignal\")\n", + " ax.set_ylabel(\"Amplitude\")\n", + " ax.set_xlabel(\"Time\")\n", "\n", - " self.interp_f = interp.interp1d(\n", - " np.arange(0, self.sample_length),\n", - " self.raw,\n", - " **interp1d_kw\n", - " )\n", - " \n", - " def __len__(self):\n", - " return self.sample_length\n", - " \n", - " def __call__(self, t_f = None, x_f = None, **kwargs):\n", - " \"\"\"\n", - " Allow this class to be used as a function.\n", - " \"\"\"\n", - " return self._translate(t_f, x_f, **kwargs)[0]\n", - " \n", - " def _translate(self, t_f = None, x_f = None, t_0 = None, x_0 = None, velocity = None):\n", - " \"\"\"\n", - " Translate the signal from (t_0, x_0) to (t_f, x_f) with optional velocity.\n", - " \n", - " Returns the signal at (t_f, x_f)\n", - " \"\"\"\n", - " \n", - " if t_0 is None:\n", - " t_0 = self.t_0\n", - " \n", - " if velocity is None:\n", - " velocity = 1\n", + " ax.plot(time, signal, label='Raw signal')\n", + " ax.plot(time2, mysignal(time2)+0.5, '.-', label='TravelSignal(periodic)+0.5')\n", + " ax.plot(time2, mysignal2(time2)-0.5, '.-', label='TravelSignal-0.5')\n", "\n", + " ax.legend()\n", "\n", - " ## spatial offset\n", - " if x_f is None:\n", - " spatial_time_offset = 0\n", - " else:\n", - " x_f = np.asarray(x_f)\n", - " if x_0 is None:\n", - " x_0 = self.x_0\n", - "\n", - " spatial_time_offset = np.sum(np.sqrt( (x_f - x_0)**2 )/velocity)\n", - "\n", - " ## temporal offset\n", - " if t_f is None:\n", - " temporal_time_offset = 0\n", - " else:\n", - " t_f = np.asarray(t_f)\n", - " \n", - " if t_0 is None:\n", - " t_0 = self.t_0\n", - " \n", - " temporal_time_offset = t_f - t_0\n", - "\n", - " # total offset\n", - " total_time_offset = spatial_time_offset + temporal_time_offset\n", - " n_offset = (total_time_offset * sample_rate )\n", - "\n", - " # periodic signal\n", - " if self.periodic:\n", - " n_offset = n_offset % self.sample_length\n", - "\n", - " # non-periodic and outside the bounds\n", - " else:\n", - " mask_idx = np.nonzero( (0 > n_offset) | (n_offset >= self.sample_length) )\n", - " n_offset[mask_idx] = 0\n", - "\n", - " # offload to scipy interpolation\n", - " if self.interp_f:\n", - " amplitude = self.interp_f(n_offset)\n", - " \n", - " # self written linear interpolation\n", - " else:\n", - " n_offset_eps, n_offset_int = np.modf(n_offset)\n", - " n_offset_int = n_offset.astype(int)\n", - "\n", - " if True:\n", - " amplitude = (1-n_offset_eps) * self.raw[n_offset_int] \\\n", - " + n_offset_eps * self.raw[(n_offset_int + 1) % self.sample_length]\n", - "\n", - " # use nearest value instead of interpolation\n", - " else:\n", - " amplitude = self.raw[n_offset_int]\n", - "\n", - " if not self.periodic:\n", - " amplitude[mask_idx] = 0\n", - " \n", - " return amplitude, total_time_offset" + " plt.show();" ] }, { @@ -194,7 +82,7 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -217,19 +105,19 @@ " self.x[key] = val\n", "\n", " def __add__(self, other):\n", - " if isinstance(other, self.__class__):\n", + " if isinstance(other, Location):\n", " other = other.x\n", "\n", " return self.__class__(self.x + other)\n", "\n", " def __sub__(self, other):\n", - " if isinstance(other, self.__class__):\n", + " if isinstance(other, Location):\n", " other = other.x\n", "\n", " return self.__class__(self.x - other)\n", " \n", " def __eq__(self, other):\n", - " if isinstance(other, self.__class__):\n", + " if isinstance(other, Location):\n", " other = other.x\n", "\n", " return np.all(self.x == other)" @@ -237,30 +125,149 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "### Receiver\n", "class Receiver(Location):\n", " \"\"\"\n", - " Receive a signal at position x and time t\n", + " A location able to trace a signal over time.\n", + " \n", + " Optionally applies a transformation to the traced signal.\n", " \"\"\"\n", - " pass" + " def __repr__(self):\n", + " return \"Receiver({})\".format(repr(self.x))\n", + " \n", + " def recv(self, travel_signal: TravelSignal) -> TravelSignal:\n", + " \"\"\"\n", + " Return a function that traces the signal as a function of time\n", + " at the receiver's location\n", + " \"\"\"\n", + " return partial(travel_signal, x_f=self.x) " ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 5, "metadata": {}, "outputs": [], "source": [ " ### Emitter\n", - "class Emitter(TravelSignal):\n", + "class Emitter(Location):\n", " \"\"\"\n", " Emit a signal from position x_0 (and time t_0)\n", " \"\"\"\n", - " pass" + " def emit(self, travel_signal: TravelSignal) -> TravelSignal:\n", + " return partial(travel_signal, x_0=self.x)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Testing" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPIAAACqCAYAAACXtRI+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAQuklEQVR4nO3dfbRVdZ3H8fcHEVAkynhIVCBdpqIktUgrZ1KzHNEpnFU2miU5FqmZzZpaiY3NOGOm0zRrKpXMykkzdFyrYaLGHkwzKx9hYkxERsYQCIGLSmIqBnznj9/vxu54z72Xyz133/Pj81rrrnP2w9n7ux8+e//2Pufco4jAzNrbkLoLMLOd5yCbFcBBNiuAg2xWAAfZrAAOslkBHOSdIOlgSb+UtEnSBS2e16ckfa2V88jzGS7pYUmv6mG810q6u9X1VOb3AUk/H6j5tZteBVnSaZLuk/Q7Sevz8/MkqdUF7ghJKyS9bQBn+UngzogYFRFf6qKeOyW9IOnZyt93+zKjiPhsRHwwT3eypJA0tDKv/trRZwN3RcTaHup5ENgo6R39MM9+IWlkXse37uDrjpW0ulV1DYQegyzp48AXgX8GXgWMB84BjgaGtbS6flbd8fvJJGBJD+OcHxF7Vf4GzY5fVVk3Hwa+2cuXfSuPP1i8G9gMnCBpn7qLGVAR0fQPGA38DnhXD+MNBz4PrATWAdcAe1SGfwhYDjwFLAAmVIYFcB7wKLAJuBQ4ELgHeAa4BRhWGf/PgcXARuBu4LW5/zeBbcDzwLOks+XkPP2zc213Af8FfLSh/geBU5os2ztJYd0I3AkcmvvfAWwFXsjze00Xr70T+GCT6R4LrM51rgeeAE4BTgL+N6+rT1XGvwS4MT9fmZfr2fz3plzH1ty9saftUpn/hcDavP4m5vU3tDLfk4CH87b5DfCJyrB98/jDmyzjWcDS/NrHgA93sfwfryz/WZXhr8z7yjPA/Xm/+HkP++EdwGXAf1frzMNWAJ/I2/q3wL8DI4CReRm2VdbnBNJJbg7wf8CTpP1w7zytzv1qVl63G4C/bdhWtwA35GVfAkyvDO+c7qa8bv+ioda/yuvtaeCHwKTuljsiegzyicCW6oZtMt4X8krfGxgFfBe4PA97a17Q1+cd60pS060a5AXAy4DDSEfU24EDSAeSh4FZedzX541+FLBbXpErOnek/PxtlWl3rvAb8gbbA3gPcF9lnCPyhhrWxXK9hnQgezuwOyl0yzvHpZug9jLIW4C/y9P+ENABzMvr8DBSOA/oIsidy1UN3Ado2NF72C6d8/+nvF32AE4GljRM4wngT/PzVwCvbxj+DPlg2sUynkw6KAs4Bniu8/WV+f9jXv6T8vBX5OE3k8IwEjicdBBpGmTSQWgbMIV0cHiwiyDfTwrp3qSgnFM9qDSM/9fAvcB+ef18BbipYf1/Na+3I0j7bedB/pK87U4i7aeXA/dWpn0q2w8Wf0nax/bJw04h7WOHAkOBi4G7dzbI7wPWNvS7m3R2eh54S95IvwMOrIzzJuDX+fnXgc9Vhu0F/B6YXAny0ZXhi4ALK93/AnwhP/8ycGlDPcuAY3oI8gGVfsNJZ7uDcvfngblNlv/TwC2V7iF5hzp2B4L8XF5fnX+XVnae54HdcveoXOtRDevilL4EuRfb5VjgRWBEZfgZ1R0u91tJaj6/rMky/gZ4S087Wh73P4GPNSx/dRnWA28k7fy/Bw6pDPss3Qf5YmBxfj6B1Dp5XUOQ31fp/hxwTTdBXgocX+neJ9c0tLL+96sMvx84rbKtflwZNgV4vpvaFwMz8/PvA2c37HPP0cNZuadr5CeBMdVry4h4c0S8PA8bAowF9gQWSdooaSPwg9y/c6U+Xnn9s/m1+1bms67y/PkuuvfKzycBH++cT57X/nke3VlVmf9m0pH+fZKGAKfT/JqwsfZteVr7Nhm/KxdExMsrf5+uDHsyIrbm58/nx2bLvqN62i4AHRHxQqX7adIBpepdpDPL45J+KulNDcNHkQ5QLyFphqR7JT2V538SMKYyypMRsaXS/RxpeceSArOqMuxxuncm6ZqdiFgD/JTUYquq3sDrnFczk4D5lXW3lHRwGN/L6TUOG9GZI0lnSlpcmfbhbF8vk4AvVoY9RTood7vP9RTke0hNhpndjLOBtMMdVtlZR0dE50KtycWRF2Ik6frnNz3MuyurgMsagrFnRNyUh0eT1zX2v5509jkeeC4i7mnyusbaRTpw9KX2/tTVcjb262m7dPWaB4EDGg7cD0TETGAc6Yx6S+cwSRNINzyXNRYjaTjwbVKLZ3w++N9K2il70kFqdu9f6Tex2ciS3gwcBFwkaa2ktaTLr9N7eYOzq/W5CpjRsK+NiIid2vaSJpGa5OcDr8zr5SG2r5dVpHsJ1fnuERHdvtXXbZAjYiPwD8BcSe+WtJekIZKmka5dOs9SXwX+VdK4XOy+kv4sT2YecJakaXnjfpZ0jbpih9dCms85ko5SMlLSyZI6zyLrSNfW3crB3UZqtnd3h/YW4GRJx0vanXTttZl0eVGnDlL91WVdB+wnaRj0aru8RESsJt10PDKPP0zSGZJGR8TvSdfDWysvORa4I7dyGg0jXcZ0AFskzQBO6M3C5VbKfwCXSNpT0hReenatmgXcRmrCTst/h5NaJDN6Mct1wCslja70uwa4LAcPSWMldXdC662RpANHR57uWbnW6nwvknRYHj5a0qk9TbTHt58i4nPA37D97uo60oX/hWzfoS8kXaDfK+kZ4MfAwfn1t5OuNb9NunFyIHBaT/NtUstC0k2hq0jNwOWka8NOlwMX52bJJ3qY3A3AVODGbua3jHSf4ErSGe4dwDsi4sUdKPuqhveRF+3Aa5vV9Rzp7uwv8rK+kXTHdgmwVtKGPGrT7dKNrwDvr3S/H1iRX38OaX10OoO043VV4ybgAtLB8GngvaQbb711Pqmpuhb4BvBvXY0kaQTpBuaVEbG28vdr0kG6uwNAZ62PADcBj+X1OYH0lusC4EeSNpFufB21A/U3m9fDpBPIPaQsTQV+URk+n3QD8ua8zh+iFwcj5QvqXY6kM4HZEfEnddcymORW0y9JN3qe6Ga8qcC1EdF4zWw12CWDLGlP0hlsbkTcUHc9Zjtrl/usdb5G7CA1a+bVXI5Zv9glz8hmpdnlzshmJXKQzQrQ398GaktjxoyJyZMn112G9dGiRYs2RMTYnscsV1sHWdJ1pG9DrY+Iw3O/vUnfbJlM+nzteyLi6e6mM3nyZBYuXNjaYq1lJPX08c3itXvT+hukb2hVzQFuj4iDSN+imjPQRVnrzF08t+4SBqW2DnJE3EX6UHnVTNJnqcmPpwxoUdZSX/6fL9ddwqDU1kFuYnznJ5Ly47ia6zFruba+Rt4ZkmaT/j8VEyc2/WKNDQJzF8/9ozPx1OunAnDuEedy3rTz6iprUGn7D4RImgx8r3Kzaxnpi/9P5P/bdGdEdPtFgenTp4dvdrWHqddP5VezfvVH/SQtiojpNZU0KJTYtF7A9m+8zAK+U2MtZgOirYMs6SbS18EOlrRa0tnAFcDbJT1K+l9bV9RZo/Wvc484t+4SBqW2vkaOiNObDDp+QAuxAeNr4q619RnZzBIH2awADrJZARxkswI4yGYFcJDNCuAgmxXAQTYrgINsVgAH2awADrJZARxkswI4yGYFcJDNCuAgmxXAQTYrgINsVgAH2awADrJZARxkswI4yGYFcJDNCuAgmxXAQTYrgINsVgAH2awADrJZARxkswI4yGYFaOtfY+yOpBXAJmArsGVX/yFsK1uxQc6Oi4gNdRdh1mpuWpsVoOQgB/AjSYskza67GLNWKrlpfXRErJE0DrhN0iMRcVfnwBzu2QATJ06sq0azflHsGTki1uTH9cB84MiG4ddGxPSImD527Ng6SjTrN0UGWdJISaM6nwMnAA/VW5VZ65TatB4PzJcEaRnnRcQP6i3JrHWKDHJEPAYcUXcdZgOlyKa12a7GQTYrgINsVgAH2awAfQqypEslfazSfZmkC/qvLDPbEX09I38dmAUgaQhwGvCt/irKzHZMn95+iogVkp6U9DrSe7a/jIgn+7c0M+utnblG/hrwAeAs4Lp+qcZe6ieX112BtYGdCfJ84ETgDcAP+6cce4mfXlF3BdYG+vzJroh4UdJPgI0RsbUfazKzHdTnIOebXG8ETu2/cgxIzenqmfiS0enxmDlw3EX11GSDWp+CLGkK8D1gfkQ82r8lGcddtD2wl4yGS35bbz026PX1rvXDwAH9XIuZ9ZE/2TXYHTOn7gqsDTjIg52via0XHGSzAjjIZgVwkM0K4CCbFcBBNiuAg2xWAAfZrAAOco3WP/MC7/nKPazf9ELdpVibc5Br9KXbH+WBFU/xpR/74+q2c4r8B/WD3cEXf5/NW7b9ofvG+1Zy430rGT50CMs+M6PGyqxd+Yxcg5998jjeOW0CI3ZPq3/E7kOYOW0CP7vwuJors3blINdg3MtGMGr4UDZv2cbwoUPYvGUbo4YPZdyoEXWXZm3KTeuabHh2M2ccNYn3HjmRefevpMM3vGwnKCLqrqF206dPj4ULF9ZdhvWRpEURMb3uOupUbNNa0omSlklaLslf6rWiFRlkSbsBVwMzgCnA6fnfE5kVqcggA0cCyyPisYh4EbgZmFlzTWYtU2qQ9wVWVbpX535mRSo1yOqi3x/d1ZM0W9JCSQs7OjoGqCyz1ig1yKuB/Svd+wFrqiNExLURMT0ipo8dO3ZAizPrb6UG+QHgIEmvljSM9GuRC2quyaxlivxASERskXQ+6TepdgOui4glNZdl1jJFBhkgIm4Fbq27DrOBUGrT2myX4iCbFcBBNiuAg2xWAAfZrAAOslkBHGSzAjjIZgVwkM0K4CCbFcBBNiuAg2xWAAfZrAAOslkBHGSzAjjITXRceVXdJZj1moPcxIarr667BLNec5DNCuAgV3RceRVLDzmUpYccCvCH525m22DnH3Gj6x9xW3rIoRz6yNKaKrId4R9x8xnZrAgOchNjPvKRuksw6zUHuYmxHz2/7hLMes3XyICkDuDxFkx6DLChBdNthXaptas6J0XELv27Pw5yC0la2C43Ydql1napc6C5aW1WAAfZrAAOcmtdW3cBO6Bdam2XOgeUr5HNCuAzslkBHOQWkXSipGWSlkuaU3c9zUhaIelXkhZLWtjzKwaOpOskrZf0UKXf3pJuk/RofnxFnTUOFg5yC0jaDbgamAFMAU6XNKXeqrp1XERMG4Rv63wDOLGh3xzg9og4CLg9d+/yHOTWOBJYHhGPRcSLwM3AzJprajsRcRfwVEPvmcD1+fn1wCkDWtQg5SC3xr7Aqkr36txvMArgR5IWSZpddzG9MD4ingDIj+NqrmdQGFp3AYVSF/0G69sDR0fEGknjgNskPZLPhNZGfEZujdXA/pXu/YA1NdXSrYhYkx/XA/NJlwWD2TpJ+wDkx/U11zMoOMit8QBwkKRXSxoGnAYsqLmml5A0UtKozufACcBD3b+qdguAWfn5LOA7NdYyaLhp3QIRsUXS+cAPgd2A6yJiSc1ldWU8MF8SpH1hXkT8oN6StpN0E3AsMEbSauDvgSuAWySdDawETq2vwsHDn+wyK4Cb1mYFcJDNCuAgmxXAQTYrgINsVgAH2awADrJZARzkQkl6g6QHJY3In+BaIunwuuuy1vAHQgom6TPACGAPYHVEXF5zSdYiDnLB8ue8HwBeAN4cEVtrLslaxE3rsu0N7AWMIp2ZrVA+IxdM0gLSfyd5NbBPRPgHrQrlbz8VStKZwJaImJf/h9jdkt4aEXfUXZv1P5+RzQrga2SzAjjIZgVwkM0K4CCbFcBBNiuAg2xWAAfZrAAOslkB/h/IuE4qizw3dAAAAABJRU5ErkJggg==\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAvEAAAF1CAYAAABh8bWyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeZhcZ3nn/e+t7la3pO6WbUleZVk2NjE2IQYLE15lsQ1J2MIyIcRmCUzCZTIZxhjDEAjJYCZAhkmGEPslMBrgzcJiglkCOMCw2MOYgEHyho1ZjBcsb1otdUvqbrX6fv84p6RSq6rX6q5q1fdzXbq66pznnLqrTkv69dP3OScyE0mSJEkLx6JmFyBJkiRpegzxkiRJ0gJjiJckSZIWGEO8JEmStMAY4iVJkqQFxhAvSZIkLTCGeEk6CkTE0oi4PiJ2RcQn5/i1LoyIu+byNeZaRNwUEa9pdh2SNFOGeEltKSIGq/6MRcS+quevaHZ9M/B7wHHAisy8dPzKiHhXROwf9763zeSFMvPGzDy3at+bI+LCqudnRsSCvwlJRLw2IjIi/t00t/tYRFw1R2VJEmCIl9SmMrO38gf4OfDbVcs+Pn58RHTOf5XTchrw48wcnWDMx6vfd2aunK/ipqOFPutXAzvKr5LUUgzxklRDOXP9qYj4ZEQMAK+MiGdGxHcj4vGIeCQiro6IrqptfjEivh4ROyLi0Yh4S7l8UUT8aUT8LCK2RcS1EXFsuW5pRHwiIraX+/1eRNQM1xFxbkT8n3LcDyLi+eXydwN/CryinGGfVuiMiM5yxvk/lDUORMQ7IuKs8v3uLj+HrnL8syPi/vLxJ4GTgS+Xr30l8K1yXWXG/+nl89dGxI8iYmdEfDkiTh33+n8cEfcAP6pR46KIuK78XB+PiBsj4klV6z9WHo8vl/V/JyJOr1r/nIj4cdlu9LdATPKZnAGsB14HPDciVlWte3ZE3B8Rb4mIrRHxcET8frnujyl+K/Kn5Xv/XLl8dUR8rhx/X0T8x6r9vav8fD9W1n5nRDytav3miLiyPOa7yrHd5boVEfGv5X53RsQXI+KUSQ+6pAXPEC9J9b0E+ASwHPgUMAq8AVhJEfCeQxHyiIjlwNeBLwInAU8Ebiz3cyXwfODXgNXAHuDqct2/B5aWy1cAfwwMjS8kIhYDXwKuB1YBbwQ+FRFnZubbgf/OoZn2f5jh+/0N4Lzyvb0d+DvgEopZ/qcCLxu/Qdm68zDw3PK131e+z+rfdnw/Il4K/GfgRWX9N1N8ttVeCDwd+MU69X0JOAs4EbgT+Kdx618O/DlFW9HPgb8AiIjjgeuAt1Icu83AMyb5LF4NfDczrwN+BoxvUVoNLKH4AeaPgA9GRH9m/h3F98p7yvf+kojoKGv/PnAKxef8nyPiWVX7e3H5fo4Bvsyh74+Kl5XbnQGcD7yqXL4I+F/AGorjtB/420nem6SjgCFekuq7KTO/mJljmbkvM7+fmTdn5mhm3gtsAH69HPtC4MHM/NvMHM7M3Zn5vXLd64A/zcyHMnMIuAp4WUQsoghdK4EzM/NAZm7MzMEatawHFgN/lZn7M/PrFGHvkmm8n5eXs9iVP18bt/69mTmQmXcAdwNfycz7M3Mn8FWKID9Tr6MItpWWn3cBF4ybNX5PZu7MzH3jNy6Pwd+X9VU+w/MjYlnVsOvKz28/8HGKH0gAXgDclpmfK9f9D2BrvUIjIihCcuWHjE9wZEvNEPCu8lh8ARim+MGtll8G+jPzPZk5kpn3AB/h8GP3fzLzq5l5gCLMnzduH+/PzEczczvFDwTnlZ/L1vJ97cvM3cB7OPQ9Keko1ip9h5LUih6sfhIRZ1MEwPMpZs87KWaUAU4F7qmznzXAFyNirGpZAscDf08xm/vPEdFPEeD+rEZv+8nAzzOz+oTRByhmdqfqE5n5mgnWP1b1eF+N58dM47XGOw34QNnKUjFGMaNdeZ0Hj9iqVM5m/yXwUoofeiqf5UqK32wAPFq1yV6gt3x8cvW+M3MsIjZPUOuvURzPfy6ffwL4rxHx5My8s1y2rQzctV5vvNOANRHxeNWyDg79pqZW7dU/nNRafxxA+UPM3wK/yaHj01enDklHEWfiJam+8VdY+Z8UbRxnZmY/8F841Fv9IPCEOvvZDPxGZh5T9aennFkdycyrMvNJwK9QtPDUujrOw8Cp5SxxxRrgoZm9tYYa/znVujLNg8AfjvsMlmTmzZNsV/H7wPOAiynam84sl0/Y2156hCKUFxsUvwFZPcH4V1P8/3hHRDwKfLus7fen8Fpw5Pt4EPjpuPfel5m/PcX9TeQtwOnABeX35MUN2KekBcAQL0lT1wfsAvaUJ1W+rmrdFyhmW18fEYsjoj8iLijXfQh4T0SsgaJHOyJeWD6+OCKeXAbL3RTtNdUzvBX/RtGT/6aI6IqIiylC7T/XGDvfHqPo1a7YAmR5cmjFh4C3V05GjYhjyj75qeqjaFnZTvFbkHdPY9svAedFxIuiuPLNGyn68o8QEUspZvv/kKJlpfLnjRQnN3dM4fXGfx7fAUYi4k0R0RMRHVGcBH3+NN5DPX0UM/M7I2IFxQ+WktqAIV6Spu5NFLO0AxSz8p+qrMjMXRQnHv4ORYj9CYd6k98HfAX4RhRXuvk3ihM4oWj1+CxFgL+L4uTYI27WlJnDwG9TnBi6jeLEx5dn5k+mUX/l6jXVf1ZMY/t63gO8s+yzvyIzByhaX24ul63LzE9TfA6fjojdwB3Ab03jNf4/it9GPEzxOf3bVDfMzMcorhjzVxQ/BKzhUBvUeP+O4vh+rPxNyaOZ+SjFyaNLKI7xZD4M/FJ5tZjrytao5wEXAPdTHL//CfRP9T1M4H0Uv5nYTvGZfLkB+5S0AMTh7ZWSJEmSWp0z8ZIkSdICY4iXJEmSFhhDvCRJkrTAGOIlSZKkBcYQXyov+XVrRHyp2bVIkiRJE/GOrYe8geI245Ne8mvlypW5du3aOS9IkiRJ7WvTpk3bMrPmfS0M8UBErAaeT3HzkCsnG7927Vo2btw453VJkiSpfUXEA/XW2U5TeD/FravHml2IdDTYtW8/L/vQd/jqXY82uxRJko5KbR/iI+IFwJbM3DTJuMsiYmNEbNy6des8VSctTN/52Xa+d/8O/tuXf9TsUiRJOiq1fYgH1gMvjIj7gWuBiyPiY+MHZeaGzFyXmetWrarZmiSpdOdDu4BiRl6SJDVe24f4zHxbZq7OzLXAJcA3M/OVTS5LWtAe2TUEwI49IwztP9DkaiRJOvq0fYiX1HhbBoYOPt6+Z6SJlUiSdHQyxFfJzBsz8wXNrkNa6LYODB98vH1weIKRkiRpJgzxkhpu2+AwZ5/YB8D2QWfiJUlqNEO8pIY6MJZs3zNyMMRvcyZekqSGM8RLaqjB4VEy4fSVvUBxcqskSWosQ7ykhhoYKi4reeLybjoWBQNDo02uSJKko48hXlJDDQ4Xob2vp4ve7k52D3mteEmSGs0QL6mhKjPvvd2d9C/pdCZekqQ5YIiX1FCDQ5WZ+E76ursOttdIkqTGMcRLaqiB4UMhvn9JJ7v3ORMvSVKjGeIlNVRl5r2vp4u+ni574iVJmgOGeEkNNVjVE9/XY0+8JElzwRAvqaEGhkZZFLB0cQf9zsRLkjQnDPGSGmpweJTe7k4igv6eTgaHRxkby2aXJUnSUcUQL6mhdg/tp6+nC4D+JV1kwuCILTWSJDWSIV5SQw0OjdLX0wkUffGVZZIkqXEM8ZIaamBo9GB4r8zIe3KrJEmNZYiX1FCDw4dm4itfveGTJEmNZYiX1FCDw6P0ljPwvZUQP+xMvCRJjTTrEB8RByLitqo/b53GtidHxHXl4/Mi4nlV6y6MiP9ntvVJml8DQ/sPttP0H5yJN8RLktRInQ3Yx77MPG8mG2bmw8BLy6fnAeuAfy2fXwgMAv821f1FRGdmmhakJhoYGj0Y3g/1xNtOI0lSI81ZO01E3B8R74mI70TExoh4WkR8NSJ+FhF/VI5ZGxF3RsRi4L8Cv1fO5v8J8EfAG8vnvxoRqyLiMxHx/fLP+nIfV0XEhoj438A/ztX7kTS5kdExhkfHDs7Ee3UaSZLmRiNm4pdExG1Vz/8yMz9VPn4wM58ZEX8D/D2wHugB7gI+VNkgM0ci4r8A6zLz9QARsQQYzMy/Lp9/AvibzLwpItYAXwWeVO7ifOBXMnNfA96PpBkaLHvfKye0Ll3cwaKwnUaSpEab63aaL5RffwD0ZuYAMBARQxFxzDRf59nAORFRed4fEX2V1zHAS81XaZupnNgaEfR2dx4M95IkqTEaEeInMlx+Hat6XHk+3ddeBDxzfFgvQ/2emRYoqXEqM+6VmfjicRe77YmXJKmhWukSkwNA3wTP/zfw+sqTiJjRybSS5s7BEN9dHeI7baeRJKnBGhHil4y7xOR/m+F+bqBol7ktIn4P+CLwksqJrcDlwLqIuCMifkhx4uusRcSpEXFDRNwdEXdFxBsasV+pHVXaZnp7Dg/xntgqSVJjzbqdJjM76ixfW/X47ylObB2/bhvw5HLZDuDp43bzlHHPf6/G61w1rYKPNAq8KTNvKXvsN0XE1zLzh7Pcr9R2DvbEdx/eTrNlYKhZJUmSdFRqpXaapsjMRzLzlvLxAHA3cEpzq5IWpkNXp+k6uKy323YaSZIare1DfLWIWAs8Fbi5xrrLyuvdb9y6det8lyYtCLVPbLWdRpKkRjPElyKiF/gMcEVm7h6/PjM3ZOa6zFy3atWq+S9QWgB2D+1ncccieroOddn19XQ5Ey9JUoPNaYiPiCsiYukE66+LiDPKx++OiAcjYnCC8WeXd4Adjog3Vy1fHBHfiogZ9fhHRBdFgP94Zn52JvuQVMzEV8/CQzETP3JgjOHRA02qSpKko89cz8RfAdQM8RFxLtCRmfeWi74IXDDJ/nZQXKXmr6sXZuYI8A1qnPg6mSguNP8R4O7MfN90t5d0SL0QX1knSZIaoyEhPiI+WPaL3xUR7yyXXQ6cDNwQETfU2OwVwL9UnmTmdzPzkYleJzO3ZOb3gVp3jvl8uc/pWg+8Cri46jKZz5vBfqS2Nzi0/7CTWsEQL0nSXGjUHVvfnpk7IqID+EZEPCUzr46IK4GLMnNbjW3WA59s0OsD3MmRl6icVGbeBEQD65Da1sDQ6GGXlwTo7S5CvSe3SpLUOI1qp3lZRNwC3AqcC5wzhW1OAhp2mZfMPACMlNd6l9QEE7fT1PoFmiRJmolZz8RHxOnAm4GnZ+bOiPh7oGcKm+6b4rjp6Aa8q4zUJAM12mkqM/O7nYmXJKlhGjET3w/sAXZFxAnAc6vWDQD1ZsbvBs6cbOcR8ZKI+MspjFsBbM1Mp/ukJhkYPnImvr8M9ZUbQUmSpNmbdYjPzNsp2mjuAj4KfLtq9Qbgy3VObL0euLDyJCL+e0RsBpZGxOaIuKpc9QRgdznmxHLMlcCfleP6y3EXAf862/cjaWbGxpLB4VH6baeRJGnONeTE1sx8TZ3l1wDX1NnsOoor17wjMw9k5luAt9QYdx7wxnJ/jwKr6+zv5cDbplO3pMbZMzJKJvSOC/GV557YKklS4zTtjq2ZuQ94B3DKJONemZkTngAbEYuBz2fmjxtYoqRpqFxCcnxPfFfHInq6FjFgO40kSQ3TqEtMzkhmfrVB+xkB/rER+5I0M4dC/JH/rPR2d9lOI0lSAzVtJl7S0WXn3hEAjlu6+Ih1/T2d3uxJkqQGMsRLaoide4oQf0yNEN9niJckqaEM8ZIaYkdlJn7ZkSG+t6fTdhpJkhrIEC+pIR7fW4T0Y5Z2HbGuv6eLXfsM8ZIkNYohXlJD7NgzwrLFHfR0dRyxbmVvN9vLdhtJkjR7hnhJDbFzz0jNfniAVX3dPL53P8OjB+a5KkmSjk6GeEkNsWPvSM1+eChm4gG2DzobL0lSIxjiJTXE9sERjq0T4lf1FSF+2+DwfJYkSdJRyxAvqSEe3T3ESf09NddVQvzWAUO8JEmNYIiXNGsjo2NsGxzmpGMM8ZIkzQdDvKRZe2z3EJlw0vLaIX5lb9Fms8UQL0lSQxjiJc3ao7uHADhp+ZKa67s7Ozi+r5uf79g7n2VJknTUMsRLmrUHy3B+8jG1QzzA2pXLuH/bnvkqSZKko5ohXtKs/eSxQbo6gtNWLK075vQVy7h/uzPxkiQ1giFe0qzds2WAM1b20tVR/5+UM1YtY9vgMDu9c6skSbNmiJc0K5nJbQ/u4pyT+ycc99Q1xwKw6YGd81GWJElHNUM8EBHPiYgfR8Q9EfHWZtcjLSQ/fmyAbYPDPPMJKyYc95TVy1ncsYh/+9n2eapMkqSjV9uH+IjoAD4APBc4B7g0Is5pblXSwvHJm39O56LgwieumnBcT1cHv/4Lq/jC7Q+xZ3h0nqqTJOno1NnsAlrABcA9mXkvQERcC7wI+GFTq5Ja0JaBIb50+yMMjR5gaP8YD+3cx2dv3cylF6zh+Dp3a632R79+Bi/90He4ZMN3edaTjqe3u5NFEfzSqcs5/7Tj5uEdSJJ0dDDEwynAg1XPNwPPGD8oIi4DLgNYs2bN/FQmtZhHdw3xX7906Ofbvp5OLnn6Gv78+VP75dX5px3H317yVP7H//4x7//6Tw8u/48XPcEQL0nSNBjiIWosyyMWZG4ANgCsW7fuiPVSO3jSSf3c+ue/wZLFHSzuWMSiRbX++kzshb90Mi/8pZPZf2CMvcMHSJLuzo45qFaSpKOXIb6YeT+16vlq4OEm1SK1tK6ORRy7bHHD9rV8adufliNJ0oz4Pyh8HzgrIk6PiMXAJcAXmlyTJEmSVFdk2hkSEc8D3g90AB/NzHdPMn4r8MB81DbOSmBbE15XE/O4tB6PSWvyuLQej0lr8ri0nmYdk9Mys+bl3wzxC0hEbMzMdc2uQ4fzuLQej0lr8ri0Ho9Ja/K4tJ5WPCa200iSJEkLjCFekiRJWmAM8QvLhmYXoJo8Lq3HY9KaPC6tx2PSmjwurafljok98ZIkSdIC40y8JEmStMAY4iVJkqQFxhAvSZIkLTCGeEmSJGmBMcRLkiRJC4whXpIkSVpgDPGSJEnSAmOIlyRJkhaYzmYX0CoiogPYCDyUmS+YaOzKlStz7dq181KXJEmS2tOmTZu2ZeaqWusM8Ye8Abgb6J9s4Nq1a9m4cePcVyRJkqS2FREP1FtnOw0QEauB5wMfbnYtkhamx3YP8dp/+D53Pbyr2aVIktqAIb7wfuAtwFi9ARFxWURsjIiNW7dunb/KJC0IX7nzUb5+9xY+8n/va3YpkqQ20PYhPiJeAGzJzE0TjcvMDZm5LjPXrVpVszVJUhv7wUPFDPzmnfuaXIkkqR20fYgH1gMvjIj7gWuBiyPiY80tSdJC88iuIrw/sGNPkyuRJLWDtg/xmfm2zFydmWuBS4BvZuYrm1yWpAVmy+5hALYNjjA2lk2uRpJ0tGv7EC9JjbBloAjxB8aSXfv2N7kaSdLRzhBfJTNvnOwa8ZI03vDoAXbt289Zx/cCsH3PcJMrkiQd7QzxkjRL2wZHAPiFE/sOey5J0lwxxEvSLO0u22fOWFXMxO/YY4iXJM0tQ7wkzdLg8CgAq49ZAhwK9ZIkzRVDvCTN0sBQEdpPLkP8wNBoM8uRJLUBQ7wkzVIltJ+4vJtFAbuHnImXJM0tQ7wkzVIlxPf3dNHb3elMvCRpzhniJWmWKj3xvT2d9PV0ORMvSZpzhnhJmqXBoVE6FgVLujroX9LF7n3OxEuS5pYhXpJmaWBoP73dnUQEfT2dB090lSRprhjiJWmWBoZH6e3uBIq++N32xEuS5pghXpJmaWBolL6eSoh3Jl6SNPcM8ZI0S4NVIb6vp9ObPUmS5pwhXpJmaXB4lL6eLgD6eroYHB4lM5tclSTpaGaIl6RZqpzYCsVM/FjC3pEDTa5KknQ0M8RL0iwNDo/SW7bTVL5Wrh0vSdJcMMRL0iztPqwnvmir8eRWSdJcmnKIj4iXRERGxNlTHH9FRCydeWmzExGnR8TNEfHTiPhURCxuVi2Sjl7DowcYGR2jr6qdBvAyk5KkOTWdmfhLgZuAS6Y4/gqgaSEeeC/wN5l5FrAT+MMm1iLpKLVnuOh9P3hiaxnmBw3xkqQ5NKUQHxG9wHqKIHxJ1fILI+LGiLguIn4UER+PwuXAycANEXFDOfY3I+I7EXFLRHy63CcRcX9EvLNc/oPKTH9EXBURHy33f2+5z8rrfj4iNkXEXRFxWY16A7gYuK5c9A/Ai2fw+UjShCptM4dObK200xjiJUlzZ6oz8S8GvpKZPwF2RMTTqtY9lWLW/RzgDGB9Zl4NPAxclJkXRcRK4M+AZ2fm04CNwJVV+9hWLv8g8Oaq5WcDvwVcALwjIrrK5X+QmecD64DLI2LFuHpXAI9nZuV/0c3AKVN8r5I0ZZWwXn2d+GK5PfGSpLkz1RB/KXBt+fja8nnF9zJzc2aOAbcBa2ts/8sUIf/bEXEb8GrgtKr1ny2/bhq3/fWZOZyZ24AtwAnl8ssj4nbgu8CpwFnjXi9q1OBFmyU1XCXEe3UaSdJ86pxsQDnLfTHw5IhIoAPIiHhLOWS4aviBOvsM4GuZeWmNddX7GL/9EfuOiAuBZwPPzMy9EXEj0DNuf9uAYyKis5yNX03xmwFJaqhKWO/rLn5R2LvYE1slSXNvKjPxLwX+MTNPy8y1mXkqcB/wK5NsNwD0lY+/C6yPiDMBImJpRDxxhjUvB3aWAf5siln+w2Rxq8QbytqhmPn/lxm+niTVNThctM1U2mgWLQp6uzs9sVWSNKemEuIvBT43btlngJdPst0G4MsRcUNmbgVeA3wyIu6gCPVTulRlDV+hmJG/A/iLcl+1/AlwZUTcQ9Ej/5FagyLi1Ii4ISLuLk+UfcMM65LUhsa300AR6O2JlyTNpUnbaTLzwhrLrq56emPV8tdXPb4GuKbq+TeBp9fY19qqxxuBC8vHV40b9+Sqp8+dQt33UpwQO5lR4E2ZeUtE9AGbIuJrmfnDKWwrqc0dDPHd40O8M/GSpLnT9ndszcxHMvOW8vEAcDdeyUbSFA0MjbK4YxE9XR0Hl/V2d3piqyRpTrV9iK8WEWspLpl5c411l0XExojYuHXr1vkuTVKLGhjaf7AfvqKvp8t2GknSnDLEl8qbT30GuCIzd49fn5kbMnNdZq5btWrV/BcoqSUNDI0eEeJ7baeRJM2xOQ/xEXFFRCydYP11EXFGecWa68s7v94VEf+tzvgXRcQdEXFbOTP+K+XyVRHxlRnW2EUR4D+emZ+dbLwkVRQz8V2HLevv6WTAdhpJ0hyaj5n4K4CaIT4izgU6ypNQAf46M8+maGlZHxG1TmD9BvBLmXke8AfAhwHKK+A8EhHrp1NcRATFlWvuzsz3TWdbSRoYGj3spFawnUaSNPcaFuIj4oPlzPhdEfHOctnlwMnADRFxQ43NXkF5/fbM3JuZN5SPR4BbKG7SdJjMHCyvAw+wjMPvxPr5cp/TsR54FXBxObt/W0Q8b5r7kNSmBodrtNN0dzK0f4z9B8aaVJUk6Wg36SUmp+HtmbkjIjqAb0TEUzLz6oi4ErgoM7fV2GY98MnxCyPiGOC3gb+t9UIR8RLgL4HjgedXrdoIvGs6RWfmTRR3lJWkaSt64g9vp6mE+sGhUY5dtrgZZUmSjnKNbKd5WUTcAtwKnAucM4VtTgIOu9RLRHRSBPurq9psDpOZnyvbbl5MccOnii0UM/+SNC9217g6TaW9xpNbJUlzpSEhPiJOB94MPCsznwJcD/RMYdN9NcZtAH6ame+fbOPM/BbwhIhYWS7qKfcpSXNubCwZHB6lv8YlJgEGhu2LlyTNjUbNxPcDe4BdEXECh99RdQDoq7Pd3cCZlScR8S5gOcXJsFQtf0lE/GX5+MzyZFQi4mnAYmB7OfSJwJ2zfjeSNAV7RkbJpObVacCZeEnS3GlIiM/M2ynaaO4CPgp8u2r1BuDLdU5svR64ECAiVgNvp2jDuaU8wfS15bgnAJVrt/8OcGdE3AZ8APi9qhNdLyr3KUlzrhLSe2tcJ756vSRJjdawE1sz8zV1ll8DXFNns+sorlzzjszcTP0TTM8D3lju773Ae+uMeyHwoqnWLEmzUQnpte7YCjBoO40kaY409Y6tmbkPeAdwyiTjXlleB76uiFgFvC8zdzawREmqqxLS612dxpl4SdJcaeQlJmckM7/aoP1spbhOvCTNi911ZuK9Oo0kaa41dSZekhayx/eOAHDs0sOvBd/T1cHijkWGeEnSnDHES9IM7dhTtNMct/TIGzr19nQyMGRPvCRpbhjiJWmGdu4ZYVEc2U4DxbLBYWfiJUlzwxAvSTO0Y+8Ixy5dzKJFR15Yq7+ni137nImXJM0NQ7wkzdDOPSMcu+zIVhqAFb2L2T44Ms8VSZLahSFekmZox56Rmv3wAKt6u9k6MDzPFUmS2oUhXpJm6PG9+zl2WVfNdSv7utk2OMzYWNZcL0nSbBjiJWmGtu8Z5rg67TSrersZHUv74iVJc8IQL0kzMDI6xrbBEU7o76m5flVfNwBbB22pkSQ1niFekmbgsd1DAJy0fJIQb1+8JGkOGOIlaQYe2VUJ8Utqrl/ZW4T4LQND81aTJKl9GOIlaQYe2bUPqD8Tf8oxRbh/cMe+eatJktQ+DPGSNAObd5Yh/pjaM/FLFndw0vIe7t+2Zz7LkiS1CUO8JM3APVsGOWl5D73dnXXHrF2xjPu3G+IlSY1niJekGfjJYwOcdULfhGNOX7WMn23d47XiJUkNZ4iXpGnaN3KAnz42yJNOnDjEn3fqMezat597tw3OU2WSpHZhiAci4jkR8eOIuCci3trseiS1to0P7GDkwBi//IQVE457+trjAPjOvTvmoyxJUhtp+xAfER3AB4DnAucAl0bEOc2tSlIr+/TGzfR2d/KM04+bcNzaFUt5wqplXLfxQVtqJEkNVf+MrPZxAXBPZt4LEBHXAi8CftjUqsa586Fd/Gzr7H4ln5NkiGTykDHpPiZ9jcnlJDuZUhSa5Xud7H1MpY6p7WN2dUzts5j95zn5cZ+H15jSPmYXlKey+f3b9/CF2x/mP1z4BJYunvif0Ijgsl87gz/5zJmtJGgAACAASURBVA/4T5+8lfVnrmTp4g4i6o+XJLWek5b3HPztaqswxMMpwINVzzcDzxg/KCIuAy4DWLNmzfxUVuXztz7Eh2+6b95fV9LhujqCl61bzRXPPmtK41+27lQ279zHR266j+t/8MgcVydJmgvPOffElgvxMduZq4UuIn4X+K3MfG35/FXABZn5n+pts27duty4ceN8lQjAtsFhdu3bP+GYqczhTTbTN7V9TLJ+kr00YrJxKvuY7Xud0ms04L1OOmSWn/dU6piP751GfJ5TKXS273Wy99nVEXR3dkxeyDgjo2Ps2DPC3pHRg8uq//Vt83+KJamlLevuqHuH7rkUEZsyc12tdc7EFzPvp1Y9Xw083KRa6lrZ233wNu6SFp7FnYs4sc7dXSVJmq62P7EV+D5wVkScHhGLgUuALzS5JkmSJKmutm+nAYiI5wHvBzqAj2bmuycZvxV4YD5qG2clsK0Jr6uJeVxaj8ekNXlcWo/HpDV5XFpPs47JaZm5qtYKQ/wCEhEb6/VFqXk8Lq3HY9KaPC6tx2PSmjwuracVj4ntNJIkSdICY4iXJEmSFhhD/MKyodkFqCaPS+vxmLQmj0vr8Zi0Jo9L62m5Y2JPvCRJkrTAOBMvSZIkLTCGeEmSJGmBMcRLkiRJC4whXpIkSVpgDPGSJEnSAmOIlyRJkhYYQ7wkSZK0wBjiJUmSpAWms9kFtIqI6AA2Ag9l5gsmGrty5cpcu3btvNQlSZKk9rRp06Ztmbmq1jpD/CFvAO4G+icbuHbtWjZu3Dj3FUmSJKltRcQD9dbZTgNExGrg+cCHm12LJEnSfPvITffxye/9vNllaBqciS+8H3gL0FdvQERcBlwGsGbNmnkqS5IkaW4N7T/AX3zphwC89PzVdHU4x7sQtP1RiogXAFsyc9NE4zJzQ2auy8x1q1bVbE2SJElacH706MDBxw/t3NfESjQdbR/igfXACyPifuBa4OKI+FhzS5IkSZof1cH9gR17m1iJpqPtQ3xmvi0zV2fmWuAS4JuZ+comlyVJkjQvtg4MHXz82O6hCUaqlbR9iJckSWpnWwaGDz7ePjjSxEo0HZ7YWiUzbwRubHIZkiRJ82brwDAn9vewa99+tg8OT76BWoIhXpIkqY1tHRxmVV83nR3Bjj3OxC8UhnhJkqQ2tmvffo5Z2sWigG2G+AXDEC9JktTGBodGObG/BygCvRYGT2yVJElqYwNDo/T1dNK/pIuBIUP8QuFMvCRJUhsbHB6lt7uLjkXBwNBos8vRFBniJUmS2tSBsSxCfE8nXR3BbttpFgzbaSRJktrUnpFi5r2/p5O+nk6GR8cYHj3Q5Ko0FYZ4SZKkNjVYts/0dnfS19MFYEvNAmGIlyRJalOVwN7X00X/ks7Dlqm12RMvSZLUpgaHix743p5O9o8Wc7teoWZhMMRLkiS1qd1V7TQHurJYts+Z+IXAEC9JktSmKj3x/T2djI4VId6Z+IXBEC9JktSmKv3vvT2dHDgY4p2JXwgM8ZIkSW2q0hPf19N1MMTvdiZ+QTDES5IktamBoVEiYGlXB1kuGxx2Jn4hMMRLkiS1qYGhUXq7O1m0KABYtrjDdpoFYsrXiY+Il0RERsTZUxx/RUQsnXlpsxMRr4+Ie8qaVzarDkmSpFY1ODxKX/ehOd3ens6DJ7uqtU3nZk+XAjcBl0xx/BVA00I88G3g2cADTaxBkiSpZQ0M7ae351CI7+vpYmDYnviFYEohPiJ6gfXAH1IV4iPiwoi4MSKui4gfRcTHo3A5cDJwQ0TcUI79zYj4TkTcEhGfLvdJRNwfEe8sl/+gMtMfEVdFxEfL/d9b7rPyup+PiE0RcVdEXFar5sy8NTPvn9nHIkmSdPQbHB6lr6fr4PO+nk7baRaIqc7Evxj4Smb+BNgREU+rWvdUiln3c4AzgPWZeTXwMHBRZl5UtrP8GfDszHwasBG4smof28rlHwTeXLX8bOC3gAuAd0RE5bvsDzLzfGAdcHlErJj6W5YkSRIc6omv6O02xC8UUw3xlwLXlo+vLZ9XfC8zN2fmGHAbsLbG9r9MEfK/HRG3Aa8GTqta/9ny66Zx21+fmcOZuQ3YApxQLr88Im4HvgucCpw1xfchSZKk0uDQKH1V7TT9PV3e7GmBmPTqNOUs98XAkyMigQ4gI+It5ZDhquEH6uwzgK9l5qU11lXvY/z2R+w7Ii6k6HV/ZmbujYgbgZ7J3ockSZIOt3tciLedZuGYykz8S4F/zMzTMnNtZp4K3Af8yiTbDQB95ePvAusj4kyAiFgaEU+cYc3LgZ1lgD+bYpZfkiRJ0zQ4vP+wnvje7k6vE79ATCXEXwp8btyyzwAvn2S7DcCXI+KGzNwKvAb4ZETcQRHqp3Spyhq+QjEjfwfwF+W+jhARl0fEZmA1cEdEfLjOuFMj4oaIuLs8UfYNM6xLkiRpwdh/YIyh/WOH9cT39XSxd+QAowfGmliZpmLSdprMvLDGsqurnt5Ytfz1VY+vAa6pev5N4Ok19rW26vFG4MLy8VXjxj256ulzp1D31cDVk40DRoE3ZeYtEdEHbIqIr2XmD6ewrSRJ0oJUuR5877jrxAPsGT7A8qXTuRK55lvbH53MfCQzbykfDwB3A6c0typJkqS5Vel9H98TD7Dbk1tbXtuH+GoRsZbikpk311h3WURsjIiNW7dune/SJEmSGqoS1Kt74vvLEO/Jra3PEF8qbz71GeCKzNw9fn1mbsjMdZm5btWqVfNfoCRJUgNVTmDt76m+TnzXYevUuhoS4iPiiohYOsH66yLijPLxuyPiwYgYHDemOyI+FRH3RMTN5ax4rX19JSIej4gvjVt+erndT8v9LC6Xvz4i/v0k9XdRBPiPZ+ZnJxorSZJ0NKjMtvfWaKfxWvGtr1Ez8VcANUN8RJwLdGTmveWiL1LcgXW8P6S4dOSZwN8A763zWn8FvKrG8vcCf5OZZwE7y/0BfBS4vF7hERHAR4C7M/N99cZJkiQdTQZqtNP02k6zYEwrxEfEB8u+8Lsi4p3lssuBk4EbIuKGGpu9AviXypPM/G5mPlJj3IuAfygfXwc8qwzYh8nMb1Bcg766rqC4IdV15aJ/AF5cjt8L3B8RtX5wAFhP8UPBxRFxW/nneXXGSpIkHRUmOrF1wHaaljfpJSbHeXtm7oiIDuAbEfGUzLw6Iq4ELsrMbTW2WQ98cgr7PgV4ECAzRyNiF7ACqLXP8VYAj2dm5TtuM4dfYWYj8KvA98ZvmJk3UdxRVpIkqW0cmok/FAf7y1l522la33TbaV4WEbcAtwLnAudMYZuTgKlczqVWkM4p1jXZtlsoflsgSZIkitn2xZ2L6O7sOLisu3MRnYvi4DXk1bqmHOIj4nTgzcCzMvMpwPVAzxQ23TfFcZuBU8vX6gSWAzumWN424JhyOyju0vpw1fqesg5JkiRRtNP0dR/elBER9PV02hO/AExnJr4f2APsiogTOPyuqQNAX53t7gbOnML+vwC8unz8UuCbmZkRcUpEfGOiDTMzgRvK7Sj38y9VQ54I3DmFGiRJktrCwNDoYa00Fb09nbbTLABTDvGZeTtFG81dFFd8+XbV6g3Al+uc2Ho9cGHlSUT894jYDCyNiM0RcVW56iPAioi4B7gSeGu5/CRgtGr7/wt8muLE180R8Vvlqj8Briy3X1Hur2I98PWpvldJkqSj3cDQ/sOuTFPR193ldeIXgGmd2JqZr6mz/BrgmjqbXUdx5Zp3ZOaBzHwL8JYa+xgCfrfG9r8MfKBq3K/WqeFealy6MiKeCtxV56RbSZKktjRYZya+r6eT3bbTtLw5v2NrZu4D3sHhV4uZzvb/b2Z+YRYlrAT+fBbbS5IkHXXqtdPYE78wTPcSkzOSmV+dj9ep89pfa9ZrS5IktaqBof30dtdop+npYnB4oMYWaiVzPhMvSZKk1vP4vv0cu7RWiHcmfiEwxEuSJLWZof0H2DtygGOXLT5iXW93J4NDoxQX/1OrMsRLkiS1mZ17RwA4rkaI7+vpYnQsGdo/Nt9laRoM8ZIkSW1mx54ixB+7tFaIL06Z3O214luaIV6SJKnN7NxTBPRaM/Ere4tl2waH57UmTY8hXpIkqc3sONhOc+SJrav6ugHYOmCIb2WGeEmSpDazc4J2mpW9RYjfNjgyrzVpegzxkiRJbWb7nhEiYPmSI2fiKyHemfjWZoiXJElqM4/tGmJlbzedHUdGwWXdnSxb3GFPfIszxEuSJLWZR3YPcfLynrrrV/Z1OxPf4gzxkiRJbebRXfs4cYIQv6q3my0DQ/NYkabLEC9JktRmHnl8iJOWL6m7/pRjl7B55755rEjTZYiXJElqI7v27WdgeJSTJpiJP23FMh5+fB/DowfmsTJNhyFekiSpjdyzZRCAM4/vrTvm9JVLGUt4cMfe+SpL02SIlyRJaiM/fWwAgLOO76s75vSVRcC/Z8ueealJ02eIlyRJaiM/eGgXvd2drD62fk/82Sf2sbhjEbf+fOc8VqbpMMQDEfGciPhxRNwTEW9tdj2SJElz5d9+tp0LTj+ORYui7pierg5+cfVybr5vxzxWpulo+xAfER3AB4DnAucAl0bEOc2tSpIkqfFue/Bx7tu2hwt/YdWkYy8++3hue/Dxgz30ai2dzS6gBVwA3JOZ9wJExLXAi4AfNrWqce7btoeHH2/cpZ4yG7YrkgbujMbWBjS0umxwcQ1+qw3dYcsf14a+18Zq5e+TRh+HRn96rX1cG7y/BlbYyv9uQuP/TjRSKx9XaGx9w6NjbPjWvaxYtpgXP/WUScf/7rrVfPDGn/H6T9zC6379DI5b1k1HHJq9j3ET+fXn9Re+43oXc/aJ/c0u4zCGeDgFeLDq+WbgGeMHRcRlwGUAa9asmZ/Kqnz8uw/w4Zvum/fXlSRJR48T+rv5u1c8jf6erknHHt/XwzUvfypv/ufbeeOnbp+H6lrXc849kQ+96vxml3EYQ3ztHxyP+Lk3MzcAGwDWrVs371MKv//MtfzmuSc2dJ/jf4Ke1b4at6tifw3/cb5xO2x0bY3/7Br4Xhu2p3J/Df/sWve4Nlpj/7429s02/Li203tt5L5a+N9NaPXPrsHvtaF7a9xnFwSnHLuEjgl64ce76BeO5+Y/fRb3bdvD7qH9jJUJaPxvCFr5ty2NcNyyxc0u4QiG+GLm/dSq56uBh5tUS11rVixlzYqlzS5DkiS1mc6ORZx1Qv3LUao52v7EVuD7wFkRcXpELAYuAb7Q5JokSZKkuuJo//XHVETE84D3Ax3ARzPz3ZOM3wo8MB+1jbMS2NaE19XEPC6tx2PSmjwurcdj0po8Lq2nWcfktMyseSkhQ/wCEhEbM3Nds+vQ4Twurcdj0po8Lq3HY9KaPC6tpxWPie00kiRJ0gJjiJckSZIWGEP8wrKh2QWoJo9L6/GYtCaPS+vxmLQmj0vrabljYk+8JEmStMA4Ey9JkiQtMIZ4SZIkaYExxEuSJEkLjCFekiRJWmAM8ZIkSdICY4iXJEmSFhhDvCRJkrTAGOIlSZKkBaaz2QW0iojoADYCD2XmCyYau3Llyly7du281CVJkqT2tGnTpm2ZuarWOkP8IW8A7gb6Jxu4du1aNm7cOPcVSZIkqW1FxAP11tlOA0TEauD5wIebXYsktauxseTqb/yUOx/a1exSJKnlGeIL7wfeAozVGxARl0XExojYuHXr1vmrTJLaxB0P7eJ9X/sJb/rn25tdiiS1vLYP8RHxAmBLZm6aaFxmbsjMdZm5btWqmq1JkqRZuGPz4wD8+LGBJlciSa2v7UM8sB54YUTcD1wLXBwRH2tuSZLUfh7csffg4z3Do02sRJJaX9uH+Mx8W2auzsy1wCXANzPzlU0uS5LazpaB4YOPH9s91MRKJKn1tX2IlyS1hq1VIX77npEmViJJrc8QXyUzb5zsGvGSpLmxZWCYNcctBWD7oCFekiZiiJcktYStA8OcfWIfANv3DE8yWpLamyFektR0Y2PJ7qH9nLGqF3AmXpImY4iXJDXdnpFRMmHFssUsXdzB7n37m12SJLU0Q7wkqekGy0tK9vZ00t/Txe4hQ7wkTcQQL0lquoGhIsT39XTS19N58LkkqTZDvCSp6SqhvbfbEC9JU2GIlyQ13UDZPtPX00X/EttpJGkyhnhJUtNVeuKLdpouZ+IlaRKGeElS0w1WtdP093R6dRpJmoQhXpLUdIef2OpMvCRNxhAvSWq6geFRImDZ4uLE1pEDYwztP9DssiSpZRniJUlNNzC0n97FnSxaFPQv6QLw5FZJmoAhXpLUdINDo/T2dALQX361pUaS6jPES5KabnB4lL4yvPd2F18HDfGSVJchXpLUdANDowfDe19P18FlkqTaDPGSpKYbGB6ltwzvfQfbaeyJl6R6phziI+IlEZERcfYUx18REUtnXtrsRMTHI+LHEXFnRHw0IrqaVYskaWIDQ/uPaKcZGHYmXpLqmc5M/KXATcAlUxx/BdC0EA98HDgb+EVgCfDaJtYiSZrA4NAofd2VE1ttp5GkyUwpxEdEL7Ae+EOqQnxEXBgRN0bEdRHxo3L2OyLicuBk4IaIuKEc+5sR8Z2IuCUiPl3uk4i4PyLeWS7/QWWmPyKuKmfQb4yIe8t9Vl738xGxKSLuiojLatWcmf+aJeB7wOoZfUKSpDl32ImtttNI0qSmOhP/YuArmfkTYEdEPK1q3VMpZt3PAc4A1mfm1cDDwEWZeVFErAT+DHh2Zj4N2AhcWbWPbeXyDwJvrlp+NvBbwAXAO6paYv4gM88H1gGXR8SKeoWX27wK+MoU36skaR6NHhhj78gBeruLf+I7FgVLF3d4dRpJmsBUQ/ylwLXl42vL5xXfy8zNmTkG3AasrbH9L1OE/G9HxG3Aq4HTqtZ/tvy6adz212fmcGZuA7YAJ5TLL4+I24HvAqcCZ01Q+98B38rM/zvhO5QkNcWe4eLOrJWZ+Mpj22kkqb7OyQaUs9wXA0+OiAQ6gIyIt5RDhquGH6izzwC+lpmX1lhXvY/x2x+x74i4EHg28MzM3BsRNwI9dWp/B7AKeF2d15UkNVnlzqy9h4X4LgaGbaeRpHqmMhP/UuAfM/O0zFybmacC9wG/Msl2A0Bf+fi7wPqIOBMgIpZGxBNnWPNyYGcZ4M+mmOU/QkS8lqIV59LytwSSpBY0WF6FpnJiKxRXqHEmXpLqm0qIvxT43LhlnwFePsl2G4AvR8QNmbkVeA3wyYi4gyLUT+lSlTV8hWJG/g7gL8p91fIhivab70TEbRHxX2oNiohTI+KGiLi7PFH2DTOsS5I0AwdDfM+hKwHbTiNJE5u0nSYzL6yx7OqqpzdWLX991eNrgGuqnn8TeHqNfa2terwRuLB8fNW4cU+uevrcKdQ96XsrjQJvysxbIqIP2BQRX8vMH05xe0nSLAzUbKfp5OHH9zWrJElqeW1/x9bMfCQzbykfDwB3A6c0typJah+VGffDTmzt7jo4Qy9JOlLbh/hqEbGW4pKZN9dYd1lEbIyIjVu3bp3v0iTpqLW7Voi3nUaSJmSIL5U3n/oMcEVm7h6/PjM3ZOa6zFy3atWq+S9Qko5SlevB93Uf6onv7elk78gBRg94XQJJqqXhIT4iroiIpROsvy4izqix/DfKu7D+oPx68RRe67Ry7G3lSal/VLXu6xFx7BRr7qII8B/PzM9ONl6S1DgDQ/vpXBT0dB36L6lykmvlGvKSpMPNxUz8FUDNEB8R5wIdmXlvjdXbgN/OzF+kuBnUP03htR4B/p/MPA94BvDWiDi5XPdPwB9PtoOICOAjwN2Z+b4pvKYkqYEGhkbp6+mk+Oe4UGmtqVxDXpJ0uBmH+Ij4YNkjfldEvLNcdjlwMnBDRNxQY7NXAP9Sa3+ZeWtmPlw+vQvoiYjuiWrIzJHMrNwQqpvD388XOPzOsvWsB14FXFzO6N8WEc+bwnaSpAYYGNp/2OUl4dA14+2Ll6TapnoZxlrenpk7IqID+EZEPCUzr46IK4GLMnNbjW3WA5+cwr5/B7i1KqDXFRGnAtcDZwL/ufKDQGbujIjuiFiRmdvrbZ+ZN1HcUVaS1ASVmfhqlVDvFWokqbbZtNO8LCJuAW4FzgXOmcI2JwETXtqlbLl5L/C6qRSRmQ9m5lMoQvyrI+KEqtVbKH4zIElqUQNDo/R2jw/xlZl422kkqZYZhfiIOB14M/CsMkBfD/RMYdN9lXER8ZKq9pV15bLVFHeH/f3M/FmN131G1TYvrF5XzsDfBfxq1eKe8jUlSS1qYHj0iHaayo2fnImXpNpmOhPfD+wBdpUz39V3UB0A+upsdzfFjDmZ+bnMPK/8szEijqH4YeBtmfnt6o0i4h8j4oLMvLlqmy9ExOqIWFKOOZaiXefH5fMATgTun+F7lCTNg4Gh/fQf0U5TObHVEC9JtcwoxGfm7RRtNHcBHwWqQ/cG4Mt1Tmy9Hriwzm5fTxHw/7xqtv34ct1TKK5EM96TgJsj4nbg/wB/nZk/KNedD3w3M/0fQJJaWK2e+P5yZt52GkmqbcYntmbma+osvwa4ps5m11FcueYdmXnYxX8z813Au8ZvEBH9wE8z88Ear/U1ioBfy6uAv6v7BiRJTZeZDNZop+nuXETnojh4IyhJ0uHm9Y6tmbkPeAdwyjS22Z2ZvzuDl7szM78xg+0kSfNk78gBDozlwR74ioigr6fTS0xKUh2zucTkjGTmV+fpdf7XfLyOJGnmKiF9fDtNsazLdhpJqmNeZ+IlSar2+L4RAI5ZsviIdb3dnV6dRpLqMMRLkppm555ipv3YZV1HrOvr6fTqNJJUhyFektQ0O/cWM/HHLTtyJt6eeEmqzxAvSWqaHXvKEL+0VojvYvc+e+IlqRZDvCSpaXaWIf6YGiF+xbLF7NgzQmbOd1mS1PIM8ZKkptm5dz+93Z0s7jzyv6OVfd3s23+APSMHamwpSe3NEC9Japqde0dqntQKsKq3G4BtA8PzWZIkLQiGeElS02zfM8KxNVppoJiJB9g6aIiXpPEM8ZKkptmye4gT+ntqrnMmXpLqM8RLkprm4cf3cdLyOiHemXhJqssQL0lqij3Do+weGuWk5Utqrj9u2WIWBWx1Jl6SjmCIlyQ1xSO7hgA4+ZjaM/Edi4IT+nvYvHPffJYlSQuCIV6S1BQPP16E8xPr9MQDrF2xjPu27ZmvkiRpwTDES5Ka4qdbBgF4wvG9dcecvmoZ9283xEvSeIZ4SVJT3LNlgOOWLWZleRWaWk5fsYzH9+5nR3lnV0lSwRAvSWqKOx/azRNPqD8LD/CLq5cDcOvPd85HSZK0YBjigYh4TkT8OCLuiYi3NrseSTra7dq7n7se3sUzTl8x4bhfWn0MXR3B9+7fMU+VSdLC0PYhPiI6gA8AzwXOAS6NiHOaW5UkHd3+9c5HGEv4tSeumnDcksUdXHD6cXzp9kc4MJbzVJ0ktb7OZhfQAi4A7snMewEi4lrgRcAPm1rVOFt2D7Fj79R6QnMa/89NayxTHzyd/U5HK9Q7nbeW09jx9PY7jcEt8TlMZ2zzP7PpfO/M0dCj+nt9594R/uqrP+Ypq5fztDXHTDr+1c9cy2X/tIk/+cwdvPT81fT1dBIEAFF8OfgVOLhOkhqlt6eTU46pfU+LZjHEwynAg1XPNwPPaFItdW341r18+Kb7ml2GJDXEqcct4X0vO4+IyQP3b5xzAq/7tTP4n9+6l+s2bZ6H6iTpcM8590Q+9Krzm13GYQzx1JyyOWJKKSIuAy4DWLNmzVzXdITfOX8155927JTHT+H/xerRc7Lf6ZQwlf/IZ7bfORo7nSqO5s+sBb53pvf5Nr/eaR23OathGoPn4DNb3LGIs0/so7Njah2dEcHbnvck/v3607lnyyCDw6PlmuKf6upfAthwI2kunDDB/SyaxRBfzLyfWvV8NfDw+EGZuQHYALBu3bp5/3/iSSf186ST+v//9u4t1NI5jOP499dmoiGHHHMccqEkJCmShHCDQpTiiguK3JAbh1IScjciUxSGHOeSC8INMxjHCUPjOM2QxNwQHhfrPxpj741k/d/X+n5qt9f6r3e3nvr19D579X/fNe23laTB2G+3ndhvt+GdSCWph5m/sBVYDRyRZFmSJcDFwKrONUmSJEkLyj+5GOn/Ksk5wD3AHLCiqm77i+O/Bj6dRm3b2Qv4psP7anHmMjxmMkzmMjxmMkzmMjy9Mjmkqua9jZdD/IgkWVNVx/euQ39kLsNjJsNkLsNjJsNkLsMzxEzcTiNJkiSNjEO8JEmSNDIO8eNyX+8CNC9zGR4zGSZzGR4zGSZzGZ7BZeKeeEmSJGlk/CRekiRJGhmH+JFIclaSD5KsT3JD73oESTYkeSfJ2iRretczq5KsSLI5ybvbrO2Z5PkkH7Xff//rjvWvLZDJzUm+bP2ytt3aV1OS5KAkLyRZl+S9JNe0dXulo0VysV86SbJTkteSvNUyuaWtL0vyauuVx9p3C/Wt1e00w5dkDvgQOIPJN8yuBi6pqve7FjbjkmwAjq8q7+XbUZJTgC3AQ1V1VFu7A/i2qm5v//TuUVXX96xzliyQyc3Alqq6s2dtsyrJ/sD+VfVGkl2B14HzgMuxV7pZJJeLsF+6SBJgaVVtSbIj8ApwDXAd8FRVrUxyL/BWVS3vWaufxI/DCcD6qvqkqn4CVgLndq5JGoSqegn4drvlc4EH2+MHmZwUNSULZKKOqmpjVb3RHv8ArAMOwF7papFc1ElNbGlPd2w/BZwGPNHWB9ErDvHjcADw+TbPv8AmH4ICnkvyepIrehejP9i3qjbC5CQJ7NO5Hk1cneTttt3GbRudJDkUOBZ4FXtlMLbLBeyXbpLMJVkLbAaeBz4Gvquqn9shg5jDHOLHIfOsuQ+qv5Oq6jjgbOCqtoVA0vyWA4cDxwAbgbv6ljObkuwCPAlcW1Xf965HE/PkYr90VFW/b5yuwAAAAY5JREFUVNUxwIFMdkMcOd9h063qzxzix+EL4KBtnh8IfNWpFjVV9VX7vRl4mkmjaxg2tb2mW/ecbu5cz8yrqk3txPgrcD/2y9S1/b1PAg9X1VNt2V7pbL5c7JdhqKrvgBeBE4Hdk+zQXhrEHOYQPw6rgSPaldFLgIuBVZ1rmmlJlraLkEiyFDgTeHfxv9IUrQIua48vA57tWIv4fUDc6nzsl6lqF+s9AKyrqru3ecle6WihXOyXfpLsnWT39nhn4HQm1yq8AFzQDhtEr3h3mpFot5e6B5gDVlTVbZ1LmmlJDmPy6TvADsAjZtJHkkeBU4G9gE3ATcAzwOPAwcBnwIVV5YWWU7JAJqcy2RpQwAbgyq17sfXfS3Iy8DLwDvBrW76Ryf5re6WTRXK5BPuliyRHM7lwdY7Jh92PV9Wt7by/EtgTeBO4tKp+7FepQ7wkSZI0Om6nkSRJkkbGIV6SJEkaGYd4SZIkaWQc4iVJkqSRcYiXJEmSRsYhXpIkSRoZh3hJkiRpZBziJUmSpJH5DTwQYA+zMeLoAAAAAElFTkSuQmCC\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "if True:\n", + " sample_rate = 3e2 # Hz\n", + " periodic = False\n", + " \n", + " t_offset = 8\n", + " t_start = 0\n", + " t_end = 1\n", + " time = t_offset + np.arange(t_start, t_end, 1/sample_rate) #s\n", + " \n", + " t_longstart = 0\n", + " t_longend = 30*t_end\n", + " longtime = np.arange(t_longstart, t_longend, 1/sample_rate) #s\n", + "\n", + "if False:\n", + " if True:\n", + " freq = sample_rate/8\n", + " signal = np.cos(2*np.pi*freq*time)\n", + " else: \n", + " from scipy.stats import norm\n", + " signal = norm.pdf(time, time[len(time)//2], (time[-1] - time[0])/10)\n", + "\n", + "\n", + "#####\n", + "# Setup Signal, Emitter and Antennae\n", + "\n", + "mysignal = TravelSignal(signal, sample_rate, t_0 = t_offset, periodic=periodic)\n", + "\n", + "source = Emitter([1,1])\n", + "emitted = source.emit(mysignal)\n", + "\n", + "antennae = [\n", + " Receiver([2,3]),\n", + " Receiver([10,10]),\n", + " Receiver([-2,-3]),\n", + "]\n", + " \n", + "#####\n", + "# Follow traces, and show geometry\n", + "ylabel_kw = {\"rotation\": \"horizontal\", \"va\":\"center\", \"ha\":\"center\", \"labelpad\": 30}\n", + "\n", + "fig, axs = plt.subplots(1,1, figsize=(2,2))\n", + "axs = [ axs ]\n", + "\n", + "### Geometry Plot\n", + "i = 0\n", + "axs[i].set_title(\"Geometry of Emitter(s) and Antennae\")\n", + "axs[i].set_ylabel(\"y\", **ylabel_kw)\n", + "axs[i].set_xlabel(\"x\")\n", + "axs[i].plot(*source.x, '*', label=\"Emitter\")\n", + "\n", + "for j, ant in enumerate(antennae):\n", + " axs[i].plot(*ant.x, '+', label=\"Antenna {}\".format(j))\n", + "\n", + "### Plot Traces\n", + "fig, axs = plt.subplots(1+len(antennae),1, sharex=True, figsize=(12,6))\n", + "axs[0].set_title(\"Traces of Emitter and Antenna\")\n", + "\n", + "# Emitter\n", + "i = 0\n", + "axs[i].set_ylabel(\"Emitter\\n at ({},{})\".format(*source.x), **ylabel_kw)\n", + "axs[i].plot(time, emitted(time))\n", + "\n", + "# Antenna\n", + "for j, ant in enumerate(antennae):\n", + " i +=1\n", + " axs[i].set_ylabel(\"Antenna {}\\n at ({},{})\".format(j, *ant.x), **ylabel_kw)\n", + " axs[i].plot(longtime, ant.recv(emitted)(longtime), label=\"Antenna {}\".format(j))\n", + " " ] } ],