random_signals/distributions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "nbsphinx": "hidden"
   },
   "source": [
    "# Random Signals\n",
    "\n",
    "*This jupyter notebook is part of a [collection of notebooks](../index.ipynb) on various topics of Digital Signal Processing. Please direct questions and suggestions to [Sascha.Spors@uni-rostock.de](mailto:Sascha.Spors@uni-rostock.de).*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Cumulative Distribution Functions\n",
    "\n",
    "A random process can be characterized by the statistical properties of its amplitude values. [Cumulative distribution functions](https://en.wikipedia.org/wiki/Cumulative_distribution_function) (CDFs) are one possibility to do so."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Univariate Cumulative Distribution Function\n",
    "\n",
    "The univariate CDF $P_x(\\theta, k)$ of a continuous-amplitude real-valued random signal $x[k]$ is defined as\n",
    "\n",
    "\\begin{equation}\n",
    "P_x(\\theta, k) := \\Pr \\{ x[k] \\leq \\theta\\}\n",
    "\\end{equation}\n",
    "\n",
    "where $\\Pr \\{ \\cdot \\}$ denotes the probability that the given condition holds. The univariate CDF quantifies the probability that for the entire ensemble and for a fixed time index $k$ the amplitude $x[k]$ is smaller or equal to $\\theta$. The term '*univariate*' reflects the fact that only one random process is considered.\n",
    "\n",
    "The CDF shows the following properties which can be concluded directly from its definition\n",
    "\n",
    "\\begin{equation}\n",
    "\\lim_{\\theta \\to -\\infty} P_x(\\theta, k) = 0\n",
    "\\end{equation}\n",
    "\n",
    "and\n",
    "\n",
    "\\begin{equation}\n",
    "\\lim_{\\theta \\to \\infty} P_x(\\theta, k) = 1\n",
    "\\end{equation}\n",
    "\n",
    "The former property results from the fact that all amplitude values $x[k]$ are larger than $- \\infty$, the latter from the fact that all amplitude values lie within $- \\infty$ and $\\infty$. The univariate CDF $P_x(\\theta, k)$ is furthermore a non-decreasing function\n",
    "\n",
    "\\begin{equation}\n",
    "P_x(\\theta_1, k) \\leq  P_x(\\theta_2, k) \\quad \\text{for } \\theta_1 \\leq \\theta_2\n",
    "\\end{equation}\n",
    "\n",
    "\n",
    "The probability that $\\theta_1 < x[k] \\leq \\theta_2$ is given as\n",
    "\n",
    "\\begin{equation}\n",
    "\\Pr \\{\\theta_1 < x[k] \\leq \\theta_2\\} = P_x(\\theta_2, k) - P_x(\\theta_1, k)\n",
    "\\end{equation}\n",
    "\n",
    "Hence, the probability that a continuous-amplitude random signal takes a specific value $x[k]=\\theta$ is zero when calculated by means of the CDF. This motivates the definition of probability density functions introduced later."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Bivariate Cumulative Distribution Function\n",
    "\n",
    "The statistical dependencies between two signals are frequently of interest in statistical signal processing. The bivariate or joint CDF $P_{xy}(\\theta_x, \\theta_y, k_x, k_y)$ of two continuous-amplitude real-valued random signals $x[k]$ and $y[k]$ is defined as\n",
    "\n",
    "\\begin{equation}\n",
    "P_{xy}(\\theta_x, \\theta_y, k_x, k_y) := \\Pr \\{ x[k_x] \\leq \\theta_x  \\wedge y[k_y] \\leq \\theta_y  \\}\n",
    "\\end{equation}\n",
    "\n",
    "The joint CDF quantifies the probability for the entire ensemble of sample functions that for a fixed $k_x$ the amplitude value $x[k_x]$ is smaller or equal to $\\theta_x$ and that for a fixed $k_y$ the amplitude value $y[k_y]$ is smaller or equal to $\\theta_y$. The term '*bivariate*' reflects the fact that two random processes are considered. The following properties can be concluded from its definition\n",
    "\n",
    "\\begin{align}\n",
    "\\lim_{\\theta_x \\to -\\infty} P_{xy}(\\theta_x, \\theta_y, k_x, k_y) &= 0 \\\\\n",
    "\\lim_{\\theta_y \\to -\\infty} P_{xy}(\\theta_x, \\theta_y, k_x, k_y) &= 0\n",
    "\\end{align}\n",
    "\n",
    "and\n",
    "\n",
    "\\begin{equation}\n",
    "\\lim_{\\substack{\\theta_x \\to \\infty \\\\ \\theta_y \\to \\infty}} P_{xy}(\\theta_x, \\theta_y, k_x, k_y) = 1\n",
    "\\end{equation}\n",
    "\n",
    "The bivariate CDF can also be used to characterize the statistical properties of one random signal $x[k]$ at two different time-instants $k_x$ and $k_y$ by setting $y[k] = x[k]$\n",
    "\n",
    "\\begin{equation}\n",
    "P_{xx}(\\theta_1, \\theta_2, k_1, k_2) := \\Pr \\{ x[k_1] \\leq \\theta_1  \\wedge y[k_2] \\leq \\theta_2  \\}\n",
    "\\end{equation}\n",
    "\n",
    "The definition of the bivariate CDF can be extended straightforward to the case of more than two random variables. The resulting CDF is termed as multivariate CDF."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Probability Density Functions\n",
    "\n",
    "[Probability density functions](https://en.wikipedia.org/wiki/Probability_density_function) (PDFs) describe the probability for one or multiple random signals to take on a specific value. Again the univariate case is discussed first."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Univariate Probability Density Function\n",
    "\n",
    "The univariate PDF $p_x(\\theta, k)$ of a continuous-amplitude real-valued random signal $x[k]$ is defined as the derivative of the univariate CDF\n",
    "\n",
    "\\begin{equation}\n",
    "p_x(\\theta, k) = \\frac{\\partial}{\\partial \\theta} P_x(\\theta, k)\n",
    "\\end{equation}\n",
    "\n",
    "This can be seen as the differential equivalent of the limit case $\\lim_{\\theta_2 \\to \\theta_1} \\Pr \\{\\theta_1 < x[k] \\leq \\theta_2\\}$. As a consequence of above definition, the CDF can be computed from the PDF by integration\n",
    "\n",
    "\\begin{equation}\n",
    "P_x(\\theta, k) = \\int\\limits_{-\\infty}^{\\theta} p_x(\\alpha, k) \\, \\mathrm{d}\\alpha\n",
    "\\end{equation}\n",
    "\n",
    "\n",
    "Due to the properties of the CDF and its definition, the PDF shows the following properties\n",
    "\n",
    "\\begin{equation}\n",
    "p_x(\\theta, k) \\geq 0\n",
    "\\end{equation}\n",
    "\n",
    "and\n",
    "\n",
    "\\begin{equation}\n",
    "\\int\\limits_{-\\infty}^{\\infty} p_x(\\theta, k) \\, \\mathrm{d}\\theta = \\lim_{\\theta \\to \\infty} P_x(\\theta, k) = 1\n",
    "\\end{equation}\n",
    "\n",
    "The univariate PDF has only positive values and the area below the PDF is equal to one. The latter property may be used for normalization."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Example - Estimate of an univariate PDF by the histogram\n",
    "\n",
    "In the process of calculating a [histogram](https://en.wikipedia.org/wiki/Histogram), the entire range of amplitude values of a random signal is split into a series of intervals (bins). For a given random signal the number of samples is counted which fall into one of these intervals. This is repeated for all intervals. The counts are finally normalized with respect to the total number of samples. This process constitutes a numerical estimation of the PDF of a random process.\n",
    "\n",
    "In the following example the histogram of an ensemble of random signals is computed for each time index $k$. The CDF is computed by taking the cumulative sum over the histogram bins. This constitutes a numerical approximation of above integral\n",
    "\n",
    "\\begin{equation}\n",
    "\\int\\limits_{-\\infty}^{\\theta} p_x(\\alpha, k) \\, \\mathrm{d}\\alpha \\approx \\sum_{i=0}^{N} p_x(\\theta_i, k) \\, \\Delta\\theta_i\n",
    "\\end{equation}\n",
    "\n",
    "where $p_x(\\theta_i, k)$ denotes the $i$-th bin of the PDF and $\\Delta\\theta_i$ its width."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAGJCAYAAABy9cILAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvFvnyVgAAIABJREFUeJzt3XmcbXdZ5/vPs/euuU6dOSfDCZwAYUyAYKTFAVGgjaikHbpF2vF6O+1t0+pVr6J4HfBqO6J9u3lpR6GdsAGVIWq4EZwARUjCEEhC4CRkOCfJGXKmmmsPz/1j7wN1phq+6+xataq+79erXklV7WevVbvWXvU7v+f3PL/ITMzMzMyqrlb2CZiZmZldDB7UmJmZ2YbgQY2ZmZltCB7UmJmZ2YbgQY2ZmZltCB7UmJmZ2YbgQY2ZmZltCB7UmNm6FRGviohXlX0eZlYN4eZ7ZrYeRcQu4G96n74yM58s83zMbP3zoMbM1qWIeBPwLqAOvDozf7DkUzKzdc6DGrM1FBH3AD+Ymf9Q9rmcLSL+ADiQmT9T9rkoIuK/AIcy87cLPs+jdAdRHz/r6x8Fvi8z7yny/GbWP15TY7YCEfFQRMxGxNSij/++gphXLP5aZj6vHwOa8x2rD89/+uc/FBF/EBHjZ31vMiJORMQ/R8QPRETtAvGnPy6/iOe3G/hu4H+c9fXrIuKfImImIj4aEU9Z5nm2A5cB953n278BvOFinbOZXXwe1Jit3Ddl5viij5vLPqE19k2ZOQ68CLge+JmzvrcFeCrwK8BPAm8+X/yij8cu4rl9L3BbZs6e/kJE7AVuA34V2Ak8eNY5n8+1wP7MnDvP924FviYiLr0oZ2xmF50HNWYFRcRPRsTB3kzF/RHx8oj4Y+ApwF/2ZiV+ovfYL8yo9P7//4qIuyNiOiLeHBF7IuK9ved6f2/m4PRxXhcRD/S+d29EfHPv6xc61uUR8RcRcSQiPh8RP3TWeV8XER/rPd/bgeGV/LyZeRB4L3DNeb53MjNvBb4d+J6IOOcxK3g9/31EfDgi3h4Rj0fEoxHx9cuEfT3wj2d97TeB38vMW3uDnbcBX7rM8zwf+HTvPEYj4k8j4p0RMd4b6NwFfN1qfyYzWxse1JgVEBHPAm4GvrQ3U/F1wEOZ+V3AI3xxduLXLvAU3wq8Engm8E10Bws/Deym+/5cPBB5APgqYCvwC8CfRMRl5ztWL/Xzl8AngSuAlwM/EhFf1zvvQeDdwB8DO4A/653LSn7mK4FXAR+/0GMy86PAgd75rta1wAuBt/fO/b8Cv7uCmPsXneMEcCPw+4seUwPONwNz9vN8KiKuAv6p95zfmplTve/fB7xgZT+Gma01D2rMVu7dvTUjpz/+A9AGhoDnRsRAZj6UmQ+s4jn/W2Ye6s1+fBD4SGZ+vDcr8C7gutMPzMw/y8zHMrOTmW8HPge8+ALP+6XA7sx8Q2YuZOaDwO8Br+l9/8uAAeC3M7OZmX8O3LGSnx/4EN1ZkV9e5vGP0R0wnRHf+3j3EnHXAr+Vme/MzA7wR8BTImKpmaRtwOSiz19O9+e7+/QxgbcCDy9zzs+nu6bm74Gfz8xfyDOrKSZ7xzKzdahR9gmYVci/ycz3n/3FiPgR4OeB50XE7cCPrmK9yKFF/z97ns/HFx3nu4EfBfb1vjQO7LrA8z4VuLz3x/y0Ot2BE8DlwMGz/mAv9wf/vD//Eq4Ajgnx1wL/96LPLwGmLrDO5bTjwJZFn+8Dbs3Mbz79hYi4Dfj/LvQEERF0U2pPA96Yme85z8O2ACfO83UzWwc8U2NWUGb+aWZ+Jd2BRNJdmErv/y+KiHgq3ZmWm4GdmbmN7tqPuMCxHgU+n5nbFn1syczT3XkfB67o/SE/bcnKoFWe75fSHdR8aJVx24ArgSOLvvxtdNNyRMRgr4ppOCK+OSJ+p/eYu+mm8E4bAmYWPe9VdBc337rE4a/q/fcVwI9FxPXnecxz6Kb0zGwd8qDGrICIeFZEfG1EDNFdrzELdHrfPkT3X/0XwxjdgcuR3nG/jzMX6p59rI8Ck71FzCMRUY+Ia3qDDYAPAy3ghyJiICK+hQunslYsIiYi4hvpLsr9k8z81Cqf4lq6Kb3XRkQjIr4B+E90Z8LIzAW6VVX/L/BdwH/uxd0GfPWi57kD+OreYukrgT8FXp+Zx3rn+QfR7cuz2POBu3vnfBPwroi4bNHPNgx8CfC+Vf5MZrZGPKgxW7nT1UWnP95Fd0bgV4CjwBN0UyU/1Xv8fwF+prem48eLHDgz76VbzfNhugOYa+kuZD3tjGNlZhv4RroLbj/fO7/fp7vI+PTg4FvolkIfo1ut9M4Cp/iXETFJd4bo9cAbge8TnudaumtfXkI3pfQLdNNW9y56zD8D3w/8p8xs9b72R8CrImKk9/nfAX8FfJbubNEfZ+bvLXqOKznz9Tt97LsBMvPdwC101wGdXsvzTcA/XORSdDO7iNxR2MzWjV466bOZ+VsX+P4uulVRfwssZOZvLPreLwOHl+so3Kv8+iTw/MxsruLcPgJ8f2Z+eqUxZra2PKgxs3UjIj4E/D+Zec6C3l6K7110G/t9hu6Mzcsyc3ptz9LM1qt1M6iJiDpwJ92KjG8s+3zMbO31qrVemJkPlX0uZlY966mk+4fpNraaKPtEzKwcvaouMzPJulgoHN09Wr6BM7t/mpmZma3YuhjUAL8N/ARfLIU1MzMzW5XS00+9nhaHM/OuiHjZEo+7iW7vCOrUv2TUWSozM9skJjl+NDN3r9Xxvu5rxvLJY+1Cz3HX3fO3Z+YNF+mUVqT0QQ3wFcCrI+JVdHcJnoiIP8nM71z8oMy8hW7fCCZiR35Z41+v2Qlmu9gv1szM1qFY22RF1GL5B13A+1pvX24bk4vq6LE2H7l9b6HnGLjsgQtt49I3paefMvOnMnNvZu6ju9ne3509oDEzMzNbznqYqVm1qAUxOLjquM7cvHa8el2KKyI7Yql9elmSmRWkzmAUuP+o99kiM+lFZk4U1Zr1T9oV/HuyrgY1mfkPwD+UfBpmZmabWgKdi7cn75pZV4MaMzMzWx86FSxILn1NjZmZmdnFUM2ZmgSE3GRteEg6nLoWB/Aal+WUkLs3s2WU8P5S15uUseZRFY0BPXjFW69eHEnSXifbKK1GNQc1ZmZm1ldeU2NmZmaVl0Dbg5q1kZnSVKU6TammrQByYUGKU8+102xJcYU4FWRmJSlSJl0TWoOA3nJjrUvIi6riTI0XCpuZmdmGUMmZGjMzM+ufBC8UXitRqxFDq08JRU2bmFJTSLD2K/OLTG+W0u3SqSszK4marq9SxVURVbw7V3JQY2ZmZv2TpBcKm5mZ2QaQ0K7emKaig5pMaAsTY2L6iQJTjfImbWLKS94IE70xVLbWuCuUmW08YiPO2sDa/xlTNlTuBhaofirQA3YzqeagxszMzPqmu6Fl9XhQY2ZmZmcJ2lSrrw54UGNmZmZnSaDAaobSVHNQU6sRQpdfuTS7SB5UJa7FiSJl2WWUV3tDSzMroEgXdbVbvLz5ZqOaf3KrxK+wmZmZncPpJzMzM6u87oaWHtSsjVqNGB9bfdzUxT+V5cgpL3V6s4ROl0XKyJ1GMjPQS7ML3X/EWDn9NDIsxZWlkx7UmJmZWcVVdabGu3SbmZnZhrC5ZmrEKqYo0rFSnKZUp1SLbEqp/pxFRsadApuFmlmflFCVKKd0inR8V9NBMzNanNIJvyRJ0K7gvMfmGtSYmZnZinhNjZmZmVVeVdfUVHNQkwktveHSqg83V8JOYjUxVVbXNqUE5Ne0UMpL3fCzSJNBM1taCVWJ8iaRRe4F4j2vtn2bdrysUoveoJ3VSz9V74zNzMzMzqOaMzVmZmbWN91duqs37+FBjZmZmZ3Da2rWSr0OE1tWHRYdMU88r5cd56TWxli9lOQOxujrVGojI/Ix29NiaaSZ9Y281q1Ad1/13lXk/pPq34RmUworo+O7KtNraszMzMxKU82ZGjMzM+urjtNPayWLlfGtljjVCBCDWol1trVJtCKXYM7OanEllHQjNoVyKbhtKmJnYDWNFGIril6wFCankIAY0O7P0RD/dArLJr7gsB6q6PapqV4yp/RBTUQMAx8Ahuiez59n5s+Ve1ZmZmabWTXX1JQ+qAHmga/NzKmIGAA+FBHvzcx/KfvEzMzMNiOXdIsyM4HTJUIDvY+l5z4jQNl8cVpLrxSaMlSrkcSqqUKpILWjZxFi51KnkcxWQO0MLKaColGgo3ldPKaYQoICaSS1M3CBVJmtTOmDGoCIqAN3Ac8A3pSZHznPY24CbgIYbkys7QmamZltMu0Kbmi5LuaWMrOdmS8E9gIvjohrzvOYWzLz+sy8frAxuvYnaWZmtkkkQZtaoY8yrIuZmtMy80RE/D1wA/Dpss/HzMxss+p4ofDqRcRuoNkb0IwArwR+tR/Hyi1jUlw0C+wI3tZyqKF2ySzQ0VPNE3emtPU/UKADqZjzL2P3YVuG+ruESv0+y9iRXj2m3GpBXBcD6OtUokCKZGRYixPXxuRQCesWRS7p1l0G/GFvXU0NeEdm/lXJ52RmZmYVU/qgJjPvBq4r+zzMzMysK4lKLhQufVAjqQWdUWEaT5ymTOVYpw85pL3EcXxSixvTF1Hn1LR2zCJlnOrUulwKrh3O+mgTpJCgnDYEcvqpjPYO4v05iqR01JTXqLY8IEcK3CtL4D41ZmZmVnmZVLKjcPXO2MzMzCovIm6IiPsjYn9EvG6Jx31rRGREXL/cc1ZypiZrQXts9VOOtQVx+rfA7HiIsblDazAYh49pBwS5cqE2rlWVAXQmtTQbrH0liSuuNhB188QSKpGKHFONVTe3LZK2KnIfkTXEdKJY/dQRlyOUI/q6S3evOOhNdCueDwB3RMStmXnvWY/bAvwwcE5T3vPxTI2ZmZmdIemmn4p8LOPFwP7MfDAzF4C3ATee53G/SLfNy9xKztuDGjMzMztHnzsKXwE8uujzA72vfUFEvAi4MjP/eqXnXKW5sC/IerCwbfXTnPV5bcqwPldgyllcXR9z4kaY6gZtAC2tyWCKcUWo0+qlVK9skuZylVKl17XA9bPWKa9Qm9kBqBtTKpsbnzY8JIV1xrWfszOk33/WWhJ0ipd074qIOxd9fktm3rKSwIioAW8Evnc1B6zkoMbMzMzWvaOZeaHFvQeBKxd9vrf3tdO2ANcA/xDdcv9LgVsj4tWZuXigdAYPaszMzOwcfd4m4Q7g6oi4iu5g5jXAa09/MzNPArtOfx4R/wD8+FIDGvCgxszMzM6S9HdDy8xsRcTNwO10y1nfkpn3RMQbgDsz81bleas5qEltfczChPbjDhRIK0ZTy903psScv1qiCMSwlifO6Rn9mOLGnfIxC6ylKKMk19aXMtZk1Yp09xXPVz5mkQ11a2rLBP2Y6tqYbIgtAepV2nYgaPexpBsgM28Dbjvraz97gce+bCXPWc1BjZmZmfVNv2dq+qV6Z2xmZmZ2HtWcqQlojax+PFbGhqOdQW3cqE6L1poFyqvbWmomBgts0iZuYsfUlHa4AptvduZX1Pvp3GNWrYxcVaUy6RKoKZ3OgtjeAaiPj0txKXbMrRW41tV2FCneK0Hv8NtpaPetuR3V2tCy3+mnfqjmoMbMzMz6JjMqmX7yoMbMzMzOUcVduis5qOnUg4WJ1U9zdsSZ0daoPmU4dFKsPljQpn9be7dLcQCDDxySY2Vi5UJNrZoqkJ5TU1fZaurHdMXVksr4najUNFKh6ieRvLlkkY7mNS3V0RnWXx91eYB8vCKltLYilRzUmJmZWf8k9HWX7n7xoMbMzMzOEk4/rZmAtjLrLA46o0hB0aB20IXt2pTq8OFZKQ6gfcVOKa5+8En5mMzN67EKcYobgBIyOpVKI6kVVwWqpuQ0kniuhSrZ1MZ0BSrZ5OtH3KQ2tmjVVoDcfK8zov8Zm98mpi/FX0lL2z+zFN0+NZ6pMTMzsw2gz3s/9UX1ztjMzMzsPDxTY2ZmZmdIwumntdJpwNzO1b/YHfGnHT2sb5g2t02bDBuY0nLhrTG9/Lwmbr7JkF5SmVPTWlyBLqtrrUgXY3W9SaU6ERexxudbZI1TbUC8ARVYBxZD2iKOGBCv2QKbSzKrdezO2lb5kMrGyAAze7Tf5cLWag0SOhVM5lRyUGNmZmb9kwltz9SYmZnZRuD00xoKYdawI2ZJmuP6L3bouDYdO79d+9U0pgt0zBVLTrNA+kktAc2mVsobYtkoAGLqobPWZetFFNmUUk0FFUkhieerlmYXKemmhGPG2KgWKG5oiZpiAxbEbujtYf31aQ9p196C+DehSHsQW5nKDmrMzMysP7oLhb2mxszMzDaAtrdJWBtZg6aQtQixcKFIWlGtuOqIxQcLYodMgOEj6pRzgSn5kWEpTN3QsiNWWwHUJgp0SxXpVV7iRphq19sCokiX51CrdMRrvciGjWKX3kLHVK+f0H4nWeBeIKe/i2SUh8VNNNX7+trvTSqrakfh6s0tmZmZmZ1HJWdqzMzMrJ+8pkYSEVcCfwTsoTvjdUtm/telYnIomXv66itKak9qc39ZYHq8NaLFNsR9KRuz+kU4MK1dDqE27QPqU+IPKlZcRYHmaXLqqkBKR00HqSkduUEcemO6IimvELMdMSjmAdQUEujVT8MFdkEUm++paeEi6afmFi2VOLurQMXVFvF9It5GZnavfXq3iI7X1EhawI9l5sciYgtwV0S8LzPvLfvEzMzMNiM33xNl5uPA473/n4yI+4ArAA9qzMzMSlLF9NO6OuOI2AdcB3zkPN+7KSLujIg725N6BYuZmZltTKXP1JwWEePAXwA/kpmnzv5+Zt4C3AIw9JQrs3Z89fnXzjatC21rXq/Dq4kp+JGjWlx7SJ8ubA9oY9y6GAd6N2L1pwyxE3ERKW7UB1Arsp5CUGSjULXzbbYLbE4qTo9nU3tjxqDeMkHuZi2WVwNySXduHdOOV6Bj98KEdv20xLJsgJlLtbi6eMlOPFCddI536S4gIgboDmjempnvLPt8zMzMNjsvFBZERABvBu7LzDeWfT5mZmabXVWb75U+qAG+Avgu4FMR8Yne1346M29bKqgjzFRmS5sabY/qJctqtercTu1ch44XSD8Na8fszOllnPWGOF19UtwkssiGlmo5b5GOueoFJB6zyOaJaml2bUgrH4YCZeRiXK1WIB2ovrZFOgqr6UvxXGcvEzfQLEDtvg4Q4ttr9nLt+pm9qsCGsb+jh24mpQ9qMvND6EskzMzMrA+qWP1U+qDGzMzM1pn0QuE1Uxtss+XKcwqkljU+rKUsWm19tHp8UpuObZ3Uqg+KXYPalHOtpXfJrE9pc8d1dVp9vkB1z8QWLXBOr35KsYtx1MXu2UU2tEytsiz1Js/yxpT1cbG6J0voCFugi3GOb5PimtvElGCB+099XtzQssBfMXWTY9riDzpXnZmPxAuFzczMbIOo4kxNdYaNZmZmZkuo5EzNyECT5+4+tOq4ux7dKx1v6xZx00Wgo1ZcXaWlLEbu0ytJGrPa9K/aNAugMa2lSeonxPG42OwPADEVlNP69UNoP6faXK4ItXJK3lwS6Mxqr22KKZ3aaIHqHvXaK7DJKE0tv9Ia1Y65MK7/O3l6j1h9WaD6qblFTCcOimnPiQKNJteYS7rNzMxsw/CgxszMzCrP2ySYmZnZhuHqpzUSQKO2+lzx1ZcekY73wOFdUhxAZ0FcZzCl/Wrmd+glp9HRLuCJh/UumSGWyDYv2yrFDRyelOIAuTtrrUAZcOfUlBQX4joMtdMu6GXkFNgksqZ221V/Jw19/ZjczbpAl+f2thE5VpEFumen+GPO7dHfX9nQYgd3aGse280C14+tSCUHNWZmZtZH6TU1ZmZmtgG4+mkNtTOYaq4+FaCmkUZH1r4Mr9kQSwYfLrA5oDg73mnoF35rRJuOHTxRQmmk2m13SN8EsTYhh0rUDsagp66ijCpXNaVTIBVEXXyDFdiENdRrVnxLdwr8RVHvP+3xIi2pReLL+pRLnpQP+aAcqavioMbN98zMzGxDqORMjZmZmfWPS7rXULNd57HJ1c/Lv2zfful4SqXVaZ87tVuKe/AJLa65Va9E6gxoE3fTC/qF3xzTKl/GxZRXbV7vtFtTu/QWSB/IVTod7TqIAqkyeePFAimdbGqbaNZG17YqCIAB7VrvTOjn2tyiVaS1RrRrtjlW4F6gFTTKFUwAT9u3+s70ANML2vvk5FwJ110B6UGNmZmZbQTuU2NmZmaVlxUt6fZCYTMzM9sQKjlTEwH12urzqB8/coV0vJFBLW8PcGpOK7EeH9d2Hz51RF8T0RoT12EcKrJmRAurzWvn2hrXd4QemNeugyiwjodY438pFegIW6jcWRTqeiX1dRW7SgN0tojtFgp0pM669nPO7dB+l2pZNkBrXF8PqHripNYz4dKtp6S4a7Y9IcUBfFyO1HlNjZmZmW0Arn4yMzOzDcIzNWukHh22Dq8+PfPcrdrU3xVDx6U4gAdntdLsI3NbtLgRLW0F8NCj2rnO79RTOvV5LW5hq3bpjk3pqaDW9lEpbuBBrWwU0DdQFFNeEXr+QN2QsFAZuZp+EkvB5a7SQLTFzVt36GXA7SGx9UFLO9eFrfofwdql2r3r8h1aKgjg6RNah9/dQ9pGs1cP6+mntVbVbRK8UNjMzMw2hErO1JiZmVkfZaE16qWp5KBmtLHAddsPrDputq119FRTSAADoXUjPjavTTlvGRTzOcDQhBbbPqSnn1RqlcXcLv1cR57Qpsdzh74rZRwTp9bF7rVyugsItaJoTr9m1c7JbNfa13ZGC1zranFPgZTX1OXaLV7dmLKmF4rSmtZe26Hdekp5qqWlPq8aPSrFHW1pywrK4uZ7ZmZmVnlJNRcKe02NmZmZbQiVnKkZr83zFeOfXXXcv0w9QzrergFtpTvAR0/uk+LmWuJGjwMLUhxAo67Nj8/t1Kd/52a1S7A5po3HJx6VwgBY2K41TyuSnKtPim9RtYpJ3ZQSQGz0VojaDE9MeXV2jGnHA1ri5q1qAz3Q00gt8cecvURPlQ1vnZPinjKuV6eO17X75VxH+13eM3W5FFcO96kxMzOzDcILhc3MzGxD8JoaMzMzq7zM7qCmyMdyIuKGiLg/IvZHxOvO8/0fiIhPRcQnIuJDEfHc5Z6zkjM1cznA/fOXrTru+aPagop2gbHfq3beLcV9auRKKe6Bab38fGJEy2lP1bROu6CXgHbWvopclvUCXXrHxW6yYvlwZAkvbIFNInNA3HhxULv11QpsTtrZof2cszv1MvvmuBbXEt/SQ/smtUBgy4i2zqkeeo7kmrHVtwYB2Dd4RIp7fO4lUtxGFBF14E3AK4EDwB0RcWtm3rvoYX+amb/be/yrgTcCNyz1vOtipiYi3hIRhyPi02Wfi5mZmXW3SSjysYwXA/sz88HMXADeBty4+AGZubhR1xjdSvMlrYtBDfAHLDP6MjMzs7WTWewD2BURdy76uGnR018BLE6fHOh97QwR8YMR8QDwa8APLXfO6yL9lJkfiIh9K338YLR46uDqOzqO1bTpzWbqL9NfH3uBFHf58Akp7rlbHpfiAE7MayXL2/foU84np7dJcUNHtfH47E59HN9paIvm6nNaV2kAJrT0U+Oo1oYgC3QUVku6o6W22oXOmJbSUY/Z2qKn59Q0UnugwGJNMXT+Ui3NNtDUr5+5hnaffeH4I/oxUyvNPtjcIcXtHdbLz8twERYKH83M64udQ74JeFNEvBb4GeB7lnr8uhjUrERvhHcTwK7LK7SgwszMrGKSlS32LeAgsHjx6N7e1y7kbcDvLPek6yX9tKzMvCUzr8/M6yd2VGYsZmZmZue6A7g6Iq6KiEHgNcCtix8QEVcv+vQbgM8t96SVHB2caI3y7qMvWnXcc8e11MzhBX0TsmeMHpbiRsVU2cm2Xok02tBKkR6eLlC9sk07Zuekdsym/vIwdFKLSzFtBYDYILqzRUsl1ia1CjhA3mVUrUQC6AyIxxzWjlmku29NLJyaXn2h5xe0xOK5kR3a5q2XbhU3YAVetmfZv1fnNdPR7z/jde16Hwjtl/mLl+i1ML8sR+r62XsvM1sRcTNwO1AH3pKZ90TEG4A7M/NW4OaIeAXQBI6zTOoJKjqoMTMzsz7K/jffy8zbgNvO+trPLvr/H17tc66L9FNE/C/gw8CzIuJARHx/2edkZma2qWXBjxKsi5mazPyO1Ty+Hh22DmjTo4rxupYKAthan5HidjS06pWjLT1V9rTx1VeUAVz2NH3K+QOPPF2Km9+qVS0MntTH8XM7tH+1ZE1f2L51v1Y5lTXt52ztEru1FVAkpRNt7c6pFjS2B/XrZ+YS7edsaLcQABYmxLg57QU6MSTmu4DPz+yS4q4a1e5bgFRFC/Dwgnaubzi6R4rr0tJzm826GNSYmZnZ+lLFvZ88qDEzM7NzeJduMzMzq7zEMzVrZqI+yyu23rPquP3zl0rH+/fbPyLFAXxifq8UVxd3JLxkQF/fMiTuLvmho8+Qj9npaG+azoi6Y2OBzSXF0Plt+o1hdo9Wrjr6qLYQIzr6P82yJnYU1hsKs7BNW6+krsVpjhfouCy+tHM79UMu7NZKjxvitX7VtmNaIDAk1rwfL9Cn4YnGViluRtxRd7Qm9mgoQwIVHNSsi+onMzMzs6IqOVNjZmZm/eU1NWtkPgd4cOGSVcd9dkYrp/ulyVdJcQDXTWibrW2paZ0u//7JZ0txALuHtY0pnzquTzkfOKVN/zZHtEt3bqd+yatpkrGldjNZRnNUm0ydfuqYFFefL5ALEhUp6VbvuuoUdXNMP9e21uSZ5lb9L4vaGfiK7dqGuvvGnpTiQE8/7SmQcu+IOeXtjWkpbpvY4qM0HtSYmZlZ9fV9Q8u+8KDGzMzMzuWZmrWRQDMLVCGs0vUTD8uxakfhO6f2SXFPHz8ixQEMhNa99lRLnFcH9m07LsXdfUKreGhtF3cVBAaPaG8XsVACgFBvKmJapjVSoDpMfEsOTOkpr/ag9i/J2V3ayXYKpMra4r6LtTn9mGpn4IlBrYv6rgGtEzrAx04+RYq77hL9/vxkW+ugfcepfVLc3mEtrWf6rwW8AAAgAElEQVQrV8lBjZmZmfXRGmxo2Q8e1JiZmdm5nH5aGwudBo/NbV913KPT2/pwNkubbmm5h2u2PCbFHZzXf0a1+mDngFYJAPDBU9qGlqNbtOnx6baeXmmNiZU2j+v/2pmf0GK3HdFSifNb9VtCp6Gda2tY/53Mb9NiG7Pa73J+txQG6Omn1iVaU0yA5+w9JMWNNrQmccdbWtUdwMt33ifFfXzmqfIxZ9vaxrhXjWhVXvdMXi7Flad6MzVuvmdmZmYbQiVnaszMzKzPnH4yMzOzDcGDmvVtuK6tGRkRN3oEaNa00lF1bcwLxh6V4gA+Ly4YeGJhQj7ma6+6Q4r7yPGnSXEfffIqKQ6gM6iVHk/u07O8Wz+n5bSPX62t5Ro6qd/F1G67nQLdGdSNKWf2aOe6oDXABqC5S7v/DE5o68cAHjulvTevvVJb0/fMkSekOICnDR6W4tTyaoDnjWs/56cmtY2K59oV+pNb0Q0tK/QKm5mZ2Vqp4t5PXihsZmZmG0IlZ2panTpPzK1+WnWwppW57p/aJcUBNNRdEEVFOi1fNaR1I368ppeR/+0RbQPO6aaWXhndrm3wBzD/yBYprtYs0IV2RIsLsXHy3E79XNVjInYFBr2L8ewe7Z+gQ08WKM9vaMccHtTT35du0TapfXhmpxT3rdvulOIAPjjzTClOLa8u4l9tfVCKe9/R517kM+mzCs7UVHJQY2ZmZn3mNTVmZma2Ech7z5WokoOaBDpCp8MnprX0QRE7hrUNLS8bPinFHW3pP2NNnGu8ZOCUfMyX7/6MFKdWH6g/I8CRfVrs7LS+o2VzWsw/iWpaI1kAxo6JlUi7CmzYKGY+o6Udc/rp+oaoo9u1e8HucX2TyJr4V0mt0vnlR79BigN42c7PSnFPNvUuxmqH31NNrT30V+58QIoDeLccKUoqmX7yQmEzMzPbECo5U2NmZmb9FF5Ts1Zm5ge5a/9TVh3XGNamjndMaNPGAI2aVv10aF5LIw0UqLb6hh2flOJOtEflYz4yr1VZ1MSfc75A86uBhlY9Nz2pp5/mdoqbaC5oN6Oh41IYACeeId4Ai9w3xenx9ogYKKatAObntc0TJ+eH5WO+9DIt3fG0Ya0R3leNfk6KA/iryefLsaqO+Ef7xJyWFn7zJ79ciuv66wKxogqmnyo5qDEzM7M+q+CgxmtqzMzMbEPwTI2ZmZmdq4IzNRUe1Kw+FzoxPicd6fikvmZk67DWwXbboLaO54VbHpHiAJ5oarv1FelivLWh/ZynWto6g+966kekOID/du/LpLjGCf31ae7SusnWjmlv7azpa0Y62pKRQr0wmlu0tVU5LK49G9DXrD1l9zEpbod4DwE4NKetzdsndhffLa4hBHh4VltfN97Q+xBsHdRe2xOD2v1nZkL/Xa45b2hpZmZmG0UVm++tizU1EXFDRNwfEfsj4nVln4+ZmdmmlwU/SlD6TE1E1IE3Aa8EDgB3RMStmXnvBYM6AVOrn9I/OailkdqT4rw68EBztxS3bZ+WKqujT/9uqWvHvG9W68oJcM3IASnuZEv7XX56+gopDqBWE9+lV03Lx4zj2jR3c4fWviBa+rXeElNBRW5+nRHtmFHXDlof1sr6AU6JpdmXj+kdu/eNaps9Pi62av6fJ18oxQFsG9BSM7c+eI18zEu2aN2a62Ka7eTRcSnOVm7ZQU1E7AN+EHg6cAz4BPCXmfnwRTqHFwP7M/PB3vHeBtwIXHhQY2ZmZnaWlaSf3gN8hi/OprwA+EBEvCkitA0wznQF8Oiizw/0vmZmZmYliSz2UYaVpJ/qmflmgIg4lpn/ISIawP8J3AJ8Tz9P8LSIuAm4CaCxdTvDh1afflpoatO/td1aWgb0ae6HTm6X4j5Uv1qKA9gzPCnFPblQoDpMrH46tqBtYrd3WG+Z+5IrHpLi9p/cJR/z0c4OKa69oFVcdQb1O9Hg5Vqabf5xfUPCwa3zUlyId9zREb3SZt9W7dp78dbPy8e8eugJKe5Ia0KKOyRWUEKBisZnflQ+5u9++GVSXG1ITEPO6pWQpahg9dNKZmreHxE39/4/ATKzlZm/DrzkIpzDQeDKRZ/v7X3tDJl5S2Zen5nX10f1m6CZmZkto+gi4XU8U/OjwE9FxJ3A5b0Zkxm6AxptFdqZ7gCujoir6A5mXgO89iI8r5mZmW0iyw5qMrMD/FJE/BbwCuCFwHbg08Dri55AZrZ6M0G3A3XgLZl5z5IxdVjYvvphYG1Om0rLJ/QN5VRHW1q1/SdbekHbl1z66PIPOo8im2jOtLVlWV82oW3U98GTenruwLRWEaJWvQCMjmrpldaQdv3MD+jVPc057dobvkyvDpuf1aq1du3QUq27xvRzfeC41lzupVoGEoA7p58mxT1r+HEpbv/MJVIcwGdOarF3p159GQ1tOqEjpndjXTRRWYUK9qlZ8V0oM2eAW3sfF1Vm3gbcdrGf18zMzDRVbL5Xep8aMzMzW4cqOKip2mSYmZmZ2XlVcqYmOjAwufr1MaNampi5XXpZ28I2bahbE0vBi7jriSuXf9B5fPnlesnp58VN7J5Y0EpO//Oev5XiAH7z8a+T4joFyiJPzI1IcTMzWoVg59SgFNcN1sJmZ/QuxpddqW0SOdjQOi4/c8thKQ7g2q2PSXEDoa9zunxQb2GgODKvd8wdqGk/5yNHCyw6EkusB8VNagv8KstRwZmaSg5qzMzMrH/KbKBXhAc1ZmZmdq4KNt+r5KAm2jB4YvVxLW0mn+EC3Xg6A9pFURtqSnGnjuvdfa/ae0SKe2BS75j74h3aFmJTYin4eyefL8UBXDEsXHRArcAc7pfs0Mrs35/P0g44oXV4BqiLpf0Lbb3LqppG2jOibWQ4XtdK7AHmU0uzbalpGz0CHGtr6aB/PK5dP9dMaCk2gAdntPvIo6F1XweoLYj3Z+32zNDaZgOLq+BMjRcKm5mZ2YZQyZkaMzMz6y+vqVkjkVATZp3VKcNpvWElza3alHzjfq26pz6gX4UHR7WOuS/dp3X3BbhCnI89OK9NOV82oKWQAPYOapU28x39bfb5aW1K/ssufUiKu/f4pVIcwJZBLTVzbE5PmZ6a07o1X7fjgBT32Ly+YeNeMX3ZTP36ubRxUop79ri2EeYHjj5DigM4Mq1XTqlqs1qyoqFlLxnQGlmXx4MaMzMzq7yKVj95TY2ZmZltCJ6pMTMzs3NVcKamkoOa2gKMH1x9a8b2kFpzX2BC6zEtVl3HUyvQVyDF2Cfnte61AHfHXinu1Ts+LsU91tTLP0+2tJz/t2y/Sz7m7y+8VIr79PHLpLhXX363FAfwT8e09RTHpvQ1NXu3a+tUbtz+MSnuw9P6Lu/bG9oO33NiKTjA52d3S3G7xMUfsy39XE+e0q6DzhGtvQPAiLjGZUDsfNDRX55yeFBjZmZmG4HX1JiZmZmVpJIzNZFJY371Q8iFCW0MVxeO9QViNqjW0gLnxRJygDyuTeM+PKandNQOpLefuFaKe9qI1jUZ4LnDWhnwHxz5CvmYe4a1+fGGuDngyZaeCmqIHYWftlNv2X3tNu36+b0nvlqK+7JtD0pxANvqWs7iUFMvIz/V0krebz/4HCmuXSD9PTY2J8W1PiW2igfa2ssjtRQBGDqh359tZSo5qDEzM7M+c/rJzMzMKi+/uFO3+rGciLghIu6PiP0R8brzfP9HI+LeiLg7Iv42Ip663HNWc6Ymobaw+mm8rGkb5w2e0oer7WFtOlYslIDUx6mzT9VSFoN1LQ7gEye06qdnTRyW4u6ZukKKA3jh8CNS3DNG9ZTXc4YPSnFHWlpH6g+d0Kt7niemEj8vbmQI8NItn5HinhzRKtk+Nr1PigPY39kjxe0e1NvQ7hzQbiRXb9eu2X+6X+8ozJR2fx4u8FcsxDSSSk1blaaPMzURUQfeBLwSOADcERG3Zua9ix72ceD6zJyJiP8D+DXg25d6Xs/UmJmZ2bmy4MfSXgzsz8wHM3MBeBtw4xmHz/z7zDy9GO1fgGX/FexBjZmZma21K4BHF31+oPe1C/l+4L3LPWk100+i7Z+ZleIWtuodk1qj2kuszjgX2PuOoce1n/OJhrYRJkDjci119cCUlrK4buujyz/oAu6YvUqKK7KJ5kBor4/asE2tmgJ4eHaHFPeiiYflY/7dqedJcVcPaxs2bm+IXdeA60YfkuKKNIx8x2NfIsUdOK69p6OuV/cMHRqU4uoL8iER99NlYKaCK2hXKbgofWp2RcSdiz6/JTNvWfW5RHwncD2wbNniphrUmJmZ2QoVH9QczczrL/C9g8CViz7f2/vaGSLiFcDrga/OzPnlDuhBjZmZmZ2p/7t03wFcHRFX0R3MvAZ47eIHRMR1wP8AbsjMFVWHeE2NmZmZranMbAE3A7cD9wHvyMx7IuINEfHq3sN+HRgH/iwiPhERty73vJWdqVEaV0Zby/d2BvWx38C0NtQdnNTOtTmqlUUCDC87sXd+C9v1Y6pdaCcGtO6jd514ihQH8GU7tG6yA2tdNwrsqE9JcV+77T75mI8saOucnj2klYIDPDSnHXOmo3XPftbw41IcQDO198nHp/Rr9tJRbXGeuqamM6ffCxZ2aPfKsUf0LsbDx7T7T2NOjJvV16yVos9LhzLzNuC2s772s4v+/xWrfc7KDmrMzMysjyq4HtqDGjMzMztHFXfpruSgptbsMHRo9aWVOaRNjY4c1Ms45y7VNlvLujalWqS8sblFixs8oU85P35S63w739Yu3Vdcdr8UB/DYnFZau6dxSj7msbbW+fazs5dKcbUCd7EDc1rKQi1bB3jR2ENS3Mm2tnHnCTEOYN+g1qX3JRMPyMf8zc+sevYegK2jWnq3vV+8iQCjWpU9gyf1a7beFGP1jFe1VHBQ44XCZmZmtiFUcqbGzMzM+mhlWx2sO9Uc1HQ61KZXPz3aqWupoGjq0+Ni41JyQJtEm9k9ph0QaDfFuKECV/6CdglOz2vdR//xsL7hXruj/U6uGTsgH3NbXduQ8JJBLeU1HOJFgL456fiEluoAaIottL98dL8UNyBW6wH8/tGvkuI+fEjrZF3E4c9r3aFrw/rr0xarTIcKbDhcm9fOt94U42aqtaNlFdfUlJp+ioh/GxH3REQnIi7UddDMzMzWWn83tOyLstfUfBr4FuADJZ+HmZmZLRJZ7KMMpaafMvM+gIhVLiXPhIXVl/lEU0tZxEmtkRlAXrFTiuuI6aehE/r074CW6ZDPFWBhQUsJzl2tpQQ7I3qq49IxrZGZ2ugN4EhLqw5TG73NdLT3CMD/fsUHpbh75y6Xj3k8tXTrew6/UIp7/tZztqZZscdmteqwy8b16rlPPaq9tjmk3UcGn9A3/x17QvsrWBNTQaCnkQaOapsjI1a12spVZk1NRNwE3AQwXNfLBs3MzGwFKrimpu+Dmoh4P3C+phmvz8z3rPR5etuV3wKwdWhPBV9qMzOzinD10/kpezeYmZlZeYJq9hisTPrpDLUajAyvOixmtB0bc0Ivk67NaCWytYa2TmXkqH4ZTu3V1n40xLU4AHO7tX8KtE5p59qe0E92sKaVY47WxJ1CgU9OXynFbW9oXbCfPaJvLvnByWdJcV86pm0UCjAtrld6zoTWa2G2wJqjObEL9qFpPd3eaYr3kYe1tTENcakJQGNG3Pz3pF4m3ZjU3pu1E9r6OuXvlq1O2SXd3xwRB4CXAH8dEbeXeT5mZmbWU8GS7rKrn94FvKvMczAzM7NzVbH5XjXTT+02nBLKrOtamWvM6mXAMSZugCeUrAO0d2kbGQJEW7uCG3N6yqsuxrbEN9tjj2gl9gA7nqmldA4saN1ZARY62lt0//QlUtxrtn9UigMYFDemfHhhl3zM8br23jy2oKWUJ1t6ef5jU1p5/olT+iaatQGtZFnsCMDo4bXfXLIxpe/iGy2xHHxYvw4qxYMaMzMz2xAqOKgpu6OwmZmZ2UVRzZmaTpLzq1+1HhPa9C+1AoVtaqw4vTkwqW9I2BnUzrU5qo+N1cKg+int0u3s1F+fJ8QqlM809JTg08ePSHFXjRyV4t47ea0UB3C0OS7FvXziHvmYH5h8thTXqGmpsmeNH5LiAO4/pqUEt4zr6e+pe7dLcXUxo1MTU0gANTGn3BkQc2VAY0p/bSVTWgq7FCVudVBENQc1ZmZm1l8e1JiZmdlG4Jma9W62QGco1RZtSp66ltIpsrnkwCmtidVIgU3amlu0qeNsaMdc2K1vfnfsSe132e6s/dI1NW11oqltMApQE++Ad848TT7mB488XYq7bscBKe6eU/rmm8eOatdPrVFgw8aW9j4ZEPfwHXtMr0SqzWspwejorw8t7ZjMiXlzsaq1NBUc1HihsJmZmW0Im2umxszMzFbE6SczMzOrPu/SvYYiiIZw6mJH4ZzWy/DkjsLimprBgye04wGt3VrJcntI3+RP3HeRha3iAY/r51q7RCv/bLb0klN1E8Thmla6/rEjz5TiAJ46cVyKG6/rG35+7Z7PSnHv2H+dFDc7rXeSjbr2FyIe1zdBHH1cixP3bmVhq/4npTGj3fOGH9TaFwAwv8ZrXDoVGyVU7HTBa2rMzMxsg6jmTI2ZmZn1TeA1NWsrV/9q56lJ6VChlmUDqJthLmgpi85ONS8DtXltznnwlFgWCd13jkAtXZ/brb9LO4e1NMD8Hv2Yh8QuxnfxFCnusvFTUhzA1kGtZcK8uGknwD0nL5PiJka0lNeuLdNSHMAjB7WNO+vzesuEjphtHTuklUnX5/VrfeigeO0NFPgzpt6fW9q9MsW40nhQY2ZmZhtBCJMHZfOgxszMzM7k6qc1lEkurL66I9tamiRP6lPyIW5MKVV3AbUnT0pxoKeuhg/rnZobM9r8+NRlWipo9LEClUi7tHd4W0wlFqF2933OxBPyMT954gop7visWCEI7BzR0kGHn9TSep1mgd/ltBbb0DNeDB/TrgN1c8mBUyV0zJ3W7z+dE9r9Moa0+3pWraNwBVVzUGNmZmZ95YXCZmZmtjF4ULOGhOZ0akpHTVsBIKTJABDPVY4DanPa1GgO6seMzoAUp06rL0zolSTz28VAsekawMkZLc2mpp+mmnpzwmZHS6+cmtUb2g02tGoSOY20oLf2qjW12LpYoAMwdEKrYlL/ha5WUHaPKR60Kd5jgRgRGxu2xU00i/wtKYFnaszMzGxjqOCgxh2FzczMbEPwTI2ZmZmdKZ1+WjPZbtMRyqxjUFwvUCQPOj4mhamdJ2NySooDoKGVdMeCnkdvjWqXYBlvttqCuB7n4RH5mPkM7feprjUpUl49MaQt/pib1dfxPCn+SgbHtPVjC9P677I+o53s6BH9Ym8PacccO6j9LqMlrjUBOKl1fFe6y39BU7x31cQkR4E1j+j7vuo8qDEzM7Oqq+reT15TY2ZmZhtCJWdqohZSKknt5ljbonUfBeSS7uyI07hiB2OAEDsn564d8jGHjsxox2xoab1OQ+8IO3hKm8pv6Rkd5g9pwY9MatfBtkvEFADw+MkJKe7qyw/Lx7z/M3uluBzQ3l8D0/q/A8cOyqFrrnFKzHUUSdWPiW+UAh3f1TSS3Bm4U7GpD+/9ZGZmZhtBFdNPHtSYmZnZmbyh5drJhBRWrddGtMqFnC3Q0lOkboRZRKHOySK1cqo2r6UPRo5KYQB0BrS3y4Ce0UFd9jY/rr2uU9P6dddoaL+Tzz66Rz5miBVpIW4yOjCpd6QWm2fTFuMAJh7QNnvM0H7OIh2FOaW9UQrdnwe1F1dNP0WjwC+zBFGgmK0sXihsZmZmG0IlZ2rMzMysz5x+Wp2I+HXgm4AF4AHg+zLzxLJxdCugVqszK07FFlixXhPTSClu0hbipoIAiOmVOFWg4Z+4odzgcW3KudbSG73N7dBe27nL9ZRFW80GzWjn2jmhT48viFPVnQk9ZdEQ009qI8XBAoU2A1PafWT0iJ4WbovNLesntHslc3qHOLXhKGKqDJCauAJEXXt/ZUvffLMMVVwoXHb66X3ANZn5fOCzwE+VfD5mZmaW9BawFvgoQamDmsz8m8w8PTz/F0BrOmFmZmYXVWSxjzKUPVOz2P8GvPdC34yImyLizoi4cyHXvhrJzMzM1re+r6mJiPcDl57nW6/PzPf0HvN6oAW89ULPk5m3ALcATMSO7AgldWoetAh1HU99fFyKy3k9p62u/ymS01Zz8LlTe32KdPQcPiaubUj9uoum9tpmQ/v3ypC6QySwsE17betH9HU8dXHpR0dfWiUbO6RdP4PH9Pd0Z0C7DmJe7JgrlkgDIK45KkJeG1NC+4tSVHBNTd8HNZn5iqW+HxHfC3wj8PLMCvZkNjMz22CquqFl2dVPNwA/AXx1ZmqbAJmZmdnFVeJi3yLK7lPz34Eh4H3RTWH8S2b+wIoiY/XTqmpptlI+/sXgNe5YKWz0+YVjqp05he7Op8VWbbPQaGrTv83dWgk5QGtEXIJW4PKpiS/t8GHtoE0xqwdAR0yVFXh9BsRuAg0xbVVf0G/yareF1pie0hl6vEANuiCfPC7HRkP7c9QpkHJX/o50iekn+Xi2UqUOajLzGWUe38zMzM7P6SczMzPbGDyoWTtKSkhNPxVZ6V4bFDu7qumnAtU9oXYULlJVNiN2ed69TYobPqK3A4gUU1cFppwHprXcjNqJOAvcEdpiR+HRQ/oxW9oetfLmkqNH9PfX4CktlzhwdFo+ZsyKVUxqKrquX+tqpWgRa97hN6u1Q6RnaszMzKz6kkItMMriVUtmZma2IXhQY2ZmZufKgh/LiIgbIuL+iNgfEa87z/dfGhEfi4hWRHzbSk65suknaZ1LCeV06nqcaIhJ/yI5W3WqsbP2eeL6YbF0dEjd9hrq49rvRGzqCsDCmBas7gjdmCnQUVirzpfjAOrikpGRJ0voXtvS3iedYb1NQ/2YWNIt7pidUwXW/4xoC6SyjLU4FUzLKPq5piYi6sCbgFcCB4A7IuLWzLx30cMeAb4X+PGVPm9lBzVmZmbWR/1tvvdiYH9mPggQEW8DbgS+MKjJzId631vxvwg8qDEzM7NzXISZml0Rceeiz2/p7eMIcAXw6KLvHQD+VdEDelDTZ3L6ae333gS1c3KthKVZYvfR9o4x+ZADx7Rp7plLJuRjDh/Xrp/2sPY7md4jhQHQEKuABwuUSaul2aotnzspx3YGtWu2NlukY26Bds3K4Qp1NF/7NJK+oaWY99x8jmbm9Wt5QA9qzMzM7EwrXOxbwEHgykWf7+19rRAPaszMzOwM3V26+zqquQO4OiKuojuYeQ3w2qJPurkGNRXq5ljG5pupbkxZZIpbjFWPWD86KUYit8wdPSy2vQWaY9r0eMxq5zo4qec9Rw9r109zXD/m8DEtPdcQX5/mNn1D1MaklkaKoyfkY6KmV8SOwupGvGVRO7fLilTgllFw1cc/mZnZioibgduBOvCWzLwnIt4A3JmZt0bElwLvArYD3xQRv5CZz1vqeTfXoMbMzMzWhcy8DbjtrK/97KL/v4NuWmrFPKgxMzOzc/Q5/dQX1R3UrGUjvQJpK3l1fQnNneRKrSK7IIryhFaFEjt36AcV009q1RRAbV6rJpm5TEuTDJ3Qr/UUU59qKgigPi+macX3V+OkXokUC2J6d8u4fExOaenWFBtqFtn8V1XomOrfkQotZZD1f6FwX1R3UGNmZmZ9kv1uvtcXHtSYmZnZOfq5TUK/eENLMzMz2xCqO1OzljnNAut31j7HXEIr4iIl3eLGebFV7NJbYMM9RsRy3hLS70MntNe109B/l9HW/lnX3KLfhtpDYkuAjhaXA2v//opJ/ZrNGW09V86La4cKlSyv/RtFbYGRnU2yFsfpJzMzM6u8hJVvI7l+eFBjZmZm5/JMzQZVxrSovNGanu6Sj6lOVQM0tB0J89hxKa62basUV4RcygvUmtpbNGvaddARNwoF6AxqU/Ijh8SdMEFOfXYa2rnW5vQOtDE5owUWeX+JaiNaF+z2tPgzlqSMEnTrLw9qzMzM7FzVm6jxoMbMzMzO5Y7CdtHI3X3FFFKRYxapeKgNiOkV8VyzQPVTiNVPRW4MtWkt9VATf8z2iFhVBtTntN9JbaFAynRe3NBSTAnGfIENEOe032UuNOVDqp2B5arEIveflvhzllFxtVk6EXtQY2ZmZpWXlNKOoig33zMzM7MNwTM1G0wpm7sV0JnVmoOp09xRZHPAWbFKp0hFmliNlGPihpaH9fRca8uQFJf1AtedVjxH7fgpLXBUqwoCyCfFKqYir8+sWHmnppHaBarDKrT572YQpNfUmJmZ2QbhQY2ZmZltCB7UmJmZWeVVdKGwBzV9VkZnYJlYbhhiV+BCxLUmnQIl3XI3YnFTQQDGxvRYQW1K717bEG+AOazfhuonxN+nuglrq8D6KHE9V0fsnl2I2jKhhI7mpZRJV600exMpdVATEb8I3Eh3PHgY+N7MfKzMczIzM7NqNt8ru6T71zPz+Zn5QuCvgJ8t+XzMzMwMumtqinyUoNSZmsxcXFc5RiV3mljaZtgwrZSfUawcDbGDMQBNsePpQIH0nFjOG1NaykstBYcCmz0W2CQStWOu2qW3QCpR7hJe0//tqd5Q1TLp2uCgeEToNMXy8zI6Cm8K5Q1Miih9TU1E/BLw3cBJ4GuWeNxNwE0Aw4yuzcmZmZlZZfQ9/RQR74+IT5/n40aAzHx9Zl4JvBW4+ULPk5m3ZOb1mXn9AFqTLzMzM1uBxOmn88nMV6zwoW8FbgN+ro+ns/GV0BW40BSueL5R06pXUp3iBlLsKBxFOp6qVTrq4Yps2NgQq1cKVBTJN87Bta/Yy+kZLVDcXLIQ8T3daZbwh8wppP6p4EtbdvXT1Zn5ud6nNwKfKfN8zMzMrKuK1U9lr6n5lYh4Ft3x4MPAD5R8PmZmZgZeKLxamfmtZR5/Q6rYVKyaRirjeGr1Ss6IaQegJtmy4EQAAAlISURBVKZJcl5rohdqg0GAySktrkB1j5y6EtNs6gasACFWwWWBPyxrvWFs1e4/tvGUPVNjZmZm600CFdwB3YMaMzMzO4v71JiZmdlG4UGN2eqonUtB7M5aoOOpvM5AjQNS7GAbY1qDyjzypBQHEFsntGMW2GRUJa85KvK7VDtSF9kkUtxsVu5+XOT12QTd1yungoOasvd+MjMzM7soPFNjZmZmZ/JCYTOBXAKqTXN35rS0A+jl4EWm5GtiiXWemtQOWCR9cOLkmh9TTbPJ0+oFUiRyyqtAyjQX9BJ06XjOIG0gWckSfQ9qzMzM7FxeU2NmZmZWDs/UWLnEDS1Lqc4Q88uhH1KufpI3bFwQK3SKHLOtT3Grr08uaB2Fi1w/8rUunqtZIV5TY2ZmZhtGBdNPHtSYmZnZuTyoMVulCq2ulzfDFNMOgF4ZJKZ0stnSjgeFKoPWnDitnm09FSQ3lyvyHlGvvQq9L61fqrlNghcKm5mZ2YbgmRozMzM7UwKd6s3YeVBjZmZm56pg+smDGqumNS4Fh4LlvKJUOyCr63/UOCAa2u0ki5SR19UMenVu1t4k0krjQY2ZmZlVX1ayT40XCpuZmdmG4JkasxWSOwoXOWZLS83ImyCGfrZqGikG9NvQWnfbLZTOEcukvUmklSIhK1ja70GNmZmZnauC6ScPaszMzOxcXihstr6VUUmippCgwPmq/8IqkusQq586aoVXAWp3aDmtByB2ay5y/ZhtNh7UmJmZ2Zky3XzPzMzMNginn8zMzGwjSM/UmK2RTVIeq5cQi+t/ClQ71MQuz/Lu5wXI66OK7GKulscW2eW9giW5tl54l24zMzOz0nimxszMzM6UuE+NmZ2rSBfamlhCrKaRagMFOgpXqPRY/Z0UagnQESfGnUKyslTw2vOgxszMzM6QFFtjV5Z1saYmIn4sIjIidpV9LmZmZpteZnempsjHMiLihoi4PyL2R8TrzvP9oYh4e+/7H4mIfcs9Z+kzNRFxJfCvgUfKPhezfiiSsuiI1TZqRVGnqf/LTO7SW0KXZ/WYhTa0VLn6yTagiKgDbwJeCRwA7oiIWzPz3kUP+37geGY+IyJeA/wq8O1LPe96mKn5LeAn6M52mZmZ2TqQnSz0sYwXA/sz88HMXADeBtx41mNuBP6w9/9/Drw8Ipb811Opg5qIuBE4mJmfLPM8zMzM7Cz9TT9dATy66PMDva+d9zGZ2QJOAjuXetK+p58i4v3Apef51uuBn6abelrJ89wE3NT7dP79+eefvjhnuCHtAo6WfRLr2Nq+PgX6tcmKZR2010c9ZhmvTzFre/1Ubw7b95+lqa/PUy/2iSxlkuO3vz//vOg61+GIuHPR57dk5i0Fn3NJfR/UZOYrzvf1iLgWuAr4ZG82aS/wsYh4cWY+cZ7nuQW4pRd7Z2Ze37+zrja/Pkvz67M0vz5L8+uzNL8+S6vK65OZN/T5EAeBKxd9vrf3tfM95kBENICtwJNLPWlp6afM/FRmXpKZ+zJzH92ppxedb0BjZmZmG8odwNURcVVEDAKvAW496zG3At/T+/9vA/4uc+m9G0qvfjIzM7PNJTNbEXEzcDtQB96SmfdExBuAOzPzVuDNwB9HxH7gGN2Bz5LWzaCmN1uzUn3NyW0Afn2W5tdnaX59lubXZ2l+fZbm16cnM28Dbjvraz+76P/ngH+7mueMZWZyzMzMzCphPfSpMTMzMyusUoOa5VoqG0TEQxHxqYj4xFmldJtSRLwlIg5HxKcXfW1HRLwvIj7X++/2Ms+xTBd4fX4+Ig72rqFPRMSryjzHskTElRHx9xFxb0TcExE/3Pu6rx+WfH18/QARMRwRH42IT/Zen1/off2qXsv//b0tALRda+28KpN+6rVU/iyLWioD33FWS+VNLyIeAq7PTPeJACLipcAU8EeZeU3va78GHMvMX+kNjrdn5k+WeZ5lucDr8/PAVGb+RpnnVraIuAy4LDM/FhFbgLuAfwN8L75+lnp9/h2+fuh1vh3LzKmIGAA+BPww8KPAOzPzbRHxu8AnM/N3yjzXjaRKMzUraalsdobM/ADdVfOLLW69/Yd0b8Sb0gVeHwMy8/HM/Fjv/yeB++h2OPX1w5KvjwHZNdX7dKD3kcDX0m35D5v4+umXKg1qVtJS2bpvmr+JiLt6XZjtXHsy8/He/z8B7CnzZNapmyPi7l56alOmVxbr7Q58HfARfP2c46zXB3z9AN0MQ0R8AjgMvA94ADjRa/kP/jt20VVpUGMr85WZ+SLg64Ef7KUX7AJ6jZyqkYNdO78DPB14IfA48Jvlnk65ImIc+AvgRzLz1OLv+fo57+vj66cnM9uZ+UK63XJfDDy75FPa8Ko0qFlJS+VNLzMP9v57GHgX3TeSnelQbz3A6XUBh0s+n3UlMw/1bsYd4PfYxNdQby3EXwBvzcx39r7s66fnfK+Pr59zZeYJ4O+BlwDbei3/wX/HLroqDWpW0lJ5U4uIsd6CPSJijO5mod7481yLW29/D/CeEs9l3Tn9B7vnm9mk11Bvoeebgfsy842LvuXrhwu/Pr5+uiJid0Rs6/3/CN0il/voDm6+rfewTXv99Etlqp8AeqWBv80XWyr/UsmntK5ExNPozs5At1v0n2721ygi/hfwMro74x4Cfg54N/AO4CnAw8C/y8xNuVj2Aq/Py+imDhJ4CPiPi9aQbBoR8ZXAB4FP8cU9yH+a7rqRTX/9LPH6fAe+foiI59NdCFynO4Hwjsx8Q+8+/TZgB/Bx4Dszc768M91YKjWoMTMzM7uQKqWfzMzMzC7IgxozMzPbEDyoMTMzsw3BgxozMzPbEDyoMTMzsw3BgxozMzPbEDyoMTMzsw3BgxozW1ZEvDwi/rjs8zAzW4oHNWa2Ei+g2/3UzGzd8qDGzFbiBcDHI2IoIv4gIn65t/ePmdm60Vj+IWZmPJ/ubtS3A7+fmX9S8vmYmZ3Dez+Z2ZIiYgA4Snfzxv+YmR8u+ZTMzM7L6SczW85zgDuAFtAu+VzMzC7IgxozW84LgH8GXgP8z4jYU/L5mJmdlwc1ZracFwCfzszPAj8JvKOXkjIzW1e8psbMzMw2BM/UmJmZ2YbgQY2ZmZltCB7UmJmZ2YbgQY2ZmZltCB7UmJmZ2YbgQY2ZmZltCB7UmJmZ2YbgQY2ZmZltCP8/Sp8WM5t9m+AAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x432 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjUAAAGLCAYAAAA/PWMAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvFvnyVgAAIABJREFUeJzt3XuY5WdV4Pvv2ruqu8kFIgQRkgiZER2VqxMz8ug5ooQxcsSIzjjg4wUPjxnnmBkdHBXFB5E5Ot515kyOM41wUEQRFbB14okwwmFULgkiAYJouCgJSISQ9LWue50/9u6hursuXevtvX/12/X9PE893bVr/y7127e33vWutSIzkSRJ6rtB1ycgSZJ0ITiokSRJc8FBjSRJmgsOaiRJ0lxwUCNJkuaCgxpJkjQXHNRI2tMi4hkR8Yyuz0PS3hfWqZG0V0XE5cAfTb59emZ+qsvzkbS3OaiRtGdFxM3A64Ah8PWZ+T0dn5KkPcxBjSRJmguuqZEkSXPBQY00QxHxvoh4atfnsZmIeEVE/J9dn0eLiPgPEfF9F2A/H42IJ29y+zsi4otb9y9pOhzUSOchIj4SEaci4viGr/98Httct/G2zPzizHzzlM7vup3v2XSMb4mI2ye/+8cj4g8j4is2HP9URByLiPsj4s8i4rsjYnDWPja7jo/a5FifFRE5+fnJiPjYToOViHg48O3Afz3r9idHxJ9O9vOOiPjcHfbzWcAjgfdv8uOfA16y3faSuuOgRjp/z8zMSzZ83dT1Cc1KRDwf+CXgJ4FHAJ8L/N/ADRvu9szMvBR4NPBTwA8BL9tkd2dfx49tcp8nAX8/+flFwL8CfjEirtzmNJ8L3JKZpzac95XALcBPAw8DPgT86A6/7uOBuzJzaZOfHQG+KiI+Z4d9SOqAgxqpUUT8UETcM5ml+EBEPC0iXsn4g//3J7MNPzi57/+cUZn8/wci4o6IOBERL4uIR0xmQI5FxBsnswanj/OCiPjg5Gd3RsSzJrdvdaxHRcTvRsTfR8SHI+LfnHXeT46IP5/s77eAQ1v8fg9hPDvxPZn52sw8kZmrmfn7mfkDZ98/Mx/IzCPAvwC+IyIeV7isTwJu2/D92yf/Hthmm68F/r+zbvt54KWZeWQy2Hk18KU7HPsJwHsBIuKiiPiNiHhtRFwyGei8E/ia8/w9JM2QgxqpQUR8AXAT8KWTWYqvAT6Smd8G/C2fmZX4mS128U3A04HPB54J/CHwI8DDGb8+Nw5EPgj8L8BDgB8Hfj0iHrnZsSZhn98H3g1cATwN+L6I+JrJeR8AXg+8Engo8NuTc9nMUxgPeF63m2uTme8A7p6c8249GXjH5FwvA36C8WDiw9ts83jgA6e/iYgHM55J+pUN9xkAm83AnL2f90TE1cCfTvb5TZl5fPLz9wNPPO/fRNLMOKiRzt/rJ+tFTn99F7AOHAS+KCIWM/MjmfnBXezz/8rMT2TmPcD/AN6eme+azAi8jvGHOwCZ+duZ+bHMHGXmbwF/DVy7xX6/FHh4Zr4kM1cy80PAS4FnT37+ZcAi8EuTWZff4cyZkY0eBnwyM9d28Xud9jHGg6aNNl7H12+x3ZOAH4iI+xgPZpLxoG27GhSXAcc2fP80xr/jHaePB7wK+JsdzvkJjNfUvAl4cWb++FnHPTY5lqQ9ZqHrE5B65Bsy841n3zhZwPpi4Isj4lbg+VusE9nMJzb8/9Qm31+y4TjfDjwfeMzkpkuAy7fY76OBR00+yE8bMh44ATwKuOesD+utPuw/BVweEQuFgc0VwH1n3bbpdTwtIg4CXwhcnZl37+JYnwYu3fD9Y4AjmfmsDfu+Bfh/tzl2AI8D/gHwC5n5e5vc7VLg/k1ul9QxZ2qkRpn5G5n5FYwHEsl4USqT/18QEfFoxjMtNwEPy8zLGK/7iC2O9VHgw5l52YavSzPzdA+ljwNXTD7ET9sqK+itwDLwDbs85y9lPKj5k91sx3hQcWKrAU1EHJhkMR2KiGdFxC9PfnQH4zDeaQeBkxu2uxq4hvFi361cPfn3OuD7I+KaTe7zhYzDepL2GAc1UoOI+IKI+OrJ7MIS49mV0eTHn2D8F/+FcDHjgcvfT477nYw//E87+1jvAI5NFjE/KCKGEfG4yUADxgOVNeDfRMRiRHwjW4SyMvMB4EXAzRHxDZPFs4sR8bURcc5aoYh4cER8HeNFub+eme/Z5e/6ZOB9W/0wM1cYZ1X9J+DbgH89+dEtwFduuOttwFdOFkxfBfwG8MLMvG9ynq+IiFectfsnAHdMzvlG4HUR8cgNv9sh4B8Db9jl7yRpBhzUSOfvdHbR6a/XMZ4N+Cngk8DfAZ8N/PDk/v8B+NHJeo5/13LgzLyTcSbPWxkPYB7PeBHraWccKzPXga9jvDblw5Pz+xXGi4xPDwy+kXEa9H2MM5Veu83xf55x6OtHGQ+sPsp41mjjmpjfj4hjk5+9EPgF4DsLv+6TmGQfbePPgOcB/8eGkNivAc+IiAdNvv9j4A+Av2I8W/TKzHzphn1cxZnXEMbX9Q6AzHw9cJjxGqDTmWHPBN68i/CipBmy95OkXolx5+7fAv47sJKZP7fhZz8J3JuZv7TDPg4wDiE9ITNXd3HstwPPy8ydBl2SOuCgRlJvTMJ8r2Nc2O8vGc/YPDUzT3R6YpL2BAc1kiRpLuyZNTWThYzviog/6PpcJEnSdEXEyyPi3ojYNJwbY/8pIu6aVF7/kp32uWcGNcD3snkDOUmSNH9eAVy/zc+/Fnjs5OtG4Je3uS+wRwY1MW46979xZjlzSZI0pzLzLZxbnHOjG4Bfy7G3AZdtLLGwmT0xqGHc/fcH+Ux9D0mStL9dwbhExGl3T27bUudtEiZFuu7NzHdGxFO3ud+NjKefuPii+Mf/6PO2a9YrSdL8eOcdy5/MzIfP6nhf81UX56fuW2/axzvvWH4fZzaQPZyZh5t2uoPOBzXAlwNfHxHPYNwJ+MER8euZ+a0b7zS5EIcBrnnioXzHrVtVdJckab4MH/nXOzVivaA+ed86b7/1yqZ9LD7yg0uZuVmrkfN1D+MimaddObltS52HnzLzhzPzysx8DOMOwn989oBGkiTtO0eAb59kQX0Z8EBmfny7DfbCTI0kSdpTkvWc7jLXiPhN4KnA5RFxN/BjwCJAZv4Xxv3cngHcxbg57Y5tV/bUoCYz3wy8uePTkCRpX0tgxHSL82bmc3b4eQLfs5t97qlBjSRJ2htGPUxI7nxNjSRJ0oXgTI0kSTpDkqz3sDekgxpJknSOaa+pmQYHNZIk6QwJrDuokSRJ86CPMzUuFJYkSXPBmRpJknSGBBcKS5Kk+dC/KjUOaiRJ0lmSdKGwJEmaAwnr/RvTuFBYkiTNB2dqJEnSGcYNLfvHQY0kSTpLsE50fRK75qBGkiSdIYFRD9fU9HJQM54Wm93E2MClR5Ik7Xm9HNRIkqTpMvwkSZJ6b9zQ0kHNXJplqOs0Q16StDd18ZnQhVE6qJEkST3X15kapwMkSdJc6OlMTbKeu5/+G4ZjOEnS2H4JI1UkwXoP5z16OqiRJEnT5JoaSZLUe31dU7OvBjWVkBV0E7aqTouaNSVJe1P1M6gbwXr27/Okf2csSZK0iX01UyNJknY2bkfUv3kPBzWSJOkcrqmZkfEIcvftQwfFB6glDmoauSSpbzJdUyNJktSZXs7USJKk6RoZfpqdWspzBxNTM04jb6mQaTq4pL5pec+bdYp1ZdlEV8Z1avr3mdD5oCYiDgFvAQ4yPp/fycwf6/asJEnaz/q5pqbzQQ2wDHx1Zh6PiEXgTyLiDzPzbV2fmCRJ+5Ep3UWZmcDxybeLk69t5+iSZD0L03jRnyq9fap+LEn7STWMZAPN6et8UAMQEUPgncDnATdn5ts3uc+NwI0AV10xnO0JSpK0z6z3sKHlnvizPjPXM/NJwJXAtRHxuE3uczgzr8nMax72sD1x2pIkzaUkWGfQ9NWFPTFTc1pm3h8RbwKuB97b9flIkrRfjVwovHsR8XBgdTKgeRDwdOCnd9quFNOsZtMV1+KM1Z4U1erHLbqI95pGLgn2x3qT0lrQjpjSXfdI4Fcn62oGwGsy8w86PidJktQznQ9qMvMO4Mldn4ckSRpLopcLhTsf1FSMp8VmOI3Xcqhy6Ko27WdemKT9pKUq8KxTs/tUURisUyNJkuZAJr2sKNy/M5YkSdpEL2dqElgtTDmOohYfbMpEqo50i2GrLmKgVjGWBP3LYKqebzWLaabLJpqFXbolSVL/Jf0MPzmokSRJ57BOzQyVJg2rhY8aZuDKq92Lmw2iPr3ZRcE/SYJ6FlNLRlE1jFQ9ZmXZRFeSYNTDlO7+DcMkSZI20duZGkmSND2GnyRJUu8lNrScmVEmS4VY6GIxPDjqoKLwYnGE3JJSOSjWI26p6FldxmMjTGl6ukjNnnV135ZjVlOzV3uW0r3ew3WWvRzUSJKk6enrTE3/zliSJGkTvZypSYLVQqrZevF4iw1ThtVRYzU1e7F4PKhPxZoKLgkaQ9HlYzakdM84NXu1T9EnMPwkSZL6LzMMP0mSpPmwnoOmr51ExPUR8YGIuCsiXrDJzz83It4UEe+KiDsi4hk77bOXMzUJrBTGY8MOKgoPyxWFZ5s1BZQztZrGxsXfc2ATTWlHfWowWT3XlorC5TBSOWzVv3DOtETEELgZeDpwN3BbRBzJzDs33O1Hgddk5i9HxBcBtwCP2W6/vRzUSJKk6UmYdpfua4G7MvNDABHxauAGYOOgJoEHT/7/EOBjO+3UQY0kSTpLTLtL9xXARzd8fzfwT866z4uBP4qIfw1cDFy30057OagZESzl7k99WJzeXG+awq1tOywOkFdbzrX4BB50MKNanaq2aJ+0sz41l6xmMEE9jLRePORS1gqcdmFcp6b5zf3yiLh9w/eHM/PwLrZ/DvCKzPz5iHgK8MqIeFzm1k/QXg5qJEnSdF2A3k+fzMxrtvjZPcBVG76/cnLbRs8DrgfIzLdGxCHgcuDerQ7on62SJGnWbgMeGxFXR8QB4NnAkbPu87fA0wAi4guBQ8Dfb7dTZ2okSdIZkrgQ4aet95+5FhE3AbcCQ+Dlmfm+iHgJcHtmHgG+H3hpRPxbxhGx52ZuH6fs5aAmCZZGhTU1xZTlljU1Q1ZL262U14w0nGvUnsBtaaO1ycL+RKal/aOL5pLVtOzxttXtau+VlbWgXRpNOZiTmbcwTtPeeNuLNvz/TuDLd7PPfl1hSZI0dZmw3sO6Og5qJEnSOaYZfpqWXg5qRhmczAO73m5QTBk8FLUQEsBwUEzpLk7FHmh4DlbTwQcNT/xZp4O3TI+bDq6+mXVqdifNJRtSupeKZSxWi9udKHxuaXd6OaiRJEnTM14o3L8/5BzUSJKkc6xPt03CVPRyUDOuKDy7abyW/hfDamPKYqbWSq6XtoN60aJR1Kd/6+Eg85+0P/SrKWUHzSXrhyyHkapZTEu5WNquCxeoovDM9W9uSZIkaRO9nKmRJEnT5Jqakoi4Cvg14BGMZ7wOZ+Z/3G6bdQYcXT+062MNq2GShsd1sRgOWixOxR5oKr5Xuz7VRqFQz5yqFkQcRv3BtImm+hQKgno4qPp7dtFcsloID+phpGoW04nRwdJ2XWlZetGVzgc1wBrw/Zn55xFxKfDOiHjDpJKgJEmaMYvvFWXmx4GPT/5/LCLeD1wBOKiRJKkjfQw/7akzjojHAE8G3r7Jz26MiNsj4vZj99WL4UmSpPnU+UzNaRFxCfC7wPdl5tGzf56Zh4HDAI9+3KV5bPSgXR+jZe1HVbXBZLWi8KHhWmk7aFj/0xBHr6aDl1NHG5rftazHkaqqVYGhYW1MsTJwF80ll7Je3qGaYl1dG3NsffefW12ZdpfuadkTg5qIWGQ8oHlVZr626/ORJGm/c6FwQUQE8DLg/Zn5C12fjyRJ+11fi+91PqgBvhz4NuA9EfEXk9t+JDNv2WqDdQYcL6R098mBqIWRqimKUE8HX22oKNwSupq1chig4X3BdHB1oZ4KXldNzW55z6umZp8shp/6VFG4rzof1GTmn9D0ti9Jki60PmY/dT6okSRJe0y6UHhmRjng2AzDT4Nic8kWh2KltN1Fxe0ADgyq1Y/rGVf17KfqY1L/y2PghKIaVMOXLU0iy1lMxdfXSvF4MPvqvlDPRjo2qn3+PNCr7CcXCkuSpDnRx5ma/gXMJEmSNtHLmZr1DI6u7X76r5NFT8UrfNFoubTdxYPadgAHisX3DhW3A1gtTldXG+cNGv7wqIYB6qXBZt9Acb9kW9WL0s0+FN3yHKiGkaphq+WGv+xXiq+UloyiWRfRO1743OqKKd2SJGluOKiRJEm9Z5sESZI0N8x+mpG1HHLfysW73m6hmLLchUNR60TesqbmULmKcT2l+2Bx29mngu+f9Saztl8ek5bU7Fkfc6m4pqa6Lgbg5KiWml1dFwP11OxqSZGja/1J6e6rXg5qJEnSFKVraiRJ0hww+2mG1jM4urr76b+FQbGiZ0MqeLUa8UXDWhjp2Kg+vXloUAt5XUS9ivGIWvipmqo6aHmRFh/L9YZjDmO24ZVZp5C3mvX5toSQ6mnkXVQU7k9zyWp6NVBujHy0mtK9Xq9+3IU+Dmr6E5CWJEnaRi9naiRJ0vSY0j1D6zng6MruV7wfGM4++6kafnrQsBYKesjwVGk7gIuKWQSXxlL5mEvFLK9DxWn1puakxTDkoJipNT7mbMMrsw53daXeXLKLisL1589qcdvV4nN9adQQfqpW9y1mMEFDFlNh+cN4u35lP6WDGkmSNA+sUyNJknove5rSvT/mmiVJ0tzr5UxNZnBqdfepcUtrtfjy2qghpbsY075kWEuTPr5Qjy9XqxGfHNTTFC/O2u+5krW1OIsN61uqaxvaKua29Pjuh5bO131aA1RPr27o0l18upfX1DR0zD5ZXFPzwNpF5WPeX1zj8kBxu8pa0C65pkaSJM0Bs58kSdKccKZmRkYZnFzZ/TTnYDD7hnLVKsaHhrXpzcsW6+Gni4rhp0sH9SnVSwe1dPDVaiXihuqs1dTsQQfHLJtxCnmremp27bq2VPethpGajln8UDpRDCOdzIbmksX06uPrDcdcKx6zsPwB4HiPwk99bZPQn4C0JEnSNno5UyNJkqYox2ndfdPLQc1oFCwVwk9RnMofjepTcNXwwaFhLbxy/2o9E6Da0PIhw5PlY1azJZaKWVOHol5VerUY6hhG/fkz6wq2XWRbtVTMHRSLg1Wva1tDy2r2U/2Y1QaT1ddltSowwAPrtfeulve8oyu18NMDxe1OLPesoaXF9yRJUt8l/Vwo7JoaSZI0F3o5U5MJa8uFafLioDMbwk/VkNeBavhpsd4w7UHFgn/HFurHPDGqha4uG9Qady41FXqrbbfeEl4pBrWrIa+W8EoXZp3F1PJYlovvNTwksw4/nRzVwyvV7Kdjaw3ZT8VspOPLte2WluvFCWfPOjWSJGlOuFBYkiTNhT6uqXFQI0mSzpDpoGZ2RsH60u5PPYoVhVfXG9bUFI95arEWmz662tDQcqFWUfh4MRYOcHJYi02fKMbuF6O+pmaRWjp4SxPNqnIV2obrsx9U0/qhnppdrQoMs0/NPrZeX193tFjd9/6V+jGra2pOFdfGrC738yO3T/ZE9lNEvDwi7o2I93Z9LpIkadwmoeWrC3tiUAO8Ari+65OQJEljmW1fXdgTc2GZ+ZaIeMz5bxDEUmE8Vhw4jg7UR5yrxfDTiYVi+Gmxnt54SbEZZjUVE+rT3Cezdn0uLlYiBlgthmbammgWqxiX6xfUNuubanPJtuq+tcdkKetVnk8UXycni6/LBxrCTw+s1rathpAATi7Vrs9qYfkDwKhSiqRDrqmZooi4EbgRYPjQyzo+G0mS5lcSvRzU7JXw044y83BmXpOZ1wwvuaTr05EkSXtMb2ZqzjCCwfLswk/ZkP00Kg4blxZr0+PHG8JPR4vhp/sX6w3ljhWPWQ1bXRT18FM1c6qaNTU+5j6JBxVVK/zWKxGXNgNgKWtvBtWqwFB/ndxfbC75wFo9/HR8dbbVfQFWillMlexboLZsokN9fPfp56BGkiRNT0/r1OyJYWNE/CbwVuALIuLuiHhe1+ckSdK+lo1fHdgTMzWZ+Zzd3D8ShkuFEWQ5+6m2HcD6sDZuXF2oPTQnikX7AI4dqE3jHj1Yz356YK02zX3ZsNYI89LBUmk7gKWsNRk91FDQrpo5Nai+o8Ts/85paaJZLTJYfURWq28iwGox/FTNYIJ6Mbzq6/JYsYAewAPLtW2rhfAA1ovZSNUw0rBn4adpi4jrgf8IDIFfycyf2uQ+3wy8mPEw6d2Z+S3b7XNPDGokSdLeMs3wU0QMgZuBpwN3A7dFxJHMvHPDfR4L/DDw5Zn56Yj47J3267BRkiSdY8rF964F7srMD2XmCvBq4Iaz7vNdwM2Z+enx+eS9O+3UQY0kSTpDMp6pafkCLo+I2zd83bjhEFcAH93w/d2T2zb6fODzI+JPI+Jtk3DVtvoZfhrBsJ6Zu2sxqk/B5aA2bhwt1GK9ywfq8eVjB4trahqaaB4txvzrlYjr6Z8HRrXU7ENRW4sD9XTwA1F7zrY0bOxC9WxXimtxmqr7FhfnVZtSApwsHvP4erGicENzyRMrs63uC8Cp2uNZKilCcS1oVxJoDz99MjOvadh+AXgs8FTgSuAtEfH4zLx/qw2cqZEkSbN2D3DVhu+vnNy20d3AkcxczcwPA3/FeJCzJQc1kiTpHFNeU3Mb8NiIuDoiDgDPBo6cdZ/XM56lISIuZxyO+tB2O+1l+CkShqcqGxaPVy8ISxbDAFlMBV8/WH9Iq83djh5qqGJcDF09sFib5q6mggMcKlYjbqkIW04HL4ZXhj2aHYd6hd9qanY1LRvqYaRqWjbAsVHt9XW0mJp9tFgVGODEcu39Z7Rcf33NOow0WC5t1p0p1prJzLWIuAm4lXFK98sz830R8RLg9sw8MvnZP42IO4F14Acy81Pb7beXgxpJkjRN029omZm3ALecdduLNvw/gedPvs6LgxpJknSuHjZ/6uWgJkYwrEzjdRF+Ks5WjxZrJ7tenE4FWCmGrk6s1Kecj60Vm9gVp8ePDeuZWovFLKaW7KfytsWw1agYturKavEvyfXim0FLKLGaeVfN9IN6ZeAHVmshr+MN7wXL1crADVV6B8UwUunzB1ioLJvQrvRyUCNJkqaopw0tHdRIkqRz9WviFujroCZhodCXsBoK6iT8tFAbIY8O1qdiq5lT1awFgKMrxSyLA7Xp8WOL9fDTodFqabulQX3OuR7uqIWthj17F6uGkapZTC2F8KphpGoGEzRkMRVfl8cb3guqzSWrGUxQDyOVt5th0dgLo38zNdapkSRJc6GfMzWSJGm6+jVxCziokSRJm3FQMxvjlO7dX+36mpqGhpbF/nejYsx2faV+ruvF1MilaiomcGK1WMW4uFbg+HrDmpqoram5bFhfZ7A0qgXhh4NiSnfDu1h1fUvLOp4Vai+w6pqalvTqY8XnXstztpqafaxYGfhUw5oaimtqWppEDgtrM6FhTU3xeJ24MA0tZ66XgxpJkjRdPStbBbhQWJIkzYleztREwsLS7oeQo2K3vsFafbiaxWOOio/M8FR9unB0oDbGXStOGwMcK1YgLTfCXKs3Bzw4qIWfjjaEDy6K2jz3MGvhp5ZQUCfhp2J8t5qafbIl/FRMzW55zh4vhpGOL9e2W12qf6REMfw9aEiTLoefymGrnk199Ox0oaeDGkmSNGWuqZEkSfMgnKmZkREMC+GnQTUTqeEqVbetVhSursoHGBWzCFYbwk+niplTR4thq2MHG0JBxXnuY8N6+ODBWZznrkWfONBQPruL8FO14nI1i6mL5pL3r9a2g/rrpPq6HC3V3wsWq80lGzKKFqqNKQufP+Pj9WiUkPQy/ORCYUmSNBf6OVMjSZKmKFxTMysxSoZLu59fz2JIp6X4XrkxZTX7qWEqdnSgdq5rDVPOy8Vp7hPFafUHVuqhoAcVs58eslBvaHliWAx3FOdgV4uNMLuyNCo+f3rUXPLYWj3kVX2drKzU3oCipblkB+GnajZSdbuFUz2L5/TsdKGngxpJkjRlPRzUuKZGkiTNBWdqJEnSuXo4U9PLQU0kDJd3v6ZmNCquqVlraWhZ27baCHP9YMP6n3KTtoYmmsvFlNxi47zjD6o33KuuiWipCPvQhYYGgQWLUV8f1YXqmppjo9pjUk3Lhvrz5+hKfR1P9XVSfV0Oi1WBoaFJZEtF4eoxi2tqBivFWgtdsKGlJEmaF30svrcn1tRExPUR8YGIuCsiXtD1+UiStO9l41cHOp+piYghcDPwdOBu4LaIOJKZd265USaD1d1XPq2mZuewPvarppGPFmrHbAkFFWfy28JPxenq5WLK6dHlhorCC7WU7pOjegjp/vVauGM9a9f1UDFtveWYLaqNKY8Xm4weX6+nVx8vNmE9VmxKCfXUbE7VwpDDk/X3gsGMq/u2bFutDFwpRaLd2fEZHxGPAb4H+IfAfcBfAL+fmX9zgc7hWuCuzPzQ5HivBm4Ath7USJIkneV8/rT6PeAv+cxsyhOBt0TEzRFR/xPiM64APrrh+7snt0mSpI5Etn114XzmJoeZ+TKAiLgvM78rIhaAfwscBr5jmid4WkTcCNwIcOjAQxisFBrvRTEUtFh/dIbL1ZDX7M91tDj7JprVCqSrS7Vp9eOHGrKfFotVaA/UQ17HF+rbVqw3LLNbrabsNVguxkyrGWnVDCaAo8Uw0omV+nO2+jqJldm/F3ST/VQLBw2Xao1fK8smOtXD7KfzeQd7Y0TcNPl/AmTmWmb+LPCUC3AO9wBXbfj+ysltZ8jMw5l5TWZes7h48QU4rCRJ2lTrIuE9PFPzfOCHI+J24FGTGZOTjAc0n7oA53Ab8NiIuJrxYObZwLdcgP1KkqR9ZMdBTWaOgJ+IiF8ErgOeBHwW8F7gha0nkJlrk5mgW4Eh8PLMfN9220QmUZnGK4afBuv1IWcOitO41e2KISSAYbFwXzVrAWB4qhh+uqgW6jhVLEYGcPzA7JtoXrpQ69ZXzUQaNWT6VQvhDaOeEXJyvVjLNRqsAAAZzElEQVSEsZjFdKwh/HS82Fzy5FJDAcZis9nq67KtuWRtu5YmkdXwU7WIXmnZRJd6WKfmvAOumXkSODL5uqAy8xbglgu9X0mSVNPH4nud16mRJEl7UA8HNXuiorAkSVKrfs7UjJJYWdv9dsU1NQzrqaqDldq4MYfFVMNiKibAQrEy8EIxbR1grRhHXyuuFVhbrj+WJw8Wq9eu1ddEVNdwjIqpmNXtoJuU7lPFNTVdNJc8tVp7/iwvFUt9A1Gs2N1FenW5Sm9xXQw0rI1ZLnz+ALHsmppp6+egRpIkTU2XBfRaOKiRJEnn6mHxvX4OajKJyvTfQnF6POvD1cFqcdlSMaV7dKC+TGp4oDr925BGXpzmHhTDbOvFCqsAS4dqYYBjxVRegGOLtXDHegdvRsuj2rUdNPw5eKzYJPL4WjU9vx5+qqZmjxpCptXX5vBU9Xi17cbbdhF+qoWDSiVFoLZsoks9nKlxobAkSZoL/ZypkSRJU+WamllJYK0w/VcMI8V6w4TWoDbdWD3isCH8tHCg2sSuIfxUzLiqbjc62NCwsRi6OtnQRPP+YvipJYupaqUafmqY466GkY4WQ4ItzSVXlotvt8VMP2h4fVWr+y7VH8vqtoOGjKJy+KkYRorK51aXHNRIkqTe62n2k2tqJEnSXHCmRpIknauHMzX9HNRkwupqYbtibDrqE1pRTM2uGqzUH9JBMTVyuNyQRl5NOS3G/NcbznVUXBPR0mX5+GI9HbxiRP35urJee321pHQfX61d2xMddMwerRSvT9Pra9bbNXTMXqq+/9TXqZTXxlS7ba+7pmba+jmokSRJU+WaGkmSpI70dKYmYVSYqqwWc4x6xUqGtXFjFJtvtqQMVsNILWmc5Wnupdp21bR1gJVDteuzcrD+Mjt+oB7uqFgb1f/OWS9uOxzUX18ni+GnY0u18NPKcr25JKdq4afhqZaSCbXtqq/plvBTtblkS5XeamXgchhppbBsQrvS00GNJEmaqh6GnxzUSJKkM/W0Tk0/BzVJbfqvErICGDQsPVqbbQOzaKmuWQyVDRsqni4sFY9ZzJpab2i4Nyie6+hB9euztFILdxSjl6wWM5ig3vd1OKi/c54sXp+VYpbgekNzyWoWU1uTyOJ2K8XtihlM422LmUjVEBLAavH9ufq+Xv0M6oqDGkmSNBd6OKgx+0mSJM2Ffs7UZJKF6b+oFtFryM4oh66qzTfX6tkZg5Xi9HgxawFgWMyyWCg26ltvqGW3UMxCWWkIP60s1h7PLDa0XF+sh0vXR7VjLgzrz59TxWJ4q6dm31xyUG7eWj4kC9Xsp1PF7KfV+mNZbi7ZFH4qZiNVs0xnvByhReCaGkmSNC8c1EiSpN7rafaTa2okSdLMRcT1EfGBiLgrIl6wzf2+KSIyIq7ZaZ/9nalZ333sNqvDzuL6BICoxlCLzTfLjdYABrVzHRQbPQIsFFNk640w649lNR28pSHh+oHa9anWLc3iuhiop3SvNvR8XaumWBe3i5WG6r4zbi4JMFwpro0prpMbtDSXrK6NqaZlw+zXxhQ+tzo1xZmaiBgCNwNPB+4GbouII5l551n3uxT4XuDt57NfZ2okSdK5svFre9cCd2XmhzJzBXg1cMMm9/v3wE8D57Xs3UGNJEk6R2TbF3B5RNy+4evGDbu/Avjohu/vntz2meNHfAlwVWb+t/M9536GnzJrFYWHxanqlpTB6jHLc/kNU7ELtXOtpmICDJZr07HV5pst4adqam1LQ8JR8TGpTnKvNoSfylPVLYcsliGoVocenqr/HbhwqrZdS/ipmppdfV22vBeUG1O2pEkXt62UFAFg1LOVt+2n+8nM3HEdzGZiXIPlF4Dn7mY7Z2okSdKs3QNcteH7Kye3nXYp8DjgzRHxEeDLgCM7LRbu50yNJEmanvNbF9PiNuCxEXE148HMs4Fv+Z+Hz3wAuPz09xHxZuDfZebt2+20t4OaLDQGa5jlLstiOCiqYauGqdhYLWaErNafRsOl2nT1QrGy6/rB+qt0/WA146p8SEYHipWBi9WzW7KfOkguJIphyEE5e6602XjbGVf3hYbsp2L4qaWhbj0TqeGY1Wyk6nY9a2g5zTo1mbkWETcBtwJD4OWZ+b6IeAlwe2Yeqey300FNRPxz4MXAFwLX7jQCkyRJMzLlJUCZeQtwy1m3vWiL+z71fPbZ9Zqa9wLfCLyl4/OQJEkbXIDsp5nrdKYmM98PELHbqeAsZQdVQlatopjFlJXsLiCK2wHl5m6xUn8aDVZqDRvr0+oNU/nVkEVDQ8LRQu2YWQ4/zf6dqCUsXA0jVZuTDosZTNBQfG+lfsxZZzGVM5igoblk/ZjV99lS9i3dfAbtN71ZUzPJb78R4FBc3PHZSJI053qWgQ4zGNRExBuBz9nkRy/MzN873/1k5mHgMMBDhg/r4aWWJKknpp/9NBVTH9Rk5nXTPoYkSbpwgm4yhlv1Jvx0hqRYmbE47Bw0PLTVmO2u1xlNtFTXHNTWYbQ00RwUY/CD5eJanJaKwsX1OC3HHNV+TbL4nO1byH+4VPs9B8V1Ki3ro8rPn6X6g1JNzR4sF0tRdNFcsqVJ5KzXxlQrxeu8dZr9FBHPioi7gacA/y0ibu3yfCRJ0sR0G1pORdfZT68DXtflOUiSpHN1lZbdop/hJ4Cc4Tz5qD6hlVGszFkMBZVTFIGoTqm2HLPYLLSacrpQrEALsF5MyV1vSumubZfFXzPWW8r7FrdreOOspjuXq/u2NJcsbrvQVIagmJpdbeJbTcuGenPJhvef8rbVMFJLyY0uOKiRJElzoYeDmq4rCkuSJF0QvZypSWqVT6OaxdQS6ip266uuro+WTIBq5tRKfco5FmtPweq0+uhgsVEo9UybhQP1vx2qDS2z+GsOGpJXqiGvpvDTjKv0NlWkLmYxVasCQz0zsRx+amku2UVGUXXbYuXtLip2l3XY6qBFLwc1kiRpyhzUSJKkeeBMzR7XSbO+6rOiOpXflAlQLL7X1ESzWuSrmDVVDFtBPYy03hSyKIafis+fcgipYdtigiAAg2Lks5z91FJ8r9yEdfbFLauvy5bqjbNuLgmUw0gzzb7tUg8HNS4UliRJc2FfzdRIkqTzY/hJkiT1n126e6CLOGgxpbtaxbhawRga0sFbYtrFdPCoNsJcqT/lByu167NQXBcDsH6gmDpabWjZ8o7Qo4rC1Sq9fWouCcCsKwO3vBcU339aKgpXPxPK6zP7thanh4Ma19RIkqS5sL9maiRJ0o4C19TMTmZvpvGq05QxqFbXbGhIWJzGzWolYiCGtdK39UaY9XMdFJthDhsqCpebGRafBuuLs29o2ZLSXU+TLm5XPB7Uw5fl6r5AVFOzi5WBW94LyqGrhlIdsw4j9aqiMPQy/NTPQY0kSZqqaGlB0REHNZIk6UxmP81WX6bxqmGkLGY/BQ0VhaMWP+ikiWZxeryaNQX1zKmWhoTDakPLDkJB1WMOGpJXBmszDj81PJbVatYtz9mZZzE1vBdUm/i2NRzuxzIGnb/eDmokSdL0uFBYkiTNBwc1umDKq+vrmTblMNKwYQq3GkYsTqvHQi3bCuqZU4OV+mNSbWhZLUHVEkksZz+1hJ9WZxtGammIWn3+xErDBSqGacvbNTS0pEeLUvuy/KGVMzWSJGk+9HBQY0VhSZI0F5ypkSRJZ0rDT9pEvaJwQ2XXquo6nqYmdsU012rMv1phFYjV2stl2LAOIxeLk6nF508XKd0tb5zVdPlqKni1KjA0pGa3vL6Ka1yqr+mm94Lie+V+Wd/SiR5eWgc1kiTpDH3t/eSaGkmSNBf6O1Mzy0qQMfux38wbYdJQxbglFbM6XV2uRFxP6a6GD2Kh/vypNtGsajnXll6qVeXU7LXidg3VfcuNKVuaRFa3rb4uW94Lqu/pDZ8FMw9d9a2CcY/S7E/r76BGkiRNTR/DTw5qJEnSmWxoOcdapgw7CF3NXMMUbjVbolz9uJo1BeXMqVhtqWJcPN9yc9KGx7J6zIYp7mo2UjV7rq255Gybt0JDNlL1MWl5L+hTFlPfwkhFLdmQXdkHn7iSJGk/cKZGkiSdq0eTZ6d1OqiJiJ8FngmsAB8EvjMz7+/ynC646jRlMWzVMoVbzpxqSXspblsOWzVlktTCSLG2WD5kNWOmOgWbw4bJ22L4qSXDohqGjGJBxHIGEzRk7DU8Z4vXpx62aqneWG3i20HG1T7Rx4XCXYef3gA8LjOfAPwV8MMdn48kSUrGf3C0fHWg00FNZv5RZp7+M+RtwJVdno8kSRqLbPvqQtczNRv978AfbvXDiLgxIm6PiNtXWZ7haUmSpD6Y+pqaiHgj8Dmb/OiFmfl7k/u8EFgDXrXVfjLzMHAY4MHx0B5G+nZpxmtxWjSt46kO56tTm9VUcKin1rY00RwW16lUnz4LDSm5xXNtSSOnWBm4vLaqpSRAdduW52x1bYzNJQUuFN5MZl633c8j4rnA1wFPy+xhTWZJkuZMXxtadp39dD3wg8BXZubJLs9FkiRNdLjYt0XXdWr+M3AQeEOM00Hflpnf3e0p7V9dNNEsp4MXp9VbUpZj1s03gVgpnm/1Icl69WPWiuGnUcPzp7htLBcfk9XV2nZQfh6U06tbtu2guaR0IXQ6qMnMz+vy+JIkaXOGnyRJ0nxwUKPO9az5Zjlbohp2aMhEKl+floyZYbHybTEW3vIeVq0r3dLQspwZVH1MWh7L6rk2hJ/6lMVUPqYhr6lxpkaSJPVf0tR1vSt7qfieJElSmTM1kiTpXP2bqHFQo56qrk8Y1CcnO+kMXu18vVB7aTf0W6/rZE1N8THpJL26P12orSg8X1xTI0mS5kMPi++5pkaSJJ1j2l26I+L6iPhARNwVES/Y5OfPj4g7I+KOiPjvEfHonfbpTI0+ozhVnaOGKr3VasTFSsTZUL22XFG4pWJuSzpvX7SESGadJt2U0l2sgt0S8upTmrSp2ftKRAyBm4GnA3cDt0XEkcy8c8Pd3gVck5knI+JfAT8D/Ivt9utMjSRJOlNegK/tXQvclZkfyswV4NXADWecQuabNvSFfBtw5U47daZGkiSdYdyle6praq4APrrh+7uBf7LN/Z8H/OFOO3VQo14qN9+shitoyG5sqGIc1eynqr4tDKyG9qrNJVtCJNVr26fqvpov7RHByyPi9g3fH87Mw7vdSUR8K3AN8JU73ddBjSRJmoZPZuY1W/zsHuCqDd9fObntDBFxHfBC4Cszc3mnAzqokSRJ55hy+Ok24LERcTXjwcyzgW854/gRTwb+K3B9Zt57Pjt1UKNOlcNIM86aatq2JYNpfVjftqIlU6ta2LDlmMVtyxlFLeHL6u/ZFPKacUaRGUzz4/wW+9Z3n7kWETcBtwJD4OWZ+b6IeAlwe2YeAX4WuAT47Uko/m8z8+u326+DGkmSdJac+hq7zLwFuOWs21604f/X7XafDmokSdI5+tgmwTo1kiRpLjhTo3YtcfSY7bi6JVU1qstbGqZw600Qa49JDBvW8LSsjSma+dqYpurQ1YrdDc8fU7PVom8lHnBQI0mSzpYQPVz37aBGkiSdy5kaaTbKqeAtGdLFYyb1lO5yReGshfWawhWD4rm2HLP6plttLtkSfjLdWZo6BzWSJOlc/ZuocVAjSZLONeWKwlPhoEbdqk7JV7OmmqqzVsMr9QyvanZP+c2opYFmMeTVlFFU/D3LWVMNyqG9LsJWhsoErqmRJElzILkQXbpnzuJ7kiRpLjhTo32lKbunWDO83HwTyiGvpiydoi7i7+Xfc78UpTOMpKIgXVMjSZLmhIMaSZI0FxzUSJKk3uvpQmEHNdL5Kq5PyIaU7iiu4+kiBaCLNOny2pgepUnblFI6f50OaiLi3wM3MB4P3gs8NzM/1uU5SZKkfhbf6zql+2cz8wmZ+STgD4AXdXw+kiQJxmtqWr460OlMTWYe3fDtxfSy04Q6MetKxHTURLP6ezaEvHqlRyEdw0jql+4GJi06X1MTET8BfDvwAPBV29zvRuBGgENcNJuTkyRJvTH1P+ci4o0R8d5Nvm4AyMwXZuZVwKuAm7baT2YezsxrMvOaRQ5O+7QlSdq/EsNPm8nM687zrq8CbgF+bIqno/2uJeulgyaa1cypahXjlhBJDKrVj3s0xW2FXu0nPXy6d5399NjM/OvJtzcAf9nl+UiSpLE+Zj91vabmpyLiCxiPB/8G+O6Oz0eSJIELhXcrM7+py+NLs9AU0ilmTpnds4MeZU1JOn9dz9RIkqS9JullN3sHNZIk6SzWqZEkSfPCQY00xzqoYtyrFOIOUt6rOlkb06fHUoJeDmr2SS11SZI075ypkSRJZ3KhsKS+a6sobHhFmh/Zy5CpgxpJknQu19RIkiR1w5kaaQ/rUwXbPp2rpB24pkaSJM2NHoafHNRIkqRzOaiRdI6WDIKWwn3aWg+zOqTZ6mebBN8xJUnSXHCmRpIknSmBUf9mNB3USJKkc/Uw/OSgRtrLZt1Es2/rf1wbI02PgxpJktR/2cs6NS4UliRJc8GZGmkedRGWMRQkzY+E7OFr2kGNJEk6Vw/DTw5qJEnSuXq4UNg1NZIkaS44UyNJks6UafE9SZI0J3oYfnJQI0mSzpHO1EiSpP6zS7ckSVJnnKmRJElnSqxTI0mS5oQVhSVJUt8lkD2cqdkTa2oi4vsjIiPi8q7PRZKkfS9zPFPT8rWDiLg+Ij4QEXdFxAs2+fnBiPityc/fHhGP2WmfnQ9qIuIq4J8Cf9v1uUiSpOmLiCFwM/C1wBcBz4mILzrrbs8DPp2Znwf8IvDTO+2380EN4xP9QcazXZIkaQ/IUTZ97eBa4K7M/FBmrgCvBm446z43AL86+f/vAE+LiNhup50OaiLiBuCezHx3l+chSZLOMt3w0xXARzd8f/fktk3vk5lrwAPAw7bb6dQXCkfEG4HP2eRHLwR+hHHo6Xz2cyNw4+Tb5Tfm77z3wpzhXLoc+GTXJ7GHeX225/XZntdne16f7VWvz6Mv9Ils5xifvvWN+Tut61wPRcTtG74/nJmHG/e5rakPajLzus1uj4jHA1cD757MJl0J/HlEXJuZf7fJfg4Dhyfb3p6Z10zvrPvN67M9r8/2vD7b8/psz+uzvb5cn8y8fsqHuAe4asP3V05u2+w+d0fEAvAQ4FPb7bSz8FNmviczPzszH5OZj2E89fQlmw1oJEnSXLkNeGxEXB0RB4BnA0fOus8R4Dsm//9nwB9nbt+7wTo1kiRppjJzLSJuAm4FhsDLM/N9EfES4PbMPAK8DHhlRNwF3Md44LOtPTOomczWnK+pxuTmgNdne16f7Xl9tuf12Z7XZ3ten4nMvAW45azbXrTh/0vAP9/NPmOHmRxJkqRe2At1aiRJkpr1alCzU0llQUR8JCLeExF/cVYq3b4UES+PiHsj4r0bbntoRLwhIv568u9ndXmOXdri+rw4Iu6ZPIf+IiKe0eU5diUiroqIN0XEnRHxvoj43sntPn/Y9vr4/AEi4lBEvCMi3j25Pj8+uf3qScn/uyYtAA50fa7zpDfhp0lJ5b8Cns44U+o24DmZeWenJ7bHRMRHgGsy0zoRQET8r8Bx4Ncy83GT234GuC8zf2oyOP6szPyhLs+zK1tcnxcDxzPz57o8t65FxCOBR2bmn0fEpcA7gW8AnovPn+2uzzfj84dJ5duLM/N4RCwCfwJ8L/B84LWZ+eqI+C/AuzPzl7s813nSp5ma8ympLJ0hM9/CeNX8RhtLb/8q4zfifWmL6yMgMz+emX8++f8x4P2MK5z6/GHb6yMgx45Pvl2cfCXw1YxL/sM+fv5MS58GNedTUlnjF80fRcQ7J1WYda5HZObHJ///O+ARXZ7MHnVTRNwxCU/ty/DKRpPuwE8G3o7Pn3OcdX3A5w8wjjBExF8A9wJvAD4I3D8p+Q9+jl1wfRrU6Px8RWZ+CePOp98zCS9oC5NCTv2Iwc7OLwP/EHgS8HHg57s9nW5FxCXA7wLfl5lHN/7M58+m18fnz0RmrmfmkxhXy70W+Ecdn9Lc69Og5nxKKu97mXnP5N97gdcxfiHpTJ+YrAc4vS7g3o7PZ0/JzE9M3oxHwEvZx8+hyVqI3wVelZmvndzs82dis+vj8+dcmXk/8CbgKcBlk5L/4OfYBdenQc35lFTe1yLi4smCPSLiYsbNQm38ea6Npbe/A/i9Ds9lzzn9gT3xLPbpc2iy0PNlwPsz8xc2/MjnD1tfH58/YxHx8Ii4bPL/BzFOcnk/48HNP5vcbd8+f6alN9lPAJPUwF/iMyWVf6LjU9pTIuIfMJ6dgXG16N/Y79coIn4TeCrjzrifAH4MeD3wGuBzgb8Bvjkz9+Vi2S2uz1MZhw4S+AjwLzesIdk3IuIrgP8BvAcYTW7+EcbrRvb982eb6/McfP4QEU9gvBB4yHgC4TWZ+ZLJ+/SrgYcC7wK+NTOXuzvT+dKrQY0kSdJW+hR+kiRJ2pKDGkmSNBcc1EiSpLngoEaSJM0FBzWSJGkuOKiRJElzwUGNJEmaCw5qJO0oIp4WEa/s+jwkaTsOaiSdjycyrn4qSXuWgxpJ5+OJwLsi4mBEvCIifnLS+0eS9oyFne8iSTyBcTfqW4Ffycxf7/h8JOkc9n6StK2IWAQ+ybh547/MzLd2fEqStCnDT5J28oXAbcAasN7xuUjSlhzUSNrJE4E/A54N/D8R8YiOz0eSNuWgRtJOngi8NzP/Cvgh4DWTkJQk7SmuqZEkSXPBmRpJkjQXHNRIkqS54KBGkiTNBQc1kiRpLjiokSRJc8FBjSRJmgsOaiRJ0lxwUCNJkubC/w/gTfXvqnPU0AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 720x432 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "K = 32  # number of temporal samples\n",
    "N = 10000  # number of sample functions\n",
    "bins = 100  # number of bins for the histogram\n",
    "\n",
    "\n",
    "# draw sample functions from a random process\n",
    "np.random.seed(2)\n",
    "x = np.random.normal(size=(N, K))\n",
    "x += np.tile(np.cos(2*np.pi/K*np.arange(K)), [N, 1])\n",
    "\n",
    "# compute the histogram\n",
    "px = np.zeros((bins, K))\n",
    "for k in range(K):\n",
    "    px[:, k], edges = np.histogram(x[:, k], bins=bins, range=(-4,4), density=True)\n",
    "    \n",
    "# compute the CDF\n",
    "Px = np.cumsum(px, axis=0) * 8/bins\n",
    "\n",
    "# plot the PDF\n",
    "plt.figure(figsize=(10,6))\n",
    "plt.pcolor(np.arange(K), edges, px)\n",
    "plt.title(r'Estimated PDF $\\hat{p}_x(\\theta, k)$')\n",
    "plt.xlabel(r'$k$')\n",
    "plt.ylabel(r'$\\theta$')\n",
    "plt.colorbar()\n",
    "plt.autoscale(tight=True)\n",
    "\n",
    "# plot the CDF\n",
    "plt.figure(figsize=(10,6))\n",
    "plt.pcolor(np.arange(K), edges, Px, vmin=0, vmax=1)\n",
    "plt.title(r'Estimated CDF $\\hat{P}_x(\\theta, k)$')\n",
    "plt.xlabel(r'$k$')\n",
    "plt.ylabel(r'$\\theta$')\n",
    "plt.colorbar()\n",
    "plt.autoscale(tight=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise**\n",
    "\n",
    "* Change the number of sample functions `N` or/and the number of `bins` and rerun the examples. What changes? Why?\n",
    "\n",
    "Solution: In numerical simulations of random processes only a finite number of sample functions and temporal samples can be considered. This holds also for the number of intervals (bins) used for the histogram. As a result, numerical approximations of the CDF/PDF will be subject to statistical uncertainties that typically will become smaller if the number of sample functions `N` is increased."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Bivariate Probability Density Function\n",
    "\n",
    "The bivariate or joint PDF $p_{xy}(\\theta_x, \\theta_y, k_x, k_y)$ of two continuous-amplitude real-valued random signals $x[k]$ and $y[k]$ is defined as\n",
    "\n",
    "\\begin{equation}\n",
    "p_{xy}(\\theta_x, \\theta_y, k_x, k_y) := \\frac{\\partial^2}{\\partial \\theta_x \\partial \\theta_y} P_{xy}(\\theta_x, \\theta_y, k_x, k_y)\n",
    "\\end{equation}\n",
    "\n",
    "This constitutes essentially the generalization of the univariate PDF. The bivariate PDF quantifies the joint probability that $x[k]$ takes the value $\\theta_x$ and that $y[k]$ takes the value $\\theta_y$ for the entire ensemble of sample functions. Analogous to the univariate case, the bivariate CDF is given by integration\n",
    "\n",
    "\\begin{equation}\n",
    "P_{xy}(\\theta_x, \\theta_y, k_x, k_y) = \\int\\limits_{-\\infty}^{\\theta_x} \\int\\limits_{-\\infty}^{\\theta_y} p_{xy}(\\alpha, \\beta, k_x, k_y) \\, \\mathrm{d}\\alpha \\mathrm{d}\\beta\n",
    "\\end{equation}\n",
    "\n",
    "Due to the properties of the bivariate CDF and its definition the bivariate PDF shows the following properties\n",
    "\n",
    "\\begin{equation}\n",
    "p_{xy}(\\theta_x, \\theta_y, k_x, k_y) \\geq 0\n",
    "\\end{equation}\n",
    "\n",
    "and\n",
    "\n",
    "\\begin{equation}\n",
    "\\int_{-\\infty}^{\\infty} \\int_{-\\infty}^{\\infty} p_{xy}(\\theta_x, \\theta_y, k_x, k_y) \\, \\mathrm{d}\\theta_x \\mathrm{d}\\theta_y = 1\n",
    "\\end{equation}\n",
    "\n",
    "For the special case of one signal only, the bivariate PDF \n",
    "\n",
    "\\begin{equation}\n",
    "p_{xx}(\\theta_1, \\theta_2, k_1, k_2) := \\frac{\\partial^2}{\\partial \\theta_1 \\partial \\theta_2} P_{xx}(\\theta_1, \\theta_2, k_1, k_2)\n",
    "\\end{equation}\n",
    "\n",
    "describes the probability that the random signal $x[k]$ takes the values $\\theta_1$ at time instance $k_1$ and $\\theta_2$ at time instance $k_2$. Hence, $p_{xx}(\\theta_1, \\theta_2, k_1, k_2)$ provides insights into the temporal dependencies of the amplitudes of the random signal $x[k]$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbsphinx": "hidden"
   },
   "source": [
    "This notebook is provided as [Open Educational Resource](https://en.wikipedia.org/wiki/Open_educational_resources). Feel free to use the notebook for your own purposes. The text is licensed under [Creative Commons Attribution 4.0](https://creativecommons.org/licenses/by/4.0/), the code of the IPython examples under the [MIT license](https://opensource.org/licenses/MIT). Please attribute the work as follows: *Sascha Spors, Digital Signal Processing - Lecture notes featuring computational examples, 2016*."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}