\n",
"\n",
"Launch in Binder [![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/esowc/UNSEEN-open/master?filepath=doc%2FNotebooks%2Fexamples%2FCalifornia_Fires.ipynb)\n",
"\n",
"\n",
"
"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# California fires\n",
"\n",
"In August 2020 in California, wildfires have burned more than [a million acres of land](https://edition.cnn.com/2020/10/06/us/gigafire-california-august-complex-trnd/index.html). \n",
"The wildfires coinciding with record high temperature anomalies, see [August 2020 temperature anomaly](../California_august_temperature_anomaly.ipynb).\n",
"![California Temperature August 2020](../../../graphs/California_anomaly.png)\n",
"\n",
"\n",
"In this example, we evaluate the UNSEEN ensemble and show that there is a clear trend in temperature extremes over the last decades. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Retrieve data"
]
},
{
"cell_type": "markdown",
"metadata": {
"nbsphinx": "hidden"
},
"source": [
"
\n",
"\n",
"Note\n",
" \n",
"In this notebook you cannot use the python functions under Retrieve and Preprocess (they are here only for documentation on the entire workflow, see [retrieve](../1.Download/1.Retrieve.ipynb) if you want to download your own dataset. The resulting preprocessed dataset is provided so you can perform statistical analysis on the dataset and rerun the evaluation and examples provided.\n",
" \n",
"
"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The main functions to retrieve all forecasts (SEAS5) and reanalysis (ERA5) are `retrieve_SEAS5` and `retrieve_ERA5`. We want to download 2m temperature for August over California. By default, the hindcast years of 1981-2016 are downloaded for SEAS5. We include the years 1981-2020. The folder indicates where the files will be stored, in this case outside of the UNSEEN-open repository, in a 'California_example' directory. For more explanation, see [retrieve](../1.Download/1.Retrieve.ipynb)."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"import sys\n",
"sys.path.insert(0, os.path.abspath('../../../'))\n",
"os.chdir(os.path.abspath('../../../'))\n",
"\n",
"import src.cdsretrieve as retrieve\n",
"import src.preprocess as preprocess\n",
"\n",
"import numpy as np\n",
"import xarray as xr"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"retrieve.retrieve_SEAS5(\n",
" variables=['2m_temperature', '2m_dewpoint_temperature'],\n",
" target_months=[8],\n",
" area=[70, -130, 20, -70],\n",
" years=np.arange(1981, 2021),\n",
" folder='../California_example/SEAS5/')"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"retrieve.retrieve_ERA5(variables=['2m_temperature', '2m_dewpoint_temperature'],\n",
" target_months=[8],\n",
" area=[70, -130, 20, -70],\n",
" folder='../California_example/ERA5/')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Preprocess"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the preprocessing step, we first merge all downloaded files into one xarray dataset, then take the spatial average over the domain and a temporal average over the MAM season. Read the docs on [preprocessing](../2.Preprocess/2.Preprocess.ipynb) for more info. "
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"scrolled": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Lead time: 07\n",
"6\n",
"5\n",
"4\n",
"3\n"
]
}
],
"source": [
"SEAS5_California = preprocess.merge_SEAS5(folder ='../California_example/SEAS5/', target_months = [8])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And for ERA5:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"ERA5_California = xr.open_mfdataset('../California_example/ERA5/ERA5_????.nc',combine='by_coords')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We calculate the standardized anomaly of the 2020 event and select the 2m temperature over the region where 2 standard deviations from the 1979-2010 average was exceeded, [see this page](../California_august_temperature_anomaly.ipynb). This is an area-weighed average, since grid cell area decreases with latitude, see [preprocess](../2.Preprocess/2.Preprocess.ipynb). "
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"ERA5_anomaly = ERA5_California['t2m'] - ERA5_California['t2m'].sel(time=slice('1979','2010')).mean('time')\n",
"ERA5_sd_anomaly = ERA5_anomaly / ERA5_California['t2m'].sel(time=slice('1979','2010')).std('time')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We use a land-sea mask to select land-only gridcells:"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"LSMask = xr.open_dataset('../California_example/ERA_landsea_mask.nc')\n",
"# convert the longitude from 0:360 to -180:180\n",
"LSMask['longitude'] = (((LSMask['longitude'] + 180) % 360) - 180)"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"area_weights = np.cos(np.deg2rad(ERA5_sd_anomaly.latitude))\n",
"\n",
"ERA5_California_events = (\n",
" ERA5_California['t2m'].sel( # Select 2 metre temperature\n",
" longitude = slice(-125,-100), # Select the longitude\n",
" latitude = slice(45,20)). # And the latitude\n",
" where(ERA5_sd_anomaly.sel(time = '2020').squeeze('time') > 2). ##Mask the region where 2020 sd >2. \n",
" where(LSMask['lsm'].sel(time = '1979').squeeze('time') > 0.5). #Select land-only gridcells\n",
" weighted(area_weights).\n",
" mean(['longitude', 'latitude']) #And take the mean\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Plot the August temperatures over the defined California domain: "
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"[]"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAESCAYAAAAG+ZUXAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3deXydZZnw8d91TvY9zdI2Sds0Kd3ohpSlBaGA4gIDuOCIyqCoiOLCjMu4vu+or6OODqMzyjAuICIiyiqICzAsUmqlrd1o2tKkLc3S5CRplnOSc5KTc79/PM8TTtMsJ8nZe30/n36anC13nyZX7nPd133dYoxBKaVUenElegBKKaWiT4O7UkqlIQ3uSimVhjS4K6VUGtLgrpRSaUiDu1JKpaGMRA8AoLy83NTW1iZ6GEoplVK2b9/eaYypGO++pAjutbW1bNu2LdHDUEqplCIiRye6T9MySimVhjS4K6VUGtLgrpRSaUiDu1JKpSEN7koplYY0uCulVBrS4K6UUgkSy5brGtyVUipB3v7fL/KFh3bH5LU1uCulVAIYYzhwvJ+cTHdMXl+Du1JKJcDxPj8DQyPUVxTE5PU1uCulVAI0dvgAqKvIj8nra3BXSqkEaOr0ArBEZ+5KKZU+Gju8FGRnUFGYHZPX1+CulFIJ0NTpo74iHxGJyetPGdxFJEdE/ioiu0TkZRH5qn37HBF5UkResf8uDXvOF0TkkIgcEJE3xWTkSimVwho7vNTFKCUDkc3cA8Clxpi1wDrgzSJyPvB54GljzBnA0/bniMhK4N3AmcCbgdtFJDa1PkoplYIGhoK09vqpj9FiKkQQ3I3Fa3+aaf8xwNXA3fbtdwPX2B9fDfzKGBMwxhwGDgHnRnXUSimVwpo8TqVMYmfuiIhbRHYCHcCTxpitwFxjTBuA/Xel/fBq4FjY05vt25RSSmHl2yF2ZZAQYXA3xowYY9YBNcC5IrJqkoePtzpwSgMFEblJRLaJyDaPxxPZaJVSKg00dngRgdqyBAd3hzGmB3gWK5feLiLzAey/O+yHNQMLwp5WA7SO81o/MsasN8asr6gY93xXpZRKS02dPmpKc2PWegAiq5apEJES++Nc4A3AfuC3wA32w24AHrU//i3wbhHJFpHFwBnAX6M9cKWUSlWNHd6YtR1wZETwmPnA3XbFiwv4tTHmcRHZAvxaRD4IvApcC2CMeVlEfg3sA4LALcaYkdgMXymlUksoZDjc6eP8urKYfp0pg7sxZjdw1ji3dwGXTfCcbwDfmPXolFIqzbT1+RkcHqG+Mnb5dtAdqkopFVdNHquyvK48tmkZDe5KKRVHjR1WcNeZu1JKpZGmTh+F2RlUFMSmYZhDg7tSSsVRo8dLXWVBzBqGOTS4K6VUHDV5fNSXxzYlAxrclVIqbnyBIG29fuorY7uYChrclVIqbg47PWV05q6UUumj0eNUyujMXSml0kajx4dLYFFZXsy/lgZ3pZSKk0aPl5rSPLIzYn9+kQZ3pZSKkyaPL6anL4XT4K6UUnFgNQyL7bmp4TS4K6VUHLT2DuIfDsW81a9Dg7tSSsVBoyf2R+uF0+CulFJx4HSD1Jm7UkqlkUaPl8KcDMoLsuLy9TS4K6VUHFiVMrFvGObQ4K6UUnHQ6PHGLd8OGtyVUirmvIEg7X2BuOXbQYO7UkrF3GuLqTpzV0qptNFkl0HqzF0ppdJIk8eLS2BhHBqGOTS4K6VUjDV6fCycE5+GYQ4N7kopFWNWpUz8UjKgwV0ppWLKahgWv26QDg3uSikVQy09gwSCIZ25K6VUOnGO1ovHuanhNLgrpVQMjZZBxuHc1HAa3JVSKoYaPV6KcjIoy49PwzCHBnellIqhJo+P+sr4NQxzaHBXSqkYavR4qSuPb0oGNLgrpVTM9PuH6egPUF8Z38VU0OCulFIx4yym6sxdKaXSSFNn/LtBOjS4K6VUjDR2+HC7hEVlSRjcRWSBiDwjIg0i8rKIfMq+fa2IbBGRPSLymIgU2bdnichd9u27RGRTjP8NSqkU8vxBD79+6ViihxEXTZ1eFpTmkpUR/3l0JF8xCHzaGLMCOB+4RURWAj8BPm+MWQ08DHzWfvyHAezb3wj8u4joOwSlFAA/eeEwX350L/3+4UQPJeaaPL64tx1wTBl0jTFtxpgd9sf9QANQDSwDnrcf9iTwDvvjlcDT9uM7gB5gfXSHrZRKVW09gwwFQ/zp5fZEDyWmnIZh8W474JjWjFpEaoGzgK3AXuAq+65rgQX2x7uAq0UkQ0QWA2eH3aeUOo0ZY2jtGQTgsd2tCR5NbCWqYZgj4uAuIgXAg8Ctxpg+4EasFM12oBAYsh96J9AMbAO+B7yIldoZ+3o3icg2Ednm8Xhm969QSqWEPn8Q39AIRTkZvPBKJ92+oamflKKaOp2j9ZJ45i4imViB/V5jzEMAxpj9xpjLjTFnA/cBjfbtQWPMPxpj1hljrgZKgFfGvqYx5kfGmPXGmPUVFRXR+vcopZJYW681a3/f+YsIhgxP7GlL8Ihip7HD7gaZrDN3sRoi/BRoMMbcFnZ7pf23C/gycIf9eZ6I5NsfvxEIGmP2xWDsSqkU09bjB+CyFZXUV+Tz2K70Tc00dXopzMmgvCC+DcMckczcLwCuBy4VkZ32n7cC14nIQWA/0ArcZT++EtghIg3AP9vPVUopWu2Ze1VJLletreavR7o53utP8Khiw6mUiXfDMEfGVA8wxrwATDS674/z+CNYlTRKKXWSth4/bpdQWZjD362dz388dZDHd7fyodfXJXpoUdfk8bGxvixhX1/rz5VScdPaO8jcwmzcLqGuooBV1UVpmZrxBoIc7/PH/YCOcBrclVJx09bjZ35J7ujnV62tYldzL0fsypJ0cXi0YVhiKmVAg7tSKo7aegeZX5wz+vmVa6oA0m727jQMS1SlDGhwV0rFiTGGtl4/VWEz96qSXM6pLU27DU2NHh8isKgsL2Fj0OCulIqLbt8QgWDopJk7WKmZg+1e9h/vS9DIoq/J46WmNJecTHfCxqDBXSkVF212yeP84tyTbn/L6vm4XcJvd6bP7L3J46M+gSkZ0OCulIoTp6dMVcnJM/fygmw21pfx2O5WjDGJGFpUhUKGps7EnJsaToO7UiouJpq5g5WaOdY9yM5jPfEeVtS19fnxD4eoS1BPGYcGd6VUXLT2DJLldlGWf+p2/MvPnEeW28Vv06BqpsnjVMpocFdKnQZae/3MK87B5Tp1w3txbiabllXwu91tjIRSOzXjHIqtOXel1GmhrWfwlEqZcFetq6KjP8DWw11xHFX0NXm8FGRnUFmYndBxaHBXSsXF2Br3sS5bPpe8LHfKb2hq9Pioq8hPWMMwhwZ3pVTMjYQMx/v8p1TKhMvNcvPGlXP5/d7jDAVDcRxddDV5vAltO+DQ4K6UijlPf4CRkBm3UibclWuq6BkYZtuR7jiNLLoGhoK09voT2nbAocFdKRVzr/Vxn3jmDrCmphiAQ3bFSao5bDdAS3SlDGhwV0rFgXMC01Qz98rCbPKy3KNBMtUkS6UMaHBXSsWBc3Zq1RTBXUSoLctP2eDe6PEiAos1566UOh209vjJy3JTlDvl4W8sLs9P2f7uTR4fVcWJbRjm0OCulIo5p497JOWBi8vzOXZikOGR1KuYaer0JkW+HTS4K6XioHWKGvdwteX5jIQMx7oHYjyq6DLGcDgJukE6NLgrpWJuqt2p4Zx8darl3dv7AviGRqjXmbtS6nQwFAzh8QamrJRxpGpwf61hmM7clVKngfY+P8ZMXePuKM3LpDg3M+WCe2OSdIN0aHBXSsWUc0hHpDN3EaG2PJ8jXakW3H3kZbmZVxTZL7FY0+CulIop55COSGfuAHXl+Rz2pFZwb+r0sbg88Q3DHBrclVIx5bQeiHTmDlBblk9rrx//8EishhV1TR5v0lTKgAZ3lUDHe/30+4cTPQwVY209fopzM8nPnnoDk2OxnbdOldSMf3iElp7BpMm3gwZ3lUDvv+uvvP32F+kd1ACfzpwNTNOxuMwO7imyqHq404cxyVMpAxrcVQI1nxjklQ4vH/3F9pTu360m19oT+QYmR215HmDlsaNlYCgYsyP8nIZhydDH3aHBXSWEf3gEbyDI2ppiXmzs4osP78GY1D47U41vJjP3wpxMyguyozZzN8ZwyXef5c4XDkfl9cZKlkOxw2lwVwnR7RsC4LpzF3LrG87gge3N/Nf/HkrwqFS0DQ6NcGJgeNozd7ArZqIU3Pv8Qdr7Auxu6Y3K643V1OmjqjiHvKzI1xViLXlGok4rTnCfk5/F35+zgFe7B7jtyYMsmJPL286qSfDoVLS0jVbKTL/2u7Y8j//d74nKODq9AQBejVG/miaPN6ny7aAzd5Ugzg9bWUE2IsK33r6GDXVlfO6B3Wxp7Erw6FS0ODXu0ymDdCwuL6DTG4hKRZWn3/p+i0UzMmPM6KHYyUSDu0qILq81cy/LzwIgK8PFHdefTW1ZPh+5ZxuHOvoTOTwVJc7u1OlsYHIsthdVj3TOPiA7k4lu31DUy289/QG8gWBSLaaCBneVIE5apqwga/S24txM7nz/OWRluHn/XS+NzrZU6nJm7vNmkJZZXG6lOZo6Z3+eamfY99Kx7sFZv164RqdSRtMySkGnL0CW20XBmI0tC+bk8dMb1tPpDfDhn2+LWemaio+23kHKC7LIzpj+yUSLyqI5cx8a/TjaeXfnl099ZYoFdxFZICLPiEiDiLwsIp+yb18rIltEZI+IPCYiRfbtmSJyt317g4h8Idb/CJV6ur1DlBVkjduHY+2CEr741hXsPNbDoY7Zz9pU4rT0+GeUbwfIyXRTXZLL4WjM3L0Bcu2j76Kdd2/y+MjJdDE/SRqGOSKZuQeBTxtjVgDnA7eIyErgJ8DnjTGrgYeBz9qPvxbItm8/G/iIiNRGe+AqtXX5hpiTnzXh/RvrywHY1dwTryGpGJjOIR3jqS3P43BXdHLui8ryKMrJiPrMvaGtj7ryAlyu5GgY5pgyuBtj2owxO+yP+4EGoBpYBjxvP+xJ4B3OU4B8EckAcoEhoC/K41Yprss3RFlB9oT315XnU5idwW4N7jEzODTC7/e0xXTzWNs0jtcbz+LyfA57vLMeo8c7REVhNovK8qMa3H2BINuOnGBjfVnUXjNappVzt2fgZwFbgb3AVfZd1wIL7I8fAHxAG/Aq8F1jTHcUxqrSSJc3MFopMx6XS1hVXczu5thsOlHw2K5WPnrvDl5ujc3cq88/jDcQnN3MvSyfPn+QEwOzq3Dp7A9QUZDNwjl5UU3LbGnsYmgkxKZllVF7zWiJOLiLSAHwIHCrMaYPuBErRbMdKMSaoQOcC4wAVcBi4NMiUjfO690kIttEZJvHE52NCip1dPuGJg3uYOXeG9r6CARTp+1rKnFmsHtitGuzrcfp4z7zmbtTOz6bvLsxhk5vgPLCbBbMyaP5xGDUFuqfPdhBXpabcxaXRuX1oimi4C4imViB/V5jzEMAxpj9xpjLjTFnA/cBjfbD3wP8wRgzbIzpADYD68e+pjHmR8aY9caY9RUVFdH4t6gUMTg0wsDQyKRpGYC1NcUMjxga2rTmPRZa7Br0WL07cvq4z6TG3eGUQx6eRcWMNxAkEAxRXpDFwjl5DI2EaO/zz/j1HMYYnj3gYWN92YyqgWItkmoZAX4KNBhjbgu7vdL+2wV8GbjDvutV4FKx5GMtwu6P9sBV6ury2btTp5i5r1lQAqB59xhpOWEF370xnrnPtFoGoKY0F7dLZjVzd8ogy+20DMDRKCzSNnp8NJ8Y5OIkTMlAZDP3C4DrsQL2TvvPW4HrROQgVuBuBe6yH/9DoAArJ/8ScJcxZnf0h65S1eju1ILJg3tVcQ7lBVnsOqZ591hwZu77j8cm9dXWO4hLoLJw8ndok8l0u1g4J29Wte7O7tTw4B6NvPuzBzoA2LQ0OTMPUzYOM8a8AExU4/P9cR7vxVpgVWpczsx9slJIsA5KXltTouWQMRAcCXG8z09deT5NnT4OHveyuqY4ql+jtcfP3KIcMtyz2ytZW5Y3q77uzk7n8oJs5pfk4HZJVCpmnjvoYUllAQvsXxjJRneoqrjrCnubPJU1NSU0erx4A8FYD+u0crzPz0jI8OZV84DYLKrOpI/7eBaXF3Ck0zfjcsjRmXthFpluF1UlObMO7gNDQbY2dSftrB00uKsE6Apr9zuVNQuKMQb2pHFJpDGGh3Y009E/+0W+SDn59vPqyijOzWRPS/TfHbX1+pk/i0oZx+LyPAaHR2jvm1mvoc7+ACIwJ8/6fls4J2/WwT2ZSyAdGtxV3HX7hsjJdJGXNXWFwdqa9F9UfWhHC//06138+qVjcfuaTr69pjSX1dXFUZ+5G2No7RmkKkozd2DGB3d4vEPMycsaTQ9Fo9b92QOepC2BdGhwV3HX6Q1Qlp89bl+ZsebkZ1FTmpu2efeOPj9ffexlwDpTNl6cmXt1SS6ra4o5cLw/qouqJwaGCQRDs6qUcTjnqc40uHd6AyelABfMyaPLNzTjVJ8xhmcPdiRtCaRDg7uKu27f0JSVMuHWLihJy4oZYwxffmQv/mCIquKc0dl0PLT0WN0aczLdrK629hMcOB69/QSz6eM+VlVxLlkZLo50zSK4F772/bZojrUxaqaz96ZOH8e6k7cE0qHBXcVdl3fq3anh1tYU09IzSJc3vfq7/25PG3/a184/vXEpZy0sHZ1Nx0NLzyDVdj58dbVVJRPN1IwT3KMxc3e5hNqyvFnN3CvCZu5OOeRM8+7PHrB21CfzYipocFcJ0O0bYk5+5LXPa0bz7ukze+/2DfF/H32ZNTXFfOjCxVSX5tLSMxjTJl7hWk4MUl1qBd6a0lxK8jKjumg9erxeFGbuYDcQm2lw7x86KS0z21r3Zw90UF+Rn7QlkA4N7iquRvt8TCMts6q6GBHYeSx98u5ffexl+vzD/Ns715DhdlFTmksgGMITh3cnxpiTZu4iEvVF1dbeQTLdQvk0folPprY8n1e7BqbdE8YXCDI4PEJ52Eaq4rzMGbf+HRwaYevh7qSuknFocFdxNTA0QiAYiqgM0lGQncEZlQVpUzHz5L52Ht3ZyscvOYPl84oARgNtPFIznd4hAsHQ6NcE6xfowfZ+/MPRWVRtsw/piFaP87ryfIZGQqPpnkiF704Nt7BsZuWQW5o6GQqG2LQsuVMyoMFdxdlrrQemN6NbU1PC7ubeuKUtYqV3cJgvPbyH5fMK+eim+tHbnRRJPBZVna9RXfpaWmFNlBdVo7WByVFb5nSHnF5q5rXgfvJkYqa17s8e8JCb6ebcxXOm/dx40+Cu4qozwqZhY62tKabLNxTXipJY+Mbv9tHlG+I771xLVsZrP37xnLmHl0E6VkV5UbW1Z3aHdIy1uGJmwd3TP/5u6AVz8mjuHiQ0jTSPMYZnDnRwwZLkLoF0aHBXcdUdYdOwsZxF1VQuiXz+oIdfb2vmIxfVndLHpTDHygPHo9a9pceasTrvFsBaVC2NwqJqcCREo8dLe58/qjP3ioJs8rPc0w/u9sy9YkzzstHWv9PYFZwqJZCOKRuHKRVNkTYNG2v5/EKy3C52N/dwxZr5sRhaTHkDQb7w0B7qK/L55GVnjPuYmtK8+KRlTgxSmJ1BcW7m6G0i1slXkc7cjTHsP97PKx1eDrX3c8jj5VCHl8OdPoZHrNlwfUVB1MYsIiyumH7FTGf/+N9v4a1/Iy3XTJUSSIcGdxVXTl+ZsmlWUWRnuFkxvzChO1VDIcO2oyc4s6qI/OzIfnT6/MPct/VV7tp8hPZ+Pw/cvJGczPHf0leX5vJqFPqMT6WlZ/CkWbtjdXUxP3q+Cf/wyIRjdNz+bCPf+eMBAFxiBcsllYVcunwuSyoLWDq3YLR+Plpqy/KnnTbq9AYozcskc0xnyvBa9/PrIjv/NFVKIB0a3FVcdXmHyM9ykxtBX5mx1tSU8NCOZkZCBncCTpp/sqGdj9yznZxMF5cur+TKNVVcsqxy3H9La88gd20+zH1/PYY3EGRjfRm3vWstZy+auBdJdUkuWxq7MMZE1JphpppPDJ6Ub3esqSkmGLJm5Ovsg1LGEwiOcOcLh7lgSRlfuXIltWX5U/4yiIa68nye2NPGUDB00nrFZMa2HnBUlViHgERa6+6UQF5//qJpjTmRNLiruOr2DTFnmvl2x5qaYu75y1GaPF7OmFsY5ZFNbU9zLxku4dqzF/D7vW08sec4eVluLlsxlyvXzOfipRUc7vTx4+eb+O2uVgxwxer53HRR3eiC5WRqSnPxBoL0DQYpzsuc8vEz1dIzOG61R/ii6mTB/fd7jtPlG+IjF9WPlnLGQ215PiEDx04MRJzy6fQOjRvcp9v6N5VKIB0a3FVcOU3DZsIJOLuaexMS3Bva+qivKODr16ziX646k61NXTy2u40/7G3jsV2t5Ga6GRweIS/LzfUbFnHjBYun9RbemU0fOzFAcV50UxqOPv8w/f7guDP36hJrUXXvFIuqP99yhMXl+Vy4pDwmY5zI4nK7Ysbjm0ZwD4wuxo81nXLIVCqBdGhwV3HV7RtiXtHMqijqKgrIz3Kzu7mHd55dE+WRTa2hrW/0h9vtEjYuKWfjknK+dvWZvNjYxZ9ePk5VSS7vO2/RjGbeNXbdeUvPYEQz/ZkYLYMcJ+cuIqyuKWH3JHntvS297Hi1h69cuTJqG5Qi5QT36TQQ6+w/ua9MuIVz8nhyX/uUr5HsB2FPREshVVx1eafXETKc22VVdOyaRRuCUMjwyN9aGBya3k7MnoEhWnv9rJh/ahoi0+3i4qUVfONtq7nlkiUzTqmMbmSKYTnkeDXu4VZXF/HKJDtV79lylNxMd0J+uZbkZVGalxnxkXuDQyP4hkZO6ggZbsGcPDq9Q/imaP3b6PHyavcAm5anRgmkQ4O7ihtjDF2+wLSaho21dkEJDW39DAVDM3r+Y7tbufX+nTyys2Vaz2tos3Zujhfco6U0L5PcTHdMyyFf2506UXAvIRgyNLT1nXJf78Awj+5q4Zqzqk4qo4ynxeX5NHm8ET12otYDjtEGYicmT808vrsNEXjTyrnTGGniaXBXcdMfCDI8YqbVNGystTUlDI2E2H/81OAzFWMMdzzXBDDt2b8T7GIZ3EXE6g4Zy5l7zyBZGa4JG3o5m6v2jpOa+c32Y/iHQ1x/fm3MxjeV+ooCmjyRzdxHNzBNEdwnKz81xvD47jbOWzyHyhmmExNFg7uKG6evzHQ3MIVbYwefXTPYSfn8K500tPWRk+ma9vMb2vooL8g+ZadjtFWX5NLcE7ta9xa7DHKifHlVcQ5z8rNOqScPhQz3/OUo59SWsrIqfhUyY9VXFtDRH6DPPzzlY50NTFPN3CdbVD3Ybm3OumJN1QxGm1ga3FXcdDt9ZabZNCxcTWkuc/KzZpR3v+PZRuYWZXPDhloOtvdPK+/ecLyPFfNjX6FTE+OZe3PP+DXuDqf979je+c+94uFo1wDXb6iN2dgiUWcvqkYye++0JxMT5dyLczMpzMmYtNb98d2tuATefOa8GYw2sTS4q7hxftim2zQsnIiwpqZ42u1/dx3rYUtTFx+8cDFnLyplJGTY1xbZ7D04EuJguzemKRlHdWkuJwaGGRia2fmeU2mZYANTuNXVxbzS4T1pUfWeLUcpL8hOeJCrr7RKIBs7ps67e/qdJnXjTyZEZNJySGMMv9vdxob6spi/Y4sFDe4qbrp9M2saNtbamhIOdXinrHIId8dzjRTmZHDduQtH6+V3RtiErKnTx1AwFJeZeyy7Q/qHR+j0BiZcTHWsqi5mJGxR9Vj3AM8c6OA95y6IeGdorCyck0eGS2jqnDq4d3oDFOdmTjrmyYL7vrY+mjp9XLE69VIyoMFdxZFzBupscu4AaxcUEzLjL/qNp8nj5Q8vH+f68xdRmJNJZVEO84pyIp79x2Mx1VFjB95YdId0DrqYauburGs4efdf/OUoLhHec17it95nul0sKsujsSOStMzUJ34tLMvj2InxW//+bncbbpfw5lWpl5IBDe4qjrp8QxRmZ8x6I4iz4/Cphqk3oAD8+M9NZLpdfOCCxWGvcWpeeSL72vrIcrui2uVwItUl1iJfcwzKIacqg3TML86hLD+LPc29+IdHuH/bMS5fOZd5UWzhOxv1FQU0RlAOOVFfmXAL5+QxFAzR0X/y8YbGGH63p42N9WWznowkigZ3FTez2cAUrrwgm7e/rpof//kwD/+tedLHdvT5eXB7C+88u+akvOnaBSUc7vTROzh11UVDWz9LKgtO6SwYC5WF2WS6JSZpmak2MDnC2/8+tquVnoFhrt+Q+Fm7o76ygCNdPoIjk+916PQOnXR26nhea/178juBvS19HO0a4MoUbC/t0OCu4qbbNxS1WdA3376aDXVlfO6B3Ww+1Dnh4+7cfIRgKMRNr6876fbR1EMEs/eGtr64pGQAXC6hqiQ3JhuZWnoGcQkRzcDX1FiLqj994TBnVBawIcK2uPFQX1HA8Ijh2BS/ACdrPeCYqBzy8T2tZLiEN6VglYxDg7uKm05vYFZlkOGyM9zccf3Z1JUXcPM928fd1NTnH+bevxzlLavmU2uX0DnWVDtNyCbPu3d6A3j6A3FZTHVUl+TSMsWuyZloOTHIvKKciN6BOIuq+4/3c/2GRTFtQTxddfaRe5NVzPiHR+gPBKescqkqycUlnFQO6VTJXHhGOSV5qZmSAQ3uKo66fEOzKoMcqzg3k7s+cA752Rm8/86XaOs9eSb3y62v0h8IcvPF9ac+Ny+TxeX5U9bLO4upK+M0cwd7I1MM0jLNExzSMR7nnU1+lpu3nVUd9bHMRn25XQ45Sd59ooOxx7Ja/+aeNHPf1dxL84lBrkzBjUvhNLiruAiFDCd80cm5h6sqyeWuD5yDNxDk/Xe+NLpzMfxAibHnlToiWVSNZ6WMo7o0l47+AIHg9JqbTSWSGnfHvKIc6srzea9dYZRMivMyKS/IniK4j38w9njGlkM+vquVTLfwxhTrJTOWBncVF33+YYIhM6umYRNZMb+I/7n+bBo9Xm6+ZztDwRAP72ihoz8w7qzdsaamhON9fjr6Jj4keX9bPxmMjGoAACAASURBVPOKciiNY8WEE4DbeiI/vHkqwZEQx/v8Ec/cRYQ//uNFfP7Ny6M2hmiqr8ifdJfqVK0HwlnB3XqnFAoZntjTxkVnVCSsOVq0aHBXceGcnTqbpmGTuWBJOf/2zjW82NjF5x7YxY+eb+LMqqJJD5RYG0Gfmn1t8Wk7EC68r3u0tPcHGAmZ0VLLSGS6XXHv2R6p+srJyyFH0zIR7Cy1Wv8GGBgK8rdjJ2jt9XPl2tStknFocFdx0eWd2cHY0/H219Xw2Tct45GdrTR1+rj54vpJFwLPrCrG7ZIJ8+5DwRCNnvi0HQhXE4O+7pMd0pGK6isKODEwPLrreSwnuEeyxjPa+rd7kMd3t5GV4eINK1I7JQMRBHcRWSAiz4hIg4i8LCKfsm9fKyJbRGSPiDwmIkX27e8VkZ1hf0Iisi7W/xCV3JymYbHeEPKxTfV8+PWLOW/xHN4yxc7C3Cw3S+cWTlgxc6jDy/CIiXtwn1ecg0ugOYoVMy12p8lIc+7Jrt6pmJlg9t7pHaIwJyOig7ud4H6ky8cTe9rYtLQi6dYZZiKSY/aCwKeNMTtEpBDYLiJPAj8BPmOMeU5EbgQ+C3zFGHMvcC+AiKwGHjXG7IzR+FWKeG2BK7bBXUT40hUrI3782ppi/vDycYwxp8zyE7GYClY6ZG5RTlR3qUa6gSlVOLuFGzu8nFN76rmmnghq3B2Lyqzg/vCOFtr7AlyRwhuXwk05czfGtBljdtgf9wMNQDWwDHjeftiTwDvGefp1wH3RGWp6eLqhne89dZDhKXbXpRvn7XM8FyYjsaamhJ6B4XGbRzW09ZGd4Ro9uzOerFr3KAb3nkHK8rPIzUqdM0AnU1WSS3aGa8KZuyeC1gMOp/XvH/cdJztNUjIwzZy7iNQCZwFbgb3AVfZd1wILxnnK36PB/ST/83wT33vqFW782Uv0R3DgQLrosjv0xWML/3SsXTDxomrD8T6WzSvEnYBFxZrS6O5SbT4ReY17KnC7hMXl+TROUDHT6Q1M2Md9LKf1rzFw6fJK8rMjSWgkv4h/0kSkAHgQuNUY0wfcCNwiItuBQmBozOPPAwaMMXsneL2bRGSbiGzzeDwz/gfEyouHOke76EWLMYYDx/s5o7KALY1dXHvHllM23qSraG9gipalcwvJznCdsqhqjKGhrZ8V8xJz6lB1aS5tvf4p+6dEqqVnkKri9AnuMHnFTGd/5DN3eC3vni4pGYgwuItIJlZgv9cY8xCAMWa/MeZyY8zZWLPzxjFPezeTzNqNMT8yxqw3xqyvqKiY2ehjZOexHt7306185ZFxfy/N2PE+P72DVhOmuz5wDi0nBrnmh5t5uXX6R8almmg1DYu2TLeLM6uKTmn/29EfoNs3FPcySEd1SR4jIUP7mG6FM2GMoXUau1NTRX1FAce6B07Z7BUIjtDnD04ruC+fV0RRTgaXLq+M9jATJpJqGQF+CjQYY24Lu73S/tsFfBm4I+w+F1aq5lfRHnCsBYIjfO6BXYQMPHOgI6qz9/3H+wFYNreQ159RwW8+ugGXCO+6YwvPHOiI2tdJRl2+QNK2Tl1TU8Lelr6TZsn7ErSY6qiOYjlkl28I/3AobRZTHfUV+YQMHB1zwLVTdjud05Nu3lTH05/eRF5WeqRkILKZ+wXA9cClYeWNbwWuE5GDwH6gFbgr7DkXAc3GmKaojzjGfvhMIwfbvXz96jMxwK9eOha11z5gB/fl9lv95fOKeOSWC1hUls+H7t7GL7e+GrWvlWy6fUNRaxoWbWsXFDM4PMKhsLf4TqXM8kQFd+dEpigclp1uNe6O8IqZcK/1lYn8+y07w52SR+lNJpJqmReMMWKMWWOMWWf/ecIY831jzFL7z+eNMSbsOc8aY86P7dCjb19rH7c/c4i3nVXN9RtquXhpBfe/9GrU8p772/qYX5xDcd5rNbRzi3L49c0beP0Z5Xzx4T18+w/7CbuUaSEUMnT7hihP0pn7Wvvwj91hx+41tPVTXZKbsC3o0TxuryXCE5hSTd0Ete6RNg1Ld8lVupBAwZEQ//zgbkryMvk/V1p10u89bxHtfQGe3h+dlMn+4/0sm3dqDrcgO4Of/MN6rjt3Af/9bCMvNnZF5esli57BYUIm9huYZqq2LJ/CnAx2huXd49nDfTy5WW7KC7Ki0h3S+QVRk2Yz97ysDKqKc06pmOnsj7xpWDrT4G778Z8Ps6ell69dvWq0FvuSZRXMK8rh3iikS4ZHrK3s4wV3gAy3iy+8dQUAf3v1xKy/XjJxzk5N1rSMyyV2h0gruPuHR2jyeFmZoMVUR3WUDu1o6RkkP8ud8o2wxjNexYzH/n5LtzTLdGlwx3pb9x9PHeQtq+bx1tWvlUJluF28+9wF/PkVD692zS732eTxWVvZJymtK8qxeozvifDg51ThNA1LxlJIx5qaEva39eMfHuFgez8hk7jFVEd1aXQ2Mjk17sl04Ea01FcU0NjhPSmV2ekNUJAdWeuBdHbaB/eRkOFzD+wmN9PNV68+85T7//6cBQhw30uzm707JwVNNHN3rKoujujot1Qy2jQsSWfuYLUhCIYMDW19CWs7MJYzc5/tGkxLT+R93FNNfUU+vqER2vteKxnt9A6d9vl20ODOz7ccYfvRE/yfK1dSWXjq2ZLzi3O5bMVcfrPtGEPBmS+sHjjeT4ZLRlf4J7K6uojWXv9oKiMdxKtp2GysXWAvqjb30tDWT36We3RjS6JUl+QSCIZG+/KMZyRkeHx36+ghJeNpOTGQdpUyDufnqSksNePp95/2+XY4zYP7se4B/u0PB9i0rIK3v27io8Tec95COr1D/Gnf8Rl/rQPH+6mvKCArY/JLvto+2zOdUjOd3iFEoDQveXO+84pyqCjMZtexHva1WW0HEt3L3OnrPll3yPv++iof/+XfeMftL550Dqij3z9Mnz84rT7uqaS+8tQj96yZuwb30za4G2P4/EO7cbuEf33b6knzkRedUUFNae6s6tAnqpQZ68xqKxWQTqmZbt8QJbmZZCRZX5lwIsLammJ2NfewP8GVMo7RjUwTLKoODo3w/adfYencAjr6A1zzw81sP9p90mNGyyDTdOZeWZhNfpb7pIqZ6fSVSWfJ+9MWQ0e7fNx0z3Y2H+ri829ZTtUU+Ui3S7ju3IW82Nh10tu/SPX5h2npGYwouKfjomqXL5DU+XbHmpoSGj0++vzB5AruEyyq3vXiYTz9Ab7xttU8/LGNFOVmct2PtvLI31pGH5NurX7HEpGTKmaGR0L0DAzrzJ3TLLj3+4f55u8beONtz7P5UCf//OblvOfchRE999r1NWS4hPv+Ov3Z+8HRnamRldatri5mbxoF907vUFLn2x1rwg7STobgXpRjtaIdb+beOzDMHc82cunySs6pnUNdRQEPf2wjr1tUwq337+S2Px0gFDKjz023GvdwTsUMvLZ4r8H9NAnuoZDh1y8d45LvPsf/PNfEVeuqePYzm/jopvqI86qVhTlcfuZcfrO9Gf/w9E6lb3CCe4QBY3V1Ma29/tGddqmu25ca1QvOTlWRyH8Rx1p1Se64G5nueL6R/kCQz75p2ehtJXlZ/PzG83jX+hr+838P8Ylf/Y3GDi9ZblfEB1ekovqKfFp7/fgCwRm1HkhX6dMlZwIvHenmq4+9zN6WPs5eVMpPb1g/WhkxXe89bxFP7DnOH/Ye55qzJl6AHevA8T4Kc6zddJFYVW3NIPe09HLJstTvUtflDVBWV5boYUypND+LhXPycAlJ09O7pjTvlIXS9j4/d20+zNVrq055h5GV4eLb71jDksoCvvn7/RhjnTSU6MXhWHIqZg53+nQDU5jk+A6OkduePMh/Pv0K84tz+P6713HV2qpZbeTYUFdGbVke9249Os3g3s+yuYURf+1V9qLq3ubUD+7BkRA9g8MpkZYB+PglSwglUW+fmtJc/tLUddIxgP/59CsERwz/9MZl4z5HRLjponoWleVz6692smSK8ttUF14x45Qrp/M7lUilbXAfCRl+tvkwlyyr4Pb3nh2V48VcLuE95y3kX5/Yz8H2fpbOnfqtuzGG/cf7uXpdVcRfpzAnk7ryfHanQd79xMAwxqROE6d3nTPegWKJU12SizcQpG8wSHFeJkc6fdz/0jGuO3chC8smL29805nzePazmxJyklQ8LSqz3m01enzk2rtStVomjXPue1t66fMHueas6qieG/nOsxeQ5XZFXBbZ2uun3x9k2TRP9FmVJouqztmpc/J1JjUTTsVMs93697YnD5LpdvGJS5dE9Py5RTlpn3/OznCzYE4ejR4vnd4AeVnutOrLPlNpG9w3N3YCsLG+PKqvOyc/i7esnseDOyJbWD1gtx2Y7gLdmppi2nr9eKJwEk8ivdY0TGdSM+GUMDafGOTl1l5+u6uVGy+spbIosvWb04VTMdM5jYOx0136BvdDnSyfVxiThZV3vK6Gfn+QZw9MffZrQ5tVKRNJCiecs6ia6rP3VGgalszCa92/88cDFOdmctNF9QkeVfKpr8jncKeP9j5/yqQAYy0tg7t/eIRtR05Efdbu2FhfRnlBFr/d1TLlYw8cn9mhD2dW2TtVUz24J3m732RXlp9FTqaLR3e18uwBDx/bVJ+WrXtnq76igEAwxN6WPp2529IyuO84eoJAMMSFZ8Sm/C7D7eKK1fN5uqGD/kkaNoFdKTODmunCnEzqKlJ/p2qXbwiXQIkGpBkREapLctl1rIe5RdncsLE20UNKSk7FjDcQpFzLIIE0De6bGzvJcAnnLo5dbfVV66oJBEP86eX2CR8zFJz8gI6prE6D9r9dPmt3ajrXWcdatd1A7FOXLT3te5RPJLzbqs7cLWkZ3F841MXaBSUUxHAjyusWllBTmstvd7VO+JhGj5dgyMx4t+Pq6mKO96X2omqXN0CZVsrMyrm1paypKeba9TWJHkrSmpOfNdp1tEJz7kAaBvfewWH2NPdwQX1sd0SKCH+3tooXDnVO2Hv9wGhPmZn1KVmd5Iuqe1t6ec+P/8K77thC6wSdC7t9qdFXJpl9/NIzePSWC8hM4q6ayaDOnr3rzN2Sdt8tW5u6CBm4YElsFlPDXb2uipGQ4Yk9bePev/94P5luGT2lfbrOrC5GxDpAIpl4+gP88wO7+bsfvMD+4/3sa+vjyv96gRcPdZ7y2C7vkJZBRkE6HpEXbfX2z5nm3C1pF9w3H+okN9PNWQtLY/61ls8rYuncgglTM/uP91FfUTDjGVdBdgZ1SdT+NxAc4X+ea+SS7z7Lgzua+eAFi3nmM5t49OMXUJafxft+upU7nms86Vi4Lt+QlkGquHDy7tp6wJJ227g2N3ZxzuI5U554FC1Xr6vmO388MO45lQeO93Pe4jmzev3V1cX8pal76gfGkDGGJ/e1840nGjjaNcClyyv50hUrRn+YinMzeeSWC/jcA7v51u/3s/PVHr5z7RpyMt30Dg5rGaSKi3ecXYPbJSyaoi3D6SKtZu7tfX4OdXhjnm8P93drrJ4xj42ZvfcODNPW659224GxVtmLqh39/kkfFxwJ8XJr9Gf4geAIH/jZS9x0z3Yy3S7uvvFc7nz/OaecBZufncEP3nMWX3rrCp5saOeaH27mpcPWLyXNuat4KC/I5kOvr9MUli2tgvuLdsuBeOTbHQvL8li3oITf7jw5uO+fYduBsdbYPcanWlT99h/2c8V/vhD10skXXunk2QMe/vENS/n9p17PxUsrJnysiPDhi+r4xQfPo2dgmOvv/CuQOk3DlEonaRXcX3ili9K8TFbG+RSdq9dVsa+tj0Md/aO3HWh3DuiYXXA/s6oIEdjT3DfhY/a19nHn5iMA3L3lyKy+3lhPNbRTkJ3BzZvqIl472FBfxuOfvHC02md+cfqeAqRUskqb4G6M4cXGTjbUl8V9w8wVa+bjEk6ave8/3k9RTgbzZtngKX90UbVn3PtDIcMXH95DSW4mf7e2it/uap2wNHO6QiHD0w0dXLS0nOyM6W2emV+cy/0fOZ/f3LxhxoejKKVmLm2C++FOH229/rimZByVhTlsrC/n0V2to5Ui+9v6WD6vKCr5vzU1JRNWzPzyr6+y81gPX75yBZ+4dAlDwRD3bzs2668JVl+bjv4Ab1gxd0bPz85wc07t7BaUlVIzkzbBfXNjFwAXxKhZ2FSuWlvF0a4Bdjf3YozhYPvM2w6Mtaq6mPa+AB19Jy+qdvT7+fYf9rOxvoxr1lWzdG4hG+rKuPcvrxIcCc366z7V0I5LSPnToJQ6HaVPcH+lk+qS3ISVQb1p1Tyy3C4e3dlK84lBvIHgrPPtjtVhZ6qG+8bvGggMh/j6NatG3yHcsHERLT2DPL2/Y9Zf96mGDtYvmkOpVrsolXLSIriPhAxbmrrYWF+WsDKo4txMNi2r4PHdrexri06ljGN0UTUsuP/5FQ+P7mzlo5vqTypLfMOKucwvzuHnW47M6ms2nxigoa2PN6zUWbtSqSgtgvu+1j56B4e58IzEpGQcV62roqM/MBpYp3tAx0TyszOorygYLYf0D4/wlUf2srg8n49uOvnghgy3i/edv4jNh7pOqt6ZrqcbrJn/TPPtSqnESovg/oLd02RDHDcvjeey5XPJz3Kz+VAXNaW5FOZEr4f56uri0R4ztz9ziCNdA3z96lXjtoD9+3Osc15/vuXojL/eUw3t1FXkjzZjUkqllrQI7i82drJ0bgGVhYk9VzI3y83lZ84DopeScayuLqajP8CWxi7++7lGrllXNeE7lfKCbK5cO58HtzdPeZjIePr9w/ylqUtn7UqlsCmDu4gsEJFnRKRBRF4WkU/Zt68VkS0iskdEHhORorDnrLHve9m+P2ZRNxAc4aUj3QkpgRzPVWutdgTRqpRxrK6xFlVv+eUOcjPdfOmKlZM+/oYNtfiGRnhwe/O0v9bzBzsZHjEa3JVKYZHM3IPAp40xK4DzgVtEZCXwE+DzxpjVwMPAZwFEJAP4BXCzMeZMYBMw/eljhHYc7cE/HEpYCeRYF55Rzvs31nLNuuqovu7K+daiardviH9+y/IpD/5eu6CEtQtK+PmWo4RCZtLHjvVUQzuleZm8bqFuPlIqVU0Z3I0xbcaYHfbH/UADUA0sA563H/Yk8A7748uB3caYXfZzuowxI9EeuGPzoU7cLuG8uuTYLJPpdvEvV53JGVFaTHXkZ2ewqqqYsxeVct05CyN6zg0bFtHU6WNz46l91icSHAnxzIEOLllWSYYeDqFUyprWT6+I1AJnAVuBvcBV9l3XAgvsj5cCRkT+KCI7RORz0Rnq+DY3drKmpjiqi5fJ6hcfPI97PnhuxO0V3rp6PmX5Wdz9YuQLq9uPnqBnYJg3rNSUjFKpLOLgLiIFwIPArcaYPuBGrBTNdqAQGLIfmgFcCLzX/vttInLZOK93k4hsE5FtHo9nRoPv9w+zu7mXC5Mk3x5rxXmZ5GVF3oI/J9PNu89dwNP72znWPRDRc55qaCfL7eKiSbo/KqWSX0TBXUQysQL7vcaYhwCMMfuNMZcbY84G7gMa7Yc3A88ZYzqNMQPAE8Drxr6mMeZHxpj1xpj1FRUzCyQHjvfjFmFjkuTbk9F7z1uEAL/YGtns/emGDs6vL4vp4eJKqdiLpFpGgJ8CDcaY28Jur7T/dgFfBu6w7/ojsEZE8uzF1YuBfdEeOMD62jns/pfLOac29kfqpaqqklwuXzmP+186hn948qWPRo+Xpk4fb1ihu1KVSnWRzNwvAK4HLhWRnfaftwLXichBYD/QCtwFYIw5AdwGvATsBHYYY34Xk9FjpR504W9yN2yspWdgmHum2NT01L52AC7TEkilUt6U772NMS8AE63gfX+C5/wCqxxSJYHz6+Zw0dIKvvFEA4HgCLdcsmTcHjxPN3Swcn7RKWfBKqVSj055TwMiwo//4WyuWVfFd/90kC88tIfhMS2Bu31DbDvarSkZpdKErpqdJrIz3PzH36+jpjSPHzxziNZeP7e/93WjC6fP7O8gZNASSKXShM7cTyMiwmfetIxvvn01mw918q47ttBuHwDy9P525hZls6qqOMGjVEpFgwb309B15y7kJzes52iXj7f9cDN7W3p57oCHy1bMjfv5s0qp2NDgfpq6ZFkl939kA8GQ4ZofbsY3NKL5dqXSiAb309iq6mIevuUC6iryKc7N1M1gSqURXVA9zVWX5PLbj19I7+DwuAd/KKVSkwZ3RU6mWwO7UmlG0zJKKZWGNLgrpVQa0uCulFJpSIO7UkqlIQ3uSimVhjS4K6VUGtLgrpRSaUiMMYkeAyLiASI/xXly5UBnlF4rWnRMkdExRUbHFJlkHBNEd1yLjDHjnlOaFME9mkRkmzFmfaLHEU7HFBkdU2R0TJFJxjFB/MalaRmllEpDGtyVUioNpWNw/1GiBzAOHVNkdEyR0TFFJhnHBHEaV9rl3JVSSqXnzF0ppU57GtzVaU9E9GzBCOh1ikyyXKeUDe7JcgGTnV6niGQmegBjiUhtoscwDr1OkUmK65RSwV1EzhSRTQAmSRYLRORcEflXEUmaa6nXKTIiskFEfgN8V0RWikjCTywRkdeJyFPA15JhPKDXKVLJdp2S5gdtMiLiEpHbgQeBL4rI10VkvXNfgsZUJCI/BH4ANBtjQomeJet1mta4Ku0xPYG1W/BTwI32fXEfn1i+BNwH/MoY8w/GmJFEjSdsXHqdIhtXUl0nSJHgDpQChcAK4L1AF/BpESkwxoQSNKYvAecDlxtjboekmCWXAAUk13X6Isl3nQDWAgeNMXcB/w48BFwtIkuNMSbeP5D2NckBXjDG/ARARM4SkYwEX69VJN91yiT5rlNSfT9BEgd3+23XUvvTYmAjkGeM8WDNTLuBW+zHxuXC2WNabn96J+ABKkXknSLyXRF5t4gsjMdYwsa0WERy7E/nkBzXabGI5Nmf/pzkuE7XichXReQq+6a/AetFpN4Y4wNeArYBH4H4/AIKG9M19k3fBqpF5N9F5CXg68DdIvLOWI8lbEwXi8h5YTftwrpOdQm8TmPH9B2s6/TdBF6na0TkiyJyhX3TThL8/TRW0gV3OzD8DvghcI+IvNEY0wS8CNxqP6wN6zfjWSJSFesLN2ZMd9tjOgBsBX4PfAw4AFwLfFZEamI5HntMtSLye+AnwL0istIYcwh4Hvgn+2Hxvk7hY7rHHtM+4M/AH0nMdRIRuRn4HHAE+I6IfAjwYv3i+ZT90B7gKSBPRObHeUzfFpEPG2O8WNfuLODTxpgrsf4/3xw20YnVmApF5CHgYeAjIlIKYIzpAu4HPmk/NJ7XaaIx+YB7gHXE/zpViMgjWD9j3cBdIvLOsMnUJ+yHxu06TSQpgvuYGeVngJ3GmA3Ao9h5K6yZ8gUistgYEwTaAT+QG+cxPQJ8yL79W8DXjDGXGmN+DHwFKy2yOE5j2mqMuQx4BviqiKwEfgacb8+0EnGdwsf0dRGpw5pp/Uu8rlM4+xfaBuBb9lvmW4BNwGVY+dElIvIGO23VBVQDvQkY0yUi8iZjzAPA240xz9sPfwqowPplFEtDwP8C7wNasX4BOx4ElovIZfG8TpONyRhzL/CuBFynemCzMeYiY8wdwKeBf7Tvu4/EXKdxJUVwx8o1OoHCBwzbtxcBDSKyBNiM9VbnuwDGmL3AIiAQ5zEVA3tFZIUxZtAYc7cT4OxZ6jzg1RiPKcP+/GX76/4AOBd4N9YPwUvAv9n3xes6jTems4GbgHxjzN3OE2J9nUTkH+y38nPsmxqw3sZnGGOeAvZirQN4gF8C37O/xy4DBMhKwJh2AZtEZIExpifsqW8EDDEIWmFjKjHGBLDeNTwFHMRKMSyzH7oL+BXxvU4TjWmp/TgxxnSHPTXW12mTnWrcjvWOD7GqYfbZfwD2YF2n78f6OkUiocFdRN4oIk9ivVV+lz2jeQE4Q0T+BrwZcGP9AF4MfBOYLyI/EJG9WD3ge8fMHuM1pp+LyOX2N5kRkatF5Gms4Nod4zEFsd4SniUia0VkLVbAWmyP7f8BNSLyX3G8ThONqQaYG/a8q2JxnexUx3wReQa4AWtB+b9EpAg4BlQCS+yH/wprobDMGPML4F7g81i/HD83JrjGa0z3A8uBMvu5l4jIDuAtwOeNMX0xHNMPRaTcGOM3xgwBW4AO4F0AxpiQMeZnWKmQLxCf6zTVmIxY1WEX2T+Xsb5O7wF+DBQbY9pFxG1X6azAmvCFX6efE4Pvp2kzxiTkD9Y39Vbgaqwc4y+Bz9j3LQMeCnvsV4Af2B/PxVo0vCoJxvQf9scbsX6jXxOHMd2HlbsutMfwONYvn/X2eG9NwHWaakwfj+V1Atz230uBX9gfZwC3A3djVVfcCVyP9cMJVvrqG2GvkZUkY/pa2DWO6v/dJGP6L+DBMY99mz3WJUA+4IrzdZpqTDlYs+L6OF+nh8Y85udY6SGAebH6fprJH+etdFyIXWttrHzUecB2Y8yj9n1PAbeJyD1YM8BjduqjASt/e6s9S27HyiMny5hcxpgXsVIQ8RrTvwO/McZ8XazcepN932ZeS7902Ncq0WN6ESvnTwyuUwbwNcAtIk9gpfFG7K8VFJGPYy0qr8T6JXMN1juJbwIhrEV67McPJcmYttqPPQQcitOYPgm0isjFxpjn7NsfFpEVwB+w1kcuARrieJ0iGdOlxkrxNSZqTFhpoMMi8jXg7SLyZmNMc7Su02zELS0jIh8AmrHKlsDKT10nr20fzgSa7Pv7scr6PikinwL+ByvvloxjimoNawRjysD6Zv4P+/PD9vNuAj4I7IDoll7Nckw3OmOKJhG5GOtdQClWEPw61rrIJSJyLoz+Ivoa8G1j5bV/BFwoIlvt5z2rYxr9Xvka8C9hz7sWay/HM8Aae0KTbGPaR5TMZExi5dxvBB7A+kVwiTGmOVpjmrV4vD3A+i37CFbZ2Q5guX3797De0m8GfgGsEjhYZAAAAyJJREFUxiotzMfKZX0C6y3s+TqmU8b0O2Cuff+tWIuo55wOY7Jf//XA9WGf3w58FHg/1rsKsCYv84DfALX2bSVAtY5p3DH9Glgc9rzX65gmHNMirJTQ94DXxWJMs/43xe0LwUL7728B99sfu7Fmwxfany/ACpxxyVel+Jh+BmTbn+edhmPKA7J5Lff5XuCb9sc7gU/YH68H7ovT/52O6fQY06/iMabZ/olbWsYY45S9fQ9YLFZN7wjQa4x5wb7vZqyywxEd05RjGgCC9nMGTsMxDRhjAvY4wCqF89gffwBYISKPY727iHpaSMd0Wo9pO6RAx9VE/EbB2pL7XNjn52JtWHqCsBVnHZOOKYLxuLHeLv8eWGLftgQr1XEhMUp36Jh0TMn+J+7H7NnVJSEReQCraiCAtTD5ijEmKqveOqbTakzOJpGfYG1TvxFrZ+AnTJRqnnVMOqZEjmnGEvEbBSu/9TxWa8xPJvo3nI4p5cd0PlYZ4QvABxM9Hh2TjikZ/sS1zj3Mx7ByaW801hbjZKBjikwyjqkZq0TuNh3TpHRMkUnGMU1b3NMy8Nrb+7h/4UnomCKTjGNSSp0qIcFdKaVUbCVLV0illFJRpMFdKaXSkAZ3pZRKQxrclVIqDWlwV6clESkRkY/ZH1fZG7OUShtaLaNOS3a74seNMasSPBSlYiJRm5iUSrRvAfUishN4BVhhjFklIu/HOkDDjXUU379jbUe/HqvdwluNMd0iUg/8EOtQ5gHgw8aY/fH/Zyg1Pk3LqNPV54FGY8w64LNj7luFdWbmucA3gAFjzFlY53j+g/2YH2H1Gzkb+AxW/2+lkobO3JU61TPGmH6gX0R6gcfs2/cAa0SkAOs82N+EdX3Njv8wlZqYBnelThXeTyQU9nkI62fGBfTYs36lkpKmZdTpqh8onMkTjdX69bB9pidiWRvNwSk1Wxrc1WnJGNMFbBaRvcB3ZvAS7wU+KCK7gJeBq6M5PqVmS0shlVIqDenMXSml0pAGd6WUSkMa3JVSKg1pcFdKqTSkwV0ppdKQBnellEpDGtyVUioNaXBXSqk09P8BOeegELnSh18AAAAASUVORK5CYII=\n",
"text/plain": [
"