Created
September 3, 2024 18:20
-
-
Save jtrive84/82376d21422fd98f15d81881669c8fbd to your computer and use it in GitHub Desktop.
Rejection sampling
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
{ | |
"cells": [ | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"### Rejection Sampling Demo\n", | |
"---\n", | |
"Author: James D. Triveri" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 8, | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"import matplotlib.pyplot as plt\n", | |
"import numpy as np\n", | |
"import pandas as pd\n", | |
"from scipy.stats import norm\n", | |
"\n", | |
"\n", | |
"# Create normal mixture model.\n", | |
"def p(x):\n", | |
" return norm.pdf(x, loc=30, scale=10) + norm.pdf(x, loc=80, scale=20) \n", | |
"\n", | |
"\n", | |
"# Proposal distribution.\n", | |
"def q(x):\n", | |
" return norm.pdf(x, loc=50, scale=30)\n" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"<br>\n", | |
"\n", | |
"Plot $p$ and $q$. " | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 2, | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAgMAAAGFCAYAAABg2vAPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABPbUlEQVR4nO3dd3ib5dn+8a+m94h3HNvZg2xICAmrzAQaKJQZyggrHdAB7cuvLS8d0ELbl9JCaSmlzJaZlhJW2BQIEAgQEjLJtJPYjuO9l6Tn98dj2c72kPRonJ/j8JHYlqVLsmWfusd12wzDMBAREZGYZbe6ABEREbGWwoCIiEiMUxgQERGJcQoDIiIiMU5hQEREJMYpDIiIiMQ4hQEREZEYpzAgIiIS4xQGREREYpzCgIiISIxTGBAREYlxCgMiIiIxTmFAREQkxikMiIiIxDin1QWIiIj05vV66ezstLqMiOByuXA4HIO+HoUBEREJC4ZhsHv3burq6qwuJaKkp6eTl5eHzWYb8HUoDIiISFjwB4GcnBwSExMH9cctFhiGQUtLC3v27AFg6NChA74uhQEREbGc1+vtDgKZmZlWlxMxEhISANizZw85OTkDnjLQAkIREbGcf41AYmKixZVEHv9jNph1FgoDIiISNjQ10H+BeMwUBkRERGKcwoCIiEiMUxgQERGJcQoDIiIiIVBdXU1OTg7FxcV9/poFCxZw1113Ba+oLtpaKBIgDW2d3PSv1XR4fJx6RC5zJ+WSkxJvdVkiEiZuv/12zjnnHEaMGNHnr7nllls48cQTufbaa0lLSwtabRoZEAmAdo+Xb/7jU15bV8F/v6zkliVrOfX377K5otHq0kQkDLS0tPDQQw9xzTXX9OvrJk+ezOjRo3n88ceDVJlJYUBkkHw+gx8+s5qPttWQ5HbwvVPGMDo7icZ2D798cR2GYVhdoogE2e7du7HZbNxzzz0ceeSRxMfHM2nSJN5//30Ali5dSlxcHLNnz97r65566ikSEhIoLy/v/thVV13F1KlTqa+vB+Dss8/m6aefDmr9CgMig/Tvlbt4eU05LoeNB66YyY/mjueRK2fhdtr5YEs1r67dbXWJIhHJMAxaOjyWvPU3xK9atQqAhx9+mLvvvptVq1ZRVFTEpZdeis/nY9myZcyYMWO/r1uwYAHjxo3jjjvuAOAXv/gFb775Jq+88kr3tMCsWbNYsWIF7e3tg3tAD0FrBkQG6fGPSgC44bRxHDcmC4CizES+feIo/vT2Fn798gZOGp9DgnvwJ4uJxJLWTi8Tf/6aJbe9/rZ5JLr7/idy9erVuFwunn/++e41Ab/+9a+ZOXMmpaWllJSUkJ+fv9/X2Ww2br/9di644ALy8vK49957WbZsGcOGDeu+TH5+Ph0dHezevZvhw4cP+r4diEYGRAZhbWk9X+yqx+2ws+Dowr0+952TxjAsPYHSulaeXLHDogpFJBRWrVrFeeedt9fiwNTU1O7/t7a2Eh9/4AXFZ511FhMnTuS2227jueeeY9KkSXt93n/+QEtLS+AL76KRAZFBeKrrj/y8yXlkJsft9bkEt4NvnjiKX7ywjhdXl3HN8SOtKFEkYiW4HKy/bZ5lt90fq1atYuHChXt9bPny5WRlZTFs2DCysrKora094Ne++uqrbNy4Ea/XS25u7n6fr6mpASA7O7tfNfWHRgZEBqi53cPzq8oAuGRW4QEvc+bkPGw2WLWzjtK61lCWJxLxbDYbiW6nJW/96fff2trK5s2b8Xq93R/z+XzcfffdLFy4ELvdzpFHHsn69ev3+9qVK1dy0UUX8dBDD3Hqqafys5/9bL/LrF27loKCArKysgb2QPaBwoDIAL24uoymdg8js5KYM+rAR67mpMZz9PAMAF5ZU37Ay4hIZFuzZg02m43HH3+c5cuXs2HDBi6++GLq6uq45ZZbAJg3bx7r1q3ba3SguLiY+fPnc/PNN3PJJZdw22238eyzz7Jy5cq9rn/ZsmXMnTs3qPdBYUBkgF7u+uN+0czCQ76K+OqUPACWKgyIRKVVq1YxYcIEbr75Zs4//3xmzpyJ1+vl3XffJT09HYApU6Zw1FFHsXjxYsAc+j/jjDM455xz+MlPfgLAMcccw5lnnsnNN9/cfd1tbW0sWbKERYsWBfU+aM2AyAC0e7x8UmzO450yIeeQlz1zylBufWk9K3fUUVbXSn56QihKFJEQWb16NVOmTOHyyy/n8ssvP+jlfv7zn3PTTTexaNEiMjIy2Lhx436Xefnll/d6/5FHHmHWrFn79ScINI0MiAzA5zvqaOv0kZUcx7jc5ENeNjc1npnDhwDwinoOiESdVatWMXXq1MNebv78+Xzzm9+ktLS0z9ftcrm49957B1NenygMiAzAh1urATh2dGafFhp9dcpQAN5cXxHUukQktAzDYM2aNX0KAwA33HADhYUHXnB8INdeey3jx48faHl9pmkCkQH4cEsVYIaBvji+qxnR5ztr6fD4cDuVw0Wigc1mo6GhweoyBk2/kUT6qbndw6qddQDdHQcPZ3R2MumJLto6fawrqw9idSIi/acwINJPK7bX4PEZFGYkUJiR2KevsdttzOzaYuhfeCgiEi4UBkT66cOt5hTBcaP71wDk6BHmIsJPig/chUxExCoKAyL95F88OKeP6wX8jh5pjgx8WlyDz6djjUUkfCgMiPRDa4eXDeXmYqFjRvYvDEzOTyPeZae2pZNtVU3BKE9EZEAUBkT6YcPuBnwGZCXHkZsad/gv6MXttDO9MB2AFds1VSAi4UNhQKQf1paaOwGmDEvt10EmfrNG9EwViIiEC4UBkX7wh4HJw9IG9PUzu8LACoUBEQkjCgMi/bC21FwvMNAwcGRROgC7alupae4IVFkiIoOiMCDSR22dXjZVNAIDDwMp8S6GZ5q9CfwLEUVErKYwINJHmyoa8fgMhiS6yE+LH/D1TByaCsD6MoUBkVhSXV1NTk4OxcXFff6aBQsWcNdddwWvqC4KAyJ9tKbXeoGBLB70m5RvhgG1JRaJLbfffjvnnHMOI0aM6PPX3HLLLdx+++3U1wf394XCgEgfDXa9gN/ErjCwXtMEIjGjpaWFhx56iGuuuaZfXzd58mRGjx7N448/HqTKTAoDIn3Us61wkGFgqPn1Wyubaev0DrouEQkPGzdu5OSTTyY+Pp5x48axdOlSbDYbq1atYunSpcTFxTF79uy9vuapp54iISGB8vLy7o9dddVVTJ06tXs04Oyzz+bpp58Oau0KAyJ90OHx8eXursWD+YMLA7mpcWQmufH6jO7rFJEDMAzoaLbmzehfy/CNGzdyzDHHcMIJJ7Bu3Tp+97vfccUVV+ByuZg4cSLLli1jxowZ+33dggULGDduHHfccQcAv/jFL3jzzTd55ZVXSEszf9fMmjWLFStW0N7ePvjH9CCcQbtmkSiytbKJDq+PlHgnhRkJg7oum83GxPxUlm2uYn15A9O6uhKKyD46W+COfGtu++YycCf1+eLXX389559/PrfddhtA99D+5s2bcbvdlJSUkJ+//32x2WzcfvvtXHDBBeTl5XHvvfeybNkyhg0b1n2Z/Px8Ojo62L17N8OHDx/8fTsAhQGRPti8xzxLYFxuyqAWD/pNHNoVBrSjQCTilZSU8Pbbb7N69eq9Pu52u5k2bRoAra2txMcfeBfSWWedxcSJE7ntttt4/fXXmTRp0l6fT0gwX4C0tLQEoXqTwoBIH2zpCgNjc5IDcn0TtaNA5PBcieYrdKtuu49WrVrVPR3Q29q1a1m4cCEAWVlZ1NYe+EySV199lY0bN+L1esnNzd3v8zU1ZsfS7OzsPtfUX1ozINIHW/aYc/tjAhQG/NsLN+5uxKvjjEUOzGYzh+qteOvHCKDdbsfr9eL19iwIfvXVV1m7dm33yMCRRx7J+vXr9/valStXctFFF/HQQw9x6qmn8rOf/Wy/y6xdu5aCggKysrIG8CD28T4E7ZpFosjmCnNkIFBhYGRWMvEuOy0dXoqrmwNynSJijRkzZuByubj55pvZtm0bzz77LNdffz1AdxiYN28e69at22t0oLi4mPnz53PzzTdzySWXcNttt/Hss8+ycuXKva5/2bJlzJ07N6j3QWFA5DA6vb7uP9iBCgMOu41xuSkAbK7QjgKRSJafn8+DDz7I4sWLmTZtGs888wyLFi0iLy+PnJwcAKZMmcJRRx3F4sWLAXPo/4wzzuCcc87hJz/5CQDHHHMMZ555JjfffHP3dbe1tbFkyRIWLVoU1PugMCByGCXVLXR6DRLdDvLTBreToDd/sPCvRxCRyHXZZZexc+dOGhsbWbx4MXv27OkeFfD7+c9/zj333IPP5yMjI4ONGzdy//3373WZl19+mVdffbX7/UceeYRZs2bt158g0LSAUOQw/H+sR2cnY7cPfieBn8KASPT64osv9usrMH/+fDZv3kxpaSmFhYV9uh6Xy8W9994bjBL3ojAgchj+xYOB2kngNybbvL7NCgMiUWfNmjVcddVV+338hhtu6Nf1XHvttQGq6NAUBkQOo3tkINBhoOv6tlY24fMZAR11EBFrVVRUWF1Cv2jNgMhhbA5wjwG/ooxE3A47bZ0+SutaA3rdIiL9oTAgcgg+n8HWysBuK/RzOuyMzDLbnW6p1FSBiFhHYUDkEErrWmnr9OF22CnK6HtHsr7qXkRYoTAgItZRGBA5hM1diwdHZiXhdAT+6TJaOwpEJAwoDIgcwrZKs9nQ6Jy+n17WH/51CJomEDH5fD6rS4g4gXjMtJtA5BD8nQf9c/uB1rvXgGEYATkRUSQSud1u7HY7ZWVlZGdn43a79Xw4DMMw6OjooLKyErvdjtvtHvB1KQyIHEJJtXlk6PDM4ISBkVlJ2G1Q39pJZVM7OSkHPuJUJNrZ7XZGjhxJeXk5ZWUWnVQYoRITEykqKsJuH/hgv8KAyCH4RwZGBCkMxLscFGYkUlLdwpY9TQoDEtPcbjdFRUV4PJ69TgCUg3M4HDidzkGPoigMiBxEh8dHaa25/39EZuB3EviNzUmmpLqFrXuaOHZ08I4oFYkENpsNl8uFy+WyupSYogWEIgexq7YFnwEJLgfZKXFBu51R2f5OhDrKWESsoTAgchA96wUSg7qQyT8F4Z+SEBEJNYUBkYMI9noBvxFZ5hREcZXCgIhYQ2FA5CC6RwaygrdeAHq2Le6sbaXTqz3WIhJ6CgMiB7G9KjQjA7kp8cS77Hh9BrtqdWCRiISewoDIQZSEaJrAbrf1rBvQVIGIWEBhQOQAOr2+7lfpI4I8TQA9gWO7woCIWEBhQOQAyupa8fgM4px2ckPQCGhElnYUiIh1FAZEDqC417ZCuz34/dFHdo0+aGRARKygMCByAP71AsE6k2Bf6jUgIlZSGBA5gOIqc2QgmG2Ie/NvLyytbaXDo+2FIhJaCgMiB7CjxnyFXhSikYHslDiS3A58BuyoaQnJbYqI+CkMiBzAzhpzJ0HhkISQ3J7NZuuektD2QhEJNYUBkX0YhsGuWvPVeWFGaKYJoGeqQOsGRCTUFAZE9lHb0klzh3mW+rD00IwMQE8/A+0oEJFQUxgQ2Yd/VCAnJY54lyNkt6sdBSJiFYUBkX10rxcI4RQB9Gxj9N++iEioKAyI7MM/MlAQosWDfkVd4aO0rhWPTi8UkRBSGBDZx07/4sEhoR0ZyEmJw+00Ty8sr28L6W2LSGxTGBDZh3+YPtQjA3a7rfs2d6rXgIiEkMKAyD6s2Fbo558qUOMhEQklhQGRXsweA/6GQ9aFAf9UhYhIKCgMiPRS2dROu8eH3QZD04N/dPG+/AFkh3YUiEgIKQyI9OJfLzA0LQGXI/RPj0JNE4iIBRQGRHrxrxcYFuLFg37d0wQKAyISQgoDIr1YuV4AoDDDDCE1zR00tXssqUFEYo/CgEgvPTsJrBkZSIl3MSTRBWh0QERCR2FApJeeHgPWjAyAtheKSOgpDIj00j0yYNGaAehZRKiRAREJFYUBkS4+n0FZndkG2KoFhKAwICKhpzAg0qWqqZ0Or9ljIC819D0G/DRNICKhpjAg0mVXnbleIC81HqcFPQb8FAZEJNQUBkS6lHWFASunCKBnW+Ou2lZ8PsPSWkQkNigMiHTxh4H8dGvDwND0eBx2G+0eH5VN7ZbWIiKxQWFApEtpV8OhYRaHAZfDTn7XuQiaKhCRUFAYEOlS2rWTwOqRAeh1YFG1woCIBJ/CgEiX0jBZMwA6ylhEQkthQKRL9wLCcBgZ0I4CEQkhhQERoKndQ31rJxAe0wQ6vVBEQklhQISeUYG0BBfJcU6Lq9HIgIiElsKACD3rBcJhVAB6RgYqGtpp6/RaXI2IRDuFARHCZ1uh35DEnhGKXV21iYgEi8KACL0XD1p3JkFvNpuNgq5dDVo3ICLBpjAgQvi0Iu5N2wtFJFQUBkQIvzUD0OvAIjUeEpEgUxgQAcrCqPugX1GmdhSISGgoDEjM83h97G4ww0BBGIWB7pbECgMiEmQKAxLzKhrb8foM3A47WclxVpfTzd9rYFdtK4aho4xFJHgUBiTm+RcP5qXFY7fbLK6mh383QVO7h9qWTourEZFopjAgMa+83pwiGJoWHtsK/eJdDnJTzZEKTRWISDApDEjMKw/DnQR+OqNAREJBYUBinn9kIC/MRgagZxGheg2ISDApDEjM868ZyA/DMFCgkQERCQGFAYl5PWsGwm+aoLC7JbHOJxCR4FEYkJhXXm/+oR0aJucS9FaolsQiEgIKAxLT2j1eqpo6AMgPx5GBrjBQVteK16deAyISHAoDEtMq6tsBiHPaSU90WVzN/vJS43E5bHR6je4uiSIigaYwIDGtrL5nW6HNFj4Nh/wcdlv3lkctIhSRYFEYkJjWvV4gDHcS+HVvL1QYEJEgURiQmOY/rTAcdxL4FWZ0jQzUakeBiASHwoDEtN31/qOLw3dkoKBrZGCXRgZEJEgUBiSm+acJwrH7oJ+2F4pIsCkMSEzzTxOE47ZCPzUeEpFgUxiQmBbODYf8/IcVVTS20dbptbgaEYlGCgMSs1o7vNS2dALhvYAwI8lNotuBYUBpnUYHRCTwFAYkZvlHBZLcDlLjnRZXc3A2m03bC0UkqBQGJGbt7nV0cTg2HOpN2wtFJJgUBiRmlXVvKwzfKQI/bS8UkWBSGJCYVV4X/t0H/bS9UESCSWFAYpZ/ZCCcFw/6aXuhiASTwoDErPLuQ4o0MiAisU1hQGJWeZ1/AWEEjAx0hYG6lk4a2zotrkZEoo3CgMSs7pGBCFgzkBznZEiiC9BUgYgEnsKAxKTmdg8NbR4AhkbAbgLQVIGIBI/CgMQk/6hASryT5LjwbTjUmxoPiUiwKAxITIqEA4r25R8Z2KXGQyISYAoDEpMi4YCiffm7EO7QyICIBJjCgMSk8u4eAxEUBjRNICJBojAgMcm/rTASGg759Z4mMAzD4mpEJJooDEhMKquPnFbEfvnp8dhs0Nrppaqpw+pyRCSKKAxITCqPoEOK/OKcDvJSzfCi7YUiEkgKAxJzDMOIqEOKetO6AREJBoUBiTkNbR6aO7xAZK0ZACjo2lGg7YUiEkgKAxJzdndNEaQnukhwOyyupn80MiAiwaAwIDGnZ/FgZI0KgFoSi0hwKAxIzCnv7j4YWesFAAqHmAFGhxWJSCApDEjMicTug37+kYGyula8PvUaEJHAUBiQmFMWgQ2H/HJT43E5bHh8RneoEREZLIUBiTnlEdhwyM9ht1HQtYhQZxSISKAoDEjM2V0fuSMDAAVd6wZ2ad2AiASIwoDEFMMwuncT5EfgmgHQjgIRCTyFAYkpdS2dtHX6AHP+PRKp14CIBJrCgMQU/6hAZpKbeFdkNRzyK+zqQrhTXQhFJEAUBiSmdB9dHKFTBKCRAREJPIUBiSnlDZG9eBB61gzsaWynrdNrcTUiEg0UBiSm+E8rjMTug35DEl0kdZ2poAOLRCQQFAYkppT7txWmR+7IgM1m044CEQkohQGJKWV1kdtwqDd/46FdWjcgIgGgMCAxpTzCGw75aUeBiASSwoDEDJ/P6NV9MLJHBvw7CnZUa2RARAZPYUBiRk1LBx1eHzYb5EV4GCjSmgERCSCFAYkZ/h4D2clxuByR/aPfvYBQawZEJAAi+zeiSD/4uw9G8k4CP/+agYY2D3UtHRZXIyKRTmFAYoa/x8DQCD2ToLdEt5PslDgASrRuQEQGSWFAYkZPj4HIDwMAIzLNqYISTRWIyCApDEjMKOsKA/kRvq3Qb3hmEgAlVc0WVyIikU5hQGLG7u41A9E1MlCsaQIRGSSFAYkZZXXR0XDIr3tkoFojAyIyOAoDEhO8PoOKhuhoOOQ3XGsGRCRAFAYkJlQ2tuPxGTjsNnKjYDcBwPAMc2SgsrGd5naPxdWISCRTGJCYUNq1rTAvNR6H3WZxNYGRluhiSKIL0PZCERkchQGJCf7TCvOjZPGgn9YNiEggKAxITCiv94eB6Fg86Kd1AyISCAoDEhP8OwmiLwxoZEBEBk9hQGJCaV10jgx09xqo0siAiAycwoDEBP+agWFaMyAish+FAYkJZVE6MuBfM1De0EZbp9fiakQkUikMSNRr7fBS29IJRE/3Qb/MJDfJcU4MA3bVaqpARAZGYUCiXlnXToLkOCep8U6Lqwksm83WPTqwrVJTBSIyMAoDEvV69xiw2aKj4VBvo7KTAdiu0wtFZIAUBiTqRet6Ab+RWeYiQoUBERkohQGJeqVR2mPAb3S2GQa2KQyIyAApDEjUK/ePDETJaYX78o8MaM2AiAxUdK2mEjmAsihtReznDwNVTe00tHWSGu+yuKLI1NLhobKxndqWTpLjHGQnx5Oa4IzKdSYi+1IYkKgXra2I/VLiXWSnxFHZ2M72ymamFaZbXVJEMAyDz0pqWbpmN+9vqWRTRdN+l8lJieP4sVmcOiGXuZNycTk0mCrRSWFAopphGN2tiIdFaRgAGJWVZIaBKoWBwzEMg9fW7eb+d7examfdXp9LcDnISHLT2NZJQ5uHPY3t/GdlKf9ZWUp+WjxXHTeSy+cMJ97lsKZ4kSBRGJCoVt3cQYfHh80GuanRuWYAYFR2Eh9vr2Fb5f6vbqVHSXUztyxZy7LNVQC4nXbOmjqUUyfkMmd0JhlJ7u7LtnZ4Wbmjlvc2VfLsylLK6tu4fekGnvi4hF+fO4Xjx2ZZdTdEAk5hQKJaedcUQU5KHG5n9A7xjsoyew1oR8HBPfPJDn7+/DraPT7cTjuLThjJlceOJDsl7oCXT3A7OG5MFseNyeLG08ex5PNS7n5zM8XVLVz20MdcNruIn581Kap/riR2KAxIVPNPEURbG+J9qdfAwXV4fNz64jqe+HgHAMeOzuT2r0/pfsz6It7lYMGsIuZPHcpdr2/iseXFPP7RDjaUN/LXS48iJ4pHnSQ2KNJKVCuLgfUCYE4TgBkGDMOwuJrw0dLh4epHP+GJj3dgs8GPTh/H49cc068g0FtKvItffm0SDy88mpR4J5+V1HLuXz5QCJOIpzAgUa13K+JoVpiRiNNuo6XDS0VDu9XlhIXGtk6ufPgT3t9SRaLbwUMLZ/K9U8ditw9+q+DJE3J44bvHMyo7ibL6Ni7+23I2VzQGoGoRaygMSFSL9h4Dfi6HnaIM/4FFWkTY0uFh4cMrWFFcQ0q8k8evPYZTJuQG9DZGZiWx+FtzmJCXwp7GdhY88BFb9dhLhFIYkKgW7a2Ie+vuRBjjQ9Yer4/vPfk5K3fUkZ7o4qlFszmqaEhQbisrOY6nFs1m8rBUqps7WPjwCvY0tgXltkSCSWFAolpPK+LoDwP+dQNb9sTuq1PDMLhlyVre2riHOKedhxYezeRhaUG9zSFJbh69ahYjMhPZVdvKVY98QlO7J6i3KRJoCgMStdo9XvY0mvPn0b5mAGBsTgpATA9VP/xBMU9/shO7Df78jaOYMTw4IwL7ykqO47GrZ5GZ5GZdWQM/WrxKCzkloigMSNSqqDeDQJzTvlczmWg1NtfsNbApRheyfbStmjuWbgDgZ2dN5PSJgV0jcDjDM5P4+8KZuB12XltXwf3vbgvp7YsMhsKARK3ebYhj4bCZMTlmGKhoaKe+tdPiakJrd30b331yJV6fwdePHMaVx46wpI6jiobwi69NBODO1zbyflenQ5FwpzAgUatnW2H0rxcAcw/80K5jmrfsiZ3RAa/P4IZnPqeqqYMjhqZyx9enWBr+vjGriAtnFOAz4MbFq6hu0lZPCX8KAxK1yutjo8dAb2NzzXUDmw9wAl+0euC9bXy0rYZEt4P7Lj2KBLe1hwjZbDZ+de5kxuYkU9nYzk/+s0brByTsKQxI1PJvK4z2VsS9je2aKtgcIzsK1uyq567XvwTgl2dPGnBnwUCLdzm4e8F0XA4bb6yv4OlPdlpdksghKQxI1IqVVsS9+cNALCwibPd4uXHxKjw+gzMn53HhzAKrS9rLpPw0bpo3HoBfvbSenTUtFlckcnAKAxK1Ym3NAPRME8RCr4E/v72FLXuayEp2W75O4GCuPX4Us0Zk0NLh5X+XrNV0gYQthQGJSoZhxMy5BL35dxSU17fR2Ba9OwrWlzXw13e2AnDbOZMZEqZbR+12G789fwpup533NlXyn5WlVpckckAKAxKVGto8NHd4gdgaGUhLcJGbGgdE77oBr8/gx89+0T098NUpQ60u6ZBGZSdzw2ljAbjtpfVUaXeBhCGFAYlK/lGBjCQ38S5rV5eH2jj/VEGU7ih4asUO1pTWkxrv5NZzJlldTp8sOmEUE4emUt/aye9e2Wh1OSL7URiQqBSLUwR+Y6J4EWFtcwe/79o98MPTx5GTEhnfX5fDzq/OnQzAvz7bxWclNRZXJLI3hQGJSv4wEEvbCv38ZxRsisJpgrve+JK6lk4m5KVw2ezhVpfTLzOGD+Girh0PP1uyDq9PiwklfCgMSFTa1RUGCobEXhgYn2eGgS93N1hcSWCtLa3niY93APDLr03C6Yi8X18/PmMCqfFO1pc38OSKHVaXI9It8p5NIn2wq9YfBhItriT0JuSlYLOZZxREy2I1wzD45QvrMAw4e1o+s0dlWl3SgGQmx/HD08cBcM+bm3TUsYQNhQGJSj1hIPZGBpLinIzINDvxbSiPjtGBJatK+bSklgSXg5u/OsHqcgbl0tnDGZmVRFVTBw+8u9XqckQAhQGJUru6ur3FYhgAmDg0FYiOMNDY1skdS80V+N87dUzErwNxOez8v67OhH9ftp2KhjaLKxJRGJAo1NLhobq5A4jNaQKAI4aa6wbWl0V+GPjbu9uobGxnZFYS1xw/0upyAuKMyXkcVZROa6eXu9/cZHU5IgoDEn1Ku6YIUuKdpCW4LK7GGhPz/SMDkb29sKqpnYc/2A6Yi+/inNHRM8Jms3HzV48A4JlPdrI5CreBSmRRGJCoE8uLB/2O6Jom2FLZRFun1+JqBu6+/26lpcPLtII05k3KtbqcgJo5IoN5k3LxGfBbNSISiykMSNTZVRvb6wUA8lLjGZLowusz2ByhnQjL6lp5/KMSAP5n3viwPIhosH58xgQcdhtvbdzD8q3VVpcjMUxhQKKOf2SgMIZHBmw2W/foQKQuIrz37c10eH3MHpXB8WOyrC4nKEZlJ/ONWUUA/PaVDTrVUCyjMCBRJ5a3Ffbm31GwPgLDwPaqZhZ/uguAm6J0VMDvB6eNJdHtYPWuet7asMfqciRGKQxI1NE0gcm/iDASw8Ddb27C6zM4ZUIOM4ZnWF1OUGUlx7Hw2BEA3P3WJo0OiCUUBiTqaAGhqXuaoKwhov7AbNzdwAurywD40dxxFlcTGotOGEWS28Ha0gbe1OiAWEBhQKJK7x4Dw2J8ZGBMTjJup53Gdg/bq5qtLqfP7np9E4YB86cOZVJ+mtXlhERGkrtndOBNjQ5I6CkMSFTx9xhIjeEeA34uh53JXVMFq3bWWVtMH32+o5Y31ldgt8GNp8XGqIDftV2jA+vKNDogoacwIFFFUwR7m144BIDVERIG7nrd7MZ3/lEFjMlJtria0NLogFhJYUCiihYP7m16UToQGSMDH26t4v0tVbgcNr5/6liry7HEol6jA2+sr7C6HIkhCgMSVTQysLcjC9MBc0dBOHciNAyD37/2JQDfmFVEYUZsfv+GJLm58rgRANz95maNDkjIKAxIVNnZNTIQ64sH/QqGJJCZ5KbTa4T1FsO3N+5h5Y464l12rj9ljNXlWOra483RgfXlDbyu0QEJEYUBiSo7a8yRgaIYfWW5L5vNxvSu0YFVO+osreVgfD6D33etFbjy2JHkpMRbXJG1hvRaO3Dff7dodEBCQmFAokpJtbmFTmGgR3cYCNN1Ay+vKWdDeQMpcU6+/ZVRVpcTFq4+fiTxLjurd9Xzoc4skBBQGJCoUd/SSUObB1AY6G1aVxhYvavO0joOxOP18cc3zFGBRSeOIj3RbXFF4SErOY4FR5tnFtz3zhaLq5FYoDAgUaOkxhwVyE6JI8EdHefeB4I/DJRUt1DT1ZApXPxnZSnbqprJSHJz9fEjrS4nrCw6cRROu40PtlRHzNZQiVwKAxI1dtSYiwc1KrC3tAQXo7OTAFhZUmtxNT3aPV7ueWszANedNJrkOKfFFYWXYekJnDN9GKDRAQk+hQGJGv4wMFxhYD+zRpqH/Xy0LXzmn5/6eAelda3kpsZx2ezhVpcTlr5z0ihsNnhtXQVb9jRaXY5EMYUBiRo7qs0wEKt71A9l9qhMAJaHSRho6fDw5/9uBeD7p44l3qVpnQMZk5PC3Im5APz1nW0WVyPRTGFAokb3yECmwsC+5nSFgfXlDdS1WL9u4LEPS6hqaqcoI5GLZhZaXU5Yu+4ks+/C86tKKa1rtbgaiVYKAxI1tGbg4HJS4xmdnYRhwMfbayytpb61k/vfNUcFbjx9LC6Hfg0dyrTCdI4bk4nHZ/D39zQ6IMGhZ6FEhQ6Pj7I6NRw6lDmju6YKLN63/tCybdS3djI2J5mvTRtmaS2Rwj868PQnO6huare4GolGCgMSFcrqWvEZEO+yk50SZ3U5YWnOqCzA2kWE1U3tPPT+dgB+NHc8DrvNsloiybGjM5kyLI22Th+Pf7TD6nIkCikMSFQo6TVFYLPpD8yBzB5l7ijYuLvRsleXf31nK80dXqYWpDFvUq4lNUQim83GohPN7oz/WF4c1odOSWRSGJCo0LNeIMniSsJXZnIc43NTAGvWDZTWtfKPj0oA+J+54xXa+umrk/MYlp5AdXMHz67cZXU5EmXU5UOiwk4tHuyTU4a7cO4ppvyzMjDSoKMZvJ3mJx0ucCdBYgYkZUNaofn/ALnr9S/p8PiYPSqDE8ZmBex6Y4XTYeea40dy20vreXDZdhYcXaRpFgkYhQGJCj0HFOnoYgAMA2q2wY7lUPY57F4DlRv5cVs9P44DirveDic+DbInQN4UyD8SiuZAxijo56v6dWX1PPd5KQA3f/UIjQoM0MVHF3L3m5vYXtXMmxsqmDcpz+qSJEooDEhU2NF1dPHwzBieJvB0QPEy2PgSfPkKNJYf8GJVpLHHl07BsAJSU9PB6e75+o4maKmBpgpo3gNt9bDzY/PNLyUfxp8BR5wNI04wRxQO47evbMQw4GvT8plakD74+xqjkuKcXDZ7OPe9s5UH3tumMCABozAgEc8wDHZ0jQzEXPdBnw+2vg1rFsOmV80/3n4ON+QfBYVHQ940yJ0IQ0Zy27Nf8sLqMr5VNIqfnnnEwa+7owVqt0PFeti9GnZ+AmUrobEMPn3YfItPg3FnwJSLYPQpYN9/GdJ7mypZtrkKl8PGTfPGB+FBiC1XHjuCB5dt57OSWj4rqWHG8MBN5UjsUhiQiFfZ1E5zhxe7DQpjZZqgpQZWPQGfPGT+wfZLzIIJ87tetR8Prv0fj7mTcnlhdRlvrKs4dBhwJ0LuJPNt6oXmxzpbofh92PAibHwZWqrgi2fMtyEj4ehrYPql3WsNvD6D37yyEYAr5oyIvbAWBDmp8Zx7ZD6LP93FA+9t42+XKwzI4CkMSMQrrjIXDw4bkkCcM8p73Fdvhff/AGv+DZ4282NxaTDtYpj0dSg8BuyHfgy+Mi4bt8POtqpmtuxpYkxOct9v35UAY0833876ozl9sO45WP2MGUpevwXe/jVMuRCOv5ElxXFsKG8gJd7Jd08eM4g7Lr0tOmEUiz/dxevrK9he1czIrBieHpOA0NZCiXjbq5oAGBHN6wVqi2HJ9fDno+Hzx80gkDsFzr4HfrQBvnonDD/2sEEAICXe1d2N8PX1uwdek91h3uZX7zRrOPsesyZPG3z+T4w/H4375e9TYKvkupPGMCTJPfDbkr2MzU3hlAk5GAY8uEwtimXwFAYk4m3vGhmIyldH9aXw4g1w7wxY9TgYXhg7D65+Db69DGZcaW4H7Ke5XQ1/lq458CLDfnMnmbV8e5lZ29i52AwvZ/ve4p24H7Ko/k/QUBaY2xLAHB0A+Pdnu9SiWAZNYUAiXnGVuXgwqsKApwOW/QH+PBM+ewR8Hhh1MlzzJly6GIpm93t7X29nTh6K22FnbWkDa0vrD/8FfWWzQdFsSuf/gwXe23jPOwUnXpyfPwr3zoT3/2jeNxm02aMymFqQRrvHxz+7mjmJDJTCgES84q6dBCOiJQxseQv+OgfeuhU6W8x1AFe9AlcsMXcGBEBGkpt5k81taU+tCHyv+9teXMdHnWP4c8GdGFe+bN6HzmZ485fw12PNHRAyKDabrXt04B/LS2jtUItiGTiFAYloPp/Bdv/IQKSvGWjcDc9cDo+fB9VbICkHvv43c9h9+LEBv7lLZhUC8PyqMprbPQG73v9+uYfX1lXgsNv41TmTsY043rwPX/+beZ+qN8M/v27e18ZBrFkQzpycR8GQBGrUolgGSWFAItruhjbaPT6cdhsFQyJ0W6FhmLsD/nIMbHgBbA6YfR1871OYtmBQ0wGHMmdUJiMyE2lq9/DSF4GZz29q9/CzJWsBuPq4EYzPM89CwGYz78v3PoVjvmPexw0vwH2zzftuGAG5/Vjjb1EM5kJCr0+PowyMwoBENP96gaKMRJyOCPxxbq6Cfy2EZ6+BtjoYOg2+9R6c8RuzoU8Q2Ww2LplVBMCTHwdmquD2lzewq7aVgiEJ/OC0cftfID4NzvyteR+HToPWWvO+/2uh+VhIv100s5C0BBfF1S28sb7C6nIkQkXgb0+RHtuqIni9wMal5mjA+ufB7oSTboZr34K8ySEr4fwZBbgddlbvque9TZWDuq53N1V2rz+484JpJMcdoo1J3mTzvp50s3nf1z9vjhJsXDqoGmKR2aLYDHUPvLfV4mokUikMSESLyJ0Eng549afw9CVmB7+cibDobTjpx33q8x9IWclxXD5nOGCeH+Ab4DBzVVM7P/73F4DZLtffx+CQHC7zPl/7lvkYNFeaj8mrN2vHQT8tnDMCt8POyh11fFYS+uOpJfIpDEhEi7idBLXF8PA8+Og+8/0534VvvmMOmVvk+pPHkBLnZH15Ay+s7v/agU6vj+ufWMnuhjZGZSfx/87o5/kD+dPNx2DOd833P/oLPHIG1Gq7XF/lpMbz9SOHAfC3d9WESPpPYUAi2rZI2kmw4UW4/0TzsJ/4dLjkaZh3OzjjLC0rI8nNt08aDcDvX/+Sts7+bVH7zdKNfLy9hiS3gwcun0GiewBdzp1x5mNxydPmY1P6GfztBNjwUv+vK0YtOtFcSPjGhgq2VTZZXI1EGoUBiVger4+dNWb3wRFZYXwAjs9r7q9/5jJor4eCo+Hb78P4M62urNvVx40kNzWOXbWt3PriOow+ru5/cNk2Hv7APCjproumMyYnZXCFjD/T7GJYcLR5AuMzl5qPnU976A9nTE4Kp/pbFL+//fBfINKLwoBErLK6Njq9Bm6nnfy0MN1W2FoHT15sdt4Dcyj8qlcgvdDSsvaV4Hbw2/OnYrPBUyt28ngfOto98N5Wfv3yBgB+ePo4zuhqYjRo6UXmY+SfNnj/j/DUgr2PZ5YDWnSi2YTo2c92UaUWxdIPCgMSsbZ2DYWOzEzCbg/OXvxBqdwED54KW94AZwKc/5A5FB7iRYJ9dfL4HH58xgQAbn1xPS8eZP1Ac1cvgTuWmkcT/+DUsXz/1LGBLcbhMh+r8x8CZzxsfh3+fipUbQ7s7USZY0ZmMM3foni51lxI3ykMSMTasscMA/06gjdUNr1uBoHqLZBaAFe/ClMusLqqw/rWiaM4Z3o+Hp/B9576nOufWMmG8gbaPV72NLax+NOdnHHPe9298H94+jhuPP0A/QQCZcoFZvfC1AKzc+HfTzEfWzkgm83WPTrwj+XFalEsfTaAlT4i4cEfBkaHWxj4+AF49cdg+KDoWLjoH5CcbXVVfWKz2bjzgmkUZSRy3ztbeXlNOS+vKcdht+3V3W5YegL/d8FUjhuTFfyi/LsNFl8OO5bDUxfDmf8HsxYF/7Yj0BmT8ijMSGBnTSv/+mwnV8wZYXVJEgE0MiARa/OeRgDGhksY8HnN/gGv3GQGgSMvhyuej5gg4Od22vnR3PEsue44jhuTSUqcszsITC1I48bTxvHajSeGJgj4JWfDFS+Yj6nhg6X/Y/Yj0MLC/Tgd9u4DjB54bxser8/iiiQS2Iy+LhsWCSOGYTDt1tdpaPPwyg9O4IihqdYW1NEC/1kEG7u2wp36Czj+xqCdKxBKhmGwp7Edt8POkCS31cXA+3+At24z359wFpz3d3CH8W4SC7R2eDnud29T09zBPQumc870YVaXJGFOIwMSkSqb2mlo82C3hUH3waZKeHS+GQQccXDBw3DCD6MiCIA5dZCbGm99EDCLgRN+ZC4sdLjNx/yxs8zvgXRLcDu48tgRANz/7rY+bxWV2KUwIBHJv16gMCOReJfDukJqS+DhuWYjoYQMWPgCTD7funpixZQLzGmDhCFmg6KH56lj4T6umDOcRLeDDeUNvDvIcyck+ikMSETa6t9JkG3heoGK9fDQXKjZBmlFcM0bUDTbunpizfA5cM2b5mNfs9UMBBXrra4qbKQnurtPpfzrOzrASA5NYUAi0mZ/GMi1KAzs+Mjsn9+0G7KPgGteg6wx1tQSy7LGmI999hHQWG5+T3Z8bHVVYeOa40fitNv4eHsNK3fUWl2OhDGFAYlIW6wcGdj0GvzjXLMjXuExcNVSSM0PfR1iSs03vwcFs8zvyT/OUS+CLvnpCZzbdYDR/RodkENQGJCIZFnDodXPwFOXgKcVxs6Fy5dAYkZoa5D9JWaY2zjHnG5+b55aYH6vhG9/xdxm+MaGiu7njci+FAYk4jS0dbKn0ey7HtKGQ8vvg+e+CYYXpl4MC57UlrZw4k6ES56CKReZ36Pnvml+z2LcmJwUTp+Yi2GY50mIHIjCgEQc/6ub3NQ4UuND1Of/3f+D135q/n/2dXDu/WF7xkBMc7jg63+DY75jvv/aT+Gd35n9CWLYt79iHlH93OellNe3WlyNhCOFAYk4WyrMMDB2sMfl9oVhwNu/hv/ebr5/8i0w7w6w66kTtux2OOM35vcK4J074K1bYzoQzBg+hFkjMuj0Gjys443lAPQbTSLOpgqzDXHQ1wsYBrzxc3jvTvP9038FX7kpapoJRTWbzfxeze0Kce//0WwVHcOB4DsnmaMDT368g/qWTourkXCjMCARZ+NuMwwcMTSIIwOGAa/+BD78k/n+mf8Hx30/eLcnwXHsd2H+Xeb/P/4rvHQD+GKzV/9J47OZkJdCc4eXf35UbHU5EmYUBiTibNzdAMCEvCCdR+DzwUs3wsf3m++fdTcc863g3JYE39HXwjl/AWzw2aOw5Dvg9VhdVcjZbDa+1bWz4JEPdLyx7E1hQCJKZWM7VU0d2GwwLjcIIwM+L7zwPfjsEcAG59wHM68K/O1IaB15GZz/INgc8MXT8Ow14I29ofKzpuZTMCSB6uYOnlyxw+pyJIwoDEhE8Y8KjMxMIsEd4DMJvB547luw6nHzj8Z5f4cjLw3sbYh1plwAFz0GdhesXwKLrwBPu9VVhZTLYee6k8xOmX97dyttnRodEJPCgESUjeXmeoEJgV4v4O00Xy2u+RfYnebJg1MvDOxtiPWOONvsReCMhy+Xms2JOlqsriqkzp8xjPy0ePY0trP4051WlyNhQmFAIsqG8iCsF/C0w+KF5qtFhxsu+idMOjdw1y/hZezp8I3F4EqErW/DExdCe+x05otzOrp3Fvz1na20ezQ6IAoDEmE2dO0kmJAXoJGBzlZ4+lL48mVwxMGCp2DCVwNz3RK+Rn0FLn8O3ClQ8j48fp55rkGMuHBmIbmpcZTXt/Hvz3ZZXY6EAYUBiRidXh9b9vi3FQZgZKCjxRwm3vIGOBPg0sUw9rTBX69EhqLZsPB5iE+DnR+bh0+1xsbJfvEuR3dXwvv+u5VOb2xut5QeCgMSMbZVNtPpNUiOczIsPWFwV9beCE9cANveAXcyXPYsjDopEGVKJBk2Axa+CAkZULYSHvsaNFdbXVVIXDKriKzkOErrWnluZanV5YjFFAYkYvT0F0jBbh9EF8C2evjneVDyAcSlmsPFI44LUJUScYZOgytfgqRs2P0FPHY2NFVaXVXQxbscfOtEs+/An/+7BY9GB2KawoBEjA2B2EnQWmsOB+9aYQ4PX7EECmcFpD6JYLmT4MqXITkP9qyDR+dD426rqwq6S2cXkZHkZkdNC0tWlVldjlhIYUAihn8nwfiB7iRorjZf9ZWtNIeFF75kDhOLAGSPh6uWQuowqPoSHvkq1Ef38Hmi28miE8zRgXve2kSHR6MDsUphQCKCYRisLTVXe08Zltb/K2jaY77a270GknLMV4FDpwa4Sol4maPNQJBeBDVb4ZEzobbE6qqCauGxw8lKjmNnTav6DsQwhQGJCKV1rVQ3d+C02/q/rbChzHyVV7nBHAa+8mXInRicQiXyDRkBVy6FISOhrsQMkTXbrK4qaBLdTr53itmV8E9vbVZXwhilMCARYc0uc1RgfF4K8a5+tCGu22kGgerNkFpgvurLHhekKiVqpBfCVa9A5lio7/oZqtpsdVVBs2BWIcPSE9jT2M4/lhdbXY5YQGFAIsIXXVMEUwv6MUVQs838JV67HdKHm0Egc3SQKpSokzq0KzweAY3l5s/Sng1WVxUUcU4HPzhtLGB2JWxsi71DnGKdwoBEBP/IwJRh6X37gspNXQvAdkDmGPNV3pDhwStQolNyjrntMHcKNPdadxKFzjtyGKOyk6ht6eSh97dbXY6EmMKAhD3DMPhiVx3Qx5GBinXw6FfNV3PZR5jzv2nDglukRK+kLFj4AuQfCS3V8OhZUPa51VUFnNNh50enjwfgwWXbqW3usLgiCSWFAQl7O2paaGjz4HbaGZd7mMWDZavMV2/NlZA3xVwsmJIbkjoliiVmwBXPQ8HR0FYHj50DOz+xuqqAO3NyHhOHptLU7uH+d7daXY6EkMKAhL3VXVMERwxNxe08xI/szk/MdrKttT1tZpMyQ1SlRL34NLNbZdGx0F4P/zwXSj60uqqAsttt3DTPHB149MNiyupaLa5IQkVhQMLeGv8UwaH6CxR/YP5ybq+Hojlw+RJIGBKK8iSWxKXAZf+GkSdCRxM8fj5se9fqqgLqpPHZzBqZQbvHx52vfWl1ORIiCgMS9r7wLx482HqBrW+bv5Q7mmDkV8xDh+IDcKqhyIG4k+Abi2H0qdDZAk9eBFvetLqqgLHZbPxsvtmH47nPS1m9s87agiQkFAYkrHl9PZ0HD7h48MtX4ckF4GmFsXPhG8+Yv6xFgsmVAAuehHFngKcNnrrE/FmMElMK0jjvKHPR7a9eWo9hGBZXJMGmMCBhbePuBpo7vCTHORmbs8/iwbXPwjOXgrcdJpwFFz9u/pIWCQVXPFz0TzjibPB2wDOXmT+TUeKmeeOJd9n5tKSWV9ZG/6FNsU5hQMLap8W1ABw1fAiO3scWf/IQ/Psa8Hlg8gVw4aPgjLOmSIldTjdc8Ij5M+jrNH8mP33Y6qoCYmhaAt860WzS9ZtXNtDuUZviaKYwIGHtk+IaAI4e3rUY0DDgvd/Dyz8EDJh5DZz3ADhc1hUpsc3hMn8GZ14NGPDSjbDsLvNnNcJ96yujyEkxDzF69INiq8uRIFIYkLBlGEb3yMDMERnmL9fXb4G3f2Ve4MSbYP5dYO/HWQUiwWB3wPw/mD+TAG/dZv6sRnggSHQ7u7ca/vntLVQ3tVtckQSLwoCErdK6VnY3tOG025g+LBme/y4s/7P5yXm/gVNuAZvt0FciEio2m/kzOe8O8/3lf4YXvgtej7V1DdL5RxUwKT+VxnYPv39dWw2jlcKAhC3/qMD0/AQSllwNqx4HmwPO/SvMuc7i6kQOYs71cM59YLPD54/DvxZCZ5vVVQ2Y3W7jF2dPAuCpFTv5rKTG4ookGBQGJGx9UlxDCi38vv1XsPElcMTBxf+E6d+wujSRQzvyUnOngcNt/uw+cQG01Vtd1YDNGpnBhTMKAPjf59bi8fosrkgCTWFAwlbJtk38y30rIxo/A3ey2fltwnyryxLpmyPOMhtguZOheBk8fAbU77K6qgH76VePID3RxcbdjTz6YbHV5UiAKQxIWGos/pzfN/yICfadeJNyzSOIR55odVki/TPyRLhqKSTnwZ718ODpsHut1VUNSEaSm5+eOQGAP7yxSecWRBmFAQk/W98m4fH55NlqKbYX4lj0FgydanVVIgMzdBpc+wZkjYfGMnOEYOvbVlc1IBfOKGTm8CG0dHi59cV1VpcjAaQwIOHl8yfgiQtxeppZ7p3IExMfgPRCq6sSGZz0IrjmNRh+PHQ0whMXwqonra6q3+x2G7/++mScdhuvravgzfUVVpckAaIwIOHBMOCd38Hz14HPw5vOE1nY+WNmHjHK6spEAiNhCFz+n65uhR5Y8h149/8irhfBhLxUrjlhJAC/eGEdTe2RvXVSTAoDYj1POyy5Dt4x92c3zPgui5q+idfuZs7oTIuLEwkgZxyc93c47gbz/f/eDs9fbz4HIsgPTh1LwZAESutauf3lDVaXIwGgMCDWaqyAR8+C1U+a+7Ln38XLud/CwM70wnRS49VmWKKM3Q6n3wpf/b35M7/qCfM50Bg5Q+6Jbid3XjANgKdW7OC/G/dYXJEMlsKAWKdsFfz9ZNi1AuLT4NJ/w9HX8v7mKgBOGJtlbX0iwTRrEVz6L4hLM58Dfz/ZfE5EiDmjM7n6OHO64MfPfkFdS4fFFclgKAyINdY+a66qbiiFzLFw7dsw5lS8PoP3tygMSIwYcxosett8DjSUms+Jtf+xuqo++39njGd0dhJ7Gtv52fPaXRDJFAYktHw+eOtX8O+rwdMKY06HRW9B1hgA1pbWU9/aSUqck2kF6dbWKhIKWWPg2jfNYOBphX9fBW//2nyuhLl4l4M/XDQdh93Gi6vLeOmLMqtLkgFSGJDQaW+EZy6DZb833z/2+/CNZ8wpgi7vbqoEzCFIp0M/nhIjEtLhG4vh2O+Z7793p/lcaW+0tKy+mFaYzvUnjQbgliVr2dMQuecwxDL9tpXQqFgPD5wMX75snjHw9b/B3F/td/zwa+t2A3DaEblWVCliHbsD5v4azr3fPNPgy5fN58ye8F+t/91TxjJ5WCp1LZ3cuHgVXl9kbZcUhQEJhdVPw99PgerNkDrMbM86bcF+F9tZ08K6sgbsNjj1iBwLChUJA9MvMdtvp+Sbz5m/nwKrn7G6qkNyO+3cffF0Et0OPthSzR/f2GR1SdJPCgMSPJ1t8OIN8Ny3zLnQUSfDt96DgpkHvLh/VGDWyAwyk+NCWKhImCmYCd9eZj5nOlvguW+az6UwPgp5TE4KvzlvCgB//u8W3t4YOVslRWFAgqVqCzw8Fz57BLDBV35inuCWdPAdAv4wcMakvBAVKRLGkrLM58xXfgLYzOfSw3PN51aYOmf6MK6YMxyAG59Zzc6aFosrkr5SGJDAMgxY+U/42wlQvhoSMsyjh0/+6X7rA3rb09jGpyW1AMxVGBAx2R3mc+eyf5vPpfLV8LcT4fPHw7aN8f/OP4LphenUt3Zy3RMraev0Wl2S9IHCgAROa525LeqF75pDmyNPhO98YG6ZOow31ldgGDCtII389ITg1yoSScacZj6XRpwAnc1mC+N/X2U+58JMnNPBXy49iiGJLtaU1nPbS+utLkn6QGFAAqP4fbj/eFj3HNidcOov4PIlkJrfpy9fuqYcgHmTNSogckCp+XDF8+Zzy+40n2v3H28+98LMsPQE7llwJDYbPPnxDv65vNjqkuQwFAZkcDqa4ZUfw6PzoX4nDBkJV78OJ/zwkNMCvZXWtfLh1moAzprSt/AgEpPsDvO5dfXrMGSE+Zx7dD688hPoCK/5+RPHZXPTvPEA/PLF9bzzpc4vCGcKAzJwJcvNVyYf32++P+NKcwV0wYx+Xc2zn+3CMGDOqEyKMhMDX6dItCmYAd9+H45aaL7/8V/h/uNgx0fW1rWP73xlNBfOKMDrM/juk5+zobzB6pLkIBQGpP/am+DVm+GRM6Fmm7kf+rJn4ex7IC6lX1fl8xn8+7NdAFw4syAY1YpEp7gU+NqfzOdeSr75XHz4DHjtf80RuzBgs9m4/etTOGZkBk3tHhY+vEI7DMKUwoD0z5evwn2z4aO/AAZMvxSuW96nRYIH8vH2GnbUtJAc5+TMyUMDW6tILBhzmvkcnH4pYMDyP8NfjoFNr1ldGWA2JHrg8pmMz01hT2M7Cx9eQXVTu9VlyT4UBqRvGsph8RXw1MXmPGVaEXzjX3DufWZf9QH612c7ATh72lAS3H1bYyAi+0hIN5+L31gMaYXmc/TJi8znbEO51dWRlujisatnMSw9gW1VzVz+0AodeRxmFAbk0DwdsPwv8OejYf3zYHOYh6lc/xGMmzuoq65t7uCVNWajoQtmFAaiWpHYNm4eXP+x+Ry1Oczn7F9mmc9hj7V/fPPS4vnHNbPISo5jfXkDlz+0gvrWTktrkh42wwjTzhVivc1vwKs/NfujAwybAWfdDUOnBuTq731rM3e9sYmJQ1N5+fvHY7PZAnK9IgKUfwEv3QCln5nvZ46FM34DY0+3tKxNFY0seOAjapo7mFaQxmNXzyI90W1pTaIwIAdStRleuxk2v26+n5gFp/4cjrysz9sFD6et08vxv3ubqqYO7lkwnXOmDwvI9YpILz6v2a3wrdugpcr82Ni5MO8OyBprWVnryxq49MGPqG3pZEJeCv+85hiyU3QeiZUUBqRHfSm8+7uuVqdesLtg9rfhxJsgPi2gN/XExyX873NrGZaewDs3nYTLoRkrkaBpq4d3/8/cBuzzmFMIR10OX/lxnxuDBdqmikYuffBjKhvbGZmVxGNXzdLWYgspDAi01MD7f4AVfwdP16lo4840z1bPGhPwm/P6DE77w7tsr2rm52dN5OrjRwb8NkTkAKq2wOu3wKZXzPed8TDrm3D8jZCYEfJyiquaufTBjymtayUr2c1DC49mWmF6yOsQhYHY1lJjvlL46K/Q3tUMpGgOnPZLKJodtJtd8nkpNzyzirQEFx/+5BSS4pxBuy0ROYCS5fDWrbBjufl+XCrM/g4c8+2Qh4KKhjaueuQT1pc3kOBy8MeLp3GGthmHnMJALGqsMPcif/owdDSZH8udYq4LGHs6BHEhX1unl1N+/w5l9W3cNG88158c+JEHEekDwzAXCb91K1SsNT/mToajr4E534XknJCV0tTu4bonVvLepkoAvn/KGG44bRx2uxYVh4rCQCyp3mqOAnz+z57pgLwpcML/wBFfA3vw5+3/8t8t3PnalwxLT+CtH32FeJd6C4hYyueDDc/De3dBxRrzY854OPJyc7Qgc3RIyvB4fdyxdCMPf7AdgJPGZ/OHi6aTkaSdBqGgMBDtDAO2v2eGgE2vAl3f7oKjzYWBY+cGdSSgt6qmdk668x2a2j3cffF0zj1SOwhEwoZhmF0L37sTSj/t+qANxp1hhoKRJ4bkd8W/P9vFzc+tocPjIy81nrsXTGf2qMyg326sUxiIVm0NsOZf5lSAfwgQup7Y14Xsid3bD57+nOdXlTG1II0l1x2nIUCRcOR/AbH8L7C5V0vj3Mkw82qYciHEpwa1hA3lDVz/5Eq2VZpnLFwxZzj/74wJJGt9UdAoDEQTw4Bdn8Bnj8G6/0Bn14EgrkSzb/kx3w7K7oC+eOmLMr775Oc47Dae/c6xTNeKYZHwV7XFPBFx1ZN7/z6ZfB4cdSUUzAzai4rmdg+/emk9T39itiwflp7Ab86bwonjsoNye7FOYSAa1O0wRwG+WAyVG3s+njUeZiyE6d+AhCGWlbenoY25d79HXUsn3z9lDD+cO96yWkRkAFprzUDw2WNQ9WXPx7OPgKkXmqMF6UVBuen3N1fxk/98wa7aVgAunFHAj8+cQFaymhQFksJApGqqhA0vmCHAvz0IzIU/k74OM66EwmNCPhWwr06vj4UPr+DDrdVMHpbKc9cdpwZDIpHKMGDHR7DyMVj3XM9CZDC3JU+50FyMnBzYV+/N7R7ufO1LHltejGFAktvBN08czbUnjNTW5ABRGIgktSWw8SXY8FJXAPB/62ww4niYepH5RBzEKYKBZBgGP3l2Dc98upMkt4Ml1x/H2NwUq8sSkUBorTNfkHyxGIrfp/v3kc0OhbPhiLNgwlkwZHjAbvKzkhpufXE9X+yqByArOY4fnDaWBUcX6kXGICkMhDNPB+xaAVvehM1v9mz78Rs6HSafb76lhd/KfP82QrsNHlp4NCdPCN2+ZREJofpSWPus+Va+au/P5U6G0aeYb0VzwBU/qJvy+QyWri3nzte+pKTaXMdQMCSBq44byUUzC0iJdw3q+mOVwkA48XnNlf8ly6F4GWx7Fzoaez5vs8Pw48y0PWE+pIfnsb+GYfCnt7bwxzc3AXDr1yax8NgR1hYlIqFRtwM2vtw1gvkhGL6ezznjzd9h/nCQc8SApzI7PD6eWrGDP721mepm83jmlDgnC2YVsvDYERQM0TkH/aEwYCVPB5R9bj5hSj405+L8bYH9ErNgzKkw5jTzyZOUZU2tfeT1GfzqpfU8+mExADecNpYbThtnbVEiYo3mKtj2Dmx923xrLN/784mZUDALCmeZa5zyjwR3//6It3V6+c/KUh56fxtbu7Yi2m1w3JgszjtqGHMn5mldQR8oDISKpwMqN0DZKihfbQ6lVazbewEOgDsFio6B4ceaf/zzpoWkM2AglNW1csMzq1ixvQaAX549kSuP0yFEIoK5+LByY08wKP4APK17X8buhLypZjgYOh3yJpu7opyH70Lo8xm8u7mSB5dt44Mt1d0fT3A5mDcpl/lT8zluTCaJbgWDA1EYCDSfD+p3QOUmcwtO5Zewew3sWQ/ejv0vn5hp/uEvOtb8N3cyOCLrh9XrM1j86U5++8pG6ls7SXI7+N0FUzlrqjVHo4pIBPC0Q/kX5rqonR/Djo+haff+l7O7IHuCGQxyJ0PuJMgaBylDD/pCqaS6mSWfl/Hc57so7lpXAOB22DlmVAYnj8/hpPHZjMxKwmbxjqtwoTAwED4vNJRBXYm5wr+uxOz7X/Wl2aRj37TrF59mpt2h0yB/uvn/jFGWb/8bKK/P4M0NFfzxjU1s3G2ubZhWmM6fFkxneGaSxdWJSEQxDKjfCTtXmM3Tdq81X0i11x/48q5EyBhtnp2QOcZ8yxhlrqVKzgW7A8MwWLWzjudXlfHmhoruXgV+WclxzBiezszhGRw1fAiT8lNj9rwUhYF9dbSY6bSxouffxnJoqugKADugfhf4Og9+HQ63+YOZNQ6yx5uLZIZOhyEjIvYPf2/bq5pZuqacp1bs6H5ypSW4+P6pY7liznBt8RGRwPAHBH8wqFgDezZAbTH4PAf/OrsLUvMhrdAMB2mFGCl57Pam8kmlg/+Wwn932qjz7j394LDbGJWVxIShqUzIS2F8bgojs5MoGJJAnDO6Q0J0hgGfFzqae97aG6Ctzuyi1f3W6/2WGmipNv/g77uA72DsLkgrMPfQpg+HjJHm3Fb2ePP9CBvqPxjDMNhV28qa0npWbK9h+dZqvqzo2eGQnujiG7OK+OaJo0hP1OliIhIC3k7zhVn1lp63qs3mSG1DKRjevl2NM4FmZwaVRio72pOp8CRRTxINhv/fRBpIooFE4pIzScvIJmNIBhlpaeSkxZOTEk9uahw5qfFkJ8fhdkbuCyFrw0BtMdRsN+eOvO3mIjtvu7morvv/XW/ejq7/t5n/72zt+mPftM+/zT09tAfKlWgOM6Xk7fPvULPl5pDhXfNVkZUUDcOgw+uj02vQ6fHR2umlsc1DY1snDW2dNLR62N3QRnldK6V1bZTXt7KjuoXG9r0TuMNu49jRmZw9NZ+vTc+P2WE1EQlDPq85mlu/C+p2miML9TvNUd7mPdC0B5orB/V3wmfYaMVNC3G0GnE0E08rcXTYEuh0xONxJOJzJoArAYcrHoc7DpszDrsrHpszDofLfHO643C443G543G5E0gtGEdWgTXt2i0NA6v/cRPTtj0QtOv3YKeVeFpsCTSSQoMt2XwjmUZbCg3430+hnhSqbelU2TJoIWGv4fwDPUT7fmTfixj7XOJAj/LhrmPfS+x/Gwe4zn0u5DPMc8L9IWAg3A47Y3KSmTliCLNGZnDc6CyG6IxxEYlk7U3maHBzZVdA2NMzatxW3zWaXIfRVo+vpQ6jrQ5HRyO2A/7mDYzlBdcy59q7gnb9h2LpWHaVLZMNvkI6cNGOiw7D2fP/A7zfjpMOw3y/HRfNRjzNJNBMPC1d6ayZeJqNeFqIpx0XMJA5+rbDXyQKuBw2UuNdpMQ7Sen6NycljqHpCeSnJ5CfFk/BkERGZSdpHYCIRJe4ZPMtc/QhL2YDusc+fV5zRKGjBTqbu/5twehoprmxnpbmRtpbGmlvbaSztQlPewuejja8nT2j3TZvB3ZvB3af/60Th68Dp9GJkZIX7Ht98Ptp5cjArtoWyur2/sN7oPV1+37owGvwbIe9zP7Xs/+FDndbtgOEi76sCRzI9Rz4PvTvftps4HLYcTvt3f+6Heb/HfbIX8woIiKDF50LCEVERKTPNPYrIiIS4xQGREREYpzCgIiISIxTGBAREYlxCgMiIiIxTmFAREQkxikMiIiIxDiFARERkRinMCAiIhLjFAZERERinMKAiIhIjFMYEBERiXEKAyIiIjFOYUBERCTGKQyIiIjEOIUBERGRGKcwICIiEuMUBkRERGKcwoCIiEiMUxgQERGJcQoDIiIiMU5hQEREJMYpDIiIiMQ4hQEREZEY9/8Bvhbv/JJxRBYAAAAASUVORK5CYII=", | |
"text/plain": [ | |
"<Figure size 640x480 with 1 Axes>" | |
] | |
}, | |
"metadata": {}, | |
"output_type": "display_data" | |
} | |
], | |
"source": [ | |
"\n", | |
"xx = np.arange(-50, 151)\n", | |
"fig, ax = plt.subplots()\n", | |
"ax.plot(xx, p(xx), label=r\"$p(x)$\")\n", | |
"ax.plot(xx, q(xx), label=r\"$q(x)$\")\n", | |
"ax.set_axis_off()\n", | |
"plt.legend()\n", | |
"plt.show()\n" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"It is necessary to scale $q$ so that it is always greater or equal to $p$ across the shared domain. \n" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 4, | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
"Scaling factor: 3.9146256816520806\n" | |
] | |
} | |
], | |
"source": [ | |
"\n", | |
"# Scaling factor. \n", | |
"k = max(p(xx) / q(xx))\n", | |
"\n", | |
"print(f\"Scaling factor: {k}\")\n" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"Visualizing target distribution with scaled proposal distribution." | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 5, | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAgMAAAGFCAYAAABg2vAPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABXdUlEQVR4nO3dd3hb5d3G8a/kve3E287eew9G2DsJhB1oC6WMDmhL+75dvFBaCnTQBWkLFEIpM0ApM+ydhIRAEidx9t527HhvSzrvH49lZ8dD8tG4P9eVSweNc34WHree6bAsy0JERETCltPuAkRERMReCgMiIiJhTmFAREQkzCkMiIiIhDmFARERkTCnMCAiIhLmFAZERETCnMKAiIhImFMYEBERCXMKAyIiImFOYUBERCTMKQyIiIiEOYUBERGRMKcwICIiEuYi7S5ARERCg8fjoampye4ywkpUVBQRERFdPo/CgIiIdFlTUxPbtm3D4/HYXUrYSU1NJTs7G4fD0elzKAyIiEiXWJbFvn37iIiIoFevXjid6oHuDpZlUVdXx/79+wHIycnp9LkUBkREpEtcLhd1dXXk5uYSHx9vdzlhJS4uDoD9+/eTmZnZ6S4DxTcREekSt9sNQHR0tM2VhCdvAGtubu70ORQGRETEJ7rSZy2d54v3XWFAREQkzCkMiIiIhDmFARERkTCnMCAiInIcBw4cIDMzk+3bt7f7NbNnz+ZPf/qT/4ryMU0tFAk3Zdtg+0LYtQQq90Bdqbk/Ph1S8qDXVOh7KvToZ2+dIgHivvvu45JLLqFv377tfs2dd97Jaaedxk033URKSor/ivMRhQGRcOBqhDWvwpePw+6lx3/uimfMbf5kmHQTjJgFkTH+rlAkINXV1TF37lzefffdDr1u5MiRDBgwgGeeeYZbb73VT9X5jroJREKZZcGaV2DORHjlFhMEnJHm0/+0/4FZD8O1L5l/sx429/Waap6ze6l5zZyJ5hyWZfdXI+JzRUVFOBwOHnzwQcaNG0dsbCwjRoxg4cKFALz11lvExMQwderUQ173/PPPExcXx759+1rvu+GGGxg9ejSVlZUAzJw5k3nz5nXfF9MFahkQCVWVe+DV78C2z8x/J+XApBth3HWQlHX811YXw4qn4Mu5ULkTXvom9DsNZj1iuhJEjsOyLOqb3bZcOy4qokPz7gsKCgB44okneOihh8jKyuJHP/oRX/va19i2bRsLFixgwoQJR7xu9uzZ/O53v+P+++9nzpw53H333XzwwQcsWbKktVtg8uTJ3HfffTQ2NhITE9itawoDIqFow9vw6nehvhwi4+CUH5p/0e1cKjYpC077CUz9Hix6CBY9aELFI6eaFoQhF/i3fglq9c1uhv+yY83qvrL2nvOJj27/n7aVK1cSFRXFa6+91jom4N5772XixIns2bOHHTt2kJube8TrHA4H9913H1dccQXZ2dnMmTOHBQsWkJfXFpZzc3NpamqiqKiIPn36dPlr8yd1E4iEEsuCBX+G52ebIJAzFr67CM78RfuDwMGiE8xrv7sIcsZAfRk8f7W5hroNJAQUFBRw2WWXHTI4MDk5ufW4vr6e2NjYo752xowZDB8+nHvuuYdXXnmFESNGHPK4d9+Auro63xfuY2oZEAkVHje883NY+k/z35NvgfPu9c3gv54D4Mb34b07zfk//DVUF8EFvwVn1/dSl9ASFxXB2nvOt+3aHVFQUMD1119/yH2LFy8mPT2dvLw80tPTKS8vP+pr33nnHdavX4/b7SYr68iut7KyMgAyMjI6VJMdFAZEQoHHA6/dBiufM/99/m/hpO/59hqRMXDRA5DWF969A5Y+Co3VcMnfQVvWykEcDkeHmurtUl9fz6ZNm1o3WgLweDz89a9/5frrr8fpdDJu3DieeeaZI167fPlyrrrqKubOncuTTz7JXXfdxUsvvXTIcwoLC8nPzyc9Pd3vX0tX6SdYJNhZFrz9ExMEHBFw+VzfB4GDnXSruYYjwlzz7Z+qy0CC0urVq3E4HDzzzDMsXryYdevWcfXVV1NRUcGdd94JwPnnn8+aNWsOaR3Yvn0706dP54477uCaa67hnnvu4eWXX2b58uWHnH/BggWcd9553fo1dZbCgEiw+/Aes34ADrj0URh1hf+vOeoKuPQRc80vH4OPfuP/a4r4WEFBAUOHDuWOO+7g8ssvZ+LEibjdbj799FNSU1MBGDVqFOPHj+fFF18ETNP/BRdcwCWXXMLPf/5zAKZMmcKFF17IHXfc0XruhoYGXn31VW6++eZu/7o6w2FZivQiQWvFs/BaSyvAzAdhwje79/rLnoQ3fmiOZz0MY6/t3utLQGhoaGDbtm3069fvmIPtAtGtt95KeXk5zz333HGfN3/+fH7yk59QWFiIs51dYg8//DCvvPIK7733ni9KPS5fvP9qGRAJVjuXwJu3m+PTftr9QQDMNU/7iTl+44ew84vur0GkkwoKChg9evQJnzd9+nRuueUW9uzZ0+5zR0VFMWfOnK6U163UMiASjKr2waPToLYEhl0MV/7bvkF8Hg+8+A1Y/yYkZMB3FkJStj21iC2CsWXAsixSUlKYN28eF110kd3ldIlaBkTCkccN/73ZBIGskabv3s7R/E4nXPZPU0ttianNY8/qcyLt5XA4qKqqCvog4CsKAyLB5rM/wvYFEJUAVz5pFgayW3QCXPEviIo3KxUuCJ6tW0VEYUAkuOxcAp/+zhzP+DOkD7K3noNlDIbpfzbHn/zW1CoiQUFhQCRYNNWZ/QYsD4yeDWNm213RkcZeA6OvNjW++j1Ts4gEPIUBkWDx0W+gbCsk5cKFv7e7mmO78A9mh8SyLfDRvXZXIyLtoDAgEgx2LIYlD5vjix+CuFRbyzmuuFS4uGVK1ZJ/mNpFJKApDIgEOldTy8I+Foz7Ogw61+6KTmzQuTD264Bl1kJwNdldkYgch8KASKBb/Dco3WDm8J8XRM3u598L8elQsh6W/N3uakTkOBQGRAJZ+Q749A/m+Lx7IS7N3no6Ii6tLbx8+geo2GlvPSJyTAoDIoHsnV+Aqx76nGpG6QebMbOhzynQXAdv/9zuakTkGBQGRALV1k9gw3yzVfD0P4LDYXdFHedwwPQ/ma9hw3zY+qndFYnIUSgMiAQijxve/T9zPOlGyBxmbz1dkTkMJn7LHL/7f1qqWALO//7v/zJr1iy/X+fAgQNkZmayffv2dr9m9uzZ/OlP/l/RU2FAJBAVPAvFhRCbAmf8wu5quu6MX0BMChSvhoLjbxcr0t0KCgoYO3as369z3333cckll9C3b992v+bOO+/kvvvuo7Ky0n+FoTAgEngaa9oW6zn9ZxDfw956fCGhJ5z+U3P80W+gqdbeekQOsnLlSr+Hgbq6OubOncuNN97YodeNHDmSAQMG8Mwzz/ipMkNhQCTQfPEw1BRDWj+YdLPd1fjO5FvM11RT3LaAkojNdu/eTWlpaWsYqKioYObMmZx66qkUFRW1+zzr16/nzDPPJDY2lsGDB/PWW2/hcDgoKCgA4K233iImJoapU6ce8rrnn3+euLg49u3b13rfDTfcwOjRo1tbA2bOnMm8efO69oWegMKASCCpK4NFLav3nXUnREbbW48vRUbDmS3jIBY9ZL5WCU2WZVp/7PhnWR0qtaCggNTUVPr27cvq1auZNGkSeXl5fPzxx2RnZ7frHOvXr2fKlClMmzaNNWvW8Pvf/57rrruOqKgohg8fDsCCBQuYMGHCEa+dPXs2gwcP5v777wfg7rvv5oMPPuDtt98mJSUFgMmTJ7N06VIaGxs79LV1RKTfziwiHbfor9BYCVmjYMRldlfjeyMvN19jcSEsehDO/bXdFYk/NNfB/bn2XPuOvR3a1rugoIAxY8bw3HPPcdttt/H73/+em28+skXuwQcfZN68eSxefOTy2rfeeiuXX34599xzD0Brs/6mTZuIjjaBfseOHeTmHvmeOBwO7rvvPq644gqys7OZM2cOCxYsIC8vr/U5ubm5NDU1UVRURJ8+fdr9tXWEwoBIoKjaB188ao7PvgucIdhw53TCWXfB81fDF4/AlO9Aco7dVUkYKygoYNWqVdx2223Mnz+fk0466ajPKykpYdOmTUfcv2PHDj766CNWrlx5yP3R0dGMGTOm9b/r6+uJjY096rlnzJjB8OHDueeee3jvvfcYMWLEIY/HxcUBZtyBvygMiASKhX8BVwP0mgqDzrO7Gv8ZfD70mgK7vjCtBIG8A6N0TlS8+YRu17U7oKCggMsuu4znnnuOioqKYz7v3nvv5d57j1wOvKCg4JDuAK/CwkKuv/761v9OT0+nvLz8qOd+5513WL9+PW63m6ysrCMeLyszXWoZGRnt+ZI6JQQ/eogEoeoiWPakOT7zF8G5wFB7ORxt0yWXPWm+dgktDodpqrfjXwd+dqqrq9m6dSu33norf/vb35g9ezZr1qzp0JfqdDpxu9243W3rZ7zzzjsUFhYe0jIwbtw41q5de8Trly9fzlVXXcXcuXM5++yzueuuu454TmFhIfn5+aSnp3eoto5QGBAJBIseBHej+cTc73S7q/G//mdA/mTTErLoIburkTC1cuVKIiIiGD58ON/61re48cYbmTlzJqWlpe0+x4QJE4iKiuKOO+5g69atvPzyy9x6660Ah4SB888/nzVr1hzSOrB9+3amT5/OHXfcwTXXXMM999zDyy+/zPLlyw+5xoIFCzjvPP+2FioMiNituhi+esIcn/6z0G4V8HI44IyfmeOvnoCa/fbWI2GpoKCAoUOHEhMTA8ADDzzAkCFDuOyyy2hqat+227m5uTz++OO8+OKLjBkzhhdeeIGbb76Z7OxsMjMzW583atQoxo8fz4svvgiYpv8LLriASy65hJ//3OzbMWXKFC688ELuuOOO1tc1NDTw6quvHnVQoy85LKuD8zBExLfeuxM+nwN5E+GmD8IjDICZAvb42bBnGZz8AzjvN3ZXJJ3U0NDAtm3b6Nev3zEHyYWTH//4x6xdu5Z33nnnkPvnz5/PT37yEwoLC3G2c4Dwww8/zCuvvMJ77713zOf44v1Xy4CInRoq4asnzfFpPwmfIADmaz2tZVXCr/5l3guRELBq1apDugi8pk+fzi233MKePXvafa6oqCjmzJnjy/KOSmFAxE7LnoSmasgYGtozCI5l0HmQPsS8B8v+bXc1Ij6xevVqRo8efdTHbr/9dnr16tXuc910000MGTLEV6Udk8KAiF1cTbDkEXN88vdDc12BE3E6zdcOZoliV/v6aUUCWXFxMV/72tfsLqNDwvC3j0iAKPwPVO+FxGwYdaXd1dhn9FWQmGXei8KX7a5GJCwpDIjYwbLMoEGAKd+GyBh767FTZIx5D8C8JxrTLNLtFAZE7LD5A9i/FqITYeK37K7GfhO/BVEJsH8NbP7Q7mpEwo7CgIgdPm9ZaGf89RCXamspASEuDSa0LN36uRYhEuluCgMi3W1vAWz7DBwRMPW7dlcTOKZ+17wn2z4175EEHS1bYw+Px9Plc2ijIpHu5v3kO/JySG3/FKOQl9obRl4Gq18yYweumGt3RdJOUVFROBwOSkpKyMjIwBFO62XYyLIsmpqaKCkpwel0tm6X3BlagVCkO1XsggfHgOWG7yyE7FF2VxRY9q2CR6eZFoLbV0FKvt0VSTvV1NSwe/dutQ7YID4+npycnC6FAbUMiHSnZf8yQaDvNAWBo8kZbd6b7QvMqoRnH7mDmwSmxMREBg0aRHNzs92lhJWIiAgiIyO73BqjMCDSXVyNbavsTfbvpiNBbdJNJgws/zec/tPwnnYZZCIiIoiIiLC7DOkEDSAU6S5rX4O6UkjKhSHT7a4mcA2dDkk5UFsCa1+3uxqRsKAwINJdlj5mbifeABFqlDumiCiYcIM5/vIxe2sRCRMKAyLdYd9K2L0UnFFmbQE5vgnXgzMSdn1hBhWKiF8pDIh0hy8fN7fDL4akLHtrCQZJ2TDsYnPsfe9ExG8UBkT8rb4cVr1kjidp4GC7eQdZrn4J6itsLUUk1CkMiPhbwXPgqoeskdB7qt3VBI/eJ0HmCGiuM++hiPiNwoCIP3k8bc3ck24CrczWfg4HTL7JHH/5uHkvRcQvFAZE/Gnrx1C2FWKSYdSVdlcTfEZdZd67si2w7RO7qxEJWQoDIv701RPmduy1EJNoby3BKCYRxlxjjr/UXgUi/qIwIOIv1cWw8R1z7J03Lx03seW92/gO1Oy3txaREKUwIOIvK58HjwvyJ0PmULurCV6ZwyB/knkvVz5vdzUiIUlhQMQfLAtWPG2Ox3/D3lpCwbiW93D50+a9FRGfUhgQ8Yedi+HAZohOhBGX2V1N8Bt5GUQlwIFNsHOJ3dWIhByFARF/WN7SKjDiUg0c9IWYJBh5qTn2triIiM8oDIj4WkMlrHnFHGsfAt/xvpdrXoGGKntrEQkxCgMivlb4sllxMGMo5E+0u5rQkT8J0oeYFQkLX7a7GpGQojAg4mvLnzK3476hFQd9yeFoG4zpfY9FxCcUBkR8qWg17F1htioeM9vuakLP6Nnmvd27HIoK7a5GJGQoDIj4knfg4NCLICHd3lpCUWIGDLnQHGsgoYjPKAyI+EpzA6x6wRyPv87eWkKZdyDhqhfMey4iXaYwIOIrG9+GhgpIzof+Z9pdTegacCYk50F9edtyzyLSJQoDIr6ycp65HXM1OCPsrSWUOSNg9NXm2Puei0iXKAyI+EJNCWx63xyP1sBBv/MOztz8vnnvRaRLFAZEfKHwP2C5IW8CZAy2u5rQlzEEcsebzYu05oBIlykMiPiCdze9MdfYW0c48b7X2slQpMsUBkS6qngt7FuJ5Yxi1qdZnPzbD/nhvBXMX7UPSzvs+c/Iy8EZCfsKYP86u6sRCWoKAyJd1fLJ9CP3WAoORLC3soHXCvZy63PLeW7pTpuLC2EJPWHQ+eZYrQMiXaIwINIVHjdNK8yI9hebT2Vyvx78+1uTuWZyLwB+/cZaCvdU2llhaPMOJFz1Injc9tYiEsQUBkS6wNr6CdH1+ym3Eokceh5P3ziZ0wdncP+lozhnWCZNLg+3PbecqoZmu0sNTYPPh7g0qN4H2z61uxqRoKUwINIF+xc8CcBb1sncPWs8MZFmfQGHw8EfrxxDXmoc2w/U8einW2ysMoRFxpixA6A1B0S6QGFApJOshipSd7wLQNOIq8hMjj3k8dT4aP5v+jAAnl+6i4ZmNWP7hXdWwbo3oLHa3lpEgpTCgEgnrfnwWWJoZKuVy8yLZh71OecNzyInJZay2ibmr9rXzRWGibwJ0HMgNNeZQCAiHaYwINJJjlWmWXpXr4tJT4o96nMiI5x8fWofAJ5avL27SgsvDkfbQEJ1FYh0isKASCeUF+9kWMNKAAac9c3jPvfqSb2IjnCycnclK3aWd0N1YWjUleZ2+wKoLra3FpEgpDAg0gnbPnkKp8NibeQw8vsPO+5z0xNjmDEmB4BnlmjdAb9I6wv5k8DywJpX7K5GJOgoDIh0Qurm1wAo6Xtxu54/e1JvAD5YV4zL7fFbXWHN2zqw+iV76xAJQgoDIh1UsWsd/Zs34rKc9Dv96+16zYQ+aaTFR1FZ38yyHeoq8IsRl4LDCXu+grJtdlcjElQUBkQ6aNen/wZgRfQ4evfq3a7XRDgdnDkkEzCtA+IHiZnQ7zRzXPgfe2sRCTIKAyIdYVmkb38dgIoBszr00nOGZwHw4br9vq5KvFq7Cv4D2iRKpN0UBkQ6oGHnMnJce6i3ohkw7aoOvXbaoHSiIhxsLa1lS0mNnyoMc0NnQEQ0lKyH4jV2VyMSNBQGRDrgwOJnAVgQMZl+uVkdem1SbBRT+/cE4EN1FfhHXCoMOs8cq6tApN0UBkTay+MmZYvpItiZexEOh6PDpzhnmAkQH6irwH9GXWFuV7+srgKRdlIYEGmv7QtJbC6lwkogddSFnTrFWUPNIMJlO8qpaXT5sjrxGnwBRCdC5U7YtdTuakSCgsKASDu5Vr4AwFvuKUwemN2pc/TqEU9+Whxuj8VyTTH0j6g4M3YA1FUg0k4KAyLt4WqEtaaLYGH8WfTqEdfpU03u2wOApdvKfFKaHIV3VsGaV8CtFhiRE1EYEGmPzR8S2VzNPqsHcQNO6dR4Aa/J/VrCwHaFAb/pfzrE94TaEtj2qd3ViAQ8hQGR9mhZ736+ewpTB2R06VTeMFCwq4KGZneXS5OjiIgyKxKCWXNARI5LYUDkRJrrsTbMB2C+e2rr9MDO6peeQHpiNE0uD6t2V/qiQjmakS2zCta9Ac319tYiEuAUBkROZPOHOJpq2W2lU5w0kl494rt0OofD0do68KW6Cvyn1xRI6QVN1bDpPburEQloCgMiJ9LSRfC2ezJjeqf55JTeQYRfaBCh/zidMPIyc6yuApHjUhgQOZ7metj4DmC6CEbnp/rktJNaWgaW7yjXlsb+5J1VsPFdaFCXjMixKAyIHM/mD6CphiIyKLAGMCY/xSenHZqdTFJsJDWNLtYXVfvknHIUWSMhfQi4G2H9fLurEQlYCgMix9PSRfCaazLgYKSPwkCE08HYXqmAmVUgfuJwwMjLzXHL/0sROZLCgMixNNXBhrYugv4ZCSTHRvns9KNbgsWq3RU+O6ccxYhZ5nbLR1CvVR9FjkZhQORYNr8PzbVUxuSyyurPGB+NF/Dyjj/Q9EI/yxgCmSPA41JXgcgxKAyIHEtLs/KimGmAo/WTvK94w8Wm/TXUN2nxIb/yLkCkrgKRo1IYEDmapjozAh14vnYcgM9mEnhlJceQkRSD22Oxdp9aB/zK21Ww9ROo03ROkcMpDIgczab3oLkOV3JvFtT2ItLpYERusk8v4XA4WmcnrNylMOBX6YMga1RLV8GbdlcjEnAUBkSOpqU5eWfO+YCDwVlJxEZF+PwybeMGKnx+bjmMt3VAXQUiR1AYEDlcU21rF8HnMdMAGJnn21YBr1HeGQV71DLgd95xA1s/hdoD9tYiEmAUBkQOt/FdcNVDWj8+rc4FYFiOf8KAdxDh1pJaqhqa/XINadFzAGSPBssN69+wuxqRgKIwIHI4bzPyiEtZX2xWBxya7Z8w0CMhmvy0OAAKNcXQ/zSrQOSoFAZEDtZY07rDXe2gmewqM1vfDs1O8tslva0D6iroBt5xA9s+g9pSW0sRCSQKAyIH2/gOuBqgR3/WW30AyE6OJS0h2m+XHN4yS2Hdviq/XUNa9OgPOWPB8sC61+2uRiRgKAyIHGztq+Z2xKWsK6oBYGiO/1oFAIa3jEdYu1dhoFuoq0DkCAoDIl6N1bDpfXM84lLWF5k/zv4aL+DlHZy4tbSWhmatROh33q6C7QuhZr+tpYgECoUBEa+N75ougp4DIWsk6/eZwYPD/NwykJUcQ1p8FG6PxabiGr9eS4C0vpA7Xl0FIgdRGBDxOmgWgQWsL/LvTAIvh8OhcQPdrbWr4FVbyxAJFAoDIgANVYd0Eewur6em0UVUhIP+GQl+v/ywlsCxVmGgexzcVVBdbGspIoFAYUAEzCwCdyP0HASZw1s/oQ/MTCIqwv8/Jt5xAwoD3SS1N+RNBCx1FYigMCBiHNRFgMPR2kUwzI/rCxzMGwbW7avCsqxuuWbY06wCkVYKAyINVbD5A3Pc8gdiQ0sYGNJNYWBgZiJREQ6qG1zsqajvlmuGveGXmNsdn0PVPntrEbGZwoDIhrfB3QTpQyBzGACb9pswMLibwkB0pJOBmeZa61pmMYifpfaC/Mmoq0BEYUDkiC4Cl9vDttJaAAZmJHZbGd4pjFp8qBupq0AEUBiQcFdfAVs+NMctI8x3lNXR7LaIi4ogLzWu20oZnqPphd3O21WwczFU7bW3FhEbKQxIePN2EWQMa+siaFn4Z2BmIk6no9tK8Y5P2FisboJuk5IHvaaa47Wv2VuLiI0UBiS8tXYRzGq9a0tJWxjoToOzTBjYfkDLEncrdRWIKAxIGKuvgC0fmePhs1rv3tTyyby7w0BmUgzJsZF4LNhaUtut1w5rwy8GHLDrC6jcbXc1IrZQGJDwteEt8DRD5nDIHNp696b9pmVgUDeHAYfD0do64J3NIN0gORd6n2SOtTyxhCmFAQlfB88iaOHxWLZ1E0DbVEaNG+hm3u8B7xbWImFGYUDCU335UbsI9lTU09DsITrCSe8e8d1e1uCWALJRuxd2L29Xwe4voWKX3dWIdDuFAQlP6+eDxwVZIyFjcOvdm1u6CPpnJBDZDXsSHK61m0AtA90rKRv6nGyONatAwpDCgISno8wigLa++gE2dBEADGoJAzvK6jSjoLupq0DCmMKAhJ+6Mtj6iTkefukhD3nXGOjuwYNe6YnRpMVHYVltrRTSTYYd3FWw0+5qRLqVwoCEn9YuglGQPvCQhzbbOHgQzIyCQZpRYI+kLOhzijlWV4GEGYUBCT/H6CKwLKv107hdYQBgcJYGEdrG+z2hKYYSZhQGJLwc3EUw4tAugtKaJqobXDgc0LdnQvfX1sI7iHBjkVoGup23q2DPV+oqkLCiMCDhZd0bYLkhezT0HHDIQ1tbugjy0+KIjYqwozoABrVsZbxR3QTdT10FEqYUBiS8HGWhIa+tLdsW90+3r4sA2roJdpXVU9fksrWWsNTaVaC9CiR8KAxI+KgthW2fmePDxgtAW8tA/wz7uggAeibG0DMhGtCMAlsMuxgcTtizDMp32F2NSLdQGJDw4e0iyBkDPfof8bB3c6D+6faGAYBBGkRoH3UVSBhSGJDw4V1M5ihdBADbvN0EGfZ2E4BWIrSdugokzCgMSHg4uIvgoL0IvJrdHnaW1QH2dxNA20qE2rDIJt6ugr3L1VUgYUFhQMLDutfB8kDuOOjR74iHd5bV4fJYxEdHkJ0ca0OBh9KGRTZLzDyoq+BVW0sR6Q4KAxIejjOLANrGC/RLT8DhcHRXVcfk7SbYU1FPbaNmFNjC+72irgIJAwoDEvpq9sP2heb4KF0EcPBMAvvHCwCkJUSTnhgDwCbNKLBHa1fBCijfbnc1In6lMCChr7WLYDyk9TnqU1oHDwbATAKvtmWJNW7AFokZ0PdUc6zliSXEKQxI6PP+Ij9GFwEcNK0wAAYPemlGQQDwtiRp3ICEOIUBCW3VxbBjkTk+ykJDXltLW7oJbF598GBaayAAHNxVULbN7mpE/EZhQEKbt4sgbyKk9j7qUyrrmymtaQKgn1oG5GCJGdB3mjnWAkQSwhQGJLS1o4vAO14gMymGxJjIbiiqfQa3bFi0t7KB6oZmm6sJY1qASMKAwoCEruqiti6C4Zcc82nbS9umFQaSlPgospI1o8B23q6CfQXqKpCQpTAgoWvt64AF+ZMhtdcxn7YtQMMAtG1nvFnjBuyTkH5QV8GrtpYi4i8KAxK6WhcamnXcp+04YMJAn54BGAY0vTAwaAEiCXEKAxKaqvbBzsXm+DhdBADbDpg9Cfqlx/u7qg7ztgyom8Bmw2aCIwL2rYSyrXZXI+JzCgMSmta+BljQawqk5B/3qd4xA30DsZugpWVgs8KAvRLSoV9LV4EWIJIQpDAgoWnNf83tcWYRAJTXNlFZb0bq9+kRgGGgZcOiPRX11GiPAnupq0BCmMKAhJ7K3bDrC8BxzL0IvLa1jBfITo4lLjrC/7V1UGp8NBlJZkaBWgdsNrSlq6BoFRzYYnc1Ij6lMCChx9uM2+dkSM457lO9gwf7BuB4AS9v64AWH7JZQk/od5o51qwCCTEKAxJ62tlFALCt1Dt4MPC6CLy8YUAtAwGgdQGiV+2sQsTnFAYktJRvhz3LzCIxJ5hFAAcNHgzAaYVeg7I0oyBgqKtAQpTCgIQW7+CuvqdCYuYJn779QODOJPDytgxorYEAkNAT+p9ujjWQUEKIwoCElkJvF8FlJ3yqZVmtqw8GQ8vA7vJ66po0o8B22tZYQpDCgISOA1tM860jwqwnfwLldc1UN5g/rn16Bu4Awh4J0fRMiAZgy/5am6uR1gWIilarq0BChsKAhA5vq0D/M0xz7gl4WwVyU2KJjQq8aYUH07LEASS+h7oKJOQoDEjo8M4iGHniLgII7JUHD6dliQNM6wJEr9pahoivKAxIaNi/HvavBWcUDJ3erpcEw+BBr7ZlidUyEBCGzgBnJBSvhtLNdlcj0mUKAxIavM21A86CuLR2vaRt8GDgjhfwGuhdeEgtA4Ehvgf0a+kqWKuuAgl+CgMS/Cyrw10EADtadisM5JkEXoNbZhTsLKujvsltczUCqKtAQorCgAS/4jVQuhEiYmDIRe16iWVZrWMGAnn1Qa+eCdGkxUdhWbClRK0DAWHo9JaugkIo2WB3NSJdojAgwc/bKjDoXIhNbtdLDtQ2Ud3owuGAXj0Cv5vA4XC0DiLUssQBIr4HDDjbHK/+j721iHSRwoAEN8s6aKGhE+9F4LW9dVphXMBPK/QamOUdN6BBhAFj1BXmtvA/5ntRJEgpDEhw21cA5dsgMg4GX9Dul7UOHgzg3QoPN7h1WWK1DASMIReZ772yrbB3hd3ViHSawoAEN+8sgsHnQ0xiu18WTIMHvbzLEqubIIDEJMKQlhBa+LK9tYh0gcKABC/LagsDHegiANh2IHgGD3p5NyzacaCWhmbNKAgYI1u6Cta8Ah6PvbWIdJLCgASv3V9BxU6ISoBB53XopcGwdfHhMpJiSImLwmPB1hLtURAwBp0LMSlQtQd2Lra7GpFOURiQ4LX6JXM7dDpEt7/v/+BphcE0ZsDMKNAgwoATGQPDZpjjQs0qkOCkMCDBye1qm1I46soOvbSkppHaJjfOIJlWeLC2ZYk1biCgjLzc3K55FdzNtpYi0hkKAxKctn0KtSUQ1wMGnNmhl3oHD+amxhETGRzTCr0Gejcs0oyCwNLvdEjIgPoy2PqJ3dWIdJjCgAQn7yIvIy6FiKgOvXRbEK08eLjB3q2M1U0QWCIiYfgsc6wFiCQIKQxI8Gmuh3VvmOMOdhFA2+DBPkGwQdHhvKsQ7jhQR6NLMwoCincBovXzzfeoSBBRGJDgs/FdaKqGlN7Qa0qHX966dXEQzSTwykqOISkmErfHYntpnd3lyMHyJ0NKL/O9uek9u6sR6RCFAQk+3lkEoy4HZ8e/hbe1/BENxm4Ch8OhZYkDldPZtmumugokyCgMSHCpr2j71NWJLgLLstjhbRkIwjAAMLilq0DLEgcg7wJEG9+Fhip7axHpAIUBCS7r3gB3E2QOh6wRHX55SXUjdd5phWnBN2YADp5eqJaBgJM9CnoOAnejGTsgEiQUBiS4tHYRXNGpl3tnEuSlxREdGZzf/gO9Cw+pZSDwOByH7mQoEiSC87ehhKfqItj2mTn2LvLSQcE8eNDLu2HRttJamlxaCz/geLsKtnwMtaX21iLSTgoDEjwKXwYsM4MgrW+nThHMgwe9clNiSYiOwOVpG/8gASR9IOSMAcsNa1+zuxqRdlEYkODR2kXQ8YGDXjtCoGXAzChoWYlQyxIHJm/rgPd7ViTAKQxIcCjdDHtXgCOibaW3Tgjm1QcPNkjjBgLbyMsBh9nFsHyH3dWInJDCgAQH72CsAWdCYkanTmGmFZpugmBcffBg3jCgZYkDVEoe9Jtmjle/aG8tIu2gMCCBz7J80kVQXNVIfbObCKcj6HYrPNzglm6CzWoZCFyjZ5vbVS+a72GRAKYwIIFvXwEc2AyRsTB0eqdP4+0iyE+LIyoiuL/1vdMLt5bW0OzWjIKANGym+Z4t3Wi+h0UCWHD/RpTw4F3adciFEJPU6dOEwuBBr7zUOBKiI2h2W60hRwJMbDIMucgcr1JXgQQ2hQEJbB53y5RCutRFALCtNQwEdxcBgNPpYEi2CUbrizRuIGCNvtrcrv4PuF321iJyHAoDEti2fQrV+yA2FQae06VTebcuDtY9CQ43NCcZgPX7tAZ+wBp4NsT3hNr9sPUTu6sROSaFAQlsK+eZ21FXQGRMl07l3fI3ZMKAWgYCX0RU22qZq16wtxaR41AYkMDVWG02JgIYc02XTuXxWK1LEfcLgTEDAEOzTcvABoWBwObtKlj/JjRq9ocEJoUBCVzr3oDmOug5EPImdOlU+6oaaHR5iHQ6yE+L81GB9vKOGdhTUU9lfbPN1cgx5U2AHgPM9/L6N+2uRuSoFAYkcK183tyOmW12g+uCbSWmVaB3z3gig3xaoVdKXBS5KbEAbCxW60DAcjjaWgfUVSABKjR+K0roqdgF2xaYY+8v0i7YVmqaZ/uHyHgBLw0iDBKjW2bCbP0EqvbZWorI0SgMSGBa/SJgQd9pkNq7y6fbGiJ7EhzOO4hwncYNBLYe/c1um5ZHyxNLQFIYkMBjWW2zCMbM9skpvQvz9M9I9Mn5AoV33IAGEQaBsdea24LntDyxBByFAQk8e5ebJVwj42DYxT45ZajsVni4YTltMwos/YEJbCMuNcsTl6w33+MiAURhQAKPt1Vg2AyzpGsXNbk87CozawyE2piBfukJREc4qWl0sbu83u5y5HhiU8x+BQAFz9tbi8hhFAYksLia2vYi8FEXwc6yOjwWJERHkJHUtYWLAk1UhLN106K1GkQY+LxdBatfAlejvbWIHERhQALLpvegvgwSs6HfGT45ZWsXQUYCji5OUQxEw3NN68mavQoDAa/f6ZCcBw0VsOFtu6sRaaUwIIGl4DlzO/pKiIj0ySm90wr7pYfW4EGvES1hYO3eSpsrkRNyRrS1eHm/10UCgMKABI7qYtj4jjke9w2fnTZUBw96jchNAdQyEDTGtHQVbP4AqovsrUWkhW8+eon4wqp5YLkhfzJkDPHZabe2rD4YaoMHvYblmOmF+yobKKttokdCtM0VBZfK+mZW7a5gY3EN5bVNVDU0ExcdQVp8NH16xDOmVyo5KbG+62JKH2jWHNj1Bax6EU75gW/OK9IFCgMSGCwLVjxjjsd93aenDvWWgaTYKPqlJ7CttJY1eyuZNijD7pIC3v6qBl5fuZfXV+5l1e4Td6/kpcYxY0wOl47La90gqkvGXmvCQMFzcPL3u7zctkhXKQxIYNj9pVlbICrezMf2kZpGF/urzajtUNm6+GiG5ya3hIEqhYHj2FVWxz8+2cx/lu2m2d22LkPvHvGMzEsmIzGGpNgo6pvdlNc2saG4mvVF1eypqOfRT7fy6KdbOXVgOj84exCT+/XofCEjLoW3fwYl62DvCsgb74OvTqTzFAYkMKx42twOn+WTtQW8vBsUpSdGkxIX5bPzBpoRucnMX7VP4waOoa7JxYMfbGLuwm24PCYEjOudyqXj8rhgZDaZSbHHfG19k5tPNuzn1YI9fLhuPws3l7JwcylnD83k15eMID8tvuMFedccWP0SFDyrMCC2UxgQ+zXVQuF/zbGPuwi2ts4kCN1WATh4EKFmFBxu0eZSfvqfVeypMIsyTRtkPtlP6tu+T/Zx0RFcOCqHC0flsKusjoc/3cJLX+3iw/X7+XzLAX5y/hBuOKVvx8cUjL3WhIHVL8F590JUaGytLcFJswnEfmtfg6Yas5lLn5N9eupQHy/g5Z1euK20ltpGl83VBAa3x+Iv72/k63O/YE9FPXmpccy9fiJP3zil3UHgcL16xHP/paN4+4fTmNy3B/XNbu55cy03P7WMyrrmjp2s3xmQ0hsaKmHt652qR8RXFAbEft6Bg2O/5vOBVG1hIDTXGPBKT4whKzkGy4L1ReoqqGpo5pv/WsqDH27CsmD2pF68/+PTOHtYlk/OPzAziXm3TOXXF48gOsLJB+uKmT5nAZuKO7BhlNMJ41um0C7/t0/qEukshQGx14EtsGMROJxtS7X6ULi0DEBbV0HhnvAOA/sq67nqkcUs2FRKXFQEf716LL+7fDTx0b7tFXU6HVx/cl9e/u7J9O4Rz+7yei5/+HOWbD3Q/pOM/Zr53t+xCEo3+7Q+kY5QGBB7FTxrbgeeA8m5Pj21ZVmtAwgHZIR+GBiZZ8LA6j3hO25gS0kNl/79c9YXVZORFMNL3zmJWePy/HrNUfkpvHrrKYzvnUpVg4vr5i7l3TXtXEwoJQ8GnmuO1TogNlIYEPu4XW1Lsvp44CBAaU0T1Y0uHA7o3bMTI76DzNheJgys3FVhbyE22by/htn/XEJRVQMDMxN55XsntwYkf+uREM1zN0/lghHZNLk93Prsct4pbGcgmHC9uV35vNmoS8QGCgNin03vQvU+iO8Jgy/w+em9XQT5aXHEREb4/PyBZnR+KgCbS2qobujgYLYgt6WkhmseW0JJdSNDs5N44ZapnZvy1wWxURH87dpxXDI2F5fH4rbnlrevhWDQeZCYBbUlsFGbF4k9FAbEPl/9y9yO/RpE+n5r4VDfoOhw6Ykx5KfFYVnh1VVQVNnAdXOXtgaB526eSs9Ee7aqjoxw8qcrx7QGgu8/v+LEYwgioszPAMDyp/xfpMhRKAyIPcp3mI1aACZ80y+X2Foa2nsSHM2YXqkAFIRJV0FlXTPXP7GUPRX19E9P4Nmbpti+N4M3EJw3PIsml4ebn/rqxDM8vN1kmz+Eip3+L1LkMAoDYo/lTwGW2d+95wC/XMI7eDAcZhJ4jW3pKgiHcQNNLg+3PP0VG4qryUqO4d/fmmxbi8DhIiOcPHTNOCb1TaO6wcU3n/iS/VUNx35BzwHQdxpgwYpnu61OES+FAel+7ua25Ycn3uC3y4TTtEIvb8vAyl2h3U1gWRa/fK2QL7aVkRgTyZM3TKZXj8AaJBobFcHj101iYGYiRVUN3Pz0Mhqa3cd+gbeFbMUz4DnO80T8QGFAut+Gt6GmGBIyYMh0v1zC7bHYcaAOCK8wMDIvGacDiqoaKKo8zifRIPfEou3M+3IXTgfMuWYcw3J8t5+FL6XER/H4dRNJiYti5a4KfvHf1ViWdfQnD50BcWlQtbutC02kmygMSPdb1jJwcNzXIdI//bt7K+ppcnuIjnSSmxo+a77HR0cyOCsJgJW7K+wtxk+Wbivj/rfWAXDHRcM4c2imzRUdX9/0BB7+2nginA5eWbGHJz/ffvQnRsW2DST8cm631ScCCgPS3cq2wZaPzPH46/12Ge/gwb4944lwhtde8WNbuwoqbK3DH0prGvn+88txeywuHZfHjaf2s7ukdjl5YDp3Th8GwH3z17FsR/nRnzjxW+Z203tQvr17ihNBYUC6m3eVtQFnQQ///SLfWmKmFfbtGT5dBF6hOqPA47H40QsFFFc1MiAjgXtnjez4ToE2+ubJfZk+Oqd1DYKy2qMsMNRzAPQ/E7Dapt6KdAOFAek+rqa2TYkm+G/gIJjV6AAGZobHGgMHG987DYAVOytodntsrsZ3/v7xZhZsKiU2ysk/vjaBhJjg2oHd4XDw+8tH0z8jgX2VDfz0P6uOPn5g0k3mdsXT0By64z4ksCgMSPfZMN+sspaYBUMu9OulwjkMDMpMJCUuivpmN2v2hsamRZ9vKeUvH2wE4DeXjGRIdpLNFXVOYkwkf7tmfOtOhy98uevIJw2+AJLzoO6A2d5bpBsoDEj3WfqYuR1/nVl1zY+2tHQTDMgIvzDgdDqY1Ne0Dny5rczmarqupLqRH84rwGPBFRPyuXJiL7tL6pLhucn8z3mDAbjnzbVsbxnf0ioisq3l7CsNJJTuoTAg3aOosGWr4gi/dxFU1DVRWmP6YweEYcsAwKS+PQD4cntwhwHLsvjFf1dTUt3I4KxEfnPJSLtL8ombpvVnSr8e1DW5+dGLBbgO784Zfx04I2HXF7BvlT1FSlhRGJDusfRRcztsptm21Y+8rQI5KbEkBlm/sq9MbAkDX+0oP/a89iDwWsFePlhXTFSEg4euGUdcdGhsOBXhdPCnq8aQFBPJip0V/OOTLYc+ISnL/KyAWgekWygMiP/VlcGql8zxlG/7/XLhPF7Aa1ReCrFRTspqm1rDUbDZX9XA3a+vAeAHZw1iaHZgLizUWflp8dwzawQAD3646cjZH96BhKtehIbQXlFS7KcwIP634mlw1UPWKOh9kt8v5w0D4ThewCs60tm63sDSbceY0x7ALMvijlcKqaxvZmReMt85wz/7V9ht1tg8po/Owe2x+PELBYcuV9znFMgYCs11sHKefUVKWFAYEP/yuOHLx83xlFugG+aFq2XAmBzE4wYO7h7445VjiIoIzV9VDoeD+2aNJDMphq2ltcz5aNPBD7a1Dix9DDyhM01UAk9o/oRJ4Nj4rtmSNTYVRl7RLZfcXKIwADCpnwkDS4NsRsH+6tDuHjhcanw097QMjHz0062sPXg66JjZEJMMBzbBlg9tqlDCgcKA+Jd34OD46yDa/7vKNTS72V1eDygMjOudRqTTwZ6Kena2bNoU6CzL4v/CoHvgcBeMzOaCEdm4PBY//+8q3J6WQZ8xSeZnB2Dx3+0rUEKewoD4T8kG2PoJOJxtzZ1+tqWkBsuC1Pgoeib4ZxOkYJEYE9m6GuGCzSU2V9M+r6/cy/trTffAA1eEbvfA0fz6khEkxUayancl/1q0re2ByTebn6GtH8P+dfYVKCEtfH7SpPst/ae5HXwhpPXplksePHgwmNat95dpg9IBWLCx1OZKTuzg7oHvnzUoYLcl9pes5FjuuMhsZvSn9zayq6ylNSetLwxt2ep7ycP2FCchT2FA/KO+AgqeN8dTbum2y27xDh4M45kEB5s2OAOARVtKj1zYJoBYlsWdrxRSUdfMiNxkvhsm3QOHu3piL6b060F9s5s7XlndtkbE1O+Z21UvQO0B+wqUkKUwIP6x7F/QXAuZw6Hf6d12WQ0ePNSovBRS4qKobnCxcnfgzlV/feVe3lsb+rMHTsTpdPDby0YRHelkwaZSXi3YYx7ofRLkjAVXAyx7wtYaJTSF50+c+JerCZY8Yo5P/n63TCf02lhswsCgLIUBMCvdnTqwpatgU2COGwj37oHD9c9I5AdnDQTgvvnrqaxvNj9D3taBpY+bnzERH1IYEN8r/A/UFEFSTrdNJwRodLnZ1rLpS7DuaucPreMGNgXeuAF1Dxzdzaf1p396AqU1jfzlfbNbIyMuNTt+1hTB2ldtrU9Cj8KA+JZlwed/M8dTvg2R3Teif2tJLW6PRVJsJNnJsd123UB3aksYKNhVYT5lBpA3Vu3jvbXFRDrDu3vgcDGREfz6ErNU8VOLt7Nmb6X5WZp0s3nCkn+YnzURH9FPnvjWlo9g/xqITvT77oSH21hcDcCQrCTNJDhIflo8/TMScHusgOoqKKlu5O7XCgF1DxzNtEEZTB+Vg8eCX762Bo/Hgok3QGQs7F0BOz63u0QJIQoD4lufzzG346+DuNRuvbQ3DAxWF8ERzh2eBcA7hUU2V2JYlsWdr66mvK6Z4TnJfO9MdQ8czZ0zhhEfHcGyHeX8Z/luSEg3qxICLPqrrbVJaFEYEN8pWm0WRnFEwJTvdPvlNxSZwYNDshQGDnfhyBwAPl6//9DNcGzyxqp9vLtG3QMnkpMSxw/PHgTA795eT2VdM5z8A7MI0ab3oKjQ5golVOgnUHzHO1Zg+CXdtsjQwbwtA5pJcKQx+SnkpsRS2+S2fSDh4d0Dw3PVPXA83zq1H4MyEymrbeKB99ZDzwHmZwxg0YP2FichQ2FAfKNyj5lFAGY6YTera3Kxs2XFNrUMHMnhcHD+yGwA3l69z7Y6LMvirlcL1T3QAVERztaNjJ79YieFeyrhlNvNg4UvQ/kO+4qTkKEwIL6x5B/gcUGfUyFvfLdfflPL+gLpidH0TIzp9usHA29Xwfvrimly2bMa4Zur9vHOmiJ1D3TQSQN6MnNMLpYFd7++BitnDPQ/Eyw3LP6b3eVJCNBPonRd7QH4qmVVtFNvt6WE1sGDahU4pgl90khPjKG6wcXnW7q/q2B/VQO/bOkeuO2sgeoe6KA7LhpKXJQZTPhqwR449UfmgeVPQU3gzBKR4KQwIF235B/QXAc5Y2DgObaUoDBwYhFOBxeMNLMKXivY263XtiyLn768ivKWxYVuPXNgt14/FOSkxHFby8qEv31rPTW5J0PueLNEsXercJFOUhiQrqmvaNud8LSfdOvSwwfb0NJNoJUHj++KCb0AeGv1PjMyvZs888VOPtlQQkykk79ePVbdA51007R+9O0Zz/7qRuZ8tLmtJW7pY9BYbWttEtz0Eylds/QxaKyCjGEwZLptZWwsUstAe4zJT2FodhKNLg+vrNjdLdfcUlLDffPXAvDzC4cySP+POi0mMoJfzhwOwBOLtrGl5xnQcyA0VMCyJ+0sTYKcwoB0XmMNLPm7OT7tf8Fpz7dTeW0TRVUNAAzWtMLjcjgcXDulNwDPL93VtkWunzS7PfzohQIamj2cOjCd60/q69frhYOzhmZx1tBMmt0Wv56/AevkH5gHPp8DzfX2FidBS2FAOu+rJ6C+HHoMMJuo2GTdvioA+vSMJyk2yrY6gsUlY/OIjXKyobia5Tsr/HqtOR9tZtXuSlLiovjjlWNwOrVMtC/8csZwoiOcfLaxhA+iz4SU3lBTrNYB6TSFAemc5vq2pYen/RicEbaVsrYlDAzL1uj09kiJi2L6qFwAnvtip9+us2xHOX//eDMA984aSXaKNo/ylb7pCdw0rR8A97y9meZTWmYWLPyLWgekUxQGpHOWPw21+yGlF4y+2tZS1u41YUBT1drv61NNV8FrBXvY1bJYky8dqGnktueW4/ZYzBqby8wxuT6/Rri79cyBZCfHsqusnkcrp7a1Dnz1L7tLkyCkMCAd52ps2yTl1Nshwt6meW/LwHDtetdu43qnMW1QOi6PxZyPNvn03G6PxQ/nFbCvsoH+GQn8ZtZIn55fjISYSO6YPgyAv322g7KJLWMHFv4Fmnwf8CS0KQxIx331L6jaA0k5MPbrtpbS6HKzeb+ZVjhMLQMd8qNzBwPw8vI9bC+t9dl5//z+BhZuLiUuKoJHvj5B4zj8aOboHKb060FDs4df7hgNqb1Ni513ETCRdlIYkI5pqoUFfzTHp/8UouztB95UXIPLY5ESF0Wu+qQ7ZHzvNM4ckoHbY/HQh75pHXjpq138/eMtAPzu8lGa6ulnDoeDX108AqcD3iwsZdPQ75oHFv1VrQPSIQoD0jFfPAK1JZDWF8Z9w+5qDukicNi04FEw87YOvFKwhy+2HujSuT7fXMov/rsagFvPHMAlY/O6XJ+c2LCcZL4x1ewS+v01Q7BS+5if0a/m2lyZBBOFAWm/+oq2LVPPuMP2sQLQNnhwmMYLdMro/FSumpiPZcH/vLSS6obOrUq4fGc53356GS6PxcwxufzPuUN8XKkcz4/PHUKPhGjWlzSwIOcGc+fCv5qWPJF2UBiQ9lv8N2iohIyhMOoKu6sB2tYY0EyCzrtrxnDy0+LYXV7PvW+u6/Drl+0o47q5S6ludDG1fw8euGK01hPoZinxUfzkfBPAfrh2MO6UPlBXCksetrkyCRYKA9I+NSWw+B/m+Kw7bV1XwMuyLM0k8IGk2Cj+fNVYHA544atdPLFwW7tf+05hEdfNXUpNo4uT+vfkiW9OIjbK/u+NcHTVxF6MzEumvBHmJV5n7lz0oNlVVOQEFAakfRb+BZprIXccDJ1hdzUA7C6vp7rBRVSEg4GZWoa4Kyb368HtZ5vxA/e8ubZ1saBjaXZ7+O1b6/jOM8uobXJz6sB0nvjmJOKjI7ujXDmKCKeDX19spnHeuWUIdT1Hmn1DPnvA5sokGCgMyIlV7oEvHzfHZ91l286EhyvcUwnAoMwkoiP1rdxVPzh7ILefMwiAB97dwC1PfdW6NbSX22Px6oo9nP2nT3n0s60A3HRqP/51wyTiotUiYLcJfdK4bHweFk5+1zzb3Pnl41DW/tYeCU+K8XJin9wP7kbocwoMOMvualqtagkDY3ql2FxJaHA4HNx+zmDioiL43TvreW9tMe+vK2Z4TjK90uKpbmxm5a5KahpdAKQnRvPri0cyfXSOzZXLwX5+4VDeW1PMU/v78938k8kp/Rw+uheu0OwCOTZ9nJLj27cSVjxrjs/5VcC0CgCs2l0BwKi8VFvrCDXfPn0A795+GheMyMayYM3eKt5ZU8SizQeoaXSR2jJY7bOfnqkgEIAyk2L54dmmhed/ymaZOwv/A3tX2FeUBDy1DMixWRa8+3+ABSMvh16T7a6olcdjsWq3aRkYna+WAV8bnJXEI9+YwM4DdWzaX82usjpioiIY2yuVwVlJRGi2QEC7/uS+zPtyJ5+X5LMq5zxGl78H798N170WUIFeAofCgBzb+vmwfQFExppWgQCyo6yO6gYX0ZFOhmRrlTt/6d0znt494+0uQzooOtLJry4ewTfmLuW24hl8EvsJzm2fwpYPYeA5dpcnAUjdBHJ0riZ4/y5zfNKtZs3zAOLtIhiek0xUhL6NRQ43bVAG54/IYqcnnbfjWmYAvfdLcLvsLUwCkn6LytF9+RiUbYWETDj1R3ZXcwRvF8EYdRGIHNOd04cTE+nkjgMX0BSdAvvXaBMjOSqFATlSXRl8+ntzfNadEBN4zfCtgwfzU22tQySQ9eoRz3dOH0AliTzoudrc+fG9UFtqb2EScBQG5Eif/NYsO5w1EsbZu0Xx0bjcHgr3mJUH1TIgcnzfPWMAealxPFxzGvsTBpuf7Q/vsbssCTAKA3Ko4jXwZct85PPvC4hlhw+3paSW+mY38dER9M/QyoMixxMbFcFdM4bhwckPK681dy5/CvYst7cwCSgKA9LG44E3fwSW2yw53P8Muys6qpUtXQQj81I0xU2kHc4fkc2pA9NZ7BrMkoSzAQve/qn5mRdBYUAOtuJp2PUFRCXAhb+3u5pjWrGzHICxvVLtLUQkSDgcDn518XAinQ5+cOAyXJEJsPtLWDXP7tIkQCgMiFFbCh/cbY7PvANS8u2t5ziW7TBhYEKfNJsrEQkeAzOT+ObJfdlPGnOdLVuQv3+3GUMgYU9hQIz3fwn15ZA1CqZ8x+5qjqmyrpmNxTWAwoBIR/3wnEGkJ8bwp6qzqIjrA7X74YNf212WBACFAYHti6CgZf+BGX+GiMBdmHJ5SxdB357xpCfG2FyNSHBJio3i5xcOpYkoflR7nbnzq7mw43N7CxPbKQyEO1cTzP+xOZ7wzYDaf+BovtpRBsCEPj1srkQkOF02Lo9xvVP5uGkYi1MuMne+/n1obrC3MLGVwkC4+/whKFkP8elw9t12V3NC3vECE/uqi0CkM5xOB/dcPBKHA75dPIumuAw4sBk++4PdpYmNFAbCWfEa+OR35vj8+yE+sD9tN7s9FOyqAGCixguIdNqo/BSumdybKhK517rR3LnoQShabW9hYhuFgXDlboZXvgOeZhhyEYy+yu6KTmjt3ioamj2kxEUxQIsNiXTJzy8cSlZyDE9VjGZD2hngccFrt2kjozClMBCuFvwJilZBXBrM+GtQ7HHu7SIY3zsVpxYbEumS5Ngo7p01CoDri6/CHZ0M+wpgyd/tLUxsoTAQjvYWwGcPmOOL/ghJWbaW015t4wUCuztDJFicOzyL6aNzKPKk8rfI682dH91nuhAlrCgMhBtXI7z6XdMkOPwSGHm53RW1i8djsXjrAQAm91MYEPGVX80cQUpcFH8pm8qOHqeCuxFevkmzC8KMwkC4+eR3sH+tmT0w/c9B0T0AsKG4mrLaJuKjIxijbYtFfCYjKYZfzhgOOJi9/+u443qa3xEfajGicKIwEE62L4RFfzXHM/4CCem2ltMRizab/dcn9e1BdKS+bUV86bLxeUwblM4+VzJ/iPm+uXPJP2Dzh/YWJt1Gv1XDRW2pafqzPDDmWhh+sd0VdcjnW0wXwSkDe9pciUjocTgc3H/pKBKiI3i0aDBr8q40D7z6Pag9YG9x0i0UBsKBx2OmEVbvg/TBMP2PdlfUIc1uD1+0jBc4eUDwtGaIBJNePeK5e+YIAK7ZPoPG1IFQUwRv/AAsy+bqxN8UBsLB5w/B5vchMhaufBKiE+yuqENW7a6ktslNanwUw3OS7S5HJGRdOTGfc4ZlUuWO4seu27CcUbD+TVj6mN2liZ8pDIS6nV/Ah/eY4wt/D1kj7K2nExZvMeMFTurfU+sLiPiRw+Hgt5eNpmdCNPNLM/ko/1bzwLu/ML9LJGQpDISyujJ4+Uaw3DDyChh/vd0Vdcqizd4uAo0XEPG3jKQY7r/MLEZ006ZJlPW9yExFful6qNlvc3XiLwoDocrjhv/eDJW7oEd/M3sgSKYRHqyuycWylm2LTx6o8QIi3eH8EdlcMSEfy3Iwu/jreHoONmOOXrpByxWHKIWBUPX+L2HzBxAVD1f+G2KDs6994aZSmlweevWIo396cI11EAlmd88cTl5qHBvL4f6k/8OKToQdC+GDwN/dVDpOYSAULX8aFv/NHM96GHJG21tPF3y4zjRLnj00C0cQtmyIBKuk2CgeumYsEU4Hj6+PYuGIlrFHi/8Ghf+1tzjxOYWBULNzCbz5I3N8+s9hxCxby+kKj8fiow0tYWBYps3ViISfCX168NPzhwBw45e5lI7+tnng1e/B7mU2Via+pjAQSip2wQtfN9sSD7sYTv+Z3RV1SeHeSkqqG0mIjmBKPw0eFLHDzdP6c9bQTJpcHq7ecj6u/ueAqx6evxrKt9tdnviIwkCoaKiC56+B2hLIGgWXPgLO4P7f+0FLF8FpgzO0BLGITZxOB3+6cgy5KbFsOdDAz5w/wsoeZX7XPHOFmbUkQU+/YUNBcwPMuxaKV0NCBlzzXNAtLHQ0H60vBuDsYcGxxbJIqEpLiGbOteOJdDp4ubCSl4f+BZLz4MAmeOEbZjdUCWoKA8HO7TJrCWxfANFJ8LWXILW33VV1WVFlA4V7qnA44IwhGXaXIxL2JvRJ42cXDAXgF++XsOr0xyAm2cwweO1Ws+y5BC2FgWBmWfDm7Wa50Iho0yKQO87uqnzirdX7AJjQO430xBibqxERgJum9WP66Bya3RbffKuO/Rc+Bs5IWP0SvPNz7WEQxBQGgtmHv4YVT4PDCZfPhX6n2V2Rz7y+ci8AM8fk2lyJiHg5HA7+eMUYRuYlU1bbxDc+jqfhoofMg0sfhffvUiAIUgoDwWrBn2HhX8zxjL8E3ZbEx7PzQB0FuypwOuCiUTl2lyMiB4mLjuCx6yaSkRTDhuJqblszGM/0v5oHP58DH99na33SOQoDwejTP5hWAYCzfwkTvmlrOb72xirTKnDygHQyktRFIBJoclLi+Oc3JhAd6eSDdft54MBJcOEfzIOfPQCfPmBvgdJhCgPBxLLgw9+0Je+z7oJp/2NvTX7wRksXwcXqIhAJWON6p/GHy83qpg9/soV5jgvh3N+YBz++17ReStBQGAgWlmX64xb80fz3effCaf9rb01+sLG4mvVF1URFODh/RLbd5YjIccwal8f3zxoIwB2vrObt5CvhzDvNgx/+2uyRojEEQUFhIBh4PPD2T01/HMCFD8DJ37e3Jj95efluAE4fnElKfJTN1YjIifz43MFcM7k3Hgt+OK+Ahbk3wDkt3ZiLHoTXb9NOh0FAYSDQNdebfcSX/hNwwIy/wpRb7K7KLxpdbl76yoSBqybm21yNiLSHw+Hg3lkjuWhUNk1uD7c8/RUrel8PF//NzHRa8Yz5HdbcYHepchwKA4GsZj88OQPWvW7WEbj8cZh4g91V+c07hUWU1TaRnRzLWUO1MZFIsIhwOvjL1WOZNiiduiY3Nzz5JRvzZsFVT0NEjFkL5ZnLob7c7lLlGBQGAlXJBnj8bNjzFcSlwXWvwagr7K7Kr55dshOA2ZN7ERmhb02RYBITGcEjX5/A2F6pVNQ1c+1jX7A+7TT4xn/bVip87CzYv97uUuUo9Bs3EG35COaeCxU7Ia0f3PgB9DnZ7qr8amNxNUu3lxHhdDB7UvAvpywSjhJiInnyhkkMy0mmtKaR2f9cwqrIkXDDW5DSG8q2wuPnwPq37C5VDqMwEEg8bvjkd/D0ZdBQCb2mwE0fQPpAuyvzu2eW7ADgnGGZZKfE2lyNiHRWanw0826eekgLwZcNeXDLx9DnVGiqhnnXmLUINNMgYCgMBIraUtOn9slvAQvGXwfXvQ4J6XZX5nelNY28+NUuAL4xta+9xYhIl6XER/HMTVOY0q8HNY0urpu7lIV7getehcktA6A/vhfmfU1bIAcIhYFAsHMJPDINtn4MkXEw6xG4eA5Ehccn5LkLt9HQ7GFMfgqnDOxpdzki4gOJMZE8ecNkTh+cQX2zm289+SWvrd4PFz0AMx8CZxRsmA8PnwJbP7W73LDnsCy109jG1WhaAhY9CJYHeg6Cq56CrOF2V9ZtKuqaOOV3H1Hb5Oax6yZy7vAsu0sSER9qdLm5fV4BbxcWAfD9swbyo3MG4yxaabZfP7AZcMApPzALFkVG21twmFLLgF32LIdHTzebDVkeGHWl6VMLoyAA8OTn26ltcjMsJ5lzhmk6oUioiYmM4G/Xjufbp/cHYM5Hm7n1ueXUpY+Eb3/WsreKZT4UzT0XitfYWm+4UstAd3M1mo2GFv4FLDfEp4fcroPtVVbbxBkPfExVg4u/Xzue6aO1Q6FIKHvpq13c8cpqmt0Ww3OSefz6ieSmxsG6N+D175t1CJyRcNJtcPrPIDre7pLDhsJAd9r4HrzzcyjbYv57xGVw0R8hITz7yf/vldU8+8VOhuUk8+b3TyXC6bC7JBHxs6+2l/Htp5dxoLaJtPgofn/5aM4bkQ1V++Dtn5hgAJDaB2b8GQaeY2/BYUJhoDsc2ALv/AI2vWv+OyHTDKIZMcvWsuy0Zm8lM+YsxLLgxW+fxOR+PewuSUS6ya6yOr7zzDLW7K0C4GtTenPn9OHERUeYNQje+l+o2mOePPwSOOdX0KO/fQWHAYUBf6orM90BSx4GT7MZPTv1u3DaTyA22e7qbGNZFlc/uoSl28uYOSaXOdeMs7skEelmjS43f3x3A48t2AbAwMxEHpo9juG5ydBYDR//Fr542IypckbCxBvh9J+GxXRrOygM+EN9OSz+hwkBTdXmvoHnwAW/g/RB9tYWAJ79Ygf/90ohcVERfPg/p5s+QxEJSws2lfDjF1dSUt1IdIST754xgO+eMYDYqAgoKoQPfgWb3zdPjk6CU38IU74DMUm21h1qFAZ8qb4CvngUFv8dGivNfVmj4Oy7YNB54FCf+Ob91cyYs5CGZg93Th/GTdPU9CcS7g7UNPKzl1fzwbpiAPr0jOdXF4/gzCEtM4y2fgLv/xL2rTT/HZtiFi+a8t2wHXPlawoDvnBgC3zxCKx4FpprzX2Zw+GMX8DQGeDUDE4wzYKX/v1z1u6rYtqgdP59w2ScGjQoIpjuw7cLi7jnjbUUVZntji8Ykc1dM4eTlxoHHg8Uvgyf/h4ObDIvioyDCdfDSbdCqvY06QqFgc7yeGD7Z7DkEdj4DtDyNmYOh9P+F4ZfqhBwEMuy+MV/VzPvy12kxUfxzu2nkZUcHissikj71TS6ePCDjTyxaDtuj0V0hJNrp/Tme2cMIDM51uzhsu4NWPjntpYChxMGnQ8TvwUDzwZnhL1fRBBSGOio8h2w8nkoeNbsKug16HwzOLD/GeoOOIo5H27iT+9vxOmAx6+fyFlDtdKgiBzb+qIqfvX6GpZsNXsXxEQ6ue6kPnzn9AH0TIwxmxxt+QgW/RW2fdb2wtTeZiGj0VdDSr4ttQcjhYH2qC2F9fOh8D+HftPFJMPoq8xgFg0MPKYXv9zFT19eBcBvLhnBN07qa29BIhIULMvi8y0H+NN7G1i+swKAuKgILp+Qx/Un9WVQVssgwpKNsOxf5kNaQ2XbCfqcAiMvh+GzNLbgBBQGjqVqH6x/E9a+BjsWmektADig32kw7utmPIBWyDquxxds5d756wD49un9+cWFw2yuSESCjWVZfLKxhD+/t5HVe9r+2J88oCfXn9yXc4ZlmUXLmupgzX+h4Dnze9vLGWl+bw++EIZcoPEFR6Ew4NXcADs/N81Omz+C/Yetj5092iwZPPpqfSO1g8vt4Q/vbuCfn20F4Fun9OPO6cM0YFBEOs2yLBZvOcC/F2/n/bXFeLxDtZJimDE6l4vH5jImPwWHwwGVu82Aw9UvQdHqQ0+UOQIGn2e6dXtNgShNbw7fMNBYDbu/NNsH71wMu74EV/1BT3BA/kQYdjEMmwk9+tlWarDZVVbH7S8UsGxHOQA/u2Ao3zm9v/kBFRHxgd3ldTz7xU7mLd1JeV1z6/19esYzc3QuZw3LZEx+qmkxKNkIG96Cje/CriUHtfQCEdGQPxn6TYPeJ0He+LBcwyA8wkBDlUmG3n/7VppP/gd/QwAk5cCAs8y//meqj6mDGl1unvp8Bw99uInqRhdJMZHcf9koZo7Jtbs0EQlRTS4Pn20s4fWVe3l/bTH1ze7Wx1Liopg2KJ3TB2dwysB0s8BZXRls/gA2vQ/bF0D1vsPO6IDMYZA3wfzLGgmZQ0M+IIRWGKgvhwNbzf7YZVtg/zooWgXl24/+/NTeJgn2nmpuM4ZqJkAnNDS7ea1gD//4ZAs7DtQBMK53Kg/NHkevHhpTISLdo67Jxftri3l3TRELN5VS1eA65PGclFjG90ljQu80xvdJY2hWIrFV28008W0LTGtx5a6jnzy1t5k67v2XMQTS+obM0vLBEwYsy/yxr9prNrCo2mOOK3aZP/wHtkB92bFfn5wPOaMhe5Tp/88dByl53Vd/iLEsi5W7K5m/ai//Xb6HA7VNAKQnxvCT8wdzxYRe2oVQRGzjcnso2FXBpxtL+HRjCWv2VuH2HPrnzumAvukJDMtOZmh2EgMzExkQX0uf+rXEFK2AvSvMh8qaomNfKK6HCQUH/0vJNy3NSdkQlxYUHzLtDQONNVBTbJpt6g5AXWnLrfdfy/21JWZ0/yF9+seQlAM9BkDP/pA+uO2Pf7x2xeuKJpeHTfurKdxTyRdby/h8y4HWVcIA8lLjuP7kPlw7pQ+JMZE2VioicqTaRhcrd1ewfEc5y3dWULCrgrKWDzFHk5kUQ+8e8WSnxNI3roEhzt30ce8gq2ELKdWbia7YgvN4H0C9ImIgKQsSs004SEg3ASK+x5G3iVkQk+jDr7r97A0D799tFozoiPh0SM6F5Dxzm5LX8sd/gNniMjrBL6WGEpfbQ5PbQ2OzuW1odlNV76KyvpmqhmYq65sprmpgT3k9u8vr2VNRz96KelyHper46AjOHpbFjNE5nD00k8gIrbgoIsHBsixKahpZv6+a9UVVrN9XzdbSWnYcqD1kQOLxpDrrGRJbzpDoUvpFlNLLsZ9cq4ge7gMku0qJd1We+CQHaZz6Q2IuuKczX06X2foR7qsSJ8OIo9KRTJUjmUpHMpWOpJbbQ/+VkEapoydNjmioBqsa2NN2LssqBUrb/ptD/3AdHHlOFH8Oz0eHP/2Qcx3nOid67eHPOPFrO1DXUZ7b7PbQ5PLg6WT8S46NZERuChP6pHHSgJ5M6JNmdhYTEQkyDoeDzKRYMpNiOW1wxiGPVdY1s6Oslp1ldRRXNbK/qoHiqgaKqxoprm6guLKB2iY3FZ44vqiL44u6ow+SjqaZDCrIcpST6agg01FOD0c1qdSQ5qghjWpSHTWkUUOao5pl+yxO744v/ihsDQMLM2ZzxaqJHXiFq+Wf+IrTAbFRESTFRpISF0VKXBTJsVGkJ8aQnxZHfo848lLjyU+LIyclVtMDRSTkpcRHMTo+ldH5qcd8Tn2Tm8p605JaUddkbuubqaxrprbJRX2Tm/pmN3Utt/VNbjY1ualrdlPf5KLR5aHZ1dJK6/LQ7Pbw/T4DbQsDtnYT7C6vY29FwyH3Hf63xnGcxw599PivNY87jvPY4a89/rm78trjPf/EX8MRVz/m44c/NSrCSUykk2jvvwinmvZFRCSIZhOIiIiIX+hjoYiISJhTGBAREQlzCgMiIiJhTmFAREQkzCkMiIiIhDmFARERkTCnMCAiIhLmFAZERETCnMKAiIhImFMYEBERCXMKAyIiImFOYUBERCTMKQyIiIiEOYUBERGRMKcwICIiEuYUBkRERMKcwoCIiEiYUxgQEREJcwoDIiIiYU5hQEREJMwpDIiIiIQ5hQEREZEwpzAgIiIS5hQGREREwtz/A8qKLACbhkLyAAAAAElFTkSuQmCC", | |
"text/plain": [ | |
"<Figure size 640x480 with 1 Axes>" | |
] | |
}, | |
"metadata": {}, | |
"output_type": "display_data" | |
} | |
], | |
"source": [ | |
"\n", | |
"fig, ax = plt.subplots()\n", | |
"ax.plot(xx, p(xx), label=r\"$p(x)$\")\n", | |
"ax.plot(xx, k * q(xx), label=r\"$k \\cdot q(x)$\")\n", | |
"ax.set_axis_off()\n", | |
"plt.legend()\n", | |
"plt.show()\n" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 9, | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
"Acceptance rate: 0.502.\n" | |
] | |
} | |
], | |
"source": [ | |
"\n", | |
"# Generate samples via RS. \n", | |
"n_samples = 10000\n", | |
"\n", | |
"# Candidate samples.\n", | |
"proposals = np.random.normal(50, 30, size=n_samples)\n", | |
"uu = np.random.uniform(0, 1, size=n_samples)\n", | |
"# mask = p(proposals) / (k * q(proposals)) > uu\n", | |
"# keep = proposals[mask] \n", | |
"keep = []\n", | |
"\n", | |
"for x, u in zip(proposals, uu):\n", | |
"\n", | |
" if p(x) / (k * q(x)) > u:\n", | |
"\n", | |
" keep.append(x)\n", | |
"\n", | |
"\n", | |
"print(f\"Acceptance rate: {len(keep) / n_samples:.3f}.\")\n" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"Plot accepted samples vs. $p$. " | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 10, | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAgMAAAGFCAYAAABg2vAPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA8eElEQVR4nO3de1zN9+MH8Nepc3RVlEuJWkVFhFwrasj9MtlG2HeY+W4zzNj2nV3NsDGX72zY5rIrX0ZuG5v5uuSS2ESI7hJCFKWc7uf3x775zTZddM55fy6v5+Ph8RhOn8+rofPq/XlfNAaDwQAiIiJSLQvRAYiIiEgslgEiIiKVYxkgIiJSOZYBIiIilWMZICIiUjmWASIiIpVjGSAiIlI5lgEiIiKVYxkgIiJSOZYBIiIilWMZICIiUjmWASIiIpVjGSAiIlI5lgEiIiKVYxkgIiJSOZYBIiIilWMZICIiUjmWASIiIpVjGSAiIlI5lgEiIiKVYxkgIiJSOZYBIiIilWMZICIiUjmWASIiIpVjGSAiIlI5lgGqsbi4OEycOBFZWVmioxARkRGxDFCN7N27F2FhYdiyZQucnJxQXl6OTz75BHl5eaKjERFRHbEMULU2bdqEQYMGISQkBJcvX4a1tTWuXr2Kt956Cy+++KLoeEREVEcsA1Sl5ORkREZG4vHHH8eOHTtgZ2cHAGjevDlWrFiBdevWYd26dYJTEhFRXWgMBoNBdAiStri4OAQEBECr1f7l95566ins2LED8fHx8PT0FJCOiIjqimWAHig/Px/169eHRqN54Gvy8vLQoUMHDB06FMuWLTNjOiIiMhaWAXqgyMhIFBUVYdu2bVW+Lj8/Hw4ODuYJRURERsc5A/S3kpKS8P3332PgwIHVvtbBwQEGgwFXrlwxQzIiIjI2lgH6Wx988AFcXV0xfvz4Gr1+8eLFCAgIQFFRkWmDERGR0bEM0F/cuHED3333HV555RVYWVnV6GOGDRuG3NxcREVFmTgdEREZG+cM0F8cO3YMY8eORUxMDJo0aVLjj+vduzfKyspw8OBBE6YjIiJjYxkgo9m4cSMiIyORkJCANm3aiI5DREQ1xMcEdJ+7d+/i4sWLD/WxERERaNu2LdLT042cioiITIkjA3SfTZs2YeTIkbh8+TLc3Nxq/fEGg6HKfQmIiEh6ODJA99m6dSs6dOjwUEUAADQaDbKyspCbm2vkZEREZCosA3RPSUkJdu7cieHDhz/0NfR6Pby8vPDtt98aLxgREZkUywDds2/fPuTn5yMiIuKhr2FjY4PQ0FDs2rXLiMmIiMiUWAbonqysLHTq1Ant2rWr03UGDhyI6Oho3L1710jJiIjIlDiBkO5jjAmAiYmJaN26NX788UcMHjzYSMmIiMhUODJAAIDCwkLcuXPHKCsBfH19ERAQgKysLCMkIyIiU+PIAAEA1q5di+effx65ubmwt7ev8/W4xJCISD44MkAAgEOHDsHf398oRQD4fYlh5WgDERFJG8sAAfi9DPTs2dNo1yspKYGbmxvWrl1rtGsSEZFpsAwQsrKykJaWZtQyUK9ePbRv356HFhERyQDLACEhIQE6nc6oZQAAQkJCEBMTA05LISKSNpYBQt++fZGfnw8XFxejXjc4OBjXrl3DhQsXjHpdIiIyLpYBQkVFBaytrY1+3aCgIGi1WiQkJBj92kREZDwsAyp3+/ZtNGrUCL/88ovRr+3s7Iy8vDwMHTrU6NcmIiLjYRlQuSNHjuDWrVvw9vY2yfVtbW05Z4CISOJYBlTu0KFDcHV1hZeXl0muv2fPHri4uCA/P98k1yciorpjGVC5I0eOoEePHibbLdDd3R3Z2dmIjY01yfWJiKjuWAZUrLy8HCdPnkTXrl1Ndg8fHx84OzsjJibGZPcgIqK60YoOQOJYWFggKSkJ9erVM9k9NBoNgoODWQaIiCSMIwMqptFo4ObmhsaNG5v0PsHBwYiNjUV5eblJ70NERA+Hpxaq2JIlS3D9+nUsWLDApPfJzs6GXq+Hu7s7TzIkIpIgPiZQsR9++AHOzs4mv0+TJk1Mfg8iInp4fEygUhUVFYiLi0OnTp3Mcr85c+Zg/vz5ZrkXERHVDsuASqWnpyM/Px+BgYFmud+FCxewdetWs9yLiIhqh2VApeLi4gDAbGUgMDAQp0+fRmlpqVnuR0RENccyoFJdu3bF2rVrTb6SoFJgYCBKSkpw7tw5s9yPiIhqjqsJyCwKCgrg4OCANWvWYMKECaLjEBHRH3BkQIUMBgNmz55t1u/S7e3tsXbtWvTo0cNs9yQioprhyIAKXbx4EY888gh++OEHDBkyRHQcIiISjCMDKnTixAkAMNuywkopKSl47733uBMhEZHEsAyoUFxcHFxcXODq6mrW+16+fBmzZ89GcnKyWe9LRERVYxlQofj4eHTo0MHs9+3YsSOA/1/WSERE0sAyoEK9evXCqFGjzH7fBg0awNPTk2WAiEhiOIGQzOqJJ55ATk4O9u/fLzoKERH9D0cGVOb69euIjY1FWVmZkPs//vjj6NWrl5B7ExHR3+PIgMp8+eWXmDhxIu7cuQM7OzvRcYiISAI4MqAyZ8+ehZeXl7AiYDAY8Ntvv+HixYtC7k9ERH/FMqAyZ8+ehb+/v9AMAwYMwNdffy00AxER/T+WAZVJSEhA27Zthd1fo9GgXbt2OHPmjLAMRER0P5YBFSkqKkLDhg2F7DHwRywDRETSohUdgMzH2tpaEm/C7dq1w/Lly1FUVARra2vRcYiIVI8jAypSUlICKSwead++PXx8fHDt2jXRUYiICFxaqCrTpk3Db7/9hpiYGNFRiIhIQjgyoCJnz56Fm5ub6Bj3FBUViY5ARERgGVAV0SsJ/ugf//gHhg8fLjoGERGBZUA1srOzkZ2dLZky4O7uLonJjERExDKgGgkJCQAgfMOhSu3atUNWVhZyc3NFRyEiUj2WAZV49NFHkZGRgVatWomOAuD3MgD8Po+BiIjEYhlQCY1GAw8PD1haWoqOAgDw8fGBTqdjGSAikgAuLVSJKVOmwNfXF1OnThUd5Z7MzEw0b94cFhbspEREIvGrsEps3boV169fFx3jPu7u7iwCREQSwK/EKnDnzh1kZWXBz89PdJT7/PTTT+jdu7ckdkUkIlIzlgEVSE5OBgD4+voKTnK/8vJy7N+/H5cuXRIdhYhI1VgGVCAxMRGA9MpAmzZtAADnz58XnISISN1YBlSge/fuWL16NRwcHERHuY+HhwdsbGxw7tw50VGIiFSNRxirgLe3N7y9vUXH+AtLS0v4+vpyZICISDCWARVYuXIlunTpgs6dO4uO8hcffPABnJ2dRccgIlI17jOgcBUVFbCzs8MHH3yA6dOni45DREQSxDkDCpeZmYmioiLJTR6slJWVhblz5/KMAiIigVgGFC4pKQkAJLfHQKX8/Hy8/fbbOH36tOgoRESqxTKgcImJibCysoK7u7voKH/L29sbOp2OKwqIiARiGVC45s2b45lnnpHMAUV/ptPp4OPjwzJARCQQJxCScE8++SRyc3Oxd+9e0VGIiFSJIwMKd+rUKRQWFoqOUaVhw4YhNDRUdAwiItXiyICC5efnw9HREd9++y2eeuop0XGIiEiiODKgYJUHFEl1JUGl8vJynDp1CtnZ2aKjEBGpEsuAglUuK/Tx8RGcpGqlpaXo1KkTduzYIToKEZEqsQwoWEpKCpo2bSq5A4r+zNraGt7e3lxRQEQkCMuAguXn5987JljqWrduzTJARCQIJxAqnMFggEajER2jWrNmzcK6deuQmZkpOgoRkepwZEDh5FAEAKBt27awsrJCSUmJ6ChERKrDkQGFysnJgb+/P9avX4/evXuLjkNERBLGkQGFSklJwfXr1+Hk5CQ6Sq1UVFSIjkBEpDosAwqVkpICAGjZsqXgJDXXqVMnzJ49W3QMIiLVYRlQqJSUFLi6usLe3l50lBpzcHBAYmKi6BhERKrDMqBQqampshoVAH5fXsgyQERkflrRAcg0PvjgA+Tn54uOUSt+fn5Yu3YtysvLJXvkMhGRErEMKJSHh4foCLXm5+eH4uJiXLx4EV5eXqLjEBGpBh8TKFBubi6eeeaZe2cTyEVISAgSExPh7u4uOgoRkaqwDChQYmIivvzySxQXF4uOUit2dnbw9fWFVssBKyIic2IZUKDU1FQAgLe3t+AktffBBx9g8eLFomMQEakKy4ACpaSkoFmzZrCzsxMdpdbOnDmD7du3i45BRKQqLAMKlJKSglatWomO8VD8/Py4vJCIyMz4cFaBhg8fLtttff38/HDjxg3k5OTA2dlZdBwiIlXgQUUkKadPn0b79u1x5MgRBAcHi45DRKQKfEygMHfu3MGWLVtw+/Zt0VEeio+PDxYvXowWLVqIjkJEpBocGVCYmJgYhISEID4+HgEBAaLjEBGRDHBkQGHkvKywUmxsLKKiokTHICJSDZYBhUlPT4eLi4sslxVW2rRpE2bNmiU6BhGRarAMKExaWpqsRwWA31cUpKeny24HRSIiuWIZUBhnZ2fZz8L38/NDeXk50tLSREchIlIFTiAkycnOzkbTpk0RFRWFESNGiI5DRKR4HBlQkPLycuj1etEx6qxx48YYOnQobG1tRUchIlIFjgwoSHx8PDp06IBff/0VnTt3Fh2HiIhkgiMDClL5jN3d3V1wkrorLi5Genq66BhERKrAMqAgaWlpsLe3R+PGjUVHqbNly5ahQ4cO4MAVEZHpsQwoSOWyQo1GIzpKnfn5+eHOnTvIysoSHYWISPFYBhRECXsMVPLz8wMAHmdMRGQGPMJYQXbs2IGCggLRMYzC09MTOp0OiYmJ6NOnj+g4RESKxjKgIDY2NrCxsREdwyi0Wi18fHz4mICIyAy4tFAhMjIy8Oyzz+KTTz5B69atRccxitLSUuh0OtExiIgUj3MGFCI5ORl79+6FlZWV6ChGwyJARGQeLAMKkZaWBq1Wq4g9Bir997//haenJ+7cuSM6ChGRorEMKERaWho8PDyg1SpnGoijoyMyMjKQnJwsOgoRkaKxDCiEkpYVVvL19QXA5YVERKamnG8jVW7q1KkoLy8XHcOoHBwc0KxZM5YBIiITYxlQiN69e4uOYBJ+fn4sA0REJsalhQpw69YtfPXVV4iMjISrq6voOEZ18uRJ2NnZwcfHR3QUIiLFYhlQgCNHjqBHjx44ffo02rVrJzoOERHJDCcQKkDl0cVeXl6Ckxjf5cuXMWXKFFy+fFl0FCIixWIZUID09HS4uLjAzs5OdBSjMxgMWL58OeLj40VHISJSLJYBBUhLS1PkqAAAuLm5wc7OjpMIiYhMiKsJFMDf3x9t2rQRHcMkLCws4OvryzJARGRCnEBIkjd27FhkZmbi0KFDoqMQESkSRwZkrqSkBBkZGfD09FTswT6jRo1CZmam6BhERIrFkQGZO336NNq3b4+YmBgEBQWJjkNERDLECYQyV7msUGnnEvxRSUkJfv75Zy4vJCIyEZYBmUtLS4O9vT0aN24sOorJVFRUYPDgwdi9e7foKEREisQyIHOVpxVqNBrRUUzG2toanp6eXFFARGQiLAMyl5OTg5YtW4qOYXI8sIiIyHS4mkDmvv/+e8UdXfx3/Pz8sH37dtExiIgUiSMDCmBpaSk6gsl16tQJzZs3R0VFhegoRESKw6WFMnbhwgX069cPGzduRGBgoOg4REQkUxwZkLGUlBSkpqbC0dFRdBSzKCsrQ2FhoegYRESKwzIgY+np6bC0tIS7u7voKGbh7u6OpUuXio5BRKQ4LAMylpaWBg8PD8VuQ/xnXl5eXFFARGQCLAMyVrnHgFpweSERkWlwaaGMzZs3D3q9XnQMs/Hz88OGDRtgMBgUvckSEZG5sQzIWOvWrUVHMCs/Pz8UFhYiKysLbm5uouMQESkGHxPIVE5ODqZPn4709HTRUcwmPDwct2/fZhEgIjIylgGZSkpKwscff4yCggLRUczG2tpaNcsoiYjMiWVApiqPLvby8hKcxLymTZuGuXPnio5BRKQoLAMylZaWhqZNm8Le3l50FLPKysrCgQMHRMcgIlIUlgGZUtuywkpcXkhEZHxcTSBTffr0Qc+ePUXHMDs/Pz9cuXIFd+7cQf369UXHISJSBB5URLLy22+/oUuXLjh+/Di6dOkiOg4RkSLwMYEMFRUVYc+ePcjLyxMdxezatGmD9evXw9PTU3QUIiLF4MiADJ0+fRrt27fHkSNHEBwcLDoOERHJHEcGZKhyWaEaJxACwM6dO7F27VrRMYiIFINlQIbS0tJgZ2eHJk2aiI4ixE8//YTFixeLjkFEpBgsAzKUnp4Ob29v1R7W4+fnh5SUFJSVlYmOQkSkCCwDMqTT6dCpUyfRMYTx8/NDaWkpMjIyREchIlIETiAk2bl8+TJatGiBHTt2YOjQoaLjEBHJHkcGZMZgMKCiokJ0DKHc3Nwwbtw4NGrUSHQUIiJF4MiAzKSnp6N169bYu3cvevToIToOEREpAEcGZCYtLQ0lJSVwdXUVHUWonJwcHD16VHQMIiJFYBmQmbS0NFhaWsLd3V10FKG++eYb9O7dG+Xl5aKjEBHJHsuAzKSlpcHDwwM6nU50FKH8/f1RVFSECxcuiI5CRCR7LAMyo9aji/+sbdu2AICEhATBSYiI5I9HGMvM6tWrkZ+fLzqGcK6urmjQoAESEhLw2GOPiY5DRCRrLAMy4+TkBCcnJ9ExhNNoNAgODkZxcbHoKEREsselhTJy/fp1PP/885g3bx7atGkjOg4RESkE5wzISHJyMrZt2yY6hqSUl5erfhMmIqK6YhmQkdTUVGg0Gnh5eYmOIgnHjx+Hvb09UlJSREchIpI1lgEZSU1NRfPmzWFtbS06iiQ88sgjKCoqwtmzZ0VHISKSNZYBGUlNTUXLli1Fx5CMJk2aoFGjRlxeSERUR1xNICPjxo1DaWmp6BiS4u/vzzJARFRHLAMyMmjQINERJMff3x/R0dGiYxARyRqXFspEXl4eNm3ahMceewyNGzcWHUcyrl27Bp1OB2dnZ9FRiIhki2VAJo4ePYrg4GDEx8cjICBAdBwiIlIQTiCUidTUVADguQR/Ulpaiscffxy7du0SHYWISLZYBmQiNTUVrq6usLOzEx1FUnQ6HY4fP45Dhw6JjkJEJFssAzLBZYUP1r59e5w+fVp0DCIi2WIZkAlPT0+Eh4eLjiFJAQEBiI+PFx2DiEi2OIGQZG/jxo2IjIzEzZs3uaqAiOghcGRABkpKSnD16lWwt/290NBQfP3116hXr57oKEREssSRARk4fvw4unXrhhMnTiAwMFB0HKIq5eXlITY2FkePHsWlS5ewZs0aAECXLl2g1+vh7OyMFi1aICgoCI8//jhcXFwEJyYi7kAoA5XLCjmB8MF27tyJO3fuIDIyUnQU1crKysILL7yAXbt2oaysDM7OzggLC7v3++PHj0dycjJycnKQkpKC77//Ht27d4eLiwuio6PRoEEDtG/fXuBnQKReLAMykJqaiiZNmsDBwUF0FMnaunUrTp48yTJgZkVFRTh8+DDCw8Ph7OyMu3fvYunSpejfvz9atmwJjUZz77UvvvjifR+r1+vvPdqZM2cO9u3bh549e2Lu3LkIDQ016+dBpHacMyADXFZYvYCAACQkJKCsrEx0FFUwGAxYu3YtWrZsiaFDhyInJwdWVlbYs2cPpkyZglatWt1XBP6OjY0NLC0tAQA///wztm7disLCQoSFhaF///64cuWKOT4VIgLLgCxcvnyZZaAa7du3R3FxMZKTk0VHUbwLFy6gb9++mDhxInr27In4+Pg6r+LQ6XQYPnw4fvvtN2zevBkODg73zuDgtCYi0+MEQhkwGAwoKiqCjY2N6CiSdevWLTg5OWH9+vUYPXq06DiKNn78eOzfvx+rVq1Cv379THqv+Ph4TJo0CV9//TVat25t0nsRqRlHBmRAo9GwCFSjYcOG+Oc//4mmTZuKjqJIRUVFOHXqFABgyZIlOHv2rMmLAABYWFigoKAA3bp1w7Zt20x+PyK1YhmQuPj4eAQGBiI9PV10FMn7/PPP0bt3b9ExFOfKlSsIDQ3FoEGDUFRUBCcnJ9SvX98s927Xrh2OHTuGfv36ISIiAu+88w4qKirMcm8iNWEZkLjExEScPHkSDRs2FB1F8vLz83Hw4EHRMRQlKSkJQUFBuHbtGn744QdYW1ubPUP9+vWxadMmzJ8/H6tXr8bVq1fNnoFI6VgGJC41NRXOzs4sAzXw008/ISwsDDdv3hQdRRHi4uLQo0cP1K9fH0ePHkWnTp2EZdFoNJg1axaSkpLg5uaGO3fuQK/XC8tDpDQsAxLHZYU1V7k7Y1xcnOAkylCvXj106dIFhw4dgpubm+g4AHDv8cTo0aMxcOBAFBYWCk5EpAwsAxKXmpoKb29v0TFkwdvbGw4ODjhx4oToKLJ27tw56PV6tG3bFrt27YKTk5PoSH8xa9YsnDhxAhERESgqKhIdh0j2uLRQ4mJjY1GvXr06n0mQry9FQXH1G/LYW2nhYKOr071E6t27Nxo2bIioqCjRUWTp9OnTePTRRzFhwgQsXrxYdJwqHThwAAMHDkS/fv2wefNm6HTy/XtLJBrLgEpk3dZj46+XkKcvfeBrHG10GN21BVwcq17GWJNiYaEBdJYWKC6rfua3MQvI22+/jXPnzrEMPITU1FSEhISgefPm2Lt3Lxo0aCA6UrV27dqF4cOHY8OGDRgxYoToOESyxTIgYRcuXMDatWvx0ksvoVGjRnW6VtZtPb44mI7cwpIHvqZ5QxuMC34E5RUP/ithaaFBcVkFok5crrJYtGhogz6tm2LrySvVFpBRXVqgWQPuoyBSTk4OgoKCoNFocOTIkTr/fTOnlJQUtGrVSnQMIlnjQUUSFhcXh7lz5+Kll14yy/3qaS2gLymv8g288k0+T19aZbFoYPv7d/rVvc4UKioqUFxczI2aamHz5s24desWYmNjZVUEANwrAp9//jmsra0xbtw4wYmI5IcTCCUsNTUVjo6Odd73vbYq38D/7sedGsw7EMlgMMDDwwNLly4VHUVWnnvuOSQkJMh6suqJEycwadIkHDhwQHQUItlhGZCwymWF1Z3+Rv9Po9GgVatWXFFQQ0uWLMGqVasAAE2aNBGcpm6WL1+OsLAwPPHEE7h48aLoOESywjIgYdxj4OF06tSJZaAG9uzZg1dffRUZGRmioxiFTqfDxo0bUb9+fTz++ONcckhUCywDEta/f38MGzZMdAzZ6dSpEy5evMidCKuQkZGByMhI9O3bF3PmzBEdx2icnJwQFRUFJycnFBQUiI5DJBucQChhr7/+uugIstS5c2cAvx/y1KdPH8FppKeoqAgjRoyAo6Mj1q9fD0tLS9GRjCowMBC//PILgN8nk1pY8HseourwX4lE5ebmIjY2FiUl5p2JrwTe3t64du0ai8AD5Ofnw9nZ+d530Ep1/vx5+Pv74/z586KjEEkey4BEHTx4EEFBQbh165boKLKj0WjQtGlT0TEkyWAwoEmTJtizZw86duwoOo5JeXh4QKPRIDIykvMHiKrBMiBRqampsLOzk/0Mb1HWrVuHXr16iY4hKZmZmQgMDMTZs2dFRzELW1tbbNiwAUlJSXjttddExyGSNJYBiUpOToaPjw+XFT4kKysrHDhwAFevXhUdRRLKy8sxduxY5ObmSuYEQnMICAjAokWL8Mknn+DHH38UHYdIsjiBUKKSkpLg6+srOoZsBQUFAQCOHj3KPesBLFy4EDExMYiOjkbDhg1FxzGrF198EQkJCXB0dBQdhUiyODIgUZaWlggICBAdQ7bc3Nzg7u6OmJgY0VGEO3nyJN555x3861//Qo8ePUTHMTuNRoOVK1eiZ8+eqKio/uAsIjXiyIBE7du3T3QE2QsKCsLRo0dFxxDu7t27GDRoEGbPni06ilB3797F4MGD8dxzzyEyMlJ0HCJJYRmQIK6NNo433nhDdARJCAkJwfbt20XHEM7W1hZNmzbF5MmTERoaimbNmomORCQZfMeRoPXr16NRo0bQ6/Wio8haQECAqh+1xMbGol+/fsjJyREdRTKWL18OnU6HadOmiY5CJCksAxKUlJQEKysrHsFrBLNnz8b3338vOobZFRcX45lnnsHt27c5ce4PnJ2dsWzZMkRFRWHHjh2i4xBJBh8TSBBXEhjP/v37cfbsWYwcOVJ0FLOaO3cuUlNTERcXB62W/8z/aOTIkdiyZQtyc3NFRyGSDH6VkKCkpKR7S+OoboKDg/HNN9/AYDCoZs+G+Ph4fPjhh3jrrbfQtm1b0XEkR6PRYOPGjaJjEEkKHxNITEVFBVJSUjgyYCRBQUHIyspCZmam6Chmc+7cOQQGBmLWrFmio0iaXq/HtGnTcPz4cdFRiITjyIDEaDQanD9/HnZ2dqKjKMIfNx/y8PAQnMY8Ro8ejVGjRnFFSjV0Oh0OHz6MgwcP4tdff4VOpxMdiUgYfrWQGI1GAw8PDzRq1Eh0FEVo3LgxFi9ejHbt2omOYnJJSUl45ZVXcPfuXRaBGtBqtVi1ahXOnDmDpUuXio5DJBS/YkhMVFQUJk2aJDqGWVmY+FH+jBkz4O/vb9qbCFZRUYFnn32W+wnUUqdOnTB9+nTMnj0baWlpouMQCcMyIDEHDhzAkSNHRMcwG9t6ltBoNMi6ra/2R76+9KHucfPmTSxevFjRx0F//vnnOHz4MFavXg1bW1vRcWTlvffeQ5MmTXiQEaka5wxIjNqWFdbTWkBfUo6tJ68gr4o3e0cbHUZ1aQEHm9o/1y0sLMQrr7wCLy8vRERE1CWuJF29ehWvv/46nn32WYSFhYmOIzv29vaIj4/nfgykahwZkJjk5GRVlYFKefpS5BaWPPBHVUWhOh4eHvDy8lLseQ/btm2DlZUVFixYIDqKbDk6OqKiogLfffcd7t69KzoOkdmxDEiIXq9HZmYmfHx8REdRnD59+ii2DLzwwgs4d+4cnJycREeRtczMTEycOBELFy4UHYXI7FgGJKSiogKffvoph3pNoHfv3jh37hyuXr0qOorR6PV6bNiwAQaDgatPjOCRRx7BzJkzsWDBAmRkZIiOQ2RWLAMSYmdnh8mTJ8Pb21t0FMXp3bs3Jk2ahJKSEtFRjGbevHkYN24c37iM6I033oCTkxNmzpwpOgqRWbEMSMjBgwc5o9lEmjRpgi+++EIxGw+dP38eCxcuxOuvvw5PT0/RcRTD3t4eH330EbZs2YJDhw6JjkNkNlxNICErV67ElStXMGTIENFRFOnmzZvYv38/nnzySdFR6sRgMOCFF16Ah4cHtxw2gdGjR6NevXoIDg4WHYXIbDgyICHnzp1DmzZtRMdQrJiYGIwcORLp6emio9TJli1bEB0djZUrV8La2lp0HMXRaDR44oknYGlpyZMNSTVYBiSivLwcSUlJLAMmFBoaCgsLC9mvKnjsscewc+dOhIeHi46iaOvWrUPLli2RnZ0tOgqRybEMSMSFCxdQXFyM1q1bi46iWA0aNEDnzp2xZ88e0VEe2sWLF6HVajFo0CDRURRvwIABAIA333xTcBIi02MZkIiSkhIMGjRI8XvoizZo0CDs3r1blqsKDh06BG9vb1VtVy2Ss7Mz5syZgzVr1uD06dOi4xCZFMuARLRp0wY7d+5Es2bNREdRtIiICISFhcnuWXBJSQmee+45dOnS5d6xzGR6zz33HFq1aoVXXnkFBoNBdBwik+FqAonIyMhAo0aNYG9vLzqKogUEBMjyZL9FixYhOTkZcXFxPJ7YjHQ6HRYtWoSdO3eitLQU9erVEx2JyCRYBiTiySefRLt27bB27dpafVy+vhQFxWVVvsbSQoOyCn5XU6moqAjR0dHo168fNBoTn59sBGlpaXj//ffx8ssvIyAgQHQc1Rk6dCiGDh0qOgaRSbEMSIDBYMD58+cxatSoWn9sQXEZNv56qcqDfFo0tEGf1k3rElFRDh8+jAEDBiAuLg4dO3YUHadaDRo0wOTJkzF79mzRUVRt+fLlcHR0xFNPPSU6CpHRcbxRAi5duoTCwsKHXlZY3Yl/d6oZOVCbsLAwODo6YseOHaKjVKu8vBzOzs5YvHgx7OzsRMdRtdjYWMycORP5+fmioxAZHcuABJw/fx4AuKzQTHQ6HQYOHCj5uQPZ2dlo3bo1Dhw4IDoKAZg/fz7y8/N5VDQpEsuABGRlZcHBwUEx++bLwWOPPYaTJ08iMzNTdJQHmjFjBnJycrgRlUS0aNECM2bMwJIlS3Dp0iXRcYiMimVAAiZMmICcnBzOEjejAQMGICQkBDdu3BAd5W/t3r0b69atw5IlS9CkSRPRceh/Xn/9dTg4OGD16tWioxAZFScQSoRWyz8Kc2rQoAEOHz4sOsbfKiwsxPPPP4/w8HA8/fTTouPQH9SvXx8xMTE8KZIUh9+KCmYwGNCmTRt8++23oqOojsFgwPHjxyX3qODmzZto0aIFPvvsM1ksfVQbb29vWFhYID4+nhsRkWKwDAiWnZ2N8+fPc6a4AKWlpRg4cCBWrFghOsp9PDw8cPDgQXh7e4uOQg9w4sQJdOjQAT/88IPoKERGwTIgWEJCAgBwkpgA9erVw5gxY/DNN9+grEz88suSkhKMGTMGJ0+eFB2FqhEYGIjw8HC8+uqrKC198B4fRHLBMiDYmTNnYGVlhZYtW4qOokrjx4/H1atXJXGS4Zw5c7Bp0yY+GpABjUaDRYsWISUlBV988YXoOER1xjIg2OnTp9GmTRtOIBQkMDAQbdu2xZdffik0x/Hjx/Hhhx/inXfeQYcOHYRmoZpp3749xo8fj9mzZyMvL090HKI64TuQYLNnz5bs8jY10Gg0ePnll3Hx4kVhGfR6PcaNG4eOHTvi9ddfF5aDam/u3Llo3LgxR3NI9lgGBGvRogVatGghOoaqPfPMM0Lvn5KSgoKCAkRFRUGn0wnNQrXTrFkz7khIisDHBAJlZmbimWeekdzSNjXKzc3F/PnzUVRUZPZ7BwQEIDU1lZNIZWzlypWYMmWK6BhED41lQKATJ07gyy+/5HeDEnDjxg289dZb+Oqrr8x2zytXrmDcuHHIzc2FlZWV2e5LxqfRaLBixQquBCHZYhkQ6MyZM2jUqBFcXFxER1E9X19fjBw5Eh9++KFZloqVl5dj7Nix2Lt3LzeuUYBnn30Wvr6+ePXVV/nnSbLEMiDQmTNn0K5dO04+kog33ngDFy9exPr1601+r3fffReHDh3CunXr4OzsbPL7kWlptVosWLAAe/fuxe7du0XHIao1lgGBTp8+jXbt2omOQf8TEBCAYcOG4YMPPkBFRYXJ7vP9999j3rx5mD9/PsLCwkx2HzKvoUOHIjQ0FNu2bRMdhajWuJpAoGnTpiEgIEB0DPqDefPm4datWyY9QTIrKwtPPfUUXnvtNZPdg8xPo9Fgx44dcHBwEB2FqNZYBgR68cUXRUegP2nbti2A35/pFxQUwNHR0WjXLioqgrW1NaZPnw6DwcDHQwpU+fflv//9L4KCgnjmCMkGHxMIcurUKWzfvl10DFmxMON75+jRozFmzBijTQa7efMmOnfujJUrVwIAi4CCXb16FYMGDcKSJUtERyGqMZYBQb755hvMnDlTdAzZsK1nCY1Gg6zb+mp/5OvrvhrgH//4B3bt2oXvvvuuzte6ffs2+vXrh+zsbPTq1avO1yNpc3V1xdSpU7Fw4UJcv35ddByiGuFjAkEqVxJQzdTTWkBfUo6tJ68gr4o3e0cbHUZ1aQEHm7rt3TB06FCMGTMGL730Erp27QpfX9+Huk5ubi6GDBmCjIwMHDhwAH5+fnXKRfLw5ptvYu3atXjvvfckd0Q20d/hyIAgZ86c4eTBh5CnL0VuYckDf1RVFGpr2bJlcHV1Ra9evXD79u2HusZrr72G5ORk7N69m3/eKuLk5IS33noLX3zxBRITE0XHIaoWRwYEyM7OxvXr1zkyIHHOzs7Yt28fdu3ahQYNGtTqY/Py8uDo6IiFCxfijTfegJeXl2lCkmRNmTIFBoMBzZs3Fx2FqFocGRDg9u3bCA8PR6dOnURHoWo0bdoUEyZMAAB89NFHWLt2LcrLyx/4+qSkJAwfPhw+Pj64du0anJycWARUysrKCq+88grs7e1Num8FkTGwDAjg4+ODPXv2wNPTU3QUqiGDwYCEhARMnDgRHTt2xNdff43jx4/f27p43bp1eOqpp+Dv749Tp05h6dKlaNKkieDUJAWrVq1C7969WQhI0lgGBEhPTxdyOh49PI1Gg6+++gqxsbFo2LAhxo8fj27dukGv1wMAVq9ejbNnz2LBggVITEzEmDFjTLpxEclHq1atEB0djU2bNomOQvRAnDMgQHh4OCIiIrB48WLRUaiWunXrhujoaOTl5SEjI+PebnN79uyBVquef075+lIUFJdV+zp7K22dV3bI3aOPPoohQ4Zg1qxZGD58OE+oJElSz1cvibh16xYuXLiAwMBA0VGoDhwdHdG+fft7P1dTEQCAguIybPz1klmWeSrBggUL0K5dO6xYsQIvv/yy6DhEf6Gur2AScOrUKQBgGSDZq1zmSdVr06YNnn32WcTHx4uOQvS3WAbMLC4uDra2tvDx8REdhYjMaPny5aobQSL54AwnM7t69So6duwIS0tL0VGIyIy0Wi0MBgO2bNmC9PR00XGI7sOaamaLFi1CWVn1E6+ISHmKiorw0ksvoVu3bti8ebPoOET3cGTAjCpPwONQIamFOU+alAMbGxt8+OGHiIqKQnR0tOg4RPewDJjR0aNH4ebmhgsXLoiOQmRy5j5pUi5Gjx6Nrl274uWXX65yN0sic+K3qGYUFxeHmzdvws3NTXQUIpMz90mTcmFhYYF///vfCA4OxoYNGzB27FjRkYhYBswpLi4Obdu2Rb169URHITIbLkH8q6CgIGzfvh39+/cXHYUIAMuAWcXFxaFLly41em1NdniztNCgrMJgjGhEZGbDhg0D8P8nXBKJxDJgJnfv3kVCQgKef/75Gr2+Jju8tWhogz6tmxorIhGZ2a5duxAZGYkzZ87Aw8NDdBxSMZYBM7GxscG5c+fQsGHDGn9MdcOrDWzV8YyVSKlCQ0Nhb2+Pf/3rX9iwYYPoOKRiLANmotFo0KpVK9ExiKrEx1PmZW9vj/nz52PChAmYOnUqQkJCREcilWIZMJNZs2bBxcUFL730kugoise17Q9P6o+nlHha4tNPP41PP/0U06ZNw/Hjx7k7KQnBMmAm3333HZcQmcEf17ZXR05vGOYk5cdTSjwt0cLCAp9++ik+/vhjFBYW3jsWm8icWAbM4PLly7h8+TK6desmOoricW27PNVmNEeJSxW7d++O7t27i45BKsYyYAbHjh0DAJYBM1LiG4ZS1XQ0Rw1zFTZv3oxTp05h7ty5oqOQyrAMmEFsbCzc3d3RrFkz0VGIJKemozlqWEp79epVzJ8/H8OHD0fnzp1FxyEVYRkwg7Fjx6JHjx6iYxBJmpTnKpjLCy+8gFWrVmHKlCmIiYmBhQWPjyHz4N80M+jQoQMee+wx0TGIVEOuK0q0Wi0+/fRTHDt2DF9//bXoOKQiHBkwsZSUFGzcuBFTp07llqNEZlDTOQgWGkBnaYHisooqX2fuVSehoaEYO3YsvvvuO0yYMMFs9yV1YxkwsT179uC9997DzJkzRUchUoXazkGo6nWiVp0sX74cdnZ2Zr0nqRvLgIkdO3YMHTp0gI2NjegoRKpS0zkIUlx5UjmKePLkSeh0OrRt21ZwIlI6lgETO3LkCAYOHCg6BhHJjMFgwDPPPAOtVovY2FjuTEgmxTJgQleuXEFaWhrCwsJERyGV45kD8qPRaLBixQqEhITg008/5VbmZFIsAyZkMBgwffp0hIaGio5CClXTN/nisgpEnbis+nX8chMUFITJkyfjzTffREREBNzd3UVHIoViGTCh5s2bY+nSpaJjkILV5mAhruOXp/nz52Pr1q2YNWsW1q1bJzoOKRTLgAn95z//QWhoKNzc3ERHob8h17Xof8Y3eWVzcHDAli1b4OnpKToKKRjLgIlcu3YNY8aMwX/+8x9ERkaKjkN/wtMNSU4qzzXJycmBVqvlniVkdCwDJnLw4EEA4ORBieLphiQ3ZWVl6Nq1K/r27YvPPvtMdBxSGG5HbCLR0dFo1aoVXF1dRUehKlQOsT/oR1VFgcictFotZsyYgc8//xyHDx8WHYcUhmXARKKjozkqQERG9fzzz6N79+6YNGkSioqKRMchBWEZMIGysjK0b98eQ4YMER2FiBTE0tISq1evRnp6Ot5++23RcUhBOGfABLRaLZcAEZFJ+Pv7Y9WqVWjTpo3oKKQgLAMmcObMGTRr1gzOzs6ioxCRAj399NMAgNLSUpSVlfHsE6ozPiYwgSeffBKzZs0SHYOIFMxgMKBv376YMWOG6CikACwDRnbhwgUkJSXxcCJ6oHx9KbJu66v8kc9VDJIh1c2pNBoNRo0ahc8++wy7d+8WHYdkjo8JjOznn3+GVqtFnz59REchiapuC2FHGx1Gd+XeBlIg9c2pnn/+eWzduhUTJ07EmTNn0LBhQ7Pen5SDZcDIfvrpJ4SEhMDBwUF0FJKwqrYQrukbEE8ZNL3abE4losBpNBqsWbMG7dq1w7Rp0/Dtt9+a9f6kHCwDRqbVavHYY4898Pd5lCxVp6ZvQDxl0HyqO/9B5AhCixYt8Pnnn6O4uNho1yT1YRkwsi1btlT5+7U5ZY7UjQcQyYfo7a1HjRp1778LCwthZ2dn1OuT8rEMGFFGRgbc3Nyg01X9D51f5ImUqbp/26b25ptvYt++fTh48GC1X4eI/oirCYxo6NChePHFF0XHICKVGjZsGH777Te88cYboqOQzLAMGMnly5dx9uxZriIgImG6deuGDz/8EIsWLcLOnTtFxyEZYRkwkp9//hkWFhbo27ev6ChEpGIzZszAkCFD8PTTT+PKlSui45BMsAwYybZt2xAcHAwnJyfRUcjIpLrpDNHf0Wg0+OqrrzB16lQ0btxYdBySCU4gNILy8nLcvHkTY8aMER2FjEzqm84Q/R1nZ2fMnj0bAHDt2jW4uLiIDUSSxzJgBJaWloiNjUVFRYXoKGRkUt90huTHnCNNcXFxCAkJQVRUFAYNGmS+G5PssAwYQUZGBjw8PGBhwacuSmWsTWe4oZS6mXukqUOHDujbty/GjBmD48ePw8fHp07XI+ViGaij/Px8+Pn5YdmyZfjnP/8pOg4Jwl0DqSbMvTmRhYUFvv32W3Tr1g3Dhw9HbGwst0qnv8UyUEc7d+5EcXExBgwYIDoKSQA3lKKaMOfmRI6Ojti+fTu6du2KGTNmYPXq1Wa5L8kLy0Adbd68GV27doW7u7voKEREf8vX1xc7duxAq1atREchieJD7jooKCjArl278MQTT4iOQkRUpbCwMDRr1gw3btzADz/8IDoOSQzLQB2kp6fDyckJjz/+uOgoREQ1smzZMowYMQL79+8XHYUkhGWgDgICApCZmQkvLy/RUYiIauSdd95Br169MGLECCQkJIiOQxLBMvCQbt68icuXL8PS0lJ0FCKiGtPpdNi0aRPc3d0xYMAAXLp0SXQkkgCWgYe0fPly+Pv7o7i4WHQUIqJacXR0xE8//QQrKyscO3ZMdBySAK4meAgVFRX48ssv8cQTT8DKykp0HCJSIFPvVNisWTMkJCTAysoKBoMBpaWlqFevnmlvSpLFMvAQ9u7di4sXL2LixIn3fi1fX4qC4rIqP467zxFRTZhrp8LKb2bef/99REdH48cff4SNjc1DXYvkjWXgIaxduxZ+fn4ICgq692sFxWXY+Osl7j5HRHVm7p0Ke/XqhQ8//BBPPPEEtm7dyhECFWIZqCWDwQBbW1tMnjwZGs3943jcfY6IjMlcOxX27NkTO3bswJAhQzBq1Chs3LiRhUBlOIGwljQaDdasWYOpU6eKjkJEZDTh4eHYsmULdu3ahYULF4qOQ2bGMlALRUVF+Prrr7mCgIgUadCgQYiOjsbMmTNFRyEzYxmohW+//RYTJkxAZmam6ChERCbRvXt32NjYID4+HhEREcjLyxMdicyAZaCGysvLsWjRIkRERPCwDyJSvMLCQhw4cABhYWG4evWq6DhkYiwDNbRjxw4kJyfjtddeEx2FiMjkgoODcejQIdy4cQMhISFISUkRHYlMiGWgBgwGAxYsWIDQ0FB069ZNdBwiIrNo27YtYmJiYGVlhcGDB6OsrOq9VEi+uLSwhl566SW4uLiIjkFEZFYeHh6IiYlBWloatFotysrKoNXyrUNp+CdaDYPBAI1Gg9GjR4uOQkQkRMOGDdG5c2dUVFQgIiICPj4+WLBgAUuBgvAxQTXWrFmDkSNHcniMiCTJ1GcY/JFGo0Hfvn2xbNkyhIeH4/r16+a7OZkUy0AVCgoK8Pbbb0On07EBE5Hk/PEMg+p+5FexrXFNaTQaTJs2Dfv27UNiYiICAwNx8OBBI3wmJBrf4aqwaNEi3Lp1C/PmzRMdhYjoL8x9hkGlnj17Ii4uDqNHj8bhw4cRGhpqlOuSOCwDD5CSkoKFCxdi2rRpeOSRR0THISJ6IHOdYfBHzZo1w759++79/Msvv0RoaCi8vb3NmoOMg48JHuDo0aPw8PDAu+++KzoKEZEkWVpawtLSEnfv3sW8efMQEBCA5cuXo6KiQnQ0qiWWgQd4+umnER8fDzs7O9FRiIgkzdbWFqdOncK4ceMwZcoU9OnTB4mJiaJjUS2wDPzJ+fPn8d5776G0tJRHeBKRYph61YG9vT1WrFiBvXv34tKlS1iwYIFpb0hGxTkDf3D37l2MHTsWer0er732GnQ640y2ISIS6Y+rDqpjb6Wt00TD3r1748yZM/dOd920aRMMBgOefPJJaDRmXAdJtcIy8D8VFRUYP348kpKScOTIEdjY2IiORERkFOZedWBjY3Pva+ju3buxZs0afPzxx1i8eDG6d+9ep2uTabAM/M+cOXOwadMmREVFoUOHDvd+PV9fioLiqjccsrTQoKzCYOKERER1I2LVwerVqxEZGYmZM2ciKCgIgwcPxjfffAMnJyez5qCqsQzg91GBpKQkzJs3DyNGjLjv9wqKy7Dx10tVtukWDW3Qp3VTU8ckIjI5U8wtCA8PR1xcHL7//nts3rwZDRo0APD7HK3WrVsb/4ZUa6ovAzdu3EDjxo2xfv36B76mujbdwJZzC4hI/mozt6C+lRb1a/E4wdLSEqNHj753zsvp06fRvn17hIWFYfr06RgyZAh3ehVI1f/n16xZg+nTp+P48eNsp0SkejWdW9DUwQojApvjTjWPUIEHT0j09/dHVFQUPvroI0RERMDNzQ0zZszAjBkz6vQ50MNRZRmoqKjAwoULMWvWLLzwwgvw9fUVHYmISDJqMhpa1wmJlpaWGDFiBEaMGIGTJ0/i888/R05ODgAgNzcXx44dQ9++fTlaYCaq+79869YtjBs3Dj/88APeeustzJkzh8tdiIgegrEmJHbs2BGfffbZvZ9v27YNEydOhLOzMyIiIvDkk0+iV69eXO5tQqrbdKiwsBDJycn48ccf8f7777MIEBFJzIQJE3DixAlMmjQJ+/fvR//+/TFp0iQAQFFREe7cuSM4ofKoYmQgKysLs2fPxvz589G8eXMkJCSgsKSi2kkyXDJIRFQ3D7M6QaPRIDAwEIGBgZg/fz7i4+NhYfH7964//vgjRo8ejaCgIISHh6Nv377o0qULHyfUkaL/7129ehX//ve/sWLFCtjY2GDcuHFo1KgRLC0tUVBcwiWDREQmZIydDzUazX17vwQFBeHjjz/Gf//7XyxZsgTvvvsuhg0bhu3bt0Ov12P//v3o0qULGjdubMxPRfEUWwZWrlyJ6dOnw9raGpMnT8brr7+Ohg0b3vcaLhkkIjIdU+x86ObmhsmTJ2Py5MkoKyvDiRMn7p2SePLkSQwePBgA4Onpia5du6JLly6YOXMmAMBgMPDR8AMorgxU7hjo7uOPV2a9jX9MeBYOjo7QA9D/r51y+J+IyHxqMtHwYR4naLVadOvW7d7Pg4KCkJ6ejmPHjuHYsWP49ddfce3aNcycORMGgwHNmzdH06ZN0bp1a3h5ecHT0xMjRoy4twmSmimuDPz/joHO0AVGYEN8DoCc+17D4X8iIumo6eOE6g5R0mg08PT0hKenJyIjI+/7vfLycrz66qtISEhAYmIioqOjceXKFYSHh6NBgwaYMGEC9u/fD1dXV7i4uMDFxQUjR45Er169cOPGDaSnp8PFxQVOTk6wt7dX3AiD4soAwOF/IiI5qcnjhLoeoqTVajF9+vT7fq2oqOjeUfXDhg2Dq6srrl27hmvXruHYsWMICgoCAPzyyy946qmn7n2cpaUlQkNDsW/fPhgMBvTv3x/169eHo6MjHBwcYGtri5kzZ8LZ2Rn79u1DZmYmbG1tYWNjA1tbW7Rs2RIeHh64e/cusrOzYWVlBUdHR9ja2j7U52YMsikDtTkwyLGavyz1rX7/tPm62r9OytmU8jopZ1Pb66ScTYmvq05NvxmvyfuFhQbQWVogt+j3o5a79RqAbr0G/OV11/L0COndH78cjMWN7OvIy8tDft5t1K/vgKzbepSWlsK2fgPcyruNS5evoLDgDvR6PSZPngwAWLVqFTZs2HDfNd9//3289dZb2L9/P4YMGQIAWLhwIV599dWafYImoDEYDHx4TkREpGKq23SIiIiI7scyQEREpHIsA0RERCrHMkBERKRyLANEREQqxzJARESkciwDREREKscyQEREpHIsA0RERCr3f/XNqIyIy+nqAAAAAElFTkSuQmCC", | |
"text/plain": [ | |
"<Figure size 640x480 with 1 Axes>" | |
] | |
}, | |
"metadata": {}, | |
"output_type": "display_data" | |
} | |
], | |
"source": [ | |
"\n", | |
"fig, ax = plt.subplots() \n", | |
"xvals = np.arange(0, 151)\n", | |
"pvals = p(xvals)\n", | |
"ax.hist(keep, 45, density=True, alpha=.65, edgecolor=\"#FFFFFF\", linewidth=1.0)\n", | |
"ax.plot(xvals, pvals, linewidth=1., color=\"#000000\", linestyle=\"--\")\n", | |
"ax.set_axis_off()\n", | |
"plt.show()" | |
] | |
} | |
], | |
"metadata": { | |
"kernelspec": { | |
"display_name": "py312", | |
"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.12.3" | |
} | |
}, | |
"nbformat": 4, | |
"nbformat_minor": 2 | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment