{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Kálmán-szűrés 1.\n",
"\n",
"## Rövid áttekintés\n",
"\n",
"* Rudolf Emil Kálmán (1930-[2016](http://www.ematlap.hu/index.php/tudomany-tortenet-2016-09/347-kalman-rudolf-es-a-kalman-szuro)) magyar származású amerikai villamosmérnök, matematikus publikálja az elméletet (1960) \n",
"* Első alkalmazása: ember nélküli amerikai Hold-szonda (1963)\n",
"* Navigációs és mérnöki alkalmazások\n",
"* Szenzor adatok fúziója (radar, lézerszkenner, kamera, INS,…)\n",
"* GNSS vevők, okostelefonok, számítógépek, játékok, stb. elmaradhatatlan része\n",
"* Tőzsde előrejelzés, idősor elemzés, függvény approximáció\n",
"\n",
"Az eredeti [publikáció](https://www.cs.unc.edu/~welch/kalman/media/pdf/Kalman1960.pdf) közel 25 ezer hivatkozást kapott. Pillantsunk bele:\n",
"\n",
"\n",
"## Meg lehet érteni?\n",
"\n",
"A Kálmán-szűrés olyan algoritmus, amely valamely lineáris dinamikus rendszerben egzakt következetést tesz lehetővé, amely a rejtett Markov-modellhez hasonló Bayes-féle modell, azonban a rejtett változók állapottere folytonos, és minden rejtett és megfigyelhető változó (gyakran többváltozós) normális eloszlású.\n",
"\n",
"A fenti definíció korrekt, de a legtöbbünk számára nem teszi érthetővé, hogy **valójában** mi is a Kálmán-szűrő. Ezért a megértés érdekében egyesével sorra vesszük a Kálmán-szűrő legfontosabb sajátosságait:\n",
"\n",
"* rekurzív\n",
"* optimális\n",
"* közvetett (rejtett) állapotot becsül\n",
"\n",
"Ezek után egy geodéziai szempontból talán könnyebben érthető módon levezetjük a Kálmán-szűrő egyenleteit és egy egyszerű navigációs példán keresztül mutatjuk be a Kálmán-szűrő működését."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Rekurzív becslés\n",
"\n",
"A legegyszerűbb példa az, amikor valamely egyetlen $x$ változóval jellemezhető rendszer állapotát mérjük egymás utáni $t_1, t_2, ... , t_n$ időpontokban. \n",
"\n",
"A méréseink legyenek $x_1, x_2, ... , x_n$.\n",
"\n",
"Az állapot legjobb becslése ez esetben (Gauss-eloszlást feltételezve) a jól ismert számtani közép vagy *átlag*:\n",
"\n",
"$$ \\mu_n = \\frac{1}{n} \\sum_{i=1}^n x_i.$$\n",
"\n",
"Tegyük fel, hogy van egy új mérésünk, $x_{n+1}$, és a célunk az állapot becslése ezen új mérés felhasználásával. Ez az $n+1$ mérésből számított átlag lesz:\n",
"\n",
"$$ \\mu_{n+1} = \\frac{1}{n+1} \\sum_{i=1}^{n+1} x_i.$$\n",
"\n",
"Viszont sokkal hatékonyabb az az eljárás, amikor felhasználjuk a már kiszámított $\\mu_n$ átlagot és csak ezt \"frissítjük\" az új mérés felhasználásával:\n",
"\n",
"$$ \\mu_{n+1} = \\frac{n}{n+1} \\left( \\frac{1}{n} \\sum_{i=1}^n x_i + \\frac{1}{n} x_{n+1} \\right) = \\frac{n}{n+1}\\mu_n + \\frac{1}{n+1} x_{n+1}.$$\n",
"\n",
"Igazából így működik a zsebszámológépünk `STAT` funkciója is: az egymás után bevitt adatok segítségével folyamatosan \"frissíti\" az átlag értékét. Ez a *rekurzív becslés* lényege: a már meglévő becslést az újabb adatok beérkezése után mindig módosítjuk, \"aktualizáljuk\". Példánkban a rekurzív állapot becslés így írható fel kicsit más alakban:\n",
"\n",
"$$\\mu_{n+1} = \\frac{n}{n+1}\\mu_n + \\frac{1}{n+1} x_{n+1} = \\mu_n + K(x_{n+1} - \\mu_n),$$\n",
"\n",
"ahol \n",
"\n",
"$$K = \\frac{1}{n+1}$$\n",
"\n",
"az ún. *erősítés*, amelyik megmondja, hogy mekkora lesz a régi átlag, azaz a $\\mu_n$ korrekciója, illetve\n",
"\n",
"$$x_{n+1} - \\mu_n,$$\n",
"\n",
"amelyik az *új információ* (idegen szóval: *innováció*), hiszen ettől függ, hogy egyáltalán változik-e az előző becslés. Ugyanis ha az új $x_{n+1}$ mérés megegyezik a régi $\\mu_n$ átlaggal, a becslésünk egyáltalán nem fog változni.\n",
"\n",
"Nézzük meg mindezt az alábbi számpéldán keresztül:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAEKCAYAAAA7LB+5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAAsTAAALEwEAmpwYAAA+N0lEQVR4nO2deZhVxbXof+t0n55AaCYRG5kUQZQWsHG4oDgxqAGJPo0JRqPX6M3okGtC5kGjRH3vJg5RuYkxTtEoEVE0GJRgJIqCIII0AgJCMw/N2DQ91Ptjn24Ph73PuIfap+v3fefrPvvsYe0aVq1atapKlFIYDAaDQV8iQQtgMBgMhuQYRW0wGAyaYxS1wWAwaI5R1AaDwaA5RlEbDAaD5hhFbTAYDJpjFHUbQUSuF5H+QcthMBgyxyjqNoCIXAmMUEqtDFoWg8GQOUZRtwGUUn9VSv1nqvNEpI+IKBEp9EOuFLI8LiJ3+vCc40RkjYj0zfC6f4rIDV7JZfO8tSJyoc3x74vIEyJi6nIeYzK3DRBTKrtEpDjhuG3lb2P8L/AdpdSaoAXJFBG5CBgGXKeUag5aHoN3GEWd54hIH+BsQAETgpXGQix8KXvJegci0gt4Qin1ih+yuI1S6jWl1FVKqaagZTF4i1HU+c81wLvA48C1LQdF5EmgF/CyiOwTke8nXigi14nIchHZKyKfishNCb9/X0Q2ichGEbkh5jY5wU6ImFX/axGZBxwA+onIQBH5h4jsFJEVMV+63bVHicgcEblfRPomumfi3RAi8jURmSci/yMiO4FfiMiHsXds+SgRORer/D8tIoUicpWILEh47q0iMiNJ2h4vIu+JyG4ReUlEOsdde6aI/FtEamPPPzfut84i8qdYuu0Skemx411F5JXYNTtF5F92DZqIRERksoisFpEdIvLXlmeLSImIPBU7Xisi74tI99hvHUXkj7E8qxGRO0WkIMn7GXRBKWU+efwBVgHfBE4DGoDucb+tBS6M+94Hy/IujH2/BDgeEGAUloIdFvttHLAZOBkoA56MXXuCgxz/BD6LnV8IdATWA9fFvg8DtgMnx85/HLgT6AK8B9xpJ2PcvW+I/f81oBH4Tuy+pQly3AhUAx3i7xV7h71A/7hz3weuSvI+NcApQDtgGvBU7LcKYAdwMVZjMDr2vVvs95nAc0AnIAqMih2/G3gkdiyK1ROSxLwCbsFqfHsCxcCjwF9iv90EvBx7n4JYvneI/TY9dm474OhYut4UdBk1nzTqcdACmI+HmQsjsZRz19j3auDWuN9bK3/s+xFKMOF+04GbY/8/Btwd99sJpFbUv4r7/iXgXwnnPAr8PPb/47FnLAVuTyYjRyrqz5Kkx1bgRLt7AU8BP4v93x9LcZcleZ8pcd8HAYdiyvEHwJMJ58/C6tH0AJqBTjb3/BXwkl0aJijq5cAFcb/1iOVzIXA98G+gMuH67kA9cQ0X8GVgTtDl1HxSf4zrI7+5FnhdKbU99v0Z4twfqRCRi0Tk3Vg3vBbLQuwa+/lYLIu4hfWJ19sQf05v4IxY97w2dv9JwDFx51wClGJZmZlwhCwichzwV+BapdQnDtc9g6W8AL4CTFdKHUjzOeuwrOCuWO92RcK7jcRSqMcBO5VSu2zudy9WD+j1mKtpssNzewMvxt17OdCEpYyfxGoUno25Vu4RkWjsmiiwKe66R7Esa4PmBB6GZfAGESkFrgQKRGRz7HAxUC4ipyqlPsSyJp2uL8bqzl8DvKSUaoj5UiV2yiasrncLx6UhVvzz1gNzlVKjk5z/v1jugVdFZJxSaj+wP/ZbGbAn9v8xCdcd9l6xtJgO/FYp9VqS570OdBWRIVgK+9Yk58Lh79wLy6rdjvVuTyqlvp54gYj0ADqLSLlSqvYwoZXaC3wP+J6InAzMEZH3lVJvJNxmPXC9Umqeg1y/BH4p1kDyq8CK2N96rN5VY4r3MmiGsajzl4lYVtYgYEjscxLwLyzlC7AF6OdwfRGWYt8GNIoVCjYm7ve/AteJyEkiUgb8LEP5XgFOFJGvikg09hkuIiclnPdtLEXzioiUKqW2YfmGrxaRAhG5HsuPnozHgGql1D3JToopsBewLNvOwD9S3PdqERkUe/9fAS8oKwLjKWC8iIyNyVgiIueKSE+l1CbgNeD3ItIp9t7nAIjIF0TkBBERrEaoKfZJ5BHg1yLSO3ZdNxG5NPb/eSIyODZIuAer8WiKPfd14P+KSIfYgOTxIjIqxTsaNMAo6vzlWuBPSqnPlFKbWz7Ag8AksaIm7gZ+EusK/3f8xTHr7rtYCnkXlitgRtzvrwH3A3OwuuvvxH6qT0e42P3HAFcBG7EGJn+D1TjEn6ewBgDXAy+JSAnwdeB2rAG6k7F8ssm4CviiHB75cbbDuc8AFwLPp2F5PonlS98MlGClF0qp9cClwI+wGrr1MXlb6ttXsRRoNZbP/JbY8f7AbGAfVnr+Xin1T5vn/g4rL14Xkb1YA4tnxH47Bqux2YPlEpmL1XCA1UAXAR9j5ekLWO4Yg+a0jCgbDDkRs4SXAsWma20wuIuxqA1ZIyJfFJEiEemEZQ2/bJS0weA+RlEbcuEmrK79aixf6jeCFcdgyE+M68NgMBg0x1jUBoPBoDmexFF37dpV9enTx4tbGwwGQ16ycOHC7Uqpbna/eaKo+/Tpw4IFC1KfaDAYDAYARGSd02/G9WEwGAyaYxS1wWAwaI5R1AaDwaA5ZlEmQ0oaGhrYsGEDBw8eDFoUQ4ySkhJ69uxJNBoNWhSDDxhFbUjJhg0bOOqoo+jTpw/WekGGIFFKsWPHDjZs2EDfvhntyWsIKUZRt2GmL6rh3lkr2Fhbx7Hlpdw+dgATh1Yccd7BgweNktYIEaFLly5s27YtaFEMPmEUdRtl+qIafvi3j6hrsFbRrKmt44d/+wjAVlkbJa0XJj/aFmYwsY1y76wVrUq6hbqGJu6dtSIgidzl0UcfZdcuu01UsjvPYAgSo6jbKBtr6zI6HiZ+9atf0blzZzp16uTKeU60b98+q+sMhkwxro82yrHlpdTYKOVjy0tzvne6vm+v+NnPnDebaWxspLCwMOV5BoNOGIu6jXL72AGURgsOO1YaLeD2sQNyum+L77umtg7F577v6Ytqcrrv2rVrGThwIDfccAOnnHIKkyZNYvbs2YwYMYL+/fvz3nvvsX//fq6//nqGDx/O0KFDeemllwB4/PHHueKKKxg/fjxjxoxxPG/ZsmWcfvrpDBkyhMrKSlauXAnAU0891Xr8pptuoqnpcJfR9u3bOeuss5g5c2ZO72gwOGEs6jZKi4XrtuWbzPed671XrVrF888/z9SpUxk+fDjPPPMMb7/9NjNmzOCuu+5i0KBBnH/++Tz22GPU1tZy+umnc+GFFwLwzjvvsGTJEjp37syPfvQj2/MeeeQRbr75ZiZNmsShQ4doampi+fLlPPfcc8ybN49oNMo3v/lNnn76aa65xtp2csuWLUyYMIE777yT0aOT7dNrMGSPUdRtmIlDK1x3SXjp++7bty+DBw8G4OSTT+aCCy5ARBg8eDBr165lw4YNzJgxg/vuuw+wwgo/++wzAEaPHk3nzp0BeP31123PO+uss/j1r3/Nhg0buOyyy+jfvz9vvPEGCxcuZPjw4QDU1dVx9NFHA9ZEoAsuuICHHnqIUaPMHrEG7zCK2uAqXvq+i4s/3/c2Eom0fo9EIjQ2NlJQUMC0adMYMOBw9838+fNp165d63ellO15J510EmeccQYzZ85k7Nix/OEPf0ApxbXXXsvdd999hDyFhYWcdtppzJo1yyhqg6cYH7XBVbzyfafD2LFjeeCBB2jZtWjRokUZnffpp5/Sr18/vvvd7zJhwgSWLFnCBRdcwAsvvMDWrVsB2LlzJ+vWWatRigiPPfYY1dXVTJkyxevXM7RhjKI2uMrEoRXcfdlgKspLEaCivJS7LxvsS9THT3/6UxoaGqisrOSUU07hpz/9aUbnPffcc5xyyikMGTKE6upqrrnmGgYNGsSdd97JmDFjqKysZPTo0WzatKn1XgUFBTz77LPMmTOH3//+956/o6Ft4smeiVVVVcpsHJA/LF++nJNOOiloMQwJmHzJL0RkoVKqyu43Y1EbDAaD5hhFbTAYDJpjFLXBYDBojlHUBoPBoDlGURsMBoPmGEVtCBUzZ85k6dKlQYvhOo2NjTz44IPU19cHLYpBQ4yiNoSGtWvX8sgjjxwWknbXXXcddo5bS49+7Wtf44UXXnDlXok8/vjjfPvb3279rpTilltuobKy8rDZlwZDC0ZRG0LD8uXL+eMf/0hBweczHxMVdRgRER588EHOOeecoEUxaEpKRS0iA0Rkcdxnj4jc4oNsrUxfVMOIKW/Sd/JMRkx5M+clMw3hY+LEifzkJz/h3HPPZerUqQBMnjyZuro6hgwZwqRJkw47f9++fVxwwQUMGzaMwYMHty5lCnDHHXcwcOBARo8ezZe//OXWxZkSmT17NmeffTYnnngir7zyCgBNTU3cfvvtDB8+nMrKSh599NHW8++55x4GDx7MqaeeyuTJkwG4//77GTRoEJWVlVx11VVHPGPbtm1cfvnlDB8+nOHDhzNv3jwA5s6dy5AhQxgyZAhDhw5l7969OaSeIeykXJRJKbUCGAIgIgVADfCit2J9TqZ7+xk85rXJsPkjd+95zGC4KPlaGY899hidO3fmwIEDVFVVcfnllzNlyhQefPBBFi9efMT5JSUlvPjii3To0IHt27dz5plnMmHCBBYuXMi0adNYtGgRjY2NDBs2jNNOO832mWvXrmXu3LmsXr2a8847j1WrVvHEE0/QsWNH3n//ferr6xkxYgRjxoyhurqa6dOnM3/+fMrKyti5cycAU6ZMYc2aNRQXF1NbW3vEM26++WZuvfVWRo4cyWeffcbYsWNZvnw59913Hw899BAjRoxg3759lJSUZJyshvwh09XzLgBWK6XWeSGMHV6ub2wIDw8//DCvvvoqhYWFbNmyhZUrV9KlSxfH85VS/OhHP+Ktt94iEolQU1PDli1bePvtt7n00kspLbVW8xs/frzjPa688koikQj9+/enX79+VFdX8/rrr7NkyZJW//Xu3btZuXIls2fP5rrrrqOsrAygdUnVyspKJk2axMSJE5k4ceIRz5g9ezYff/xx6/c9e/awd+9eRowYwW233cakSZO47LLL6NmzZ8ZpZsgfMlXUVwF/8UIQJ/J5b79QksLy9YK5c+cya9Ys5syZQ1FREaNGjeLgwYNJr3n66afZtm0bCxcuJBqN0qdPHw4ePEgma9sk7vQtIiileOCBBxg7duxhv/3973+33Rl85syZvPXWW8yYMYM77riDZcuWHfZ7c3Mz77zzTmvD0cLkyZO55JJLePXVVznzzDOZPXs2AwcOTFt2Q36R9mCiiBQBE4DnHX6/UUQWiMiCbdu2uSWf4zrGbqxvbAgHu3btomPHjhQVFVFdXc38+fNbf4tGozQ0NBxxze7duzn66KOJRqPMmTOndWnSkSNH8vLLL3Pw4EH27duXdPus559/nubmZlavXs2nn37KgAEDGDt2LA8//HDrMz/55BP279/PmDFjeOyxxzhw4ABgLYfa3NzM+vXrOe+887jnnnuora1l3759hz1jzJgxPPjgg63fW9w4q1evZvDgwfzgBz+gqqqK6urq7BLPkBdkYlFfBHyglNpi96NSaiowFazV81yQDbDWN473UYN/6xsb9GDcuHE8/PDDVFZWMmDAAM4888zW32688UYqKysZNmwYTz/9dOvxSZMmMX78eKqqqhgyZEirNTp8+HAmTJjAqaeeSu/evamqqqJjx462zx0wYACjRo1iy5YtPPLII5SUlHDDDTewdu1ahg0bhlKKbt26MX36dMaNG8fixYupqqqiqKiIiy++mF/+8pdcffXV7N69G6UUt956K+Xl5Yc94/777+db3/oWlZWVNDY2cs455/DII4/w29/+ljlz5lBQUMCgQYO46KKL3E/YPCTojZW9Iu1lTkXkWWCWUupPqc51e5nTfE38sJBvy2nu27eP9u3bc+DAAc455xymTp3KsGHDghYrY/ItX3IlMfAALKPOr/XQcyXZMqdpWdQiUgaMBm5yU7B08WJvP0Pb5cYbb+Tjjz/m4MGDXHvttaFU0oYjyefAg7QUtVLqAOA8xG4whIhnnnkmaBEMHpDPgQdmZqLBYMgL8jnwwChqQ1p4sWWbIXsyyY+2MrM3yI2VvSbTOGpDG6SkpIQdO3bQpUsX21hhg78opdixY0dasxUzmdkb9kH7FlnD/A5OmM1tDSlpaGhgw4YNKSeZGPyjpKSEnj17Eo1Gk543Ysqb1Nj4aCvKS5k3+fzW72GPmMgHco76MLQdnKyqvn37Bi2aIQvSHWALKmJCByteBxlSYRS1oRWzAFb+cWx5qa1FnTjAFkTEhA7lTQcZ0sEMJhpaSWZVGcJJugNsQURM6FDedJAhHYyiNrSSz3GobZWJQyu4+7LBVJSXIli+aTu/cxAREzqUNx1kSAfj+jC0km432RAu0pnZG0TEhA7lTQcZ0sFY1IZW8jkO1ZCaiUMrmDf5fNZMuYR5k8/33EerQ3nTQYZ0MBa1oZV8jkM16IcO5U0HGdLBxFEbDIQjRMuQ35g46hBgFEVwhCVEy9B2MT5qDWhRFDW1dSg+VxT5uiaDboQlRMvQdjEWtQbk8zq6YSAsIVq6YHp//mMUtQb4pShMBbMnLCFaOmDcRMFgXB8a4MesMONecSYsIVo64NT7+8WMZW1iKdWgMIpaA/xQFMYP60y6s/cMzr282roGYwR4iHF9aIAfsZzGD5scsy9neji5iRIxYyzuYhS1JnitKHTxwxo/eW4EnX63jx1wxLrVThgjwD2M66ONoIMf1vjJc0OH9LNzE3Uqs9+8wAzGuoexqENKppZVNu4Vt603E4aYG7qkX2Lvz2l3GDMY6x5GUYeQbEOkMnGveBGGZfzkuaFr+oVlvYwwYxR1CPHDsvLiGbr4yYMi1x6KDunn9A5mMNZbjKIOIX5YVl48w24gqq10kbPpoSQqxfMGdmPawprA0s9MdjkcPwd2zWBiCPFjgowXz2jL8cqZxrHbDRxOW1jD5adVBJZ+Jhb/c/we2DUWdQjxwzL16hlttYucaQ/FSSnOqd7GvMnnuy5fOujqIw8Cvwd2jaIOIX4M3pgBInfJ1L+so1LUwUceJPGuDqdV/L3KH6OoQ4oflmlbtX69INMeio5K0YwxpJ7o41X+hF5RBz1Tq61g0jk3Mu2h6KgU23Ivy87VkYiX+RPqrbicAu0TB1iMksmNdNPZ4C6m3PpLsvTuO3mmo7tDwJX8ydutuNJx6JuQotzRZUZcW8Nt15NR/M6k0hNOrqiK8lJfBnfTCs8TkXIReUFEqkVkuYic5bVg6ZBswGX6ohpGTHmTW55bHEhIUcvz82F9Xh0HtgyZocM6ITqTKvTQbq0cwUpHP+p3uhb174C/K6X+j4gUAWUeypQ2Tq1cx9JoSse/l0om36z4ZANbxkrTD7s8Mb2i5KQyRuL98zW1dQi0ukL8qN8pLWoR6QCcA/wRQCl1SClV64k0GeK0IpwIgY3OQvoTA8JidTul83kDuxkrTTOcLGenNaRNr8ginQleE4dWMG/y+VSUlx7hr/a6l56O66MfsA34k4gsEpE/iEg7zyTKAKeZbrUHGpJel2p0NlcFmo6rIExdUad0nlO9zcxU0wwnI6FAxPZ8v8L9dDdKMlkGOAhXYDquj0JgGPAdpdR8EfkdMBn4afxJInIjcCNAr1693JbTEbsBl5buiR0VKbrnbrgt0omBDVtX1C6db31use25xkrLnlxdSU5p36QUpdGCQML9wuAKzCT0MIgY93Qs6g3ABqXU/Nj3F7AU92EopaYqpaqUUlXdunVzU8aMcWodf/ulIcybfH7SwuHGegbptM75MEDnx5ojbQk3ellOad/SCwpinZCwrBHS4tpYM+WSpHoiiE04UlrUSqnNIrJeRAYopVYAFwAfeyaRC+QSmO+GAk3n+TrOPMsUHSdlhBk3elnJ8iSomab5YJTEE8TEn3SjPr4DPB2L+PgUuM4ziVwim0I5fVENERGabCYBZapAUz0/H5Rcy/v9YsYyauuscYGSqFmQMVv8MhL8JkijxKuoJL8bvbQUtVJqMWA7YyZfaOl22ilpLxSojhUqW+obm1v/33WgQTv/Y1hwS6HptkZLUEZJrr5xnUJPQz0z0U2c5vIXiGTty0uV0bpVqGwI26CozuRDL8uOoIySXMqmbgOgRlHHcOpeNiuVtZLWKaO9It/8j0HipNAARkx5UwvLLluCMEpyKZu6GSB5p6jT6a7YneO2H023jPYKp3SLiNB38szQKpagSLXDd742+F6QS53WzQDJq5GfdMKbnM45b2A3V0NudMtor7ALVQIrblf3iTxhICyhbTqSSxidbqGneaWo0ynUTuf8Zf56V/ej0y2jvSJx1qLdDDijWLKnrTT4XpDLHp1BxEonI7SuDzv3RTqFOtnMrWkLa1ybBJCvA0N2xHfX+06eaXtOTW1dm3WF5BI9oHu8vU6REXZk6xvXLSorlIrayW9XXhZll806H/GF2qngg7s+ZDcyWvdKYEey9I13hUDb8LHm6mPOtcH3sgzlu/9cp6isULo+nNwXSpGyu+LkU23BzS5lulNS7QjTok3xpEpfaFuukFx9zLl0370uQ8Z/7h+htKidlOnuugb+50tDUsYuA3zvrx+6MgPRK8IaNZLYk/B7t2bdcGu2YTZ57nUZMv5z/wilok7mt0unULf8rrMPOcyVID4PRkx5U2sfq9ek447zCq/LkO7+83wilK4PN0Zkc+lS+oGbUSNBrgWs2+i5n0xfVMO+g41HHI8WiC/v73XkUVvOW78JpUXt1oisToMFibgVNRL0gI9uo+e5ksng3L2zVtDQfKTzp11RoS/v73XkUb7lrc6IsvHT5kpVVZVasGCB6/fVEa9H1XO9t5Prwa/dk/OJxEYPLMXn1BPrO3mmrY9egDVTLvFO0DjCGDnUVhGRhUop28XvQmlR64LX1qobFn8YfN1hUSaZDs7p4MPVudfoNWEpV+kQSh+1LmQbnuSnz1j3GZJhCkPMtNEzPtzgCFO5SgftLWqdW8VsrFW/fcY6zZC0y0svQsi8KjOZWsjGh+sO2eRnWMNbndBaUaer1IJS5tl0bf0uQLooC6e8tFsDHLJ3zXjZEGbT6LVl14MbZJufYXD5ZYLWijodpRZkVEM2FTeIAqSDsnDKywKXtj5L9Rw3GkJdGr18IF3jKtv81GF8wE20VtTpKLUguzjZVNx8K0DpkmwxrNJogWuuGS8aQp3db7kSxLtlYlylm5+J73HewG5MW1ijhcvPDbRW1OkotaC7OJlaqzr5jP3EKS8r4nzVbigLtxvCoOPQnWRyY7Gvmto6BFpDCOPfDQ7ftLhTWZSfjz/ZlXfOxLhKJz/t8mjawhouP62COdXb8qKB1VpRp6PUwmahttXuc7K8dNM143ZDqNuglBsbtsZfn+h0qmto4hczlrG/vvGwyTq7DjRw+wsfpv2cZGRiXKWTn055NKd6W97MFdBaUaej1MJooergM/YbLxqoZJalW89JpVT8dh3k2nA4beIcT4sVnUhDk3KlgcrEuEonP4PuVfuB1ooaUiu1tmqhhhE3G6hUlmX8YPO9s1Zw63OLsyobyZRKEG6RXJVSrsrLDeWXqXGVqtzo0Kv2usHWXlGnQ1u0UNs6bkYEJatkyZRKEG6RXJVSso0dwHq3kmjEdsW/TJ6TDLeNq6B71X402PkxM7F+Lyx6Clb8PWhJDD6Ra0RQC6lmsCVbZTGILneusx3trm/Z5bLl3X4+/mSikSP3vnRz1b9cN9WIn9kLpFwJ08vZwH5soKCtRZ12V2LnGph6LhyshbKucPsqsNlgNefnGLTCrYigdKxipx5bEF3uXK3RTK73KuojF5ys17svG+w4cOi1xetHg62lorZL2F/8bSG9173Aul5f5N7XV7UWssd7vkz/+r0w7Fr44M+wdxN0ODbr5wQdemVID7cignKpZEF1uXN19aVzva7uxGzcTV67qPxosLV0fdgl7Iim9xm6+Ge89uLTrd3UrbV76bRyGpu6j4Ihk6wTN1khRChluUP2bcvoOWbPt3CQzCUBViO8v/7IRfvtlLkd6VSyVDIESZCbRXhJNg2r1xavH4tvaWlR2yXgMbIDgNPVh8yiEoALIh/QVXbz3zvO4r7uJwNiKeoBF8HSafDSt+DMb8G4u9J+TrLjBr1wsvrs1o0G++57rlZxtpanncsN3Blgy+eeYjbWq9cWrx+RZ1oqaruE7S61AIyMtMycUlxX+Hc2qs68uPck7ituD137w6YlcGg/vP5T67RlL8KYO2HRk9DvXOjUO+lzWo4bwotTrHCZzc4qQYR32inS25//EMSKVW45lq1y1W2Sjptk07Amu8atMSqvXUVaKmq7hK2I7AJgQGQDR7OLoZFVnBGp5icN13FMeTvrpB6nwrp3YO49sHcjDP86vP+/8PpP4N2HYMjVMPGhpM8pjRZw+5gT4YMnoXYdnHgRVAzLaIDSECyZ9pTctIqznXRit2VXtso1n3uK2TSsTtcAoel5pKWoRWQtsBdoAhqdtotxC7uEPaPdIer3llN8qJYvFrzNlwve5JPmCqZHRnNnS2va41T46Hn49/0w9Ktw4c8tS/rdmHKufgWafgsF0SOeU1NbR4EIdQ2N7J/5Y2iabl3z1r3QdxSM+C70GArtuuT0bibKxHv86Cnl4l7IRGFmo1zzoadot8hSLut22DXGI6a8GZqeRyaDiecppYZ4raRbSIyz7Nq8g+L+51Ff1JkfRv/CcbKV3xffwJ2XDfk8UXucav3teByMuxuKj4ITx1rHRt5qhfB9OveI57QMBjSrJu4sfIxJTdP5S/NoZl70bxg3BTZ/BE9dDvceDzO/B/X7snonL3adyNdBo1zwY3Anl4HoTBRmNso17DvL2NWTp979zPXdWsLU89DS9XEESsHezdChguLx98KOVRSccjm/7dr/8PMqqmDAJXD29ywlDTD6V1D5JTjhQnjvD/Dxi9D/ws+vWfYiha88weOykQ5FBzgp8hmPNI5nSuNVVMzZwiWTvwFDr4YNC6B6Jrz/B1g7D67/O5SWZ/QabvsO83nQKBf88Dsnq+Spek23jx3Arc8ttt34Np5slWsu759Lj8+t3mI665G4YfmGqeeR1i7kIrIG2IW12NajSqmpyc53fRfyulr4TW8Y82v4j29nf59pX4eVr8Oty6CoHcy5C966h+2qA6vVsSiEWU1V/KlpHCD2u0Wvmg3PXAW9/wMmvQCFRWk/3nlXasWaX4+Df/wUGuvh4nshUmBzpkX8MpV2eLbDeGM9RKIQ0TKq01ecdncvL41S39iccqfyPpNnOt5bIBC3WKa7rGd6bbqK3KmeJJLrbu65vK8XuLEL+Qil1EYRORr4h4hUK6XeSnjIjcCNAL169cpJ4CPYu8n6e9Qxud3n9Bvho7/Cv+6zFM5b98DQq7ns44l8tvvQEafbtqwnXAgTHoDp/wV/Hg/DvmqFAp75Teg/+sjzleLV+ct4YM4qoIB4b9M5kQ+5o/BPlEaa4ImTYN3b1g/NDTDiFsvl8vFLoJqg6wA447+Y/snBpFtYQeo9G1NWlk0fQmknOLgH3nkICgotJb3sReh8PIy5w/5d02XHatiyFNp1s9xVRe2yv1cO5GIBOg1Ei5BWr6kiyfrcQS3NmUuPL9W1mfT+Uq1HEn9eLoRpQbe0FLVSamPs71YReRE4HXgr4ZypwFSwLGpXpWxR1GnOOHTkuOFw6ldg3v2W8ht6NUx4kNt6bcws5GfIl60ByZdvtmK1C4qsaJOvvQI94xrEXWvZ8pdvcfHWt7kY2FDcldlNwyigmdMj1QyIbGC1Oha69oN186weQ90uqyH54AnrHu27Q3EH+HgGvPswK5qvpq5hJABfKXiTKwrmUkwDbzefwkONl7Kb9o4FOLGybKrdT/TFG9g3u4b2o74D3U+GhX+GJc9+flFxB4gUQtMhOPUqy+3z9P+h+sRv8J/rRrN3906+1H4xX+/0AUfvWQoDL7EU8KYPoeI0KySyy/GweSlsWgw1C2HlP2hdCbmwFDr3g11rLBdXWRdLjh6VLNjfjReXbOf9fV3Z3+EEbh830LfV91LhVMlvfW6x7fmJjaedogc4cKiR6YtqAlEWGftslbKioQ7tp/PuZYwtWEFX2c0hCpnffBLVzcexp/YAfSfPJCJCkaqjA43UU2R9Ghp44bXXmXioAZrqobAECov53aADPPH+FvY0FVBPEc0qQqnUU8ZBIig+UT05WNiBH599NDQ1tAYHANBwEDYvgaL2UFgM+7fDge3W37pdlgFy1DGWcdDlBGuQccixlvGw4X3Y8Ta8AbQ/Bo4+ySr7+zaDFFj32bXWqhON9bD9E9i+EhoPWmVYBKJlcOWfXc+blK4PEWkHRJRSe2P//wP4lVLKcQUk110fi5+B6d+A7y6yEiQX9m6B358JPYfDVU+3ZnIq68r2975NsGudFb/92FirIEx8BAZeDCteg2k3cOBQI1MbL2avKuPsyEecEVlOPVFWqON4u+hs+o/7JhOq+lmx30XtQCnenj2dN95bxOr9pSwtGoKKFNC1bg13lT7N8OYPmd88kCIaGRpZxUfNfdiljmJkZCkHKGaWOpNep49neNUZVoHqUGFZxMA37nqQK+v+SqkcYpPqTJRGvlAwn7VSQR8VG5iRCIy8zVK2jXVw2tegpNyqlJEINB5i3RM30fuzv3FAFROlkag08ZnqjvSs4rhtc61C3PVE2FZtNYitiBXHPvgKGHAx7N9muZJ2roEuJ1hy7t0Mmz9CbfsEobn1yhrVhZfVOfQb8w3GjDwjtzKAs+siV4s2k/tOX1Rz2HoaLXjS/a7fB5/OsRrK2vVWeevcz1Jch/ZBv1Fc89xaTt4/n3LZSxn1lEk9HThANBrl3Asuhk59rQZ73b8tw2JbtWWkNH3eG21QBRTQTEQ+1yvbVQeaiLTOhQCoU0UU03DYeVkhEejYEzr0tMrrthXQcCD964s7QHMTNOyP3S/mclQOPdZIITQ3Ws/t1Mcq54XFlgKXiGVYfeW57F4liesjHUXdD3gx9rUQeEYp9etk17iuqN+6D968A368GaIuOPoP7rFa3DR9rWn5snathb9eA5s+ZDud6Kxq+STSj+vrbmYjXY+4p5N/zWlWnXVNM98oeJnxBe/QSIRpTefw56YxKCKcKOv5btksxsp8oo37P7+oU1/Lr1+ziOZFT7OFTnymjuZ42UhX2cPvGydwb+OXWHNLH0txdqiAbicmTY8Rd7/Bf+z7OyfKBg5SxOymYXyojqeivIx5t50JKKvRqdtlKYYdn1rWScVpUFSWRorDeXe/RmTPeopp4JTIGsZF3mdU5EMKREGfs63Y9kjUklnEWpCr+yDoPdJqZPZthn1bYP8OK9qncz84prK10bLzg3ZiD0Miq/nTD2/I2s2WsqzU1VqWY7TESssjFLuil2zl2KOiPDv5K5aLaMMCy+Jr381SEp/OhT010FBnleNOvS3XWLcTrfyOllo9mn1b4eBu2PkprH3bUmSRQijvZfVkdq62rEEpOEwxHVDF1FFEHcXso4ye7RXt96//XMSi9tDrTMtt1dwIJR15b3dHbp9fxrqGjrSjjhGRpfSQnRRziH6yiQKa+VQdy0GilHKIctlHHcXsL+nBj779DSgt59VFa/nDnI+p3bOXnkdFuP7MHpzbr4P1jGhZzJBpgq3LrYZGBHbXWL2xPRut9+7UF/qNsiztpkNWuWjXxfpbWg4Hdlplpn6vdZ9da6w06Xoi9DrL6v1JAezfaqU9wFE9QDVDSUerQWiqt5RyYXFWZcSJnBR1NriuqGf+txUfPXmde/fMgHStpBkL1vDJjHs4tnkTuziKBxq/SD3FtgMjTpab07Piid/nDhIUQWM9bFvBewvf55+LP+HSQ68wILKBxoJSXmIUP9t/BfspRWiml2xlnepORXlZRvsWOg+K5ja4k+oZPdjBFQVzua1iOexYCc2NHCzqxP76JjqqPRRKzAIvLLWUkh1SAN1P5oXN3dnXaE1iKqKB7lLLyMhSiiVm3XbpD8edYVXuQ/vhwA6rcWh3tDURqiBqVda6WqshkIhV2dt35931dcxY9Blr90cpb1/G109qYGjJZqvR2rDAUigVpwHwwaebKKGBYg5RLA10YD8dJCZ7i/UGlnJvPBhLiFMt+aKlUL/H6pFsX2n/zsUdLEXe83Q45TIrMirWSNDcZN1fKVg5Cw7u5rX6U7lz7o4jy8H+7ZYLUik4elBrgxdPfK8zHa0SX251G9gLAjcGE4Nl7yarVfOZVNEVib6738xeQ82hLxxxnp1idfJ/pxPDqbAUva1SLSxm+uYu/HB+BXUNx/AoIxgo69lQeBzjT+tL88IaaGhCEWGdOobSaAHnDeyWkb/Wi7CmRNdSeVn0iMXrN9GFaUdN4rZvng9Njby0uIbJ05dT19BEEQ2cKOs5P7qUCb2LOWFgpVVmyrpaoZpbP2b50g9YuGojx9es4PzIO6293EYi7FAdeb75XPqOvJIR7TdbFujKWZaSjpZaCm/5jJgkcTlaWGK5hpoOQd1OAM6MfSgCDgEfAtF2lsU/6vuxnsYHUBClsbAdNQ2WL7ZeRdmvSliuenFUaQk/Gi5wzGArwqjjcZZ13NxkP+mquRl2fwa1n1kydz/Fcgkkm1EbKfg8umjQpQBcBFx0ls257bpanyTETypxMjgKRGhWynbrtLBMPgmC8CjqDv4q6mQuiBYSFZOTkk2qWG3umcqiTuVHjS/0TRSwTPWBBphTvY27Lxt8hOWcaSVxe3lPu4G9aESIFkjr2hdHPKOgkHv+sbr1mkNEWar6sfRQP6bVlDLvqwn+4I0d+eHH7ZPmZ3lplMVjx1hfRnz3yBP2brH8ueW9rS64aj7cFdfSrW44YFnDdbUxf31/S9HauNo2OlmSlwyGxLRPFrcfiflMO/VxPsdHnMqIk4UcpsknQaC/om5qgK3VMHSSr49NFXRvp5iclGwmA1RO0QDJnptIskJvN5U23UiFFtwOa3Ja+6K8NEq74kJXNjVNZxLFbodNXVs5qjvQPfbFpuqUdbY+GRCmELFMyPS9dJ58osOyD/or6o2LrBHZPiP9fWySlrzCIbPcsDQTC3jH0igiUHugIe1Ckmmhz6aSuLlamFNa765rYPHPxzhel4nc6Vhmidf5VUF1XaQ/VzJ5r6D3PXRCl9m/+ivqtf+y/vb2V1FnYx27ZR3lWnEzLfRBV5JsralM5E5nU9f463SpoDrjZkNmV3fOG9gtpx3k3UAX37k2itox09f8C44+OedV6zIlW+Wlg3WUaYMR1JrM8b2GpP5oBzKR2y4/W4YE7XpIulRQXfGiIYuvO7o0lLr4zrVQ1E6ZIs2HuHT9fGvJUp8Ju+8w0wbDzwYmMb9r6xqIRoROZdGMXDyQvtyZ5qcuFVRX3GzI7Iw0LxvKTHoCuvjOtVDU985aQXNDHd8v/BsLm/vzRvNp1DU08eqs17i04QD0PTsQuXSwjvMRp8HDsqJCFv3M2SedK5nkZ7IKqsPgEgQ7yOVWQ+ZkpDkN/ObaUGZqqQftFmxBi2XQNtbWUU+Uywve4tKCf7ce77s/tlFtr/8ISDKDF4TBWnVa07kl5tzttZEzxYu1zTMhl02B43GynAsc4r9ztWQzXUd84lA9NjDWQlFbiS/Maz6FsyLLaJlMMKJ4tTUDy2f/tMFb3KrkXuJUQedUb9Ni5/pcNi5wA7c2J3BqnJuU8mTzg2yMhMRNTILoPWmhqFsy/d/NJ9NN9jBA1lMajXB64SprGq8hbcKw40tYdiCxq6CpNgzwK+2D7pW4ZWk6Nc4t93Pbkg2DkWCHFj7qlsT/89/3QP2jXNJuBSeNqqT4jV1w3OmByqaLPzIdvBgp9+L9wzxQ6+S77lga9TVKQYdBLjfGcJL5gL0YI9LF55wp+i3KdP8wa8nLkydaS5t+cz4cPdBV+dIlbAvFuL10Z9je3w+c0qQkGjlibRLwbiMAv/PGS4PFb2NIV+MrXIsy9TsXljzXunwiXZMvueklYYuldbs7HLb394NcNwzwWg6vlLSXvYXE+GmvJ7mEMZpLP0V92rXUL/kbxavfYE7Tqfzknn8G1uIF7QfMFLe7w2F7f7+wq+hOqyx66YrwS+H41WDrMslFR7QYTIxn+uaujKq7j981fpGHGycEFv4E7g48+DHQ5PYgXVADL2EYEE3ELu0FS9mE5R2c8KvBDjqSRWe0U9T3zlrB5oZS/qfxCt5TJwHBZZZbis+vmFe3Yz6DiM4IOj44W+LTHg5fgzws7+CE2w22U0NsenDOaOX6mL6oJu1F+v3ALT+gn75eN7vDbrx/pgM3YfaLt6S93aBuWN7BDjcjJZK5N3SIZNEVbRR1SwY6EVRmuTHQEWZLIRfFn43PMZ200nXUvoUw57cdbg5cJmuIwxo65wfaKOpkC7vrkFm5DHSEzVJwSxFmYx2nSqswDDiFLb/Twa2eWqpNLSDY+HpdjQBtfNTJrA0d4nZzGegIy0w8cNdHnI1lmSqtwjDgFKb89ptU/u4gp2unU/aDGujWRlEnm0oatJKG3LqzqQb5dIpycFMRZjMIlSqtnNK7JbpChzQMaiEfncqREzo3YqnKfpAD3dq4PnT3T+XanXXqOurWlXfTv+rF5gtO+dASCgfBp2HLc4OwBnUpR0744d7I1n2RquwHOdCtjUWty3KCTnhlCQTVlXeyvtwMxfIiT53ilRMXQtDNHeI1YXAJteCleyMXqzdV2Q9ykFgbixpSWyFBOvq9sgSCyPxk1pfbPZtMLctUeWyXDzqFdAaFSQOLXKzeVGU/yEFirRR1MnTo2nnRnXUz89NtyJIV5pYFhIJoENPN48R8cFqMKsxRFpkwfVGNba8C2kYaxJd7pyXm0h1LAueyH6R7NjSKOswTIZLhVuZn0pClsuKDWrQm2zzWYXwjyN7evbNW2CooAW3GeLzCbhVBO3IdS2r5DYIxYkKjqHWbROBWxQxi9qOucb7Z5nHQ8bdB9/ac0kf59PwgSTb/ogU3G+2gjJjQKGqdlIvbFdONzM9EyelggdqRSx4HuXRl0L09p3SraANuj2SNuMBhjbauk1nSQZuoj1ToFH+p4wh7JtEafkTYZBPTq1MeZ4Ifvb1k6RnWdHODZPMv4qNKwrrYVwuhsaiD7t7Go5sbBjK3kr20QLPtceiUx5ngdW8vVXqGNd3cIN1yH3SvJ1dCo6hBn50ZdHLDtKBTZc2lUuiSx5ngtSspnfQMY7q5QbrlXkfjKhPSVtQiUgAsAGqUUl/wTiT9caqY5w3sxogpbwa2r5wuldXNShEGv6LXjWTYlYzXpFPudTSuMiETi/pmYDnQwSNZQoNdxTxvYDemLazxZOQ/6KiCTEm3UqRSwm69tx/K3stGMuxKRgd0HUBPl7QGE0WkJ3AJ8AdvxQkPidNg51Rv82yAUcfBy2SkM7iVzuCOG+8d9kEkaNuDhW6h+xIVqUjXov4t8H3gKKcTRORG4EaAXr165SxY2PCyexq2rm86roB0/K5uvHfYB5EgeXpm21sIg0vJbXRxDWZDSkUtIl8AtiqlForIuU7nKaWmAlMBqqqqnGZy5i1edk/D2PVNVSnSUcJuvHc2yl5HJWaXntm6hvx0pemYlmEkHdfHCGCCiKwFngXOF5GnPJUqhHjZPc3Hrm86cd9uvHemqwGGyVWSrWvIL1ea12kZhvW33SKlolZK/VAp1VMp1Qe4CnhTKXW155J5hFeZ66UPLNW9w1hg01HCmaapXTo4LYvastFAYlrpMh6QTp5m6xryy5XmZVqGqUF1g1DFUeeK110+L31gTvd2452C6J6mG9KWbpo6pcPdlw3m7ssGc++sFdTU1h22ypxdWumwXGi6eZqta8gvV5qXDUI+jD1kQkZTyJVS/wxzDLUu1pKb5PpOQVomiZEzuVSwVBV33uTzqSgvTbrBQMtyoXb4OR6Qbp5m6xryy5Xm5iYUieg0wO5HjzY0a324gU6Z6xa5vlO+NF7ppEM6Wy3psFxounmarbvNr1A1LxsELxuBTPDL0GlTrg+3unw6jWTn+k75MoswnXRIdY4uy4VmkqfZutv8CFXzcsamLhNY/HLBtCmL2o0WXrdBjFzfyS3LJOh0SScdUp2TbCU2Pwds8ynKx033VuJ9dZjA4lcvvU1Z1G608LoNYuT6Tm5ZJkGnSzrpkOqcZGu4+DmFX6cFtnRGhwksfg3MilLuz02pqqpSCxYscP2+OtB38kxHP+aaKZf4LU5S0nVFuOGyCFO6JMMuLVoiRhKpKC9t3WPS4A86uR1b5LFr3LOx7kVkoVKqyu63NmVRu0FYZglmErbnhmXiR7oEtbjSrc8ttj03zIPQYUTHxcn86v20KR+1G4TFf+h3NIfX6RKkD1yXCIO2jq4RSl754eMxijpD7AYxLj+tgntnrdBqZqDfoYheD+4EWUnD0jjnO/kYXpsuxvWRBfHdYx27YxCMi8bLwZ0gK6kZ3NMDP8u0br5wo6hzJOhoByd0iTN1i6DHBtxqhHRTAGHCrzKto/FlFHWO+GnpZVLJvbACg1Qy+dDw6KgA0kGXxsWvno2OxpdR1Dnil6WXTSV30xURtJLJB/eDjgogFUHneyJ+xE7r6AvP68FEP2aT+TXQFPSId9DP18WqywUdFUAqgs73INAxyidvFbVf4Vx+TWUNupIH+fygp6e7hY4KIBVBlzuvSGbE6Rjlk7euDz+7mX50x4IeTAvy+WF0GdgRRj970OXOC1K5c3R0s+WtRZ1vlkDQrXyQz8+XvNRlIaFMCLrceUE67hw/JrFkQt5a1PlmCQTdygf5/HzKSx0WEsqEoMudF4Sx4c9bRR3GbmYqgq7kQT0/H/MyTARd7twmjA1/3ro+wtjNNNhj8tLgJmF055hlTg0GQ5tDx3BPs8ypwWAwxBE2d07euj4MBoMhXzCK2mAwGDTHKGqDwWDQHOOj1gwdBzkMBkOwtBlFHQYFqNtKZQZDNoShroWNNuH6CMuiPm1xpTJDfhGWuhY22oSiDosCDOPUVoMhnrDUtbDRJhR1WBRgGJfBNBjiCUtdCxttQlGHRQGGcWqrwRBPWOpa2GgTijosCtCsaWEIO2Gpa2EjZdSHiJQAbwHFsfNfUEr93GvB3CRMSzWGbWqrwRBPmOpamEi5KJOICNBOKbVPRKLA28DNSql3na4xizIZDAZDZuS0KJOyNPm+2Ndo7OP+knsGg8FgsCUtH7WIFIjIYmAr8A+l1HxPpTIYDAZDK2kpaqVUk1JqCNATOF1ETkk8R0RuFJEFIrJg27ZtLotpMBgMbZeMppArpWpF5J/AOGBpwm9Tgalg+ajdEjATzNRVg8GQj6QT9dENaIgp6VLgQuA3nkuWIWadDIMhO4yBoz/puD56AHNEZAnwPpaP+hVvxcocM3XVYMgcszZHOEgn6mMJMNQHWXLCTF01GDInmYFjrGp9yJuZiWbqqsGQOcbACQd5o6jN1FWDIXOMgRMO8kZRm3UyDIbMMQZOOMirHV7MOhkGQ2aYtTnCQV4paoPBkDnGwNGfvHF9GAwGQ75iFLXBYDBojnF9GDLGzGQzGPzFKGpDRpip+gaD/xjXhyEjzFR9g8F/jKI2ZISZyWYw+I9R1IaMMDPZDAb/MYrakBFmJpvB4D9mMNGQEWYmm8HgP0ZRGzLGzGQzGPzFuD4MBoNBc4yiNhgMBs0xitpgMBg0xyhqg8Fg0ByjqA0Gg0FzRCnl/k1FtgHrMrysK7DddWHcQVfZjFyZYeTKHF1ly0e5eiulutn94ImizgYRWaCUqgpaDjt0lc3IlRlGrszRVba2JpdxfRgMBoPmGEVtMBgMmqOTop4atABJ0FU2I1dmGLkyR1fZ2pRc2vioDQaDwWCPTha1wWAwGGzQQlGLyDgRWSEiq0RkcoByHCcic0RkuYgsE5GbY8d/ISI1IrI49rk4ANnWishHsecviB3rLCL/EJGVsb+dfJZpQFyaLBaRPSJyS1DpJSKPichWEVkad8wxjUTkh7Eyt0JExvos170iUi0iS0TkRREpjx3vIyJ1cWn3iM9yOeZdwOn1XJxMa0Vkcey4n+nlpB+8L2NKqUA/QAGwGugHFAEfAoMCkqUHMCz2/1HAJ8Ag4BfAfwecTmuBrgnH7gEmx/6fDPwm4HzcDPQOKr2Ac4BhwNJUaRTL1w+BYqBvrAwW+CjXGKAw9v9v4uTqE39eAOllm3dBp1fC7/8X+FkA6eWkHzwvYzpY1KcDq5RSnyqlDgHPApcGIYhSapNS6oPY/3uB5YDO63leCvw59v+fgYnBicIFwGqlVKYTnVxDKfUWsDPhsFMaXQo8q5SqV0qtAVZhlUVf5FJKva6Uaox9fRfo6cWzM5UrCYGmVwsiIsCVwF+8eHYykugHz8uYDoq6Algf930DGihHEekDDAXmxw59O9ZNfcxvF0MMBbwuIgtF5MbYse5KqU1gFSLg6ADkauEqDq88QadXC05ppFO5ux54Le57XxFZJCJzReTsAOSxyztd0utsYItSamXcMd/TK0E/eF7GdFDUYnMs0FAUEWkPTANuUUrtAR4GjgeGAJuwul5+M0IpNQy4CPiWiJwTgAy2iEgRMAF4PnZIh/RKhRblTkR+DDQCT8cObQJ6KaWGArcBz4hIBx9Fcso7LdIL+DKHGwS+p5eNfnA81eZYVmmmg6LeABwX970nsDEgWRCRKFYmPK2U+huAUmqLUqpJKdUM/C8edfmSoZTaGPu7FXgxJsMWEekRk7sHsNVvuWJcBHyglNoSkzHw9IrDKY0CL3cici3wBWCSijk1Y93kHbH/F2L5NU/0S6YkeadDehUClwHPtRzzO73s9AM+lDEdFPX7QH8R6RuzzK4CZgQhSMz/9UdguVLq/8Ud7xF32heBpYnXeixXOxE5quV/rIGopVjpdG3stGuBl/yUK47DrJyg0ysBpzSaAVwlIsUi0hfoD7znl1AiMg74ATBBKXUg7ng3ESmI/d8vJtenPsrllHeBpleMC4FqpdSGlgN+ppeTfsCPMubHaGkao6kXY42grgZ+HKAcI7G6JkuAxbHPxcCTwEex4zOAHj7L1Q9r9PhDYFlLGgFdgDeAlbG/nQNIszJgB9Ax7lgg6YXVWGwCGrCsmf9MlkbAj2NlbgVwkc9yrcLyX7aUs0di514ey+MPgQ+A8T7L5Zh3QaZX7PjjwH8lnOtnejnpB8/LmJmZaDAYDJqjg+vDYDAYDEkwitpgMBg0xyhqg8Fg0ByjqA0Gg0FzjKI2GAwGzTGK2mAwGDTHKGqDwWDQHKOoDQaDQXP+P4zewyNpLdWUAAAAAElFTkSuQmCC",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"\"\"\" Kálmán-szűrő: rekurzív átlag becslés\n",
"\"\"\"\n",
"\n",
"npt = 200\n",
"x = 5 + np.random.randn(npt)\n",
"n = np.arange(1,npt+1)\n",
"K = 1.0/(n+1)\n",
"mu = np.zeros((npt))\n",
"mu[0] = x[0]\n",
"\n",
"for i in range(1,npt):\n",
" mu[i] = mu[i-1] + K[i-1]*(x[i]-mu[i-1])\n",
"\n",
"plt.plot(n,x,'o',label='mérések')\n",
"plt.plot(n,mu,label='átlag becslés')\n",
"plt.legend()\n",
"plt.title('Átlag rekurzív becslése')\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"A példa jól bemutatta a Kálmán-szűrő egyik nagyon fontos jellemzőjét: **rekurzív becslést** szolgáltat a mérések függvényében a rendszer állapotára.\n",
"\n",
"Most vegyük szemügyre a Kálmán-szűrés másik fontos sajátosságát: a szűrő a különböző információkat *optimális* módon hasznosítja a becslés során."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Különböző információk optimális felhasználása\n",
"\n",
"Tegyük fel, hogy $x$-et *két különböző eszközzel* mérjük (pl. ha $x$ távolság, akkor megmérjük mérőszalaggal és kézi lézer távmérővel is). Legyen az első eszközzel végzett mérés eredménye $x_1$, a másiké pedig $x_2$. A két mérőeszköz nem azonos pontossággal dolgozik. Ezt a tényt a méréseik különböző $\\sigma_1$ illetve $\\sigma_2$ szórásaival (középhibáival) fejezhetjük ki.\n",
"\n",
"Célunk az, hogy *optimális* becslést adjunk $x$-re a két $x_1$, $x_2$ mérés és $\\sigma_1$, $\\sigma_2$ szórásaik segítségével. Legyen ez a becslés a két mérés *súlyozott átlaga*:\n",
"\n",
"$$\\hat{x} = wx_1 + (1-w)x_2.$$\n",
"\n",
"A kérdés tehát az, hogy milyen $w$ súly esetén lesz az $\\hat{x}$ becslés $\\hat\\sigma$ szórása minimális?\n",
"\n",
"Mivel feltételezzük azt, hogy mind $x_1$, mind $x_2$ *Gauss-eloszlású*, ezért a valószínűségszámítás [nevezetes eredménye](https://www.statlect.com/probability-distributions/normal-distribution-linear-combinations) szerint a lineáris kombinációjuk is normális eloszlású lesz, amelynek a szórása\n",
"\n",
"$$\\hat\\sigma^2 = w^2\\sigma_1^2 + (1-w)^2 \\sigma_2^2.$$\n",
"\n",
"Az optimális $w_{opt}$ súly ennek a kifejezésnek a minimuma lesz, vagyis ahol a $w$ szerinti derivált zérus:\n",
"\n",
"$$\\frac{\\partial \\hat\\sigma^2}{\\partial w} = 2w_{opt}\\sigma_1^2-2(1-w_{opt})\\sigma_2^2 = 0.$$\n",
"\n",
"Az optimális súly tehát\n",
"\n",
"$$ w_{opt} = \\frac{\\sigma_2^2}{\\sigma_1^2+\\sigma_2^2}.$$\n",
"\n",
"Ezek után megadhatjuk $x$ optimális $\\hat{x}$ becslését:\n",
"\n",
"$$ \\hat{x} = \\frac{\\sigma_2^2}{\\sigma_1^2+\\sigma_2^2} x_1 + \\frac{\\sigma_1^2}{\\sigma_1^2+\\sigma_2^2} x_2.$$\n",
"\n",
"Az optimális becslés szórását is kiszámíthatjuk az alábbi összefüggéssel:\n",
"\n",
"$$ \\hat\\sigma^2 = \\frac{\\sigma_1^2\\sigma_2^2}{\\sigma_1^2+\\sigma_2^2}.$$\n",
"\n",
"Számunkra, földmérők számára egyébként ez a jól ismert súlyozott legkisebb négyzetes becslés, hiszen $\\hat{x}$ minimalizálja a javítások súlyozott négyzetösszegét:\n",
"\n",
"$$\\hat{x} = \\underset{x}{\\arg \\min} = \\left[ p_1 (x_1 - x)^2 + p_2 (x_2 - x)^2 \\right],$$\n",
"\n",
"ahol $p_1$, $p_2$ a súlyok, melyek arányosak a szórásnégyzetek reciprokaival. Valóban, a fenti kifejezést $x$ szerint deriválva és a deriváltat zérussal egyenlővé téve a \n",
"\n",
"$$ p_1(x_1-\\hat{x}) + p_2(x_2-\\hat{x}) = 0$$\n",
"\n",
"egyenletet kapjuk, melynek megoldása $\\hat{x}$-re\n",
"\n",
"$$ \\hat{x} = \\frac{1}{p_1+p_2}(p_1 x_1 + p_2 x_2)$$\n",
"\n",
"lesz, ami megegyezik a korábbi megoldással, ha $p_1=1/\\sigma_1^2$ és $p_2=1/\\sigma_2^2$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Rekurzív optimális becslés\n",
"\n",
"Az előző optimális becslést fogalmazzuk át úgy, hogy rendelkezzünk először az $x_1$ méréssel, és ezt \"frissítsük\" az új $x_2$ mérés beérkezése után egy korrekcióval. Nyilván az első lépésben $\\hat{x_1} = x_1$, $\\hat\\sigma_1^2=\\sigma_1^2$, ahol az alsó index a lépésre utal. A második lépésben frissítsük a becslést\n",
"\n",
"$$\\hat{x}_2 = \\hat{x}_1 + \\frac{\\hat{\\sigma}_1^2}{\\hat{\\sigma}_1^2 + \\sigma_2^2} (x_2 - \\hat{x}_1),$$\n",
"\n",
"$$\\hat{\\sigma}_2^2 = \\left( 1 - \\frac{\\hat{\\sigma}_1^2}{\\hat{\\sigma}_1^2 + \\sigma_2^2} \\right) \\hat{\\sigma}_1^2.$$\n",
"\n",
"A zárójelben található $(x_2 - \\hat{x}_1)$ tényező az *új információ* (innováció), az előtte található szorzótényező, \n",
"\n",
"$$\\frac{\\hat{\\sigma}_1^2}{\\hat{\\sigma}_1^2 + \\sigma_2^2}$$\n",
"\n",
"pedig a $K$ *erősítési tényező*. Foglaljuk most össze a rekurzív optimális becslésre vonatkozó egyenleteinket:\n",
"\n",
"$$\\hat{x}_2 = \\hat{x}_1 + K(x_2 - \\hat{x}_1),$$\n",
"\n",
"$$\\hat{\\sigma}_2^2 = ( 1 - K ) \\hat{\\sigma}_1^2,$$\n",
"\n",
"$$K = \\frac{\\hat{\\sigma}_1^2}{\\hat{\\sigma}_1^2 + \\sigma_2^2}.$$\n",
"\n",
"Láthatjuk, hogy ez alakjára nézve teljesen hasonlít a korábbi rekurzív becslésre, amely így nézett ki\n",
"\n",
"$$\\mu_{n+1} = \\mu_n + K(x_{n+1} - \\mu_n),$$\n",
"\n",
"csak most sikerült megtalálni a $K$ erősítési tényező *optimális* értékét annak érdekében, hogy a becslés hibája (szórása) a lehető legkisebb legyen. Ezt az optimális becslést folyamatosan elő tudjuk állítani mindig amikor egy új mérés érkezik.\n",
"\n",
"Már most megjegyezzük, hogy amikor az állapotot nem egyetlen $x$ változóval, hanem több változóból álló $\\mathbf{x} = (x_1, x_2, ...)$ vektorral írhatjuk le, akkor is hasonló az eljárás, csak ez esetben varianciák helyett kovariancia mátrixokkal fogunk dolgozni. Ez többek között azt vonja majd maga után, hogy reciprok képzés helyett inverz képzésre lesz szükség."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Közvetett állapot becslés\n",
"\n",
"Idáig feltételeztük azt, hogy a rendszer $x$ állapota *közvetlenül* mérhető. A gyakorlatban sokszor erre nincs lehetőségünk, csak arra, hogy a rendszer állapotától függő valamilyen $\\mathbf{z}$ mennyiségeket megmérjünk. A geodéziában ezt a helyzetet igen jól ismerjük: ez a közvetett mérések feldolgozásának (II. kiegyenlítési csoport) esete. Ezt a feladatot oldja meg a Kálmán-szűrő is: a közvetlenül nem mérhető (vagyis \"rejtett\") rendszer állapotot szeretnénk becsülni a méréseink alapján. \n",
"\n",
"Először a *statikus* rendszerekkel foglalkozunk, vagyis olyanokkal, amelyekben a rendszer állapota időben nem változik. De mind a statikus, mind a dinamikus (időben változó állapotú) rendszerek esetében feltételezünk egy nagyon fontos dolgot: a rendszer $\\mathbf{x}$ állapota **lineáris** modellen keresztül kapcsolódik a mérhető $\\mathbf{z}$ mennyiségekhez:\n",
"\n",
"$$ \\mathbf{z} = \\mathbf{H} \\mathbf{x} + \\mathbf{r}.$$\n",
"\n",
"Ebben az összefüggésben a $\\mathbf{H}$ mátrix a mérési modell (számunkra ez a jól ismert alakmátrix), továbbá $\\mathbf{r}$ a zérus átlagú ($\\bar{\\mathbf{r}}=\\mathbf{0}$) és Gauss-eloszlású mérési hibákat (mérési zajt) jellemzi. Ezeket a hibákat (zajt) az $\\mathbf{R}$ *zaj kovariancia mátrix* írja le\n",
"\n",
"$$ \\mathbf{R} = \\mathbb{E}\\left[ (\\mathbf{r} - \\bar{\\mathbf{r}}) (\\mathbf{r} - \\bar{\\mathbf{r}})^T \\right]$$\n",
"\n",
"ahol $\\mathbb{E}[.]$ a várható értéket jelöli. Célunk optimális $\\hat{\\mathbf{x}}$ becslést adni a rendszer állapotára, amely minimalizálja a kovariancia mátrix inverzével súlyozott mérési zajt (mérési hibákat). Ez az igen jól ismert paraméteres kiegyenlítés (a paraméterek vektora most a rendszer $\\mathbf{x}$ állapota), ahol a $\\mathbf{P}$ súlymátrix arányos a zaj kovariancia mátrix $\\mathbf{R}^{-1}$ inverzével."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Rekurzív becslés**\n",
"\n",
"A hagyományos paraméteres kiegyenlítéssel szemben most tételezzük fel azt, hogy a mérések folyamatosan (nem egyszerre) érkeznek. Minden egyes $k$ mérési időpontban érkezik egy új $\\mathbf{z}_k$ mérés és szeretnénk a rendszer $\\hat{\\mathbf{x}}_k$ állapotára, illetve a rendszer állapot $\\hat{\\mathbf{P}}_k$ kovariancia mátrixára optimális becslést adni Ezt a minimalizációs feladatot kell tehát megoldani:\n",
"\n",
"$$\\hat{\\mathbf{x}}_k = \\underset{x}{\\arg \\min} = \\left[ (\\mathbf{z}_k-\\mathbf{H}_k \\mathbf{x})^T \\mathbf{R}_k^{-1} (\\mathbf{z}_k-\\mathbf{H}_k \\mathbf{x}) \\right].$$\n",
"\n",
"Ennek a minimalizációs feladatnak a megoldása a *statikus Kálmán-szűrő*. A szűrő kiszámítja az új $\\hat{\\mathbf{x}}_k$ becslést és annak $\\hat{\\mathbf{P}}_k$ kovariancia mátrixát az előző $\\hat{\\mathbf{x}}_{k-1}$, $\\hat{\\mathbf{P}}_{k-1}$ becslésekből és a $\\mathbf{z}_k$ mérésből valamint annak $\\mathbf{R}_k$ kovariancia mátrixából. \n",
"\n",
"Megemlítjük, hogy tipikus példa a statikus Kálmán-szűrésre a fix helyzetű GNSS vevő méréseinek a feldolgozása. Ez esetben a rejtett állapot ($\\mathbf{x}$ paraméter vektor) a vevő ismeretlen koordinátája és a vevő órájának igazítatlansága. A $\\mathbf{z}$ mérések pedig a műhold-vevő áltávolságok. Ha az összes mérés egyszerre rendelkezésre áll, akkor ez a hagyományos legkisebb négyzetes kiegyenlítéssel megoldható. Amennyiben a mérések epochánként folyamatosan érkeznek, akkor rekurzív kiegyenlítésre, vagyis statikus Kálmán-szűrésre van szükség.\n",
"\n",
"A következőkben levezetjük a statikus Kálmán-szűrő egyenleteit a geodéziában ismert *csoportokban történő kiegyenlítés* segítségével. Látni fogjuk, hogy a statikus Kálmán szűrő összefüggései igazából teljesen megegyeznek a közvetítő egyenletekkel törénő kiegyenlítés során alkalmazott folyamatos csoportképzés egyenleteivel.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Statikus Kálmán-szűrő egyenletei\n",
"\n",
"A levezetésben követni fogjuk Detrekői: Kiegyenlítő számítások tankönyvének 5.8.2-es pontjában alkalmazott gondolatmenenetet. A jelöléseket viszont a Kálmán-szűrő esetében szokásos jelöléshez közelítjük, így könnyebben látjuk majd a kapcsolatot a két eljárás között. Csak a lineáris esetre korlátozódunk, tehát nincs szükség a közvetítő egyenleteknek a paraméterek előzetes értékei helyén végzett linearizációjára.\n",
"\n",
"A folyamatos csoportképzés esetében a méréseket különböző időben végezzük, s a később végzett mérések kiegyenlítésekor felhasználjuk a korábban végzett mérések eredményeit is. Abból indulunk ki, hogy az $i=1$ alkalommal a $\\mathbf{z}_1$ mérési eredmények ($n$ db.) alapján meghatároztuk a kiegyenlített $\\mathbf{x}_1$ paramétereket és azok $\\mathbf{P}_{11}$ kovariancia mátrixát a mérések $\\mathbf{R}_{11}$ kovariancia mátrixát felhasználva, a \n",
"\n",
"$$ \\mathbf{r}_1 = \\mathbf{z}_1 - \\mathbf{H}_1\\mathbf{x}_1 $$\n",
"\n",
"lineáris javítási (reziduál) egyenletekből kiindulva az\n",
"\n",
"$$ \\mathbf{N}_1 \\mathbf{x}_1 - \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 = \\mathbf{0}$$\n",
"\n",
"normálegyenletek alapján. Ha az $ \\mathbf{N}_1= \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{H}_1$ normálegyenlet együtthatómátrixa nem szinguláris, akkor az ismeretlen paraméterek értéke\n",
"\n",
"$$\\mathbf{x}_1 = \\mathbf{N}_1^{-1}\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 .$$\n",
"\n",
"A kiegyenlített paraméterek kovariancia mátrixa\n",
"\n",
"$$ \\mathbf{P}_{11} = \\mathbf{N}_1^{-1} = (\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{H}_1)^{-1}.$$\n",
"\n",
"A kiegyenlítés elvégzése után az $i+1=2$ alkalommal további ($s$ db.) mérést végzünk, ezek eredménye $\\mathbf{z}_2$, kovariancia mátrixa $\\mathbf{R}_{22}$. Tételezzük fel, hogy ezek a mérések függetlenek a korábbi mérésektől. Határozzuk meg az új $\\mathbf{x}_2$ paramétervektort, a javítások (reziduálok) új $\\mathbf{r}_r$ vektorát és a kiegyenlített paraméterek új $\\mathbf{P}_{22}$ kovariancia mátrixát."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Az együttes feldolgozáskor a normálegyenlet\n",
"\n",
"$$ (\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{H}_1 + \\mathbf{H}_2^T\\mathbf{R}_{22}^{-1}\\mathbf{H}_2)\\mathbf{x}_2 - ( \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 + \\mathbf{H}_2^T\\mathbf{R}_{22}^{-1}\\mathbf{z}_2) = \\mathbf{0}.$$\n",
"\n",
"A megoldás érdekében vezessük be az új $\\mathbf{y}$ ($s$ méretű) paramétert:\n",
"\n",
"$$ \\mathbf{y}= \\mathbf{R}_{22}^{-1}\\mathbf{r}_2 = \\mathbf{R}_{22}^{-1}\\mathbf{z}_2 -\\mathbf{R}_{22}^{-1}\\mathbf{H}_2 \\mathbf{x}_2 .$$\n",
"\n",
"Ekkor az $\\mathbf{y}$ segítségével a\n",
"\n",
"$$ \\mathbf{r}_2 = \\mathbf{z}_2 - \\mathbf{H}_2\\mathbf{x}_2 $$\n",
"\n",
"javítási (reziduál) egyenlet átalakítható csak paramétereket tartalmazó kényszerfeltételi egyenletté:\n",
"\n",
"$$ \\mathbf{R}_{22}\\mathbf{y}= \\mathbf{z}_2 - \\mathbf{H}_2\\mathbf{x}_2,$$\n",
"\n",
"amiből a szokásos felépítésű kényszerfeltételi egyenlethez jutunk:\n",
"\n",
"$$ \\mathbf{z}_2 - \\mathbf{H}_2\\mathbf{x}_2 - \\mathbf{R}_{22}\\mathbf{y} = \\mathbf{0}$$\n",
"\n",
"A bevezetett $\\mathbf{y}$ vektort felhasználva a normálegyenlet átalakítható:\n",
"\n",
"$$(\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{H}_1)\\mathbf{x}_2 - \\mathbf{H}_2^T\\mathbf{y} - \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 = \\mathbf{0}$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"A kényszerfeltételi egyenleteket hozzávéve ehhez az átalakított normálegyenlethez a következő $n+s$ méretű blokkos egyenletrendszerhez jutunk:\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{H}_1 & - \\mathbf{H}_2^T \\\\ - \\mathbf{H}_2 & - \\mathbf{R}_{22}\\end{bmatrix} \\begin{bmatrix} \\mathbf{x}_2\\\\ \\mathbf{y}\\end{bmatrix} = \\begin{bmatrix} \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1\\\\ -\\mathbf{z}_2\\end{bmatrix}. $$\n",
"\n",
"Az $\\mathbf{x}_2$-re vonatkozó megoldás meghatározása érdekében invertáljuk az egyenletrendszer $\\mathbf{A}$ blokkos együttható mátrixát. Ehhez a Gauss-Jordan kiküszöbölést fogjuk alkalmazni. Az inverz meghatározása során kiindulunk a kibővített\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{A} & | & \\mathbf{I} \\end{bmatrix} $$\n",
"\n",
"mátrixból és a Gauss-Jordan eliminációval ebből kihozzuk a\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{I} & | & \\mathbf{A}^{-1} \\end{bmatrix} $$\n",
"\n",
"alakot ($\\mathbf{I}$ az egységmátrixot jelöli).\n",
"\n",
"A kiinduló kibővített mátrix (az egyszerűsítés érdekében $\\mathbf{N}=\\mathbf{N}_1$)\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{N} & - \\mathbf{H}_2^T & \\mathbf{I} & \\mathbf{0} \\\\ - \\mathbf{H}_2 & - \\mathbf{R}_{22} & \\mathbf{0} & \\mathbf{I}\\end{bmatrix} $$\n",
"\n",
"első sorát szorozzuk végig balról $\\mathbf{N}^{-1}$-el, majd adjuk hozzá a második sorhoz az $\\mathbf{H}_2$-vel balról szorzott első sort. Mivel $\\mathbf{H}_2\\mathbf{I}=\\mathbf{H}_2$, ha az $\\mathbf{I}$ $r\\times r$ méretű egységmátrix, eredményül a\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{I} & - \\mathbf{N}^{-1}\\mathbf{H}_2^T & \\mathbf{N}^{-1} & \\mathbf{0} \\\\ \\mathbf{0} & - \\mathbf{R}_{22}- \\mathbf{H}_2\\mathbf{N}^{-1}\\mathbf{H}_2^T & \\mathbf{H}_2\\mathbf{N}^{-1} & \\mathbf{I}\\end{bmatrix} $$\n",
"\n",
"kibővített mátrixot kapjuk. Vezessük be az $\\mathbf{S} = (\\mathbf{R}_{22}+ \\mathbf{H}_2\\mathbf{N}^{-1}\\mathbf{H}_2^T)^{-1}$ jelölést. Balról a második sort $-\\mathbf{S}$-el végigszorozva a\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{I} & - \\mathbf{N}^{-1}\\mathbf{H}_2^T & \\mathbf{N}^{-1} & \\mathbf{0} \\\\ \\mathbf{0} & \\mathbf{I} & -\\mathbf{S}\\mathbf{H}_2\\mathbf{N}^{-1} & -\\mathbf{S}\\end{bmatrix} $$\n",
"\n",
"mátrixot kapjuk meg. Végül adjuk hozzá az első sorhoz a balról $\\mathbf{N}^{-1}\\mathbf{H}_2^T$-vel szorzott második sort. Eredményül a\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{I} & \\mathbf{0} & \\mathbf{N}^{-1} - \\mathbf{N}^{-1}\\mathbf{H}_2^T\\mathbf{S}\\mathbf{H}_2\\mathbf{N}^{-1} & -\\mathbf{N}^{-1}\\mathbf{H}_2^T\\mathbf{S} \\\\ \\mathbf{0} & \\mathbf{I} & -\\mathbf{S}\\mathbf{H}_2\\mathbf{N}^{-1} & -\\mathbf{S}\\end{bmatrix} $$\n",
"\n",
"kibővített mátrixot kapjuk, amelynek jobb oldalán a keresett inverz mátrix áll.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Ezzel az inverz mátrixszal balról végigszorozva a megoldandó egyenletrendszert megkapjuk a keresett megoldást\n",
"\n",
"$$\\begin{bmatrix} \\mathbf{x}_2\\\\ \\mathbf{y}\\end{bmatrix} = \\begin{bmatrix} \\mathbf{N}^{-1} - \\mathbf{N}^{-1}\\mathbf{H}_2^T\\mathbf{S}\\mathbf{H}_2\\mathbf{N}^{-1} & -\\mathbf{N}^{-1}\\mathbf{H}_2^T\\mathbf{S} \\\\ -\\mathbf{S}\\mathbf{H}_2\\mathbf{N}^{-1} & -\\mathbf{S}\\end{bmatrix} \\begin{bmatrix} \\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1\\\\ -\\mathbf{z}_2\\end{bmatrix}. $$\n",
"\n",
"Ha a továbbiakban a $\\mathbf{K} =\\mathbf{N}^{-1}\\mathbf{H}_2^T\\mathbf{S} $ jelölést alkalmazzuk, akkor az $\\mathbf{x}_2$-re vonatkozó megoldás\n",
"\n",
"$$ \\mathbf{x}_2 = \\mathbf{N}^{-1}\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 - \\mathbf{K}\\mathbf{H}_2\\mathbf{N}^{-1}\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 + \\mathbf{K}\\mathbf{z}_2.$$\n",
"\n",
"Ezt egyszerűsíteni tudjuk, ha figyelembe vesszük, hogy $\\mathbf{x}_1 = \\mathbf{N}_1^{-1}\\mathbf{H}_1^T\\mathbf{R}_{11}^{-1}\\mathbf{z}_1 $. Végeredményben\n",
"\n",
"$$ \\mathbf{x}_2 = \\mathbf{x}_1 - \\mathbf{K}\\mathbf{H}_2\\mathbf{x}_1 + \\mathbf{K}\\mathbf{z}_2 = \\mathbf{x}_1 + \\mathbf{K}(\\mathbf{z}_2 - \\mathbf{H}_2\\mathbf{x}_1 ).$$\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Az $\\mathbf{S}$ mátrixot visszaírva (mivel $\\mathbf{N}^{-1}=\\mathbf{P}_{11}$)\n",
"\n",
"$$ \\mathbf{K} = \\mathbf{P}_{11}\\mathbf{H}_2^T(\\mathbf{R}_{22} + \\mathbf{H}_2 \\mathbf{P}_{11} \\mathbf{H}_2^T)^{-1}.$$\n",
"\n",
"A kiegyenlített paraméterek kovariancia mátrixa pedig az invertált hipermátrix bal felső blokkja (vagyis az $\\mathbf{x}_2$-höz tartozó blokk)\n",
"\n",
"$$ \\mathbf{P}_{22} = (\\mathbf{I} - \\mathbf{K}\\mathbf{H}_2) \\mathbf{P}_{11}.$$\n",
"\n",
"Ezek megegyeznek a statikus Kálmán-szűrő egyenleteivel. A $\\mathbf{K}$ a Kálmán-féle erősítési mátrix, $\\mathbf{H}$ a mérési mátrix, $\\mathbf{x}$ pedig a rendszer állapota. A szűrő minden új mérés után frissíti a rendszer állapot becslését és az állapot becslésének $ \\mathbf{P}$ kovariancia mátrixát.\n",
"\n",
"A *dinamikus* Kálmán-szűrő egyenletei mindössze annyiban térnek el ettől, hogy a rendszer $ \\mathbf{x}$ állapota is változik:\n",
"\n",
"$$\\mathbf{x}_{k+1} = \\mathbf{F}\\mathbf{x}_k + \\mathbf{B}\\mathbf{u}_k + \\mathbf{q}$$\n",
"\n",
"a $k$-adik időpontról a $k+1$-edik időpontra, ahol $\\mathbf{F}$ az állapot átmenet mátrixa, $\\mathbf{u}$ a rendszer gerjesztése, $\\mathbf{B}$ a gerjesztés (vezérlés) mátrixa és $\\mathbf{q}$ a zérus átlagú rendszer zaj, aminek $\\mathbf{Q}$ a kovariancia mátrixa.\n",
"\n",
"A dinamikus rendszer állapot-átmenete utáni $\\mathbf{P}^-$ kovariancia mátrixot egyszerűen a jól ismert kovariancia terjedéssel (hibaterjedés) számíthatjuk ki:\n",
"\n",
"$$\\mathbf{P}^-_k = \\mathbf{H} \\mathbf{P}_{k-1} \\mathbf{H}^T + \\mathbf{Q}.$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Érdemes megjegyezni, hogy a Kálmán-szűrő egyenleteit valószínűségszámítási megfontolások alapján közvetlenül is le lehet vezetni (lásd [Pradu](https://math.stackexchange.com/questions/840662/an-explanation-of-the-kalman-filter) válaszát)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Példa a Kálmán-szűrésre: járműnavigáció\n",
"\n",
"Példaként vegyünk egy egyenes úton haladó járművet. A modellünkben a rendszer (jármű) $x$\n",
"állapotát a $p$ helyzete és $v$ sebessége jellemzi. Az $u$ bemenetet az általunk alkalmazott ismert gyorsulások jelentik és az $y$ kimenet a mért helyzetinformációval azonosítható. Mondjuk azt, hogy gyorsulásokkal vezéreljük a rendszert és a helyzetét $T$ másodperces időközönként mérjük. Ez esetben az elemi kinematikai ismereteink szerint a $v$ sebesség a következő egyenlettel írható le:\n",
"\n",
"$$v_{k+1} = v_k + Tu_k.$$\n",
"\n",
"Viszont ez a modell még nem teljes. Ugyanis a sebesség véletlenszerűen is megváltozhat elkerülhetetlen külső tényezők miatt. A sebesség $\\tilde{v}_k$ zaja időfüggő valószínűségi változó, ezért \n",
"\n",
"$$v_{k+1} = v_k + Tu_k + \\tilde{v}_k$$\n",
"\n",
"a sebességre vonatkozó modellünk. Ehhez hasonló egyenletet írhatunk fel a $p$ helyzetre:\n",
"\n",
"$$p_{k+1} = p_k + Tv_k + \\frac{1}{2}T^2u_k + \\tilde{p}_k,$$\n",
"\n",
"ahol $\\tilde{p}_k$ a helyzet bizonytalansága (zaja). Az állapotvektor mint említettük a helyzetből és a sebességből áll:\n",
"\n",
"$$x_k = \\left[ p_k \\atop v_k \\right].$$\n",
"\n",
"Mivel tudjuk, hogy a mért $z$ mennyiség a jármű helyzete, az alábbi lineáris rendszeregyenleteket írhatjuk fel:\n",
"\n",
"$$x_{k+1} = \\begin{bmatrix}1 & T\\\\0 & 1 \\end{bmatrix} x_k + \\begin{bmatrix}T^2/2\\\\T \\end{bmatrix} u_k + w_k$$\n",
"\n",
"$$y_k =\\begin{bmatrix}1 & 0\\end{bmatrix} x_k + \\nu_k.$$\n",
"\n",
"A mérés $\\nu_k$ zaja a mérőeszköz hibájából adódik. Célunk pedig a jármű $x$ állapotának minél pontosabb megismerése.\n",
"\n",
"A Kálmán-szűrés esetében viszonylag könnyű a mérési zaj, esetünkben a jármű pozíciója hibájának ($R$ \"mátrix\", mely egyetlen elemből álló skalár) a meghatározása, mivel a felhasznált mérőeszköz pontosságát általában ismerjük. Legye a példánkban a helyzet meghatározás szórása ± 3 méter. Nehezebb feladat a rendszer zaj $Q$ kovariancia mátrixának felvétele. Erre az ún. diszkrét fehérzaj gyorsulás modellt (DWNA) alkalmazzuk. Ha ismerjük a gyorsulás $w_{acc}$ mérési zaját (legyen ez a példánkban ± 0.2 m/$\\mathrm{s}^2$) akkor a rendszer zaj kovariancia mátrixa\n",
"\n",
"$$Q = \\begin{bmatrix}¼T^4 & ½T^3\\\\½T^3 & T^2 \\end{bmatrix} w_{acc}^2.$$\n",
"\n",
"A jármű helyzetének és sebességének becslését a `kalman(duration,dt)` függvény végzi el Kálmán-szűréssel, és a jármű pontos, mért valamint becsült helyzetét egy ábrán szemléltetjük a szimuláció lefutása után. Az ábráról jól látszik, hogy miért is nevezhető ez az eljárás szűrésnek."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAe4AAAGDCAYAAADtffPSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAAsTAAALEwEAmpwYAAB8qUlEQVR4nOzdd3gVVf7H8fe5Pb1XEkjovUlRQcRGsfeuKPay6upvrbvq2l1ddde6dsUu2AsqCIJSpNfQQgLpvd2UW8/vj3sTCCQhQBII+b6eJw+5M3NnzkxCPvecOXOO0lojhBBCiM7BcKgLIIQQQojWk+AWQgghOhEJbiGEEKITkeAWQgghOhEJbiGEEKITkeAWQgghOhEJbiH2QSmVqZQ6+TAoR5xSaoFSqkop9W+l1P1KqTf961KUUlopZTrU5WytPc/nUJdHiM6i0/wnF+JgKKUygWu11nP8ry8GXgXO1lr/dijLth+uB4qBUL2PARiUUvOBD7TWb3ZEwQ5Qq89HCLGL1LhFl6OUmga8DJzWiUIboAewsSNCTillbMd9K6WUgYM4n87UsiBEW5PgFl2KUup64N/AZK31Iv+yXkqpX5VSJUqpYqXUh0qp8Gbe/7BS6nOl1Af+Jt51Sqm+Sqn7lFKFSqkspdSk3ba/WimV5t92u1Lqht3WTVRKZSul7vK/N08pdXUzx30XmAbcrZSyK6VO9pflgya2fRw4DnjJv+1L/uX9lVK/KKVKlVKblVIX7r5/pdSrSqkflFLVwAlN7He+UupJpdSfSqkKpdTXSqnI3dYfrZRapJQqV0qtUUpN3OO9jyul/gBqgPebOB+rUuoFpVSu/+sFpZR1j2t1j1IqH3inI38WSqkA/+2JHf5z/10pFbCv8xaiXWit5Uu+jvgvIBOYBRQAw/ZY1xs4BbACMcAC4IU93nuy//uHgTpgMr5bTe8DGcADgBm4DsjY7b2nAb0ABRyPL7RG+tdNBNzAI/73nupfH9HMObwLPLbb64fxNYcDpAAaMPlfz8d3a6B+2yAgC7jaX+6R+JqpB+227wpgHL4P9LYmjj8fyAEG+/c3a7fjdwNK/Odg8F/PEiBmt/fuBAb5j29u4nweAZYAsf6fwyLg0T2u1dP+n1NAR/4s8LXQzPefpxE41l+OFs9bvuSrPb6kxi26klPwBcO63RdqrbdprX/RWju01kXAc/j+sDdnodb6J621G/gcX8g8pbV2AZ8AKfU1dq3191rrdO3zG/AzvtpwPRfwiNbapbX+AbAD/drkbBs7HcjUWr+jtXZrrVfiC97zd9vma631H1prr9a6rpn9zNBar9daVwP/AC70N6tfDvygtf7B//5fgOX4Aq3eu1rrDf7ju5rY92X4rkWh/+fwT+CK3dZ7gYf8P6da/7J2/1n4m/WnA7drrXO01h6t9SKttaOV5y1Em5LgFl3JjUBf4E2llKpfqJSKVUp9opTKUUpVAh8A0S3sp2C372uBYq21Z7fXAMH+fU9VSi3xN0+X4/uDvvu+S/yhU6+m/r1trAcw1t+cW+4vy2VA/G7bZLViP7tvswNf7TTav/8L9tj/eCBhP/af6N/n7vtP3O11URMfKDriZxEN2ID0JsrcmvMWok1JcIuupBA4CV8t65Xdlj+Jr5l5qNY6FF8tSu399v3jvz87C3gWiNNahwM/tMW+W2HPDl9ZwG9a6/DdvoK11je18J6mJO/2fXd8tdRi//5n7LH/IK31U/ux/1x8Qbj7/nP3s3xNOsifRTG+JvleTaxrzXkL0aYkuEWXorXOBU4EpiilnvcvDsHXLFqulOoG/K2NDmfBdx+0CHArpaYCk1p+S5spAHru9vo7oK9S6gqllNn/NVopNWA/93u5UmqgUioQ3/3gmf4a7gfAGUqpyUopo1LK5u/wlbQf+/4Y+LtSKkYpFQ086N9vWzjgn4XW2gu8DTynlEr0n98x/g8DbXHeQuwXCW7R5Wits/CF9/lKqSfx3Usdia9z1vfAF210nCrgNuAzoAy4FPimLfbdCv/Bd35lSqn/+ssyCbgYXy02n10dvfbHDHydyvLxNR/fBg3X9CzgfnzhmIXvA9D+/I15DN/94bX4+iGs9C87aG3ws/g/f5mWAaX4rp2hjc5biP2itJZxD4QQ+6Y6x6AuQhzx5FOhEEII0YlIcAshhBCdiDSVCyGEEJ2I1LiFEEKITkSCWwghhOhEOsUMO9HR0TolJeVQF0MIIYToECtWrCjWWsc0ta5TBHdKSgrLly8/1MUQQgghOoRSakdz66SpXAghhOhEJLiFEEKITkSCWwghhOhEOsU97qa4XC6ys7Opq2tu2mCxv2w2G0lJSZjN5kNdFCGEEM3otMGdnZ1NSEgIKSkp7Da1sjhAWmtKSkrIzs4mNTX1UBdHCCFEMzptU3ldXR1RUVES2m1EKUVUVJS0YAghxGGu0wY3IKHdxuR6CiHE4a9TB7c4cMHBwYe6CEIIIQ6ABPcRxO12H+oiCCGEaGcS3AchMzOT/v37c+211zJ48GAuu+wy5syZw7hx4+jTpw9//vkn1dXVTJ8+ndGjRzNixAi+/vrrhvced9xxjBw5kpEjR7Jo0SIA8vLymDBhAsOHD2fw4MEsXLgQaFxDnjlzJldddRUAV111FXfeeScnnHAC99xzD+np6UyZMoWjjjqK4447jk2bNgGQkZHBMcccw+jRo/nHP/7RgVdJCCFEW+q0vcp3989vN7Axt7JN9zkwMZSHzhi0z+22bdvG559/zuuvv87o0aP56KOP+P333/nmm2944oknGDhwICeeeCJvv/025eXljBkzhpNPPpnY2Fh++eUXbDYbW7du5ZJLLmH58uV89NFHTJ48mQceeACPx0NNTc0+y7BlyxbmzJmD0WjkpJNO4rXXXqNPnz4sXbqUm2++mV9//ZXbb7+dm266iSuvvJKXX365LS6REEKIQ6Bdg1sp9VfgWkAD64CrgUDgUyAFyAQu1FqXtWc52lNqaipDhgwBYNCgQZx00kkopRgyZAiZmZlkZ2fzzTff8OyzzwK+3vA7d+4kMTGRW2+9ldWrV2M0GtmyZQsAo0ePZvr06bhcLs4++2yGDx++zzJccMEFGI1G7HY7ixYt4oILLmhY53A4APjjjz+YNWsWAFdccQX33HNPW14GIYTostZkleP2ejmqR2SHHK/dglsp1Q24DRiota5VSn0GXAwMBOZqrZ9SSt0L3AscVIq0pmbcXqxWa8P3BoOh4bXBYMDtdmM0Gpk1axb9+vVr9L6HH36YuLg41qxZg9frxWazATBhwgQWLFjA999/zxVXXMHf/vY3rrzyykY9vvd8ZCsoKAgAr9dLeHg4q1evbrKs0mtcCCHa3ou/biW3vI4fbj+uQ47X3ve4TUCAUsqEr6adC5wFvOdf/x5wdjuX4ZCaPHkyL774IlprAFatWgVARUUFCQkJGAwGZsyYgcfjAWDHjh3ExsZy3XXXcc0117By5UoA4uLiSEtLw+v18uWXXzZ5rNDQUFJTU/n8888B36Aqa9asAWDcuHF88sknAHz44Yftd8JCCNHFVNa6CQvouBEn2y24tdY5wLPATiAPqNBa/wzEaa3z/NvkAbFNvV8pdb1SarlSanlRUVF7FbPd/eMf/8DlcjF06FAGDx7c0DHs5ptv5r333uPoo49my5YtDbXm+fPnM3z4cEaMGMGsWbO4/fbbAXjqqac4/fTTOfHEE0lISGj2eB9++CFvvfUWw4YNY9CgQQ2d4f7zn//w8ssvM3r0aCoqKtr5rIUQouuoqHURGtBxXcZUfU2wzXesVAQwC7gIKAc+B2YCL2mtw3fbrkxrHdHSvkaNGqX3nI87LS2NAQMGtHGphVxXIYTYP8c8OZfxvaN55oJhbbZPpdQKrfWopta1Z1P5yUCG1rpIa+0CvgCOBQqUUgn+giUAhe1YBiGEEKJdVda6CD0SmsrxNZEfrZQKVL5eUScBacA3wDT/NtOAr9uxDEIIIUS7cXu8VDs9hNqOgODWWi/F1zS+Et+jYAbgdeAp4BSl1FbgFP9rIYQQos1kldYw4pGf2ZDbvn16Kut8I1aGdeA97nY9ktb6IeChPRY78NW+hRBCiHaxeHsJZTUuVu0sZ1BiWLsdp7LWBXDENJULIYQQh8S6bF9NO6t036NPHozKOn9wHwlN5UIIIcShsjbHH9xl+w7ur1blkFNee0DHqaz1NZVLjbsTKC8v55VXXjnUxRBCCLEHp9tLWp5v/oqs0pYDubTayR2fruajpTsO6FgV/qbyI2IAliNdc8FdPwKaEEKIQ2NLQRVOt5eoIAs799FUXh/wxVXOAzpWQ1N5B3ZOk+A+QPfeey/p6ekMHz6c0aNHc8IJJ3DppZc2TC4yePDghm2fffZZHn74YYBmp90UQgjRNtb5m8knD46notbVEK5NaQhuu+OAjtXQOc3Y/DHa2hExrSc/3gv569p2n/FDYGrzT6o99dRTrF+/ntWrVzN//nxOO+001q9fT2pqKpmZmc2+7/rrr29y2k0hhBBtY212OWEBZsb3juajpTvJKq1ptmf5xvrgrj6AGrfXS0zuXD6xvE3grFiY9s3BFLvVjozgPgyMGTOG1NTUFrdpadpNIYQQbWNtdgVDk8LoHhkI0GJwp+VVYQzcRlFNSusP4KyG1R/Bklc5tzSdPEM0qs9loDV0wCyMR0Zwt1Az7ij1k4QAmEwmvF5vw+v6aTj3Ne2mEEKIg1Pn8rA5v4rrJ/QkOaI+uHd1UKtxujn9xd954NQBHNcnhu218wns8RmVzmgKq48hNqjxvFcur4vt5duJD4onTJlg2Zvwx3+gpgS6HcU7CQ/xQcVQ5h57coed45ER3IdASEgIVVVVTa6Li4ujsLCQkpISgoOD+e6775gyZUqjaTcvuOACtNasXbuWYcPabmB6IYToyrYUVOH2aoZ0CyMs0EyozdTokbD0wmq2F1Xz1u8ZVOl0zHGzsLi74zDlc+3P1zEp5RRWFa6izFGG1+sl256Nw+NgkDWGGVk7MVcXQa+T4Pi7IXksv727DGf4Vzyx9E/uHXMvBtX+XcckuA9QVFQU48aNY/DgwQQEBBAXF9ewzmw28+CDDzJ27FhSU1Pp379/w7oPP/yQm266icceewyXy8XFF18swS2EEG2kvnadEu1rBU2ODGzUs3xHaTXGoC2srFnH5uUb0O4wzkx8mPdXLCbX9D5vrHuD/pH9SQ5OxqAUx1qiCdyxiP9RxDuxiVx/wofQfWzD/spqy6gM/JmS2uM7JLRBgvugfPTRR82uu+2227jtttv2Wp6amsrs2bPbs1hCCNFl5Vf6bk3Gh9oASI4IZFNRDl7tRaH4PP1NArt/jvZYcVT1w108mfHju/POwkKeP/ZLRiRHEGwOgu3zYd4TkP0nRSEDWBPTm1drN2Cr2sz6+d/SN6I/1wy9jHz1E17l5KZhN3XYOUpwCyGEOGLkV9RiNRkIDzSzIHsBW4wvURKVxpRZL9EzvCcrK//AWD2WvsZp/JlVydCkMOL8Ie8pryC4+EfS57xBL9cWdGg3Vgx7hAuW9uS4qGDCAv7JM8ufQWkLP+74mmJHNjW230g0HU3viN4ddo4S3EIIIY4Y+ZUO4sNsLMpdxC1zbyHUFIOj6CSS42pZkruEKPdkotQ5nH9Ud/7MWMvg+EC6ZX3L++bXGP/dRtAenN5k/u65mqKIC/htZSUaL0UVJj654BOqXdVc/PI2XOEz+TDtQ7RSjAy5qEPPUYJbCCHEEaOgoo64UCv/W/sa8UHx3D3oHa57bw039T+GoacEM+HphaT0DubU/uFkhS7kpu3fELg+h1QVw4qkKwkfcylTPyxifO9oft9cTGyIlYl9I1iSUUJ8UDxOt5di+yZ01Znccu4Q3vgti6QePTr0HCW4hRBCHDHyKmtJTSpgVeEq7htzH31jIwDYXFDF4G5hUJnDJZXfEPzy19zlLIOYo+Cs5zjtYzg3JpmjvBFAEX8/fQA5ZbX0iAripw35zN6QT63TQ7HdgdYABoYHX4SzZDmhto6NUgluIYQQB23Wimy6RwUyOiXykJVBa01BhYOApO+JtEVybp9zsRqtdA/W1K75Ctfmpfxu/RFjNtD/VDj6JugxDpQiOng+RXYHGcXVAKREBdE/PhSAtdm+e+D5lXUUVe0aNGvVzjKgY2cGAwluIYQQbeCx7zdid7h5/qLhnD40sWH5yp1lBFoMrK34CbPBzNm9z0a10+hipdVOPNat5LnWcPuQ27EVbIDlb/OLZybW3DqclnDe8kzluEvvY+DAIY3eGxVsocTuIMNooFt4ADazsWFdfJgvuPPKaynabUzz1VnlQMcHt0wy0gHeffddcnNzG15fe+21bNy48aD3m5mZ2eIjac256qqrmDlz5kEfXwgh6lU7PGgNf/l4FW8s2I7WmgVbirjwjZ+44vtreXTJozy46EEeX/o4bq+72f04PU5qXPueQ3t3Bf5HwLYWF2BL/JQYcwyXrpkNb5wI678ku9upXOJ8gH/2+5In3ZeRmNJvr31EB1spsTvZXlxNanRQo3UJYQEA5FXUkVvuO1ZcqJU19cFtkxr3Eefdd99l8ODBJCb6PoW++eabbbLf+uC+9NJL22R/QghxIJxuL06Pl1tP6M22QjuP/5DG79uKWZ5ZiDX5dWoN5dx91AMUO3J5Z/07pJWmccvwWxgVN4qt5VvJqMgg157L+uL1LMlbAsDtI2/nkv6XtDioicvr4uF5b/Dl1u84e8CxZBSuwWyq5MUdBQRig5MfhlHXYLQbWfzsfNauKSQswEx4oGWvfUUFW1iy3YGnso6zhndrtC4hbFdTeWFlHSE2E0OTwvllYwEAYR04pSdIcB+UDz74gP/+9784nU7Gjh3LK6+8wjXXXMPy5ctRSjF9+nSSk5NZvnw5l112GQEBASxevJipU6fy7LPPMmrUKIKDg7nllluYM2cOERERPPHEE9x9993s3LmTF154gTPPPJPMzEyuuOIKqqt9915eeukljj32WO69917S0tIYPnw406ZN47bbbuPee+9l/vz5OBwObrnlFm644Qa01vzlL3/h119/JTU1Fe3rWSGEEG2i1ukBIDLIwquXj+St3zN46sdNhCX8gdNaRM3Oq4kZdTxXHpVAn/A+vLDyBW745QYMyoBX75rXITEokTN6nkFOdQ5P/fkUv+78lRdOeIEQSwirs0r578LfOPUoTbZ9J+WOclYXria9Ip0gQxizMz7FpeC2Ejs9hlwLJ9wFQVEA9LBq4kNt5FfWMTQpuMlziAqyUlbjm5pzzxq3zWwkItBMbnkthVUOEsMC6BmzaxupcR+Ap/98mk2lbTuvdf/I/twz5p5m16elpfHpp5/yxx9/YDabufnmm3nsscfIyclh/fr1AJSXlxMeHs5LL73UENR7qq6uZuLEiTz99NOcc845/P3vf+eXX35h48aNTJs2jTPPPJPY2Fh++eUXbDYbW7du5ZJLLmH58uU89dRTPPvss3z33XcAvP7664SFhbFs2TIcDgfjxo1j0qRJrFq1is2bN7Nu3ToKCgoYOHAg06dPb9PrJYTouqqdvqbvIKsRpRTXHteToSlwy28PMzHhBOZtH8Tv24qZOiSBM3qdweSUyXy3/Tuyq7LpF9mPPuF9SAhOIMDka5LWWvPVtq94ZPEj3DjnRh48+kH+9vtD5Ls2sGwJmJSJUGsoccrK3flwee068szR/JQ4meeKJ3L1HeeDcVdNXSnFMb2i+HJVDj2igpo6BaJDrA3fp8bsvU18WAD5FXUUVNURH2ajV8yuDwDSOa2TmDt3LitWrGD06NEA1NbWMmXKFLZv385f/vIXTjvtNCZNmrTP/VgsFqZMmQLAkCFDsFqtmM1mhgwZ0jCvt8vl4tZbb2X16tUYjUa2bNnS5L5+/vln1q5d23D/uqKigq1bt7JgwQIuueQSjEYjiYmJnHjiiW1wBYQQXZXWmivf/pPLxnZnyuAEavw17gDLrkiZmfEaXu3h7jF/ozY7n9+3FTessxgtnNvn3Gb3r5TinD7nEGoN5f/m/x/nf3s+RmzU5Z9Bd8tAZh9TjXnJy1C+g23eRP5pvJk/LCcyPDiWwJAiTMa9m9eP7hnJl6tySIkKbPKY0UG7ms97Ru8d3IlhNvIq6iiorGNItzB6+cPdYjI06sjWEY6I4G6pZtxetNZMmzaNJ598stHyxx9/nJ9++omXX36Zzz77jLfffrvF/ZjN5oYelgaDAavV2vC92+37FPv8888TFxfHmjVr8Hq92Gy2Zsv04osvMnny5EbLf/jhh3brxSmE6HpqXR4Wbi2mX1yIP7j9NW6LL8CW5y/nx4wfuWHoDSSHJDOut5s5aYVkldaQHNl0cDblpO4n8cIJL/Bj5o/kbRrJ0Ko5TDPeh3l2JXQbxbuh1/Pcjl5cM643W+dswWarahijfE/j+8RgMRmanZc7Ktj3t9dsVHQLD9hrfXyYjaUZpdgdbhLCAugZ7atxd3QzOUiv8gN20kknMXPmTAoLCwEoLS1lx44deL1ezjvvPB599FFWrlwJtDwFaGtUVFSQkJCAwWBgxowZeDyeJvc7efJkXn31VVwu332aLVu2UF1dzYQJE/jkk0/weDzk5eUxb968Ay6LEELU3wuubyKvdvj+JgVaTLi9bp748wkSghK4Zsg1AIzvHQ3AH7vVulvr+LC+PFVr4pWc27jL9CnbjL14IOxpii/+nv/m9GNi/3jG9vQ9O74up6Lh0a09dQsPYOl9JzF5UFyT66ODfTXu7pGBTdbYE8Js2B3uhu8jgixEBlkI7eCOaXCE1LgPhYEDB/LYY48xadIkvF4vZrOZ5557jnPOOQev19fZor42ftVVV3HjjTc2dE7bXzfffDPnnXcen3/+OSeccAJBQb4mmqFDh2IymRg2bBhXXXUVt99+O5mZmYwcORKtNTExMXz11Vecc845/PrrrwwZMoS+ffty/PHHt92FEEJ0OWXVTgCq6nxBVrPbPe7PNn/G1rKtPDfxuYZ71r1jg4kLtfL7tmIuHtO9xX1/tSqHxeklPD0xEBb9B9Z8gva6maePIXvgdYSlHsWHX6zjw8fnAjB5UDxDk8IwGhQer262xg0QEbR3b/J69TXuPTum1at/JAwg0V8j7xsXjPcQ9PWV4D4IF110ERdd1Hhw+fpa9u7OO+88zjvvvIbX8+fPb/jebrc3fP/www83el/9uj59+rB27dqG5fUfCMxmM3Pnzm30nieeeIInnnhirzK89NJL+zgbIYRonXJ/jbu+Blrtv8dtMLh4dc2rjE0Yy8ndT27YXinF0KRwthS03PK4Kb+S12d9z63qc/T6ZSijBUZcgX3UTdz6whbuT+jPJWO6MyAhlF83FZJdVsOJ/WMJsBjpHx/ChtxK4sP2buZujVCbifBAMwMTQptcn7BbTb7++2fOH4bnECS3BLcQQohG3B4v//ppM9cel0psyN412LIaX4272h/cNf5/F+Z/T7mjnFuG37JXv5rukYEs3FqE1rrJPjd1NVWsfucOvjZ+RS1Wdgy8gZRT74LgWHLzq4AtDaE8PDmc4cnhjd4/snuEP7ite+27NZRSfPeX8UQFNf3++EbB7SvH/tyvb0tyj1sIIUQjmwuqeH3B9oYBRvZUXtO4qdxX4/bwRfpHDI8ZzojYEXu9p3tkIHUub6MhQwHmphVwz0szyP/XGC52zKIg5SxOcD7HFxHXQHAsAHkVtYCvZ3dzRvYIBxo3ae+vpIhAAixN9xCv329EoLnZbTqK1LiFEEI0Uj+RRmGlo8n1ZXs0ldc43JhC15Ffk8f9Y+9r8j3d/bXTrNLaXbV4r5es7//FY/Z3qTGFs/Totxl78nnEvLCgYThRgPwK3zCjzXU8Azh1SALVDk+7TXISYDESHmg+qA8GbaVT17hlBLC2JddTCAFQbPfVqAurmgvuxk3ldocba9QCUsNSOT656c6vyZG+wMsq9Y9DXpWP5/1zuMr+JhkR4wi7cxljT/b1BRqeHM6a7PKGv0m5FXUoBXEtdDyzmoxcfnQPjIb2e/S1d0wwfeKaHnmtI3Xa4LbZbJSUlEjYtBGtNSUlJc0+Iy6E6Dp21bjrmly/Z+e0vNptGGy5XNr/0mbHFk+K8NW4d5bWwKYf4JVjIGsJ97muoWDqmxC4q6Y8LDmc8hoXO0p8IZ9fUUtMsBVzE49pdaQ3rhzFY2cPPqRlgE7cVJ6UlER2djZFRUWHuihHDJvNRlJS0qEuhhDiECv234feV43b5dE43B7Sa+eDNjI1dWqz+7SZjSSGGBmx4UlYOAvih/Bm7APMXKF5MCWq0bbDksIBWJNdTkp0EHkVdY16dR8qLT1O1pE6bXCbzWZSU1MPdTGEEOKIUx/cBc3UuEtraghMeQl31SDKao4n17MYm3MoYdamRyUDoLqYN3iMQaXr4Oib4eSH+eG1ZQxPNuzV2atvXDABZiOrs8o5a3g38irq6B1z6JuoDxedtqlcCCFE+6hvKi+2O5p8TrnQswxjQDbW2J94ZvmTuKki3Hts8zvc+DW8eix93Zv5p+l2mPIklW4D63IqOKZn1F6bm4wGhnQLa+igll9R12LHtK5GglsIIUQj9TVur4YS+97N5dWW31HuKDy1Sfyc9Q1GHUq0ccjeO6rMg08ug8+uhOA4Ph36Nu9Wj8Xh9vDn9lK8Go7utXdwA4zoHs76nEq2Fdr944NLcNeT4BZCCNFIUZWDWP80l3ve595Sug0CthOjj6c2+3IirbEE100kyLLbwCVaw4r34OWxsG0OnPwwXPcrAd1HoDXklNWyKL0Ei8nAyO4RTZbh0rHd8WrN/V+uAyChiYk/uioJbiGEEA1cHi9lNS4GJfqG/tzzPvdHaZ+htZEBwSei3eE8MHQGxqqTCLL671OXpMN7Z8C3t0HCULhpEYz/KxjNdPdPqbmloIqvVucwoU90s1Ni9ogK4tKx3fkzoxRAaty76bSd04QQQrS9Uv8EIoMSw5i3uYjCKge17lqeWfYMG0s2srl0C+7KQfTslgBspdapqXF6CTYDv78A858EowXO+A+MuBIMu+qH9YOwPP/LVkqrndx4fK8Wy/KXE/swc0U2NU6PBPdupMYthBCiQX3HtAH+yTbyKqq5Z8E9zNwykzBrGMclnIqj8DSSI3xN19W1dUxwLOC2bdfAnIeg98lwy59w1FWNQhsgJtiKxWRgc0EVY1IiGbWPUc5iQqzcckJvIgLNLQ6+0tVIjVsIIbqgWqenyTG368cSjw+zERlkYl7RG2S65nFZ79u5bfRV/LGthG/cy0kONXKpcS6nLbibSw3ZlNIDLnwfBpwJTUwiAmAwKJIjAkgvquamE1qubde7eWIvrhmfesgHXzmcyJUQQoguZllmKUP/+RPpRfa91tXXuJ2qAJXwKpmuXziv5+W89m0C7y7KpKzGySTDMkZ9fSJPmN/CbgjlBudf+WzsLBh4VrOhXW9ItzCGJ4czsW9Mq8qqlGr2PnhXJTVuIYToYhZtK8Hl0SxKL6HXHgObFNsdKFM5ty+8Eo9JEV17JanG84GNbN62jSmON7jQMgdP4CCm268jqNsJ/FSSx3iruVXHfvaCYbi9TU/tKVpHatxCCNHFrM+tAGDVjrK91hVXOQmKXEutu5ZjAh6mpnQECzYXcYFxPo9mT6d78QKedV+E4Yb5pNmGNzwuFmhpXT3QZDRIDfogSY1bCCG6mA05vuBeuXPv4C6yOzCFrmZwzDB6GXuy2r6MaZmPMsG8mqXe/nwUexd/lEfwfyYLwVZTQ3A3PA4m2p0EtxBCdCEldge5FXXEhVrJLKmhxO4gKnjX4CnZ9u14LLmcmjqN5I1zuM78GBbcLOp7N5etHYrKMdAzxjfZRpDVxNaCKqD1NW5x8KSpXAghupD1uZUAXDqmBwCrdpY3Wp/nXgxaMWn1Vxy37n626W6c6X6aYefdQ3igFa+GiEDf/ewQm4lqpweQGndHkuAWQoguZL2/mfySMcmYDKpRc7nWGq/5d8bWuYjeNo/s0fdxgfMh4lIHEWQ1cVQP3/Ck4YG+GnewdVctW2rcHUeCWwghupANuRX0iAokNtTGgITQXTVut5N3Z55Prbma412BcMMCLMfdgVYGTuwfC8BRPXwDptTXuIN2C+4gCe4OI1daCCG6kHU5FQztFg7AyO7hfL4iG3fpDmZ9dSnPGSvpXRWG4ai3IbYvscA3t4xnQEIIQEONO6KpGrc0lXcYqXELIUQXUVHjIqu0lkHdfMOZjugewYnu3/nm3RN4zFhJX7qzKvtuesbtmmpzSFIYJv+oZUOTwogJsdI3zhfkIbbdm8oluDuK1LiFEKKL2OB/fntwYhg4qjht+yPYIr7h7shoor29WbFlGteM783RPZseQ9xmNrL0vpMwGHyDp9Q3lSsFNpMEd0eR4BZCiC5iR2kNAP08W+C1W1hTm8d9CXEE695kbLmCs4b14IFTB7Q4qll9aMOupvJAs7HRctG+JLiFEKKLKKms4RbjV8R+NouqsETuT+lPoiWYD6e+y5osB8f0jNqvAK5vKg+0SpR0pHa92kqpcOBNYDCggenAZuBTIAXIBC7UWu89fI8QQoi246rlxHV3M9D8Gww6jyeiIyjcOZf3T3qR8IBQju+7/7us70keJPe3O1R7d077DzBba90fGAakAfcCc7XWfYC5/tdCCCHaS00pzDiH/hULeNlyDYuOuZbvdvzMDcNuYGjM0APebXB9jVseBetQ7Xa1lVKhwATgKgCttRNwKqXOAib6N3sPmA/c017lEEKILq08Cz44D8oy+E/4fSy2TWDbls+ItEVy7eBrD2rX9fe4ZdS0jtWeNe6eQBHwjlJqlVLqTaVUEBCntc4D8P8b245lEEKIrqtgI7x1ClTlw+Vf8K1nLKHBtczPms+Zvc7EbGzdVJzNaeicJjXuDtWewW0CRgKvaq1HANXsR7O4Uup6pdRypdTyoqKi9iqjEEIcmfLXwbungdYw/UdIPY7iKgd281Lc2s05fc456EPUN5VLjbtjtWdwZwPZWuul/tcz8QV5gVIqAcD/b2FTb9Zav661HqW1HhUTE9OOxRRCiCNM7mp47wwwB8LVP0DcIBxuD5V1LrLd8xkZO5KeYT0P+jD1Ne4As9S4O1K7BbfWOh/IUkr18y86CdgIfANM8y+bBnzdXmUQQogjSU55LWl5lS1vlLsK3j8TLCFw9fcQ1QuAErsTY+B2Kty5nNvn3DYpj9VkwGRQUuPuYO39MekvwIdKKQuwHbga34eFz5RS1wA7gQvauQxCCHFEePTbjWzKr2T+305oeoOcFTDjHLCFwbTvIKJHw6qiqjosMb8QYo5kUsqkNimPUoorj0nh+H7SKtqR2jW4tdargVFNrDqpPY8rhBBHoq2FVWSV1eL2eBvGD2+Q7Q/tgHC46jsI795o9YKchZgCMzk39U4CTAFtVqYHzxjYZvsSrSOTjAghRCfg9njZWVqDx6vJr6xrvDJ7Bcw4GwIj4KrvG0K7xlXDzsqdlNWV8VXmm3idUZzX57yOL7xoU9KjQAghOoGc8lpcHg1AVmktSRGBvhV5a+CDcyAw0hfaYUl8m/4tr615jayqLDS6YR+OootJCAs6FMUXbUiCWwghOoHtxdUN32eX1QBRUJgGM87BYw2l8IK3CAuM5ON1b/HCyhfoHzGIm4bfREJQAnannZ83FLHS0R+bWTqSdXYS3EII0QlkFO0e3LVQkg7vn8VSq4Wnknuw7ZerG9YPCJnAumVT+d9Jk4gMsgCwdM0qooPLO7rYoh1IcAshRCeQUVxNiM1EsNVEVf5WFnx0Nx+GmlhkNdJNe7h3zL04PA7CLGEsXJnKn6485mws4MLRyQAUVzmIDrYe4rMQbUGCWwghOoGM4mpSok2YzT+yyPEjn4caiLGGc/ugaVwx8Aqsxl2h/OaPCwGYvSF/V3DbHfSKCT4kZRdtS4JbCCEOM1+uysagFGcN79awLKO4mtD4V8hUWxhU5+X2kbdx8rDpe4037vZ4SS+0YzYqft9aTFWdixCbmWK7g7E9Izv6VEQ7kMfBhBDiMPPmwgwe+HI9lXUuAOpcHsLr5pGttjC90k3Rjjs5edh1TU4SkllSg9Pj5fyjknB6vPy6qRCXx0tZjUuayo8QEtxCCHGYqapzY3e4+WjpTgAK18ymf/THBHohafDLpHsTySvf9Sx3eY2TNVnlAGwpqALgkjHdiQmx8tOGfEqrnQAS3EcICW4hhDjM2B1uAN75IwPXhu/g5+nMCbZxUtKZdOs+FICsshoAPF7Nte8t54LXFlNR42JTfhVKQd+4ECYNjGPepiKWbC8BJLiPFBLcQghxGNFaU1XnYlBiKGPt8zDMvJLnoxLxYuS6UbeQ7B94Jdsf3G/9vp3lO8pwerz8tCGfLflVpEQFYTMbuXpcKkFWI3d8uhqAmBDLoTot0YYkuIUQ4jDicHtxeTQ3Jm7necsrPBySwhybG5P9RFIjEkkIs2E0KLLLatlSUMWzP21h0sA4ekQF8u3aXLYUVNE3ztd7vHdsMF/ePI4+sb7XcaG2Q3lqoo1Ir3IhhDiMVNa5GKa2MTXtCb5O7MNX1lqiOJqz+90IgMloID7UxvLMMr5YmUOIzcQT5w7h3T8yeWX+NgBOH5bYsL/kyEBm3XQsq3aW7xomVXRqUuMWQoh2UOfyMOn53/hiZfZ+vc+ZvYb3LE+zPiiGxwPcHJ1wNL9c/gp3nNKvYZukiAAWby+hss7Fu1ePITrYyunDEvBq8GroFxfSaJ8hNjMT+srUm0cKCW4hhGgHm/Kr2FJg57Hv06iodbXuTQUbif/6YkqVlb8lJBJqDePpCU/v9dhXn7hggixG3ps+hiFJYYAvrOubxPvFy0ArRzIJbiGEaAfrcyoAKK128tKvW/f9BnsRfHAeHmXmwvBjKXDm8+RxTxJp23vQlPumDuDX/5vIyO4RDcuUUlw8pjuxIVZSomQGsCOZBLcQQhwErTWz1+fjdHsbLd+QW0F4oJkLjkri3UWZZOw2u1e9SmclFY4KtMcNs66B2lJ+HPMYdeHrmZR8LkcnHN3kMYOspiY7mk0fl8Kie0/EZJQ/7Ucy+ekKIcRBSC+yc+MHK/h+XW6j5etyKhicGMbfpvTDq2HWCt+97lx7Lo8sfoTTvjiNcR+PY/wn4xn9wSierliDnvoMX1T+CdrAFf2n73dZlFIS2l2A9CoXQoiDUFbju3+dXrirRu10e9mcX8X08anEhthIjghgfdFmnvrzSz7f/DkAxyUdxzl9zsFSksH6DR/zQVgoblcWa8rn4iofTUp4wiE5H3H4k+AWQoiDYK/zjXK2vdjesGxLQRUuj2ZItzBcHhfOmFdY4dnEmk0mTu15KrcOv5WE4AQoy4TvJ6DDkzGnjueTzZ9gwISz5HiCrMZDdEbicCfBLYQQB6HKPzzp9qJdNe76jmmDE8OYtXUWVWoT3tJJzLnxfqIConwbuergsytBg7rwAx4OS8KrDOzIC2K1IUqavEWz5DdDCCEOQpV/Bq/Mkmq8Xg3A+twKQqwmYkMNvL72dbrZBlJdcAIe9269vWffA3lr4JzXIDIVs9HMk8c9SbLhTEJsUqcSzZPgFkKIg1DfVF7n8pJf6Zuxa31OJYO6hfL51s8oqi3ivJ7XAorMYt/44qz+GFa8C+P/Cv1PbbS/KoeLYKsEt2ieBLcQQhyE+pm8ADKKq3F7vKTlVdI/IYC317/NMQnHMKnnsQBkFldD/nr47q+Qchyc8Pe99ldV5ybEtvc820LUk491QghxEKrq3CgFWsP24mqsJgMOtxdTyHpKS0uZPmQ63cIDMBkUOQUFsPg6sIXB+W+Dce8/wb7glj/Nonny2yGEEAfB7nCTEGqjrMZFRlE1ueW1mAyKbbXz6BbcjTHxYzAoAz0ibJyQ9hDU7ICrvoPg2Cb3V1Xnolt4QAefhehMJLiFEOIgVNW5CLGZCQ+0sL3YTl55HcNSvKwsXMZNw2/CoAygNX83vMXwqt9h8pPQ49gW9ueWe9yiRXKPWwghDoLd4SbYZiI1JogVO8rYXFBFRPxaAM7qdZZvo7mPcELVd7yhz8I15kamv7uM537ejNa6yf1JU7loiQS3EEIcBLu/htwrOoiqOjfgJaNuPmMTxpIYnAgL/w2/P8eWpPN43HEhL/66jV83FfLfX7fx0DcbGh4hA3B7vNQ4PdI5TbRIglsIIQ5Clb+GnBrje0Y7sdtWCmvzOK/vebDkNZj7CAy5kPzxTwCKl+dtY0xKJDdM6Mn7i3dw/5fr8PjDu76HutS4RUvkt0MIIQ6C3d8LPDU6GPCiIn6hZ3BPTinJ9w2y0v90OPtVUiucAHi8mvtPG8CwpDCsJgP//XUbDreXZ84f6q+xQ7AEt2iB/HYIIcRBqO9MNiAhhHHDcljrzOYfEedh/PYO6H1yw2NfieFGQqwmJvaPZXhyOAB3TuqH1WzkmZ82kxBm4/ShiQCESnCLFshvhxBCHCC3x0uty0Ow1YzJAPaAH+ltjGfS/P9Cj3Fw4QwwWQEwGhTf/GU8CWGN59G+5YTe/LGtmAVbizi+bwyA3OMWLZJ73EIIcYCqHR7A17T92ZbP2F6xnZuzNmOIHwyXfgKWwEbbp0YHYTPvPevXkG5hbMm3U1bja06Xe9yiJRLcQgixH7YWVDHuqV/ZUVJNlcM3wYgyVvLfFS9wjMPNyeZIuGwWWENavc+BiaE4PV5W7SwHkOe4RYskuIUQYj/MWLKDnPJa1uVUNPQCX5D/Ek5XNQ/YPajLv4TgmP3a56DEUACWZJQC0lQuWibBLYQQrVTn8vD16lwA8ivqsNe5CbVuZ0XlEqZX1dHjkpkQmbrf+02NDibAbGyYx1uaykVLJLiFEKKV5qQVUFHrax4vqKyj2l7FWeFvo7TmokkvQMLQA9qv0aDonxCCx6uxGA1N3gcXop4EtxBCtNJny7NJDLORHBlAYUU1/X6/nc1BNQwISCK63+kHte/65nJ5hlvsiwS3EEK0Ql5FLQu3FnHeUUkkhZq5eOc/oWgBaVYLx/Q856D3PygxDJBmcrFvEtxCCNGEjbmVnPfqIqrqfE3jS7eXojWcOjCau6uf5Zi6hbzd4zQATkk58aCPV1/jluAW+yLBLYQQTViaUcKKHWWsy/Z1GNtWaCfCUE2/X6Yxomo+T3ku5zezF68jmoHRvQ/6eH3jQjAalDwKJvZJglsIIZpQbHcAsLmgyvc6N52vbQ/jyVrCE30v4+1wI3nODRhqB6OUOujj2cxGhiWFkRQRuO+NRZcmH+2EEKIJxVW+Ucw251eBvZBbdt5FlKrgsxNv5+P0z7BEKZQ3iADn6DY75nvTx2AySH1KtEyCWwghmlBf487Oy8U74yaivCXMHPQC7+W9Rd+woaxYcgkBZhNxkQFtdkwZeEW0hny0E0KIJhTbHQRSx/8V/R2KtnCD669sjqymsKaQaQOuBRS1Lo+ErehwEtxCCNGEiqpq3rA8x2DSmTvoSRZ6B7G45HOGRg9lcs8JDdtJZzLR0SS4hRBdwndrc3n4mw2t2lZ7vdxZ+x/GGdZzt+t6XsofgDl8BUV1eVw/9HqsZiPRwRZABkwRHU+CWwjRJcxen8+HS3fg8ep9buuY9y/ONPzB791v4gvvBNZkFxIQO4fhMcOZkOSrbceF+ubVDpEat+hgEtxCiC6h2O7A5dEUVNa1vOHGb7AtfJJZnvEUDb+FxDAblsiFaGMld426q+HRr3h/cEtTuehorQpupVSEUmqQUqqnUkrCXgjR6RRV+XqJZ5XWNL9RSTp8fQv26GHc57qO6BAbKfFuLFG/kWQZw/DY4Q2bxoX5g1uaykUHa/Y3TikVBtwCXAJYgCLABsQppZYAr2it53VIKYUQ4iAV233PZWeV1TK2qQ1ctfDZNDAY+XPUczizC6gjj+2mp8GtOaP7tY02r69xS69y0dFa+qg4E3gfOE5rXb77CqXUUcAVSqmeWuu32rF8Qghx0Jxub8N0nHvWuLcV2pmzMZ8bSp5GFayDy2aSVRSNwbqGB5c9jjKYqNlxA+NOHtDoffFyj1scIs3+xmmtT2lh3QpgRbuUSAgh2lhJtaPh++yy2kbrXpm/jcQ1L6HMn8OJf4c+p1CcsRlr9FwAPj7tQ7bkmBmaFNbofdJULg6VVv3GKaWGAim7b6+1/qKV7zUCy4EcrfXpSqlI4FP//jKBC7XWZftVaiGE2A/1w5cCZJXtqnF7vBpL2tf8n/lzspPPIOm4/wMgs2InppD1XNRvOinhyaSE773PYUlhHNMzimHJTawUoh3ts6OZUupt4G3gPOAM/9f+zBh/O5C22+t7gbla6z7AXP9rIYRoN/XDl/aMDiJ7t6byrat+42Hviyzz9mVGzP+Bv8d4Ws13gIFLB1za7D7DAy18fP3RdAtvuyFPhWiN1tS4j9ZaDzyQnSulkoDTgMeBO/2LzwIm+r9/D5gP3HMg+xdCiNao71E+vHs4X67Kwen2YqnOJWn2dIoI5+nQf0COL9Dzq/Mp1AuJ1GOJDYw9lMUWokmtebRrsVLqgIIbeAG4G/DutixOa50H4P+3yf8ZSqnrlVLLlVLLi4qKDvDwQggBRf4a94juEWgNBTk74P2zMLhqeD7mUQb17c4G+2wu+f5SJs2chMZDP9v+NCwK0XFaE9zv4QvvzUqptUqpdUqptft6k1LqdKDQ35Ftv2mtX9daj9Jaj4qJiTmQXQghBOBrKg+2mugTG0wElUR+cQHeihyudPyNwN5Wfqq6FVPsl1TV1XHjsBtx7bydlLDeh7rYQjSpNU3lbwNXAOtoXHPel3HAmUqpU/E9/x2qlPoAKFBKJWit85RSCUDh/hZaCCH2R7HdSXSwhe4hmncs/8JalcPco15i+R+KwvL/EmmLoGjDRVx4yhmc078b//rk54axyIU43LSmxr1Ta/2N1jpDa72j/mtfb9Ja36e1TtJapwAXA79qrS8HvgGm+TebBnx9oIUXQojWKK5yEBtkIuHnmxmiMpjV81Ee2RBNXOpPVDkr+PfEfxGu+rMmq6JhoJboYOshLrUQTWtNjXuTUuoj4Fug4WHI1j4O1oSngM+UUtcAO4ELDnA/QgjRKsVVdfxDvYEq/InnLTfy1uYeuKzLCbAs5YYhNzAgagDDk+2szipr6MgmwS0OV60J7gB8gT1pt2UaaHVwa63n4+s9jta6BDip1SUUQoiDdJb9Yybo72D8nazKnIwzZxlBSZ8wIGogNwy9AYBRKZHM3VTI32auASAmRIJbHJ72Gdxa66s7oiBCCNEe3Cs/4lb9CWkxUxlw0oOEzPyKAMN7JAQl8urJr2I2+sYav3pcCh6vl0+WZRFgNpIoz2eLw1RLk4z8Hd9EIqXNrD8RCNRaf9dehRNCiIOSPg/jd3/hd88gskY+RkH2byx3PEFUQATvn/oWkbbIhk1tZiO3ntiHmyf2xu50EyqTh4jDVEs17nXAt0qpOmAlu2YH6wMMB+YAT7R3AYUQoiV//XQ1YQFmHj5zUOMV+evg0ytwhPfmpty/crl7OU/Ne5IBkQN46aSXiA6IbnJ/BoOS0BaHtWZ7lWutv9ZajwNuBDYARqAS+AAYo7X+q9ZaRkYRQhxSf2aUMietoPHCimz48AKwhrDyuNepIoA5ee/TL6Ifb09+u9nQFqIzaM097q3A1g4oixBC7NPyzFJiQ2x0jwrE69UUVNbh9mrKa5yEB1qgthw+OB+c1TB9NjlZoRgDfie/NoubRzxCoDnwUJ+CEAelNc9xCyHEYaHO5WHa23/y7182A1Ba48Tt1QBsyK30hfZHF0LJNrhoBsQNosjuwBy+jEBTEJNTJh/C0gvRNiS4hRCdxvzNhVQ7PWT5Z/gqqKxrWJe+fRu8exrkrITz3oSeEwHIqyzHFLqWqalTpLYtjggyA7wQotP4bm0eALnlvsCuD+4ESpiy7G9AOVz2GfQ6seE9a8rnoAwuzu1zbkcXV4h2sc/gVkrFANcBKbtvr7We3n7FEkKIxmqdHuamFWI0KAqq6nB5vBRUOkhSRcwMeIJAVxVc8zUFkT3IL1pDsDmYmVtmkuH9iDBDb4ZEDznUpyBEm2hNjftrYCG+x7887VscIYRo2rzNhdS6PJw1PJGvV+eSX1FHXcE2PrE8ijI6Od1yDt3XvcLqopVodMP7nGXHcO3R/4dS6hCWXoi205rgDtRa39PuJRFCiBZ8vzaP6GAL54zoxtercynZuZFz117PS5GKz8Oj8LAQsz2Zm4ffzMCogdiddrbkWHkhrY5xPRMOdfGFaDOtCe7vlFKnaq1/aPfSCCFEM/5IL2bSwDi6RwYyVqXR/4db+dmq+CQ8iJOST+GbBX248MSTGBcTQ7DVRPekQOatWENEYAG9Y4MPdfGFaDOt6VV+O77wrlNKVfm/Ktu7YEKIrunPjFIq61yNltU43ZTXuOgRGUDt2me5P/hfZJuDeSgqmkDdg2cmPkGMpSf/+mkzp/53IRf8bxE1TjfLd5RxVI9IaSYXR5R9BrfWOkRrbdBa2/zfh2itQzuicEKIrqW02snFry/mvT8yGy3Pr6gjjlKO3nQzV+38nMuT4jg31ohbORkddDNmg5lLxnRnbGoUt53Ym4JKB0/8kEZGcTVjUiMOzckI0U5a9TiYUupMYIL/5XyZWEQI0R7W51Tg1bCtyN5ouWvdF3xvvZe/mEKxmIIJrr4AjyWD3Pwkeo3uDcAdJ/dt2H5roZ0PluwEfNN1CnEk2WeNWyn1FL7m8o3+r9v9y4QQYr+U1zi58LXFbPcHs9aamSuyqaj1NY2vzS7DHL6UjWXL8Wov2Ath1nX0W/AXXgyPYb3VzD/GP0a/oJMp3XEO7soRxIfZ9jrOvVP7YzYqbGYDgxPDOvQchWhvralxnwoM11p7AZRS7wGrgHvbs2BCiCPP+pxK/sws5ddNhfSMCWZLgZ3/+3wNd57Sl1tP7MVX2S9gS1hIPnD6x+9ya2E+p1ZW8N+U05nFOib3mMqU1Cn8uX4jdocbgLhQ617H6REVxL1TB1Bid2AxyQCR4sjS2pHTwoH6ebnl46sQ4oDkVdQC/nHFgbXZ5RgsBXy/fSuZxp0U6IVEFQ/lMs9a5oQXcU9EIF/2OZolxetQtT15fPwjAHSLCGjYZ2zI3jVugGvGp7bz2QhxaLQmuJ8EViml5gEK373u+9q1VEKIw0pueS3vL97Bku0lVNW5uPzoHlw8ujsBFuN+7Se/wjdE6YbcCgC+Sv+YoF4zyAXydxg4oSSM/1Z9R46O4vgxjzAvoprX1rxGCL2w1dyEzeQL6W7hu8K6qaZyIY5krZnW82Ol1HxgNL7gvkdrnd/eBRNCHD7eWLiddxdlMrJ7BGEBZv757UZmLN7Bz3+dgMnY+qbovMo6lKmcDHsmz/y5nLV1Mwix9+SiEhdXqZUEuLNY1+taLtx4LI+GjOWGYUmc1fssrn93E2Fhu2rZieG+700GRWSgpc3PV4jDWbPBrZTqr7XepJQa6V+U7f83USmVqLVe2f7FE0IcDnLKaukbG8Ksm44F4NX56Tw9exP5lXUkRbR+xq0NFfMJ6v0GSmneT4PjqxQvFM/HpS38oI/lJdcZfHbOFTg3zSWj2NeBLT4onsKKDfSN3VWzrg/u2BArBoM8oy26lpZq3HcC1wP/bmKdBk5sYrkQ4giUW1FLwm7N00nRbkyhK8koGdYQ3E6PE4vRgld72Va+jdWFq9lRuYNyRzlTUqYQYglhO28T44rhL8W5DPYW4XbEs37YA7xVMYZvt1TTLTyAmBAr3SMDySz2Td3p8WqK7A7iQ3cdPyrIgtVkIE6ayUUX1Gxwa62v9387VWtdt/s6pZT8bxGiC8ktr2NoUjgAHq+HGdsfJaDbGv66+DvG7TyazWWbybHnEGAKwGQwUeWsAsBmtGExWvgm/RsMKGLdMCtvJSXeBF5Qd/Gtcwg/jz+e0dtL+HbLBoZ08/V9TYkKJKO4GoBiuwOPVzcKaaUUPWOC6REp82uLrqc1ndMWASNbsUwIcQSqc3korXaS6A/Ot9a/RVrZGhxFp9A3xUFaaRqDogZxVu+zsDvt1LprGRYzjFHxo0gMSsRdW8bs2bcyN38p00pdrO51N/+rOo7FmRUEmI30ignG4r9PPjTZH9zRQSzNKEVr3dChLSG0cX3hrWmjsJn3r3OcEEeClu5xxwPdgACl1Ah8HdMAQgH5mCvEYWzRtmJ6RAfRLTxg3xvvQ2657xGuxPAANhRv4NXVrzI1ZSo/bj+RIZZkHj5zUPNvzliI+csbOaMyh+OHXs34ncfwUP9jGZBbyeLMCgYmhmI0KFKig/jo2rEMTQ4HIDU6iBqnh8IqB/mVvuDes/d4YhucmxCdUUs17snAVUAS8Nxuy6uA+9uxTEKIg6C15tr3l3PeyCQePXvwQe8vt7w+OK08uuRRIm2RPHD0A6xasbqhNrwXjxsW/At++xdE9YJrfma9syeVfy4lYbcArm8aBzi2d3TD96nRQQBkFFc3HCMuVO7QCQEt3+N+D3hPKXWe1npWB5ZJCHEQiu1Oapyehprywcr1D5qy2f4bG0o28MT4JwizhhEfZmsYUKWRnUvh+zuhYD0MuxROfQasweSv9D2YEh9ma7hfPbJH0xOApET5gjuzuJr8yjrMRkVUkDz2JQS03FR+udb6AyBFKXXnnuu11s818TYhxCFWH9h5zdWGD2B/yuDg/c2vMCR6CKf1PA2AxLAANuVX7drQ64F5j8PCf+MNSeTL3k9Sm3Aao8s0/eJ3lSc+zEagxcTsO46jb2xIk8dMDA8gPNDMV6tzSAgLIDbEJo99CeHXUlN5kP9fmYFeiE4kxx/cBZVtE9x55XWEJSyiuLaIF054HoPydSSLD7NRbHfgdHuxOMth1rWQPhdGXMFzxqt56fd8WL8eg4LPbjiG/Io6wgLMBFp8f3b6xzc/O7DRoLh7cn/u/3IdAWYjAxNlJmEh6rXUVP4//7//7LjiCCEOVn2Nu6TaicPtwWpqvue13eFGa02IzdzsNjsrivCG/MYpPU5hWMywhuWJ4Ta0htLtK4n6bjqGqhzcU5/DPXwa7z85l1OHxHPPlP6c+p+FzFqZTVGVs9H97X25eHQyX67KZllmWaNnuIXo6lozree/lFKhSimzUmquUqpYKXV5RxROCLH/sst23XcurHS0uO3tH6/ijk9Wt7hNuvM7tHJw87CbGy2PD7VxoXEeMZ+djqO2mgvq/s59mUcxc0U2lXVurhnfkx5RQUwaFM8P6/LJKq3Zr+A2GBRPnjsEs1GRLM9rC9GgNc9xT9Ja362UOgffsKcXAPOAD9q1ZEKIA7J7p7T8yrpmQ09rzbLM0hZ7axfXFFNj+40elmPpHdF714qKHEb9fjPHm+dTFD6Ge7md9TVGVq7KYfaGfEZ0D+cof8ezM4Yl8OWqHCpqXc12RmtO79gQfrz9OGKlxi1Eg9bMDlDfhnYq8LHWurSljYUQh1ZOeW3D89stdVDLKq3BGfExRbZ38Gpvo3Vaa5bkLeGmObeCcjMx7rL6FbDqA3jlGALz/+RB1zRm9H2R33INTB+fynF9oqlxehpNqTm+dwzhgb4/I/tT467XOzaE0Baa8oXoaloT3N8qpTYBo4C5SqkYoG16vQgh2lxueW1DbbegheB+b/3nmMNX4AlawcebPm5YXueu4y9z/8p1P19HYU0hdbkXMiS2N1TmwkcXwte3QNwg1E1/8KXpNL5bm4/bqxnXO4qXLhnJfy4eztTBCQ37s5gMDa9lCk4hDt4+g1trfS9wDDBKa+0CqoGz2rtgQoj9V+N0U1bjol98CIEWY8OoY0u2l7Ct0N6wXX51Pl/ueBV3dSruqv48t/w5VhasZHPpZq75+Rp+y/kVd/EUzoh8EW/lMAblfQkvHw0ZC2HK03DV9xDZk/gwG9uLqzEbFaN6RBIWaOas4d0w7vHo1vlHdUMp6BMrD6kIcbD2eY9bKWUGrgAmKKUAfgNea+dyCSEOQP397aSIAOJDbeRX1OH1am6Y+REDEoL56PJpZFZmcu+Ce/FoD3V554PXQnTUK0ybPQ0Aq9FKcMXV5BX1ZdG8n/nO8g7Jf+yAHuPgzBd9I6H5JYQHsLXQzvDkcAIszfdeP6pHJMseOJnoYGv7XgAhuoDWdE57Fd997lf8r6/wL7u2vQolhDgw9T3KE8MDiA+zkV9Zx9ztq/HEvskGr5czvpxJfnU+AeYArKVXYiOWWo+HR8e8RqF7PaGWUJIC+3De82t5I+FrTiz9jAIi8Z73DobB54BqXJOun/jjmJ5R+yybhLYQbaM1wT1aaz1st9e/KqXWtFeBhBAHrn5c8W7hvhr3kowinl7xNNoTiK3qdGLitjIoehDXD7qDk59ZxXF9Ili4tZhAQzwX9hsIwKol8/jW8gB9y3LI7n0x38fdxA1Dhjd5vPo5uo/ute/gFkK0jdYEt0cp1UtrnQ6glOoJeNq3WEKIA5FTXoPRoIgNsRIXZqPU/ANVddtx5F1BXc1g3rrtfowGxaJtxQAc2yuahVuLKa9xgscFC55l2G/PUKDCqD7/U5IGT+GGFo53fN8Y1mVXMLL7/j3mJYQ4cK0J7r8B85RS2/FN7dkDuLpdSyWEOCC55XXEh9owGhQ7PN9ijp6DoXokbvsgQFNsdxAXamNjXiUAx/pryrpgAyx8GPLXsjT4FB73TuO7wVP2ebwR3SN466rR7XhGQog97TO4tdZzlVJ9gH74gnuT1rrl4ZiEEIdETlkt3SICeGb5MywsnoGrYjh1eecyLDmcNVnl5FXUERdqY0NuJfGhNlKjbNxo/IaTF86CgDC46APu+iqQUSmRh/pUhBDNaM2QpzbgFuBh4EHgJv8yIcRhJqe8lqiwWmZsnMHEhNOpy70QtInThsQDNMxtnZZXycToCkI+PJ17zZ+QHnEc3LKUoqRJ5FbUMTQprKXDCCEOodYMwPI+MAh4EXgJGAjMaM9CCSH2n9PtJb+yDodlJQBXD74aMGBQMGWQbwCUgso6PB4Px5V8zqP5N6JKtnKfup0ZyY9AUDTrcsoBGJoUfmhOQgixT625x91vj17l86RXuRCHn8ySajxeTb5nKQMiBzAsvg9Gwzb6xYWQFBGA2aioKdiO6+2beMC4iNzoCSRe/jqL/7eJoXVuANZkVWBQMEim0RTisNWaGvcqpdTR9S+UUmOBP9qvSEKI1tBa8+mynQ3zbm8tsKNMZWTVbGJSyiSMBsXgbmGcNCAWg4LrAhdw9bpLMRWs4W+u68md+i6EJhAWaKG81gXApvxKUqODCLK25jO9EOJQaE1wjwUWKaUylVKZwGLgeKXUOqXU2nYtnRCiWWuzK7hn1jreXZQJwLZCO+awdQBMTpkMwFc3H8udY4Phg/O42/Uq20z9mDnmMz73TKRXbAgAYQFmKvzBvbO0lh5RQR1/MkKIVmvNx+p9PxMihOhwX6zMBmBtdjkAWwurCIxYR/+oQSSHJAOg1s+C7+8Ej4tPY2/jf9UnMNYeSmRQLRFBFgDCA8zsLKkGILushjEp8ky2EIez1jwOtqMjCiKEaD2Xx8u3a/MAWJtVgderWVe2GE9wFmf2ug+8Xpj3GCz8NySPhbNfZetiB3lLd7Kt0E6vmF216voad0WNi6o6N0kRTc/fLYQ4PLSmqVwIcZj5bXMRpdVOJg+Ko8rhZn1eCSXWTwkzJnFBzzPhi2t9oT3ySt9MXlG9iA+zUevysC6ngt67zdIVHugL7p2lNQAkRwYcqtMSQrSCBLcQndCXq3KIDLJw20l9APjXkldQ5jLOi78S8wfnwfpZcPLDcMZ/wWgGds2FXefy0itmV3CHBZjxakjL942mJjVuIQ5v0nVUiE6m2uHml7QCLhmdTJ/YIIJi57Om6mfMVf24ueQpqM6DC96FQec0el986K5xk/YMboANORUAJEtwC3FYk+AWopPZVmjH6fYypmcYN869HkPUMqyVffi4ZAUWmwemfQvdx+71vrh9BPf63EpCbCbCAs3tfxJCiAMmTeVCdDLpRXYAttX9wrL8ZZxsPpUfi5cRhQd11Q9NhjbsCm6LyUC3iF33scMDfb3LN+ZWSjO5EJ2ABLcQh5nsshomPjOPLH9nsT2lF9kxmeqYmf4OY0N78kz6u9RoG0/GPw/xg5vdr8VkIDrYQs/oIIwG1bC8vsZd6/KQHCEd04Q43LVbcCulkpVS85RSaUqpDUqp2/3LI5VSvyiltvr/lYdGhdjNhtxKMktqWJ1V3uT69MJqorr9QbmjnDvTFuGNGci5zn8S0m3APvc9qkck43pHN1oWvlvTuNS4hTj8tec9bjdwl9Z6pVIqBFihlPoFuAqYq7V+Sil1L3AvcE87lkOITqW02glAbnltk+s3F2/HFTGH0+zVDOw5CX3u61y/OJ+TBsTuc9+vXXHUXsvqa9wgj4IJ0Rm0W3BrrfOAPP/3VUqpNKAbcBYw0b/Ze8B8JLiFaFAf3Hn+KTh353Q6CLH9i1rt4Y4ep8Opz6MMRm44vtcBH89mNmI1GXC4vdKjXIhOoEPucSulUoARwFIgzh/q9eHeZDVBKXW9Umq5Ump5UVFRRxRTiA6VVVrDA1+uw+H2NFpeYvcFd86eNW63k5mfnE1GYB2XmIcTf9p/wGBsk7LU17qTpMYtxGGv3YNbKRUMzALu0FpXtvZ9WuvXtdajtNajYmJi2q+AQhwis9fn8+HSnazNrmi0vLTaATRuKnfWVfLZx6fxomcnMTVhjD/uv6AUbaUhuKXGLcRhr12f41ZKmfGF9oda6y/8iwuUUgla6zylVAJQ2J5lEOJwtb3YN7HHxtxKRqdENiwvqXaCcpFTWcTWsq3MTv+Wrza8TyEeehrjWJt7Nb1j23a+7PBAMxGBZoJlOk8hDnvt9r9UKaWAt4A0rfVzu636BpgGPOX/9+v2KoMQh7OMYt/z2BtzGzdE5dWlEdz3RTwGJ+d+AwYNx9TV8eigq/gi7zRybcWNOpS1hd6xwQRaJLSF6Aza83/qOOAKYJ1SarV/2f34AvszpdQ1wE7ggnYsgxCHrYz6GnferuDWWlNkmQleGxSeyIMRv3N84VZiz/4fDD6PZ1/5o9HMXm3lsbOH4NW6zfcrhGh77dmr/HeguZtwJ7XXcYXoDKodbgoqHVhMBjYXVOHyeDEbDczZMQevZQfday/gqeovGOzIxnDBOzDgDLTWpBdVc/rQhDYvj9GgMDb731UIcTiRtjEhDoH62vbEvtHMy/uSs756ldEJI/gzbxk4onm9/CtiVBa/j3yBCQPOACCzpIaKWhc9dxtnXAjR9ciQp0K0k22FVXy0dGeT6zKKq8FQS0nwq9jiv8XtNvBz5s9k27O4t7SAeGc217vuYrnVN+64y+Plzs9WE2I1MXlQXEeehhDiMCPBLcQBmrOxoNlgBvjkzyzu/3IdtU7PXusyiquxxvxCun0lrsKzmRD0OL9N/YTPy42cU1vM2gmvszl4bMMjYc/8tJlVO8t56ryh8siWEF2cBLcQB+iV+dv457cbsDvcTa4vr3UBkFW292QhaUU7sUQs5azeZ9E3YDL5WelYZpxN76pCpjnvQfU8nsRwG7nltazLruD1Bdu5/OjunNYO97eFEJ2LBLcQB8Dr1WzKr8Lh9jJnY0GT21T6g3tnyd7BvcY+CwXcMPQGxkfbubfgLnR1Mb+N/R/LdH+igqwkhgeQW17LO39kEGQxcveU/u15SkKITkKCW4gDsLO0hhp/E/g3a3Kb3KbCH9w79piec2flTqpMf5BiPYHEqmJu3X4LwbqaonM+ZatlIACRwRa6hQeQU17Lt2tzuWBUMqG2tn12WwjROUlwC3EA0vzPXh/TM4oFW4oor3HutU1FQ427umFZZkUm1/x0Hdpr4uKAQfDuaZjMZs53PsRKd09Kq51YTAaCLEYSwmy4PBq3VzPt2JQOOS8hxOFPgluIA5CWX4VBwV2T+uL2amavz99rm4amcn+Ne0PJBqbNnobdVc3w7GO5eO29EJqId/rPZKpk1maXU2J3EhVkQSlFYrhvwo8T+sWSGt32g64IITonCW4hDkBaXiWp0UEc1SOC1OggvlyVs9c2lXW+Tms7SmtYkreE6bOnYzXamG4fy/ueGTjjhsPVP2KL6k7/hBDWZldQWu0gMsgCwMDEUEJtJm48iCk7hRBHHgluIQ5AWl4lAxJCUUpxyZhklmaUsnR7ScN6t8eL3eHGYKomj2+5ec7NxAcmcFlWN67LfZOt4eOwXf0NBPomFxmaFM6a7HJKqp0NwZ0UEciahyYxJjWyyTIIIbomCW4h9lNlnYvssloGJPhm6LrymBTiQ208PXsT2j/ed2WdG1PIWoJ6P4Up+meOij2avxWEM61iFtuTz6XfbV+DZdfz2MOSwqiqc7Mpr4oof3ADqDaculMIcWSQ4BZiP23KqwJgoD+4bWYjd5zch5U7y/nF/2hYcXUN1rgfCDEkUp1+B4+Umhif/yXfBZ1Lz+lvg7HxaMNDk8IBcHq8RAZZO+5khBCdjgS3EPtpU76vR3n/hJCGZecflURqTCD/+XULAD9mfI/BXM4pcVfyoHc2CZvf4033VHLHPABN1KL7xAYTYDYCEBVs2Wu9EELUk+AWYj+l5VUSHmgmPtTWsMxkNBCW9DU7Ah5gQfYCvs6cgac2kZsLZnOV6WdmqDN5zH05Jw6Ib3KfJqOBwd18NfjIIAluIUTzJLiF2E/rcyoZ6O+YVu/rbV+TXjcXlJNb5t5CUV0OV5TXErd9Fu+YLuIftReREhXU4lza9c3lEtxCiJZIcAuxHxxuD5vyKxtCFmB7xXYeX/o4qcFDsW+7l7O6X8Rx7iDuqVtD1XEP8mvCNYDixP5xLXY2G9Hdt8+43WryQgixJwluIfZDWl4VLo9mWFJYw7LHljyG1WjllkEPYtVw55alvJKVxoOuqzFPuIPukb7e4ycNiG1x31MHJ/D6FUc12rcQQuxJgluI/bA2uxyAocnhAPyZ9yfL8pdx47Ab6R8Rzf/MzxOZt4DvU+/nMyZjNRkY3zuaoUlhjE5p+Xlso0ExaVC8PAImhGiRad+bCCHqrcmqIDrYQmKYDa01L69+mdiAWM5PPR3TZ1fTw7iGef3+we/WSYQGFKKUYuqQBKYOkek4hRBtQ4JbiP2wNrucoUnhKKVYlLuIlYUreWDU3VhnToftc3mQGzEETaWiqo6wAPnvJYRoe9JULkQr2R1uthXZGeq/B/3u+neJC4zj3NXfwrY5cOaL/BEylcKqOipqXYQFyDScQoi2J8EtRCutz6lAaxiWFE5RTRFL85dytseKZetPcNq/YeSVxIXaKKx0UFHrIlSCWwjRDiS4hWilho5pSWH8lPkTXu3l1PQlcPy9MPpaAGJDrBRU1VFZ65YatxCiXchNOCH2Iau0hnmbC5m1IoekiACigq38sGEGAxxOeg65FCbe27BtfY3bZjZKcAsh2oUEtxD7cN37y9mU75u16/oJPdmx+VvW1eRylyUWTnuu0djjMSFWHG4vDrdXglsI0S4kuIVoQVm1k035Vdx+Uh/uOLkPqjCNVz+/CxVsZcrpb4CxcTjvPupZqE2CWwjR9uQetxAtWLGjDIBxvaNRpdtxzTibWYEWxsQMJT66/17bx4bsmpJTatxCiPYgwS1EC5ZllmIxGhga4YQZ5/CdRVNgVFw1/KYmt29U45bgFkK0AwluccQqsTvQWh/UPpZlljKqmw3bzMvx2At4OyGFAZEDGJc4rsntY0Olxi2EaF8S3OKIlFdRy9gn5vLrpsID3kedy8P6nDL+4X4RspczZ+JtZNbkc82Qa5odTzzQYiLY6us6Eiojpwkh2oEEtzgiZRRV4/ZqlmWWHfA+VmeVc4f6lAGlc6k96UFeKlpCSmgKJ3c/ucX31de6pcYthGgPEtziiFRQVQdAWl5lq9/zyZ87WZ1V3vC6evHb3Gz6BsewafzbXENmZSb3j70fo8HY4n7qO6hJcAsh2oMEtzgiFVQ6gP0L7n9+u5E3Fmz3vUj/lYlbn2CZaSRLR57Dp5s/5cqBV3JM4jH73E9cqA2DoqHJXAgh2pIEtzisLNxaxG0frzroTmX5Fb4ad2GVgxK7Y5/bVzvc1Lo8pOVXQsFG9GfTSCeJL/v+g4eWPkLfiL7cPvL2Vh17VEoko1MiZV5tIUS7kOAWh5X5m4v4Zk0u+ZV1rdp+U34lxz45l+yymkbLC6t2vT8tr2qf+yn2h7u7OAP94fl4TQFMq/sbWdYfKa0r5dFxj2IxWlpVpiuO7sGnN+y7Zi6EEAdCglscVspqnABsyt932AJ8uSqH3Io6Nu0RzvkVdQxICAVa11xebHeSSDEfmh/H47Cz5JjXKAyoZHX5bK4YcAUDowbu55kIIUT7kOAWh5XyGhcAm1sZ3HM2FgBQtEdzeEGlgwHxIcSFWlsV3LW5G/jE8ihhys780f/jpwojtm6fkRiUyM3Db97PsxBCiPYjvWfEYaXcX+Pe0org3l5kJ72oGoCiql3BrbWmsKqO2FAbAxJC2ZBXwY7KHSzPX05JXQl17jpOSD6BITFDfG/Y+gtj5l5FkcHIRZ6/0Keqmt8rHsdkdPPcCa8SaA5s+xMVQogDJMEtDiv1Ne7WNJXPTfMNrmIxGhoFd2m1E5dHEx9qJd+ZwTL3a5z+ZXHDeoXijXVvcHrP07nLHUT0b8+yNaQ3l4RZ0OZPyK4EvGEcG/ggg6IGte0JCiHEQZLgFoeV+nvc24rsuD1eTMbm7+b8klZA//gQvFo3Cm7fo2BellS8zR/lX6F1DNP738U5/SeSGJyIy+vizVUv837aB6xwu3huwCRudtag3fkMNF/NuiwntZUpjJ06oL1PVwgh9pvc4xaHDa9XU1Hrolt4AE63l8ySmma3Lat2sjyzlFMGxhETYm10j7ugsg5T2Cr+KPqKqd3PoybjNoIdE0gJS8FitBBUns3tSz/j/dwCXNZQLqlLo8y7kzD7NE5JPpuasiFoTwiDu4V1xGkLIcR+keAWh42qOjdeDWNTI4GWO6h9viILr8YX3MHWRjXu3IoqrNFz6B3Wn6eOf5Djeifw37lbKaiohRXvwRsnQl05gy6ZxUdnf8WY+DHEOS8i0TKK/vGhDfsZmBja1KGFEOKQkuAWh436ZvJRKZEYFGwuaDq4y6qdvPTrNo7vG8PQpHBfjbtq10xgC/K+x2Ap4/aRt2EwGHjs7MFEeIoo/N/Z8O1t0G0kXP8bpIwnITiBtya/hbdiHDHBVvrFhwCQGh0kI58JIQ5L8pdJHDbKa30d0+LDrKREB7E5v+nHuP7761bsDjf3nzqA9PJ05lX9HU94N1bkJ1LiyGNF5ecoRyrHJ48Hr5ce2z9hjvUfeKpdLOn3N8ZefB9qj/HGi+wOjgu2EBNiJS7UyrAkaSYXQhyeumxwF9sdrNhRxuRB8Ye6KMKvvsYdFmChX1xIkz3LM4qrmbF4BxeNTqZvXDA3/HIXJc4dWKLSufrn+QAYdQAxzvNRJdvgm9tg5yJMKRP4u/MaPltr5gTHCp6/aDjhgb6R0BxuD1V1bqKDfZODfHjtWMICWjdKmhBCdLQu21T+8dKd3DBjBRX+Wt6+lNgdXPHWUgpaORSn2H/1z3BHBJrpFx9CZkk1dS5Po23+NXsTFpOBv57Sl99zfmdx3mLO7nEt1dvu5Zp+9/PxaR/TreIprvVshFfHQeEGOOtlDNO+4enrzuKfZw5i4dZiXp63rWGfJXbfcaP9s3r1jg0hxv+9EEIcbrpscNePTV0/GcW+LN9RxsKtxfy+tXjfG4sDUlbt+xAVEWghNToIrSGrdFfP8mWZpfy4Pp8bj+9FZJCJfy//N91DunNen4vQ7lD6BE5ksMvDQ+UPcknZq9DrRLhlGYy4HJRCKcW0Y1MY3yea2RvyG+6J1/8uRAVJLVsIcfjrssFd6h/oI7eitlXbZ5f5tttebG+3MnV15bUulILQADM9ooIAGh4J01rz2PdpxIVaufa4VN7b8B7pFen89ai/khgWTJIqZMiCG+B/E+jvTeeXXg/AJR9DSNxex5k6OJ6s0lo25PruodcHd7TUsoUQnUCXDe6yal/zaGtr3PWzT6UXVrdbmbq68honoTYzRoMiJco3zOiOEt/1/nVTIWuyyvm/Sf1Ir0jjpVUvcUqPUzgp+UTCN33MbMu9xJUtp/KYexjv+A/F/S6GZqbVPGVgPAYFP23IB3wTjADEBEtwCyEOf102uEv9wZ3XyuDO8de404ukxt1eympcRASaAQgPtBBqM7HDX+NevqMMs1Fx0qAw7l5wNzGBMTyUcjbqnakYvr2NTYbe/LvXu6xMuZZKgkiOaH588cggC2NTo/hxfX1w+2vcEtxCiE6gywZ3fQ/mvPL9ayrfUVKD2+Ntt3J1ZeU1zoae3gAp0UHs8N/j3pRXSc+YQB5a/AB51bk87Q4l7N0zoDQdzvgPj0Q8yVZHOL9tKcJqMjAqJaLFY00dEs+2QjvbCqsornISZDESYDG2+B4hhDgcdMng1lo31LjzW9lLPLushhCrCafHS3ZZLQ63h6dnb2o0YtfBenneNl77Lb3N9tfZlNe4CPfXuAG6RwY2NJVvzq8iNOJL5mfN5+7iEkakL4IJf4PbVsFRVxEdGkCR3cFvm4sY2zMKm7nlEJ400PcY4OfLsym2O4iS2rYQopPoksFd6/LgcPtqzfVN5auzyrnwf4updXr22r6yzkVlnZtje0cBvuby3zYX8er8dOZtKmzYbn1OBQ733u9vrW9W5zJrRfYBv7+zK6txErF7jTsqiOyyWsoKsjjG+DSb+JXz7DVcMuAKuH0NnPh3sPpGOosJtrK1wM724mom9o3Z57Hiw2ycO6Ibb/+RwZrscqKDpUe5EKJzOCTBrZSaopTarJTappS6t6OPX1/bDraayCuvRWvN7PX5/JlR2uQwm/X3tyf4AyG9yM68zb7ALvHvq6zayVkv/8Gny7IOuFyFVXXsKK3B69UHvI/OrGKPGnevEBd3Gj5m4Qfj+SU2l5GE8cClc1BTn4TgxuEcE2Jt+DA2sd++gxvg/tMGEGA2sqOkRu5vCyE6jQ4PbqWUEXgZmAoMBC5RSg3syDLUPy88ICGEaqeHKoebjXm+R4Mymnjcq/7+9qDEMKKCLKQXVjNvUxEApdX+58Er6/B4NVsLDqzzmtPtpazGhdPtJa8LDvLi8nipcrh9NW6Pm9W//oPv153PwuQ/+XtUKMaaZB4+/VvMEalNvr9+wJTkyABSo4NadczoYCt3T+nv+14eBRNCdBKHosY9Btimtd6utXYCnwBnddTBtctBsb0G8DAgIQjwkFVWxYbcUsBDRtHej3vl+B8FS4oIoFdMML+kFZBfWQuG2oZRt+r/3VHa/FSULSnebVrKzOK2feSssKqOx7/fiNN9+HaqK/c/V9/HsYHq1ydwd8ZMtlltbHP3xl1xLO6im0iJDG/2/fXBPbFvLKqZx8CacumY7lx+dHdOHZxwUOUXQoiOcijGKu8G7N6enA2M7aiDv/TD1bxevo6QAfBlOYQMgIt/BpIhBHg/L5Bfv4whwhpBVEAUCUEJrM5xEpS8gctm/5fKIKjDQHC3IpTRwZqayWg9jBJ/zXtnyYGF7u6d3DJLqhnXO7rV762sc7G9qJrhyeFNrn9/0Q7eWJjBpEHxjE6JPKDytbfK0gKeNL3B1D/n8UR8EvkmM+9NfptLXi6gzuVlTEpki4Hcw//c96RBew+40hKDQfHY2UMOquxCCNGRDkVwN/XXd6+bukqp64HrAbp3795mBx+bPJHSkgrs5cWkBDpQteWY8KCAKmUlwxSMLcRMmbOGzLoyFuUuotZdizkghmGxo9hSUE55ZRFhnrEYlYNCw0889WcUlWU9MEeuJsc+GLfHi8m4f40ZhVUHXuOesXgHz/+yhTUPTSLIaqLW6WF9bgWjUyLRWvP9ujzAN0FHWwd3dlkNoQFmQm3mfW/cFI8Llr1Fj1+foIexil/6nccnzuVcOuBSRsSPpEfkAjYXVNE/IaTF3QxKDGPB306ge1Tzz28LIcSR4FAEdzaQvNvrJCB3z4201q8DrwOMGjWqzXprjRl5PQuLj+e9rdtIe3QKA/7xA/3MJQzxrOfC6J1cUrSchOKNu8qRMJR3yo+iMqgfd/SbxLzoUK5esYYLTu7DjpJqFhTb+GjTRwDY4kB75zJjfSBXDzt/v8pVX+MOtZkahvncn/e6vZqM4moGdwvjs+VZPPTNBj66dixhgWYy/B8E2roJHuDi15dwXJ8Ynjz3AGqtW+fAT/dD8WYqYo/lkqyz6B6+gYjKCG4bcRvgq0lvLqhqmCe7JRLaQoiu4FAE9zKgj1IqFcgBLgYu7cgClFb7BvqwmoxEBweQVhWDPfJUBh/dk/O/Xs/yvwwi2pkFuatRG79meu1HUAu8+hATTTbmhPUkxn06nxiH80Peqbx782W8vSCXORvs2BJm8dzqf7KlcgUPjH2AYEtwq8pUWFWHUnBUj4j9Dtj6WbXqg3uTfx7rJ3/cxLje0RgNiohAC5kH2IzfHIfbQ3ZZLat2lu1z23/N3sS6nAr+fcEwYh074ecHYOvPENkTLv6YuZWD2ZK9nJLiNzmv77kEmn0hXN8E3r8VwS2EEF1Bhwe31tqtlLoV+AkwAm9rrTd0ZBl8zwv7mnYTwgMorHIwKCHM3xtZsc0RSnTPCZA6gepRN3PsQ7P4+7EBXJBSh8pbQ++dS2DJv7gBGGdIoU/GjXxU3ZcoayjFO69j6vgN/JDxMWuK1vD8xOfpF9lvn2UqqnIQGWihd2wwf6SX4PVqDIbWdbIq909NWl+z3lZox2oysC6ngs0FVRzbKwqTQZFRvH81+TqXh9VZ5QxLCm9yVLH6cd63Ftqpc3maHfTE7fHy4dKdUFvGvP/8hwu8P4E5AH3yIxiPvhFMVsoXpGMK3ojT62Bq6tSG945JjeL7tXn0jw/dr7ILIcSR6lDUuNFa/wD8cCiODb4ad6R/CseEUBtrgIGJoaT4Z6TKKK7m6J6+wVZyymupIBhrjxEwNBGGXujbSVU+K2e/h3ndR1h/uZfnMLEiYBwzzEfR3XsG10+ZzN9++xtXz76al09+mRGxI1osU2GVg5gQKynRQQ2PhHULD2jV+dT3yN5eZEdrzdZCO+eM6MbqrHI25Vdx+tAENuVXsTSjFK31Pntde7yaB75cx3dr87A73Dx4+kCmj9/7Mazc8rqG7bcUVDE0KZxfNxUQHmhhZPddQ46uzizgHOe33BP8FRa3nY/dJ/Bc7fnoeTFc7dhJn7gQvlyViyVsDQlBCQyLGdbw3lMGxnHKwP3rcCaEEEeyLjlyWlm1q2GErvgwGwCDEkNJDA/AYjI01FyL7Q7+/tV6APrF7dFUGxJP2eCrOMP5BJvO+p4vjZMZ4lzFy8bn+OuqKYz49VlmpFxAlC2c63++nq+2fYXH2/yoakX1we3/8LCjiebyGqe7ydnMKnarcZdUOymvcdE3LoR/njmI4cnhTBmUQGp0EDVOT5NDtLo8Xp76cVPDwDTpRXY+WZbF0T0jCTAb2dnMI255u02Juj6nEpfHyx2frObZnzb7d1wLy9+mz6cTedj8PuZuw6idPo9uV7zG384dz7CkMP79yxZu/GAFZbVlmIO3MiVlCgbVJX8thRCiVQ5JjftQK61xMqJ7OOAbD1spX6/k+ukktxdVk1lczaVvLKGk2skLFw1vsnNUfa09x9aHvzuuIOPo+zBmL6Jv6W+cmbOChE3f8Y7Jwu3JPfjHH//gnbVv8OC4Rzgq7qi99lVU5aBnTBAp/sFDMkqqOXa3R8K8Xs30d5eRWVzD4vtObFRrrr/Hvb2oumEAmN6xwYztGcVXt4wDaNSaEBtqa3TsTXlVvPZbOkkRAVx+dA/SC337uOPkvuwoWU1uMxOx1A8XG2gxsiG3gtTMICrrXFCwHn75AVbOgNpS8ox9eTHyNv5+1c0EK8VE//svHtOdtLxKMkpKWG3/hI83e5iSOqXpH5oQQgigCwa31pqyaicR/tC9aHQygxJDG2reqdFBbMyr5IYZK6hxeZh547EMSQprcl9RQb5BP3aU1OB0e4kIDaQg6TjuyenBGffOQOWuIjrtW2Zs+pZfaot4weXiph+v4rVupzFy7G18mwGxIVbGpEY21LgTQm1YTIaG6Szrfbh0B0u2lwKQXlRN71hfpzevV1NR6yLEZqKqzs2S7SUA9Ilr3CmuPrgzS6oZ678NUK/K4auxb/UP91o/dWnPmCDiw2zNTsSSV1FLWICZgQmhrM+pYGj1Ir6xvMpQdwZ6kRHVdwolQ69lyowa7pswYK/5sfOr8/kh50NmbpmJ3WXnlB6nMCByQJPHEkII4dPlgrvK4cbt1UT6m8qDrKZGQZYaHcxPGwowKHhv+phmQxsgyj8xxRZ/4EUHW7GZjdS6PBTZXcQmjYKkURhOfpjJRZux/fwGzzjncHP2t/x97Ufk2cfyTdTl9Js+FafHS2SQoriuiB6RgSzYUkS38ABiQ6zYLEae/HETAxJCScur5M+M0obgrnK48WoYnhzOwq3F/LyxgGCrifg9atURwRqz0ddBLa+ilke/28jjZw8hIshCtcPjPw9fYKcXVdMtPIBAi4nEsADS8vYevx0gr7yOhDAbp4RmMSrtSYYWp7NDxfEP11VceOVfGNKvNz8t3Qms44T+sQ3vK6sr48VVL/Ll1i/RaCb1mMS0QdMYFD1of36UQgjRJXW54C7z38etr3HvqT4Q75rUj+P6tDxZRaDFiNVkaJiYJCrY2tB8vqO0ZleTtFKkqySmrT+Ro/tOQoU8x32GXAxRaZxZdQf6vRc5y9qHGVlpvJRRQlBMDBVlyTy2IAVPTQ+8zmiCrVZevWIY57+6mD8zSrh0rG9Qmgp/x7SR3SNYuLWYtLxKhiWHNzSla615dc2rvLbmNUJSe7GoaCzrZlv4M7uU83YmcNKAROz+GveW3WrcPWN8NfSEcBvFdgdOtxeLqfG9Z29pJn/3fs74TT+RTwT/57qB0DGXMWNxNsPsNoYA8zcX0i08gD6xwWitmbV1Fs+veJ5qVzUX9ruQqwZdRWJwYut+eEIIIbpecNd3wIoManqkr9OHJhAZZGZi39gm1+9OKUVUkKXhvnJ0sIVAi++SPvFDGqE2M9dP6Mm43tG8tygTgJXpXpb9/Sue+W0OP6Z9xleRq8lyFFKQnI12Grk9qA8bA4NZbs6gzLESALPBSrA5mDO/LcWYFMTvBafj8Q7FaDBSXus7n4GJoViMBpweL71jfB8+tNY8v+J53tnwDuMSx7EiJ50M/SHUgS0eFhX046QBV2CvcwO+mc6K7Q7SC+1cMMo3Rk6C/xZCQWUdyZH+AU5yVsCfb/BG5adoZaJs+A2ctGQM1QTw2/jezPgzh22FdjxezZLtJUwdnEBRbREPLXqI33N+Z3T8aO4fcz+9I3rv189OCCFEFwzuMn9Hrt3nfd6dzWzkxP6tf/woKthKrr+TVrS/xj00KYwSu5Os0lpu+Wgln99wDDNXZNMzJojtRdUs3V5BZnYsIe7p5OX/xoq4L1DeAN4w92bs1j/AVY0OjCJzwFTWJw4kzVNNlauK+KB4vt2ygBw+YfrsDN6Z+nrDo2BRQRZ6RAWytdBOr9gAZmfM5u31b5NWmsZF/S7i/rH389h3abyzdDVoI4GpL7K6ZCFwBVUOd8P5/L61mGqnh17+loeEsABAU7F9OclrF8Gm7yB/HdocxPvuSahxt3Pl5GPwrviJ/lGB9IgKIiUqiPQiO5vyK6msc5OaWMUl3/+VSkcl9425j4v7Xyw9x4UQ4gB1ueAura4PuraZxjFytyb3yCALZqOBb24dD/ieqz7tv79z3quLqHF6ePaCYVz51p/M3pDPssxSLhvbg9VZp7B6Ryx4rQy470owuiH9V9Taz0hd/RmpK5ycETsIBp0Nwf05efRpnP75q6zkB77b/h3UjgIgxKboGRPE1uISfiv/F6/uWEpKaAqPHPsIZ/c+G6UUqTHBaHc4Jw+IZWHpYNIty6h2VVO9W3D/4B/XvE+EEWfad6SvepMvg9Yw+LtSQEHSKJj6DJlJZ/DIiyv5d2x3jAbFX0/pQ1KEjc82f0ZA9Ha2Fg1nSXoExqAtvJ85E6vJwgenftCqwWiEEEI0r8sF96573Ac4KcYeovzBHR5oxrzHxCI9Y4J58IyB3PfFOkb1iGBk9wgm9I3m69W5eLyaY3tFEWAxsHJnD2xmAyFWEygz9D/N91VTChu+gNUfw7zHAegH/EIcd7hC+M8fD/N8XX+O7pbJRT/Z6eUJIznFwaaKOu7rfTEX9zwDgy0UasvAEszolAjiQ23cPaU/8/43BE/kIhZkL8Bel0y8zcUEVnLctmXcZdmJdVYxl0eHk2a1EJwYys22y7nirNsg2HffP3dbMeC7Bw5w7ugw7l94P4vzFgOgI2bxcrqVwO4VhFmTeH3S6ySHJCOEEOLgdLngLq1xYjYqgq1tc+r1PcujmunsdvHoZBwuD8f08j2TfUK/WH5Yl49BwZiekVjNBl6el05MiHXvEc0CI2H0tb6v2nLIWw05KyldPJfpJdk8EO/lDtNaiqyKqTVOtpiKMBoNvJpfwrEZ/wL+1Wh3/Y1WlliD4ZMQvvd6udEDv8y9h7sKKhgQ4mFMXS3BzmDeC0piZrwFk9HC/cNv5Ykl/+N5zy/0qpjKsf7g9j3D7eXrrJf456rF5FbnYjVaefiYh8kriuClFe/gUR5GRh/HG2dOaxh7XAghxMHpcsF9+0l9uOrYlH0O+9lakf4m9+jgppvelVJcNW7XcKET+/k6vQ3pFkaozcyoHpGYjYqYZt7fICAcek6EnhNZ4Did/8zdyumjfmZ+9jy8RWfzr/97BADtrEFVZEFlDjiq/F92cO72vaOKgoodjHTaWRBgZ3FiBHaDB4jAU5uEMSCbETEjePq4p0kITuCNOUFUhb/GDXNu4IyeZ3DnqDvJK6/FGvs9P+z4g5O6n8Q5fc7hlB6n0Cu8F+sCK3j2W9+gLWcfM1RCWwgh2lCXC26b2djsZBgHor6m3Vxw7ykmxMp1x6UyuJvv+fAAi5HThiSQ2MpxycE3mIrWcMOAB3AUncwGd1DDOmUJhJh+vq8WPJm1gPCAHBw8gZEUkuvOJiFhB0vcs0k1ncnbk/+JyeD79UgOSaK48j6mjF/P2+ve5ucdPxOuBmOJWsHlAy7nnjH3NNp3r9hd5Tl6j8FehBBCHJwuF9xtLbIhuJtuKm/KA6cNbPT6hYtbnoBkT/XDouaVa7yOeMICmh7ZrCUhNhM4Uph96WxuejediEAb1ww8nblvjmLy5H4NoQ2+R8I25FbwlxF/4axeZ/HGujf4auu32NyDuGvUXXvtO9Biolt4AFprkiJa/4FECCHEvklwH6SGe9ytrHG3hRT/HNU7Sqopr3ERHrj/He1CbGYKKuvoFtwNu2MryZEmhiaHMyY1khP6NX6GPSEsgGK7E4fbQ/fQ7jw67lGWLD+G7uFRjQJ+d1cc04MAs7HNbkkIIYTwkeA+SPVjnNcPVNIRwgMthAWYySypprzWRd89xiVvjRCbiW2FvsfA7HVuQqwmgq0mPrvhmL22bRiEpcJBZkk17y/OZFu+Zlxq88e98fhe+10mIYQQ+ybBfZASwgL49PqjGe6fbayjpEQFsqOkhvIaF2EBrW+mr+eblMT3TLvd4W6xl339I19v/r6d9xfvIDbEyg0TenFNE3N0CyGEaF8S3G1gz9m2OkKPqCBW7iyjotZ5wE3lVXVuPF5NjdNDsK2F4A7z3ad+f/EORqdEMOOasW3awU8IIUTrybiTnVRKVCDZZbW4PJrwgAMJbhNur6bE7gBoucYdZkMp6BUTxBtXjpLQFkKIQ0hq3J1Uj6hdj1wdaI0b6gdSaTm4g6wm3rxyFEO6hRHezBjvQgghOoYEdyeVEr1rUJMDuccd6m8az6vwDZTSUlM5wEkDWj/xihBCiPYjTeWdVMpB17jrg3vfNW4hhBCHDwnuTioyyOKblIS2aSoP2UeNWwghxOFBgruTUkrRw99cHn4ATeX1Nezccl9TeZDUuIUQolOQ4O7E6juoSVO5EEJ0HRLcndiI5HC6hQcc0ONZDU3l/hp3iLVt5icXQgjRvqSa1YlNH5fKlcekHNB762vYBVW+57iDrPJsthBCdAZS4+7EDAaFxXRgP0KjQRFsNeHxagLMRkxG+VUQQojOQP5ad2H197mlY5oQQnQeEtxdWH1wy6NgQgjReUhwd2H1HdSkR7kQQnQeEtxdWH1NW4JbCCE6DwnuLqyhxi1N5UII0WlIcHdhUuMWQojOR4K7C5PgFkKIzkeCuwsLlaZyIYTodCS4uzCpcQshROcjwd2F1Qe2PMcthBCdhwR3F1bfqzzIIsEthBCdhQR3F9bQVC41biGE6DQkuLuwoUlhXDImmbGpkYe6KEIIIVpJqlpdWKDFxJPnDj3UxRBCCLEfpMYthBBCdCIS3EIIIUQnIsEthBBCdCIS3EIIIUQnIsEthBBCdCIS3EIIIUQnIsEthBBCdCIS3EIIIUQnIsEthBBCdCIS3EIIIUQnIsEthBBCdCIS3EIIIUQnIsEthBBCdCJKa32oy7BPSqkiYEcb7jIaKG7D/XVVch0PnlzDtiHXsW3IdTx4bXUNe2itY5pa0SmCu60ppZZrrUcd6nJ0dnIdD55cw7Yh17FtyHU8eB1xDaWpXAghhOhEJLiFEEKITqSrBvfrh7oARwi5jgdPrmHbkOvYNuQ6Hrx2v4Zd8h63EEII0Vl11Rq3EEII0Sl1ueBWSk1RSm1WSm1TSt17qMvTGSilkpVS85RSaUqpDUqp2/3LI5VSvyiltvr/jTjUZT3cKaWMSqlVSqnv/K/lGu4npVS4UmqmUmqT/3fyGLmO+08p9Vf//+f1SqmPlVI2uY77ppR6WylVqJRav9uyZq+bUuo+f95sVkpNbosydKngVkoZgZeBqcBA4BKl1MBDW6pOwQ3cpbUeABwN3OK/bvcCc7XWfYC5/teiZbcDabu9lmu4//4DzNZa9weG4buech33g1KqG3AbMEprPRgwAhcj17E13gWm7LGsyevm/zt5MTDI/55X/Dl0ULpUcANjgG1a6+1aayfwCXDWIS7TYU9rnae1Xun/vgrfH8pu+K7de/7N3gPOPiQF7CSUUknAacCbuy2Wa7gflFKhwATgLQCttVNrXY5cxwNhAgKUUiYgEMhFruM+aa0XAKV7LG7uup0FfKK1dmitM4Bt+HLooHS14O4GZO32Otu/TLSSUioFGAEsBeK01nngC3cg9hAWrTN4Abgb8O62TK7h/ukJFAHv+G85vKmUCkKu437RWucAzwI7gTygQmv9M3IdD1Rz161dMqerBbdqYpl0q28lpVQwMAu4Q2tdeajL05kopU4HCrXWKw51WTo5EzASeFVrPQKoRppz95v/HuxZQCqQCAQppS4/tKU6IrVL5nS14M4Gknd7nYSveUjsg1LKjC+0P9Raf+FfXKCUSvCvTwAKD1X5OoFxwJlKqUx8t2hOVEp9gFzD/ZUNZGutl/pfz8QX5HId98/JQIbWukhr7QK+AI5FruOBau66tUvmdLXgXgb0UUqlKqUs+DoNfHOIy3TYU0opfPcU07TWz+226htgmv/7acDXHV22zkJrfZ/WOklrnYLv9+5XrfXlyDXcL1rrfCBLKdXPv+gkYCNyHffXTuBopVSg///3Sfj6rsh1PDDNXbdvgIuVUlalVCrQB/jzYA/W5QZgUUqdiu9eoxF4W2v9+KEt0eFPKTUeWAisY9f92fvx3ef+DOiO7w/BBVrrPTttiD0opSYC/6e1Pl0pFYVcw/2ilBqOr4OfBdgOXI2vEiLXcT8opf4JXITvqZFVwLVAMHIdW6SU+hiYiG8WsALgIeArmrluSqkHgOn4rvMdWusfD7oMXS24hRBCiM6sqzWVCyGEEJ2aBLcQQgjRiUhwCyGEEJ2IBLcQQgjRiUhwCyGEEJ2IBLcQQgjRiUhwC3GE8U97efNurxOVUjPb6VhnK6UebGH9EKXUu+1xbCG6KnmOW4gjjH8imO/80zW297EWAWdqrYtb2GYOMF1rvbO9yyNEVyA1biGOPE8BvZRSq5VSzyilUpRS6wGUUlcppb5SSn2rlMpQSt2qlLrTP9PWEqVUpH+7Xkqp2UqpFUqphUqp/nseRCnVF3DUh7ZS6gKl1Hql1Bql1ILdNv0W3zCvQog2IMEtxJHnXiBdaz1ca/23JtYP5v/bu2OQqqI4juPffwRGIM0uLkVjKg2B2tDUINjUktDS0tLW4ORWi5OLg5OzjuEmiIKEEJFUNLpHQfmgCMK/w7kPb/EyuD6M+/p+pnfuPe9w7nL//C/3vR88oOQCPwO+VUlbL4GH1ZwV4Elm3gSeAss91pkCXtfGC8DdzBwDZmvHXwG3z3A9kmou/usNSDp3W5nZAToR8ZXSEUP5L/obVXzrJLBe8icAGOqxzgglG7trF1iNiDVK2lTXR0p0pKQ+sHBL/58ftc9HtfER5Z5wAfiSmeN/Wec7cKU7yMzHEXELmAHeRMR4Zn4GLlVzJfWBj8qlwdMBhpt+OTMPgYOIuA8l1jUixnpM/QBc6w4i4mpm7mXmAvCJkxzi68C7pvuR9CsLtzRgqi53t3pRbLHhMnPAo4jYB94D93rM2QEm4uR5+mJEvK1ehNsB9qvjd4CNhvuQ9Bt/DiapsYhYAl5k5uYfzg8B28B0Zv48181JA8qOW9JZPAcun3J+FJi3aEv9Y8ctSVKL2HFLktQiFm5JklrEwi1JUotYuCVJahELtyRJLXIM0KoDTDpSKVUAAAAASUVORK5CYII=",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"def kalman(duration, dt, accelnoise = 0.2):\n",
" ## function kalman(duration, dt) - Kalman filter simulation\n",
" ## duration = length of simulation (seconds)\n",
" ## accelnoise = acceleration noise (m/sec**2)\n",
" ## dt = step size (seconds)\n",
"\n",
" measnoise = 3 ## position measurement noise (m)\n",
"\n",
" F = np.array([[1, dt],[0, 1]]) ## state transition matrix\n",
" H = np.array([[1, 0]]) ## measurement matrix\n",
" x = np.array([[0, 0]]).T ## initial state vector\n",
" xhat = x ## initial state estimate\n",
"\n",
" Q = accelnoise**2 * np.array([[dt**4/4, dt**3/2],[dt**3/2, dt**2]]) ## process noise covariance\n",
" P = Q ## initial estimation covariance\n",
" R = measnoise**2 ## measurement error covariance\n",
"\n",
" dimt = int(duration/dt) ## number of epochs\n",
" pos = np.zeros(dimt) ## true position array\n",
" poshat = np.zeros(dimt) ## estimated position array\n",
" posmeas = np.zeros(dimt) ## measured position array\n",
"\n",
" k = 0\n",
" for t in np.arange(0,duration,dt): \n",
" ## Simulate process\n",
" ProcessNoise = accelnoise * np.array([[(dt**2/2)*np.random.normal(), dt*np.random.normal()]]).T\n",
" x = np.dot(F,x) + ProcessNoise\n",
" ## Simulate measurement z at epoch t\n",
" MeasNoise = measnoise * np.random.normal() \n",
" z = np.dot(H,x) + MeasNoise \n",
" ## Innovation: z - H.xh\n",
" Inn = z - np.dot(H,xhat) \n",
" ## Covariance of Innovation: S = H.P.H' + R\n",
" S = np.dot(np.dot(H,P),H.T) + R \n",
" ## Gain matrix K = F.P.H'.inv(S)\n",
" K = np.dot(np.dot(np.dot(F,P),H.T),np.linalg.inv(S))\n",
" ## State estimate: xh = F.xh + K.Inn \n",
" xhat = np.dot(F,xhat) + np.dot(K,Inn) \n",
" ## Covariance of prediction error: P = F.P.F' + Q - K.H.P.F'\n",
" P = np.dot(np.dot(F,P),F.T) + Q - np.dot(np.dot(np.dot(K,H),P),F.T) \n",
" ## Save some parameters in vectors for plotting later\n",
" pos[k] = x[0] \n",
" posmeas[k] = z \n",
" poshat[k] = xhat[0]\n",
" k = k + 1\n",
"\n",
" return pos,posmeas,poshat \n",
"\n",
"## Simulate\n",
"\n",
"duration = 100\n",
"dt = 0.5\n",
"p,pm,phat = kalman(duration,dt)\n",
"\n",
"## Plot the results\n",
"t = np.arange(0,duration,dt)\n",
"plt.figure(figsize=(8,6))\n",
"plt.plot(t,pm,label='measured')\n",
"plt.plot(t,p,label='true')\n",
"plt.plot(t,phat,label='estimated')\n",
"plt.legend()\n",
"plt.xlabel('time (s)')\n",
"plt.ylabel('position (m)')\n",
"plt.title('Kalman filter performance')\n",
"plt.show() "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Ez az egyszerű példa jól bemutatta a Kálmán-szűrés alapvető algoritmusát. A gyakorlatban azonban általában a mérések függése az állapottól nem lineáris, vagyis nem írható le egyszerűen egy $\\mathbf{H}$ mérési mátrixszal. Ezzel az esettel a továbbiakban foglalkozunk, ami el fog vezetni a kibővített Kálmán-szűrés (Extended Kalman Filter, EKF) és a szagtalan Kálmán-szűrés (Unscented Kalman Filter, UKF) eseteihez."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.16"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
![](\"REK.jpg\")
\n",
"* Első alkalmazása: ember nélküli amerikai Hold-szonda (1963)\n",
"* Navigációs és mérnöki alkalmazások\n",
"* Szenzor adatok fúziója (radar, lézerszkenner, kamera, INS,…)\n",
"* GNSS vevők, okostelefonok, számítógépek, játékok, stb. elmaradhatatlan része\n",
"* Tőzsde előrejelzés, idősor elemzés, függvény approximáció\n",
"\n",
"Az eredeti [publikáció](https://www.cs.unc.edu/~welch/kalman/media/pdf/Kalman1960.pdf) közel 25 ezer hivatkozást kapott. Pillantsunk bele:\n",
"![](\"KalmanP.png\")
\n",
"\n",
"## Meg lehet érteni?\n",
"\n",
"A Kálmán-szűrés olyan algoritmus, amely valamely lineáris dinamikus rendszerben egzakt következetést tesz lehetővé, amely a rejtett Markov-modellhez hasonló Bayes-féle modell, azonban a rejtett változók állapottere folytonos, és minden rejtett és megfigyelhető változó (gyakran többváltozós) normális eloszlású.\n",
"\n",
"A fenti definíció korrekt, de a legtöbbünk számára nem teszi érthetővé, hogy **valójában** mi is a Kálmán-szűrő. Ezért a megértés érdekében egyesével sorra vesszük a Kálmán-szűrő legfontosabb sajátosságait:\n",
"\n",
"* rekurzív\n",
"* optimális\n",
"* közvetett (rejtett) állapotot becsül\n",
"\n",
"Ezek után egy geodéziai szempontból talán könnyebben érthető módon levezetjük a Kálmán-szűrő egyenleteit és egy egyszerű navigációs példán keresztül mutatjuk be a Kálmán-szűrő működését."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Rekurzív becslés\n",
"\n",
"A legegyszerűbb példa az, amikor valamely egyetlen $x$ változóval jellemezhető rendszer állapotát mérjük egymás utáni $t_1, t_2, ... , t_n$ időpontokban. \n",
"\n",
"A méréseink legyenek $x_1, x_2, ... , x_n$.\n",
"\n",
"Az állapot legjobb becslése ez esetben (Gauss-eloszlást feltételezve) a jól ismert számtani közép vagy *átlag*:\n",
"\n",
"$$ \\mu_n = \\frac{1}{n} \\sum_{i=1}^n x_i.$$\n",
"\n",
"Tegyük fel, hogy van egy új mérésünk, $x_{n+1}$, és a célunk az állapot becslése ezen új mérés felhasználásával. Ez az $n+1$ mérésből számított átlag lesz:\n",
"\n",
"$$ \\mu_{n+1} = \\frac{1}{n+1} \\sum_{i=1}^{n+1} x_i.$$\n",
"\n",
"Viszont sokkal hatékonyabb az az eljárás, amikor felhasználjuk a már kiszámított $\\mu_n$ átlagot és csak ezt \"frissítjük\" az új mérés felhasználásával:\n",
"\n",
"$$ \\mu_{n+1} = \\frac{n}{n+1} \\left( \\frac{1}{n} \\sum_{i=1}^n x_i + \\frac{1}{n} x_{n+1} \\right) = \\frac{n}{n+1}\\mu_n + \\frac{1}{n+1} x_{n+1}.$$\n",
"\n",
"Igazából így működik a zsebszámológépünk `STAT` funkciója is: az egymás után bevitt adatok segítségével folyamatosan \"frissíti\" az átlag értékét. Ez a *rekurzív becslés* lényege: a már meglévő becslést az újabb adatok beérkezése után mindig módosítjuk, \"aktualizáljuk\". Példánkban a rekurzív állapot becslés így írható fel kicsit más alakban:\n",
"\n",
"$$\\mu_{n+1} = \\frac{n}{n+1}\\mu_n + \\frac{1}{n+1} x_{n+1} = \\mu_n + K(x_{n+1} - \\mu_n),$$\n",
"\n",
"ahol \n",
"\n",
"$$K = \\frac{1}{n+1}$$\n",
"\n",
"az ún. *erősítés*, amelyik megmondja, hogy mekkora lesz a régi átlag, azaz a $\\mu_n$ korrekciója, illetve\n",
"\n",
"$$x_{n+1} - \\mu_n,$$\n",
"\n",
"amelyik az *új információ* (idegen szóval: *innováció*), hiszen ettől függ, hogy egyáltalán változik-e az előző becslés. Ugyanis ha az új $x_{n+1}$ mérés megegyezik a régi $\\mu_n$ átlaggal, a becslésünk egyáltalán nem fog változni.\n",
"\n",
"Nézzük meg mindezt az alábbi számpéldán keresztül:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAEKCAYAAAA7LB+5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAAsTAAALEwEAmpwYAAA+N0lEQVR4nO2deZhVxbXof+t0n55AaCYRG5kUQZQWsHG4oDgxqAGJPo0JRqPX6M3okGtC5kGjRH3vJg5RuYkxTtEoEVE0GJRgJIqCIII0AgJCMw/N2DQ91Ptjn24Ph73PuIfap+v3fefrPvvsYe0aVq1atapKlFIYDAaDQV8iQQtgMBgMhuQYRW0wGAyaYxS1wWAwaI5R1AaDwaA5RlEbDAaD5hhFbTAYDJpjFHUbQUSuF5H+QcthMBgyxyjqNoCIXAmMUEqtDFoWg8GQOUZRtwGUUn9VSv1nqvNEpI+IKBEp9EOuFLI8LiJ3+vCc40RkjYj0zfC6f4rIDV7JZfO8tSJyoc3x74vIEyJi6nIeYzK3DRBTKrtEpDjhuG3lb2P8L/AdpdSaoAXJFBG5CBgGXKeUag5aHoN3GEWd54hIH+BsQAETgpXGQix8KXvJegci0gt4Qin1ih+yuI1S6jWl1FVKqaagZTF4i1HU+c81wLvA48C1LQdF5EmgF/CyiOwTke8nXigi14nIchHZKyKfishNCb9/X0Q2ichGEbkh5jY5wU6ImFX/axGZBxwA+onIQBH5h4jsFJEVMV+63bVHicgcEblfRPomumfi3RAi8jURmSci/yMiO4FfiMiHsXds+SgRORer/D8tIoUicpWILEh47q0iMiNJ2h4vIu+JyG4ReUlEOsdde6aI/FtEamPPPzfut84i8qdYuu0Skemx411F5JXYNTtF5F92DZqIRERksoisFpEdIvLXlmeLSImIPBU7Xisi74tI99hvHUXkj7E8qxGRO0WkIMn7GXRBKWU+efwBVgHfBE4DGoDucb+tBS6M+94Hy/IujH2/BDgeEGAUloIdFvttHLAZOBkoA56MXXuCgxz/BD6LnV8IdATWA9fFvg8DtgMnx85/HLgT6AK8B9xpJ2PcvW+I/f81oBH4Tuy+pQly3AhUAx3i7xV7h71A/7hz3weuSvI+NcApQDtgGvBU7LcKYAdwMVZjMDr2vVvs95nAc0AnIAqMih2/G3gkdiyK1ROSxLwCbsFqfHsCxcCjwF9iv90EvBx7n4JYvneI/TY9dm474OhYut4UdBk1nzTqcdACmI+HmQsjsZRz19j3auDWuN9bK3/s+xFKMOF+04GbY/8/Btwd99sJpFbUv4r7/iXgXwnnPAr8PPb/47FnLAVuTyYjRyrqz5Kkx1bgRLt7AU8BP4v93x9LcZcleZ8pcd8HAYdiyvEHwJMJ58/C6tH0AJqBTjb3/BXwkl0aJijq5cAFcb/1iOVzIXA98G+gMuH67kA9cQ0X8GVgTtDl1HxSf4zrI7+5FnhdKbU99v0Z4twfqRCRi0Tk3Vg3vBbLQuwa+/lYLIu4hfWJ19sQf05v4IxY97w2dv9JwDFx51wClGJZmZlwhCwichzwV+BapdQnDtc9g6W8AL4CTFdKHUjzOeuwrOCuWO92RcK7jcRSqMcBO5VSu2zudy9WD+j1mKtpssNzewMvxt17OdCEpYyfxGoUno25Vu4RkWjsmiiwKe66R7Esa4PmBB6GZfAGESkFrgQKRGRz7HAxUC4ipyqlPsSyJp2uL8bqzl8DvKSUaoj5UiV2yiasrncLx6UhVvzz1gNzlVKjk5z/v1jugVdFZJxSaj+wP/ZbGbAn9v8xCdcd9l6xtJgO/FYp9VqS570OdBWRIVgK+9Yk58Lh79wLy6rdjvVuTyqlvp54gYj0ADqLSLlSqvYwoZXaC3wP+J6InAzMEZH3lVJvJNxmPXC9Umqeg1y/BH4p1kDyq8CK2N96rN5VY4r3MmiGsajzl4lYVtYgYEjscxLwLyzlC7AF6OdwfRGWYt8GNIoVCjYm7ve/AteJyEkiUgb8LEP5XgFOFJGvikg09hkuIiclnPdtLEXzioiUKqW2YfmGrxaRAhG5HsuPnozHgGql1D3JToopsBewLNvOwD9S3PdqERkUe/9fAS8oKwLjKWC8iIyNyVgiIueKSE+l1CbgNeD3ItIp9t7nAIjIF0TkBBERrEaoKfZJ5BHg1yLSO3ZdNxG5NPb/eSIyODZIuAer8WiKPfd14P+KSIfYgOTxIjIqxTsaNMAo6vzlWuBPSqnPlFKbWz7Ag8AksaIm7gZ+EusK/3f8xTHr7rtYCnkXlitgRtzvrwH3A3OwuuvvxH6qT0e42P3HAFcBG7EGJn+D1TjEn6ewBgDXAy+JSAnwdeB2rAG6k7F8ssm4CviiHB75cbbDuc8AFwLPp2F5PonlS98MlGClF0qp9cClwI+wGrr1MXlb6ttXsRRoNZbP/JbY8f7AbGAfVnr+Xin1T5vn/g4rL14Xkb1YA4tnxH47Bqux2YPlEpmL1XCA1UAXAR9j5ekLWO4Yg+a0jCgbDDkRs4SXAsWma20wuIuxqA1ZIyJfFJEiEemEZQ2/bJS0weA+RlEbcuEmrK79aixf6jeCFcdgyE+M68NgMBg0x1jUBoPBoDmexFF37dpV9enTx4tbGwwGQ16ycOHC7Uqpbna/eaKo+/Tpw4IFC1KfaDAYDAYARGSd02/G9WEwGAyaYxS1wWAwaI5R1AaDwaA5ZlEmQ0oaGhrYsGEDBw8eDFoUQ4ySkhJ69uxJNBoNWhSDDxhFbUjJhg0bOOqoo+jTpw/WekGGIFFKsWPHDjZs2EDfvhntyWsIKUZRt2GmL6rh3lkr2Fhbx7Hlpdw+dgATh1Yccd7BgweNktYIEaFLly5s27YtaFEMPmEUdRtl+qIafvi3j6hrsFbRrKmt44d/+wjAVlkbJa0XJj/aFmYwsY1y76wVrUq6hbqGJu6dtSIgidzl0UcfZdcuu01UsjvPYAgSo6jbKBtr6zI6HiZ+9atf0blzZzp16uTKeU60b98+q+sMhkwxro82yrHlpdTYKOVjy0tzvne6vm+v+NnPnDebaWxspLCwMOV5BoNOGIu6jXL72AGURgsOO1YaLeD2sQNyum+L77umtg7F577v6Ytqcrrv2rVrGThwIDfccAOnnHIKkyZNYvbs2YwYMYL+/fvz3nvvsX//fq6//nqGDx/O0KFDeemllwB4/PHHueKKKxg/fjxjxoxxPG/ZsmWcfvrpDBkyhMrKSlauXAnAU0891Xr8pptuoqnpcJfR9u3bOeuss5g5c2ZO72gwOGEs6jZKi4XrtuWbzPed671XrVrF888/z9SpUxk+fDjPPPMMb7/9NjNmzOCuu+5i0KBBnH/++Tz22GPU1tZy+umnc+GFFwLwzjvvsGTJEjp37syPfvQj2/MeeeQRbr75ZiZNmsShQ4doampi+fLlPPfcc8ybN49oNMo3v/lNnn76aa65xtp2csuWLUyYMIE777yT0aOT7dNrMGSPUdRtmIlDK1x3SXjp++7bty+DBw8G4OSTT+aCCy5ARBg8eDBr165lw4YNzJgxg/vuuw+wwgo/++wzAEaPHk3nzp0BeP31123PO+uss/j1r3/Nhg0buOyyy+jfvz9vvPEGCxcuZPjw4QDU1dVx9NFHA9ZEoAsuuICHHnqIUaPMHrEG7zCK2uAqXvq+i4s/3/c2Eom0fo9EIjQ2NlJQUMC0adMYMOBw9838+fNp165d63ellO15J510EmeccQYzZ85k7Nix/OEPf0ApxbXXXsvdd999hDyFhYWcdtppzJo1yyhqg6cYH7XBVbzyfafD2LFjeeCBB2jZtWjRokUZnffpp5/Sr18/vvvd7zJhwgSWLFnCBRdcwAsvvMDWrVsB2LlzJ+vWWatRigiPPfYY1dXVTJkyxevXM7RhjKI2uMrEoRXcfdlgKspLEaCivJS7LxvsS9THT3/6UxoaGqisrOSUU07hpz/9aUbnPffcc5xyyikMGTKE6upqrrnmGgYNGsSdd97JmDFjqKysZPTo0WzatKn1XgUFBTz77LPMmTOH3//+956/o6Ft4smeiVVVVcpsHJA/LF++nJNOOiloMQwJmHzJL0RkoVKqyu43Y1EbDAaD5hhFbTAYDJpjFLXBYDBojlHUBoPBoDlGURsMBoPmGEVtCBUzZ85k6dKlQYvhOo2NjTz44IPU19cHLYpBQ4yiNoSGtWvX8sgjjxwWknbXXXcddo5bS49+7Wtf44UXXnDlXok8/vjjfPvb3279rpTilltuobKy8rDZlwZDC0ZRG0LD8uXL+eMf/0hBweczHxMVdRgRER588EHOOeecoEUxaEpKRS0iA0Rkcdxnj4jc4oNsrUxfVMOIKW/Sd/JMRkx5M+clMw3hY+LEifzkJz/h3HPPZerUqQBMnjyZuro6hgwZwqRJkw47f9++fVxwwQUMGzaMwYMHty5lCnDHHXcwcOBARo8ezZe//OXWxZkSmT17NmeffTYnnngir7zyCgBNTU3cfvvtDB8+nMrKSh599NHW8++55x4GDx7MqaeeyuTJkwG4//77GTRoEJWVlVx11VVHPGPbtm1cfvnlDB8+nOHDhzNv3jwA5s6dy5AhQxgyZAhDhw5l7969OaSeIeykXJRJKbUCGAIgIgVADfCit2J9TqZ7+xk85rXJsPkjd+95zGC4KPlaGY899hidO3fmwIEDVFVVcfnllzNlyhQefPBBFi9efMT5JSUlvPjii3To0IHt27dz5plnMmHCBBYuXMi0adNYtGgRjY2NDBs2jNNOO832mWvXrmXu3LmsXr2a8847j1WrVvHEE0/QsWNH3n//ferr6xkxYgRjxoyhurqa6dOnM3/+fMrKyti5cycAU6ZMYc2aNRQXF1NbW3vEM26++WZuvfVWRo4cyWeffcbYsWNZvnw59913Hw899BAjRoxg3759lJSUZJyshvwh09XzLgBWK6XWeSGMHV6ub2wIDw8//DCvvvoqhYWFbNmyhZUrV9KlSxfH85VS/OhHP+Ktt94iEolQU1PDli1bePvtt7n00kspLbVW8xs/frzjPa688koikQj9+/enX79+VFdX8/rrr7NkyZJW//Xu3btZuXIls2fP5rrrrqOsrAygdUnVyspKJk2axMSJE5k4ceIRz5g9ezYff/xx6/c9e/awd+9eRowYwW233cakSZO47LLL6NmzZ8ZpZsgfMlXUVwF/8UIQJ/J5b79QksLy9YK5c+cya9Ys5syZQ1FREaNGjeLgwYNJr3n66afZtm0bCxcuJBqN0qdPHw4ePEgma9sk7vQtIiileOCBBxg7duxhv/3973+33Rl85syZvPXWW8yYMYM77riDZcuWHfZ7c3Mz77zzTmvD0cLkyZO55JJLePXVVznzzDOZPXs2AwcOTFt2Q36R9mCiiBQBE4DnHX6/UUQWiMiCbdu2uSWf4zrGbqxvbAgHu3btomPHjhQVFVFdXc38+fNbf4tGozQ0NBxxze7duzn66KOJRqPMmTOndWnSkSNH8vLLL3Pw4EH27duXdPus559/nubmZlavXs2nn37KgAEDGDt2LA8//HDrMz/55BP279/PmDFjeOyxxzhw4ABgLYfa3NzM+vXrOe+887jnnnuora1l3759hz1jzJgxPPjgg63fW9w4q1evZvDgwfzgBz+gqqqK6urq7BLPkBdkYlFfBHyglNpi96NSaiowFazV81yQDbDWN473UYN/6xsb9GDcuHE8/PDDVFZWMmDAAM4888zW32688UYqKysZNmwYTz/9dOvxSZMmMX78eKqqqhgyZEirNTp8+HAmTJjAqaeeSu/evamqqqJjx462zx0wYACjRo1iy5YtPPLII5SUlHDDDTewdu1ahg0bhlKKbt26MX36dMaNG8fixYupqqqiqKiIiy++mF/+8pdcffXV7N69G6UUt956K+Xl5Yc94/777+db3/oWlZWVNDY2cs455/DII4/w29/+ljlz5lBQUMCgQYO46KKL3E/YPCTojZW9Iu1lTkXkWWCWUupPqc51e5nTfE38sJBvy2nu27eP9u3bc+DAAc455xymTp3KsGHDghYrY/ItX3IlMfAALKPOr/XQcyXZMqdpWdQiUgaMBm5yU7B08WJvP0Pb5cYbb+Tjjz/m4MGDXHvttaFU0oYjyefAg7QUtVLqAOA8xG4whIhnnnkmaBEMHpDPgQdmZqLBYMgL8jnwwChqQ1p4sWWbIXsyyY+2MrM3yI2VvSbTOGpDG6SkpIQdO3bQpUsX21hhg78opdixY0dasxUzmdkb9kH7FlnD/A5OmM1tDSlpaGhgw4YNKSeZGPyjpKSEnj17Eo1Gk543Ysqb1Nj4aCvKS5k3+fzW72GPmMgHco76MLQdnKyqvn37Bi2aIQvSHWALKmJCByteBxlSYRS1oRWzAFb+cWx5qa1FnTjAFkTEhA7lTQcZ0sEMJhpaSWZVGcJJugNsQURM6FDedJAhHYyiNrSSz3GobZWJQyu4+7LBVJSXIli+aTu/cxAREzqUNx1kSAfj+jC0km432RAu0pnZG0TEhA7lTQcZ0sFY1IZW8jkO1ZCaiUMrmDf5fNZMuYR5k8/33EerQ3nTQYZ0MBa1oZV8jkM16IcO5U0HGdLBxFEbDIQjRMuQ35g46hBgFEVwhCVEy9B2MT5qDWhRFDW1dSg+VxT5uiaDboQlRMvQdjEWtQbk8zq6YSAsIVq6YHp//mMUtQb4pShMBbMnLCFaOmDcRMFgXB8a4MesMONecSYsIVo64NT7+8WMZW1iKdWgMIpaA/xQFMYP60y6s/cMzr282roGYwR4iHF9aIAfsZzGD5scsy9neji5iRIxYyzuYhS1JnitKHTxwxo/eW4EnX63jx1wxLrVThgjwD2M66ONoIMf1vjJc0OH9LNzE3Uqs9+8wAzGuoexqENKppZVNu4Vt603E4aYG7qkX2Lvz2l3GDMY6x5GUYeQbEOkMnGveBGGZfzkuaFr+oVlvYwwYxR1CPHDsvLiGbr4yYMi1x6KDunn9A5mMNZbjKIOIX5YVl48w24gqq10kbPpoSQqxfMGdmPawprA0s9MdjkcPwd2zWBiCPFjgowXz2jL8cqZxrHbDRxOW1jD5adVBJZ+Jhb/c/we2DUWdQjxwzL16hlttYucaQ/FSSnOqd7GvMnnuy5fOujqIw8Cvwd2jaIOIX4M3pgBInfJ1L+so1LUwUceJPGuDqdV/L3KH6OoQ4oflmlbtX69INMeio5K0YwxpJ7o41X+hF5RBz1Tq61g0jk3Mu2h6KgU23Ivy87VkYiX+RPqrbicAu0TB1iMksmNdNPZ4C6m3PpLsvTuO3mmo7tDwJX8ydutuNJx6JuQotzRZUZcW8Nt15NR/M6k0hNOrqiK8lJfBnfTCs8TkXIReUFEqkVkuYic5bVg6ZBswGX6ohpGTHmTW55bHEhIUcvz82F9Xh0HtgyZocM6ITqTKvTQbq0cwUpHP+p3uhb174C/K6X+j4gUAWUeypQ2Tq1cx9JoSse/l0om36z4ZANbxkrTD7s8Mb2i5KQyRuL98zW1dQi0ukL8qN8pLWoR6QCcA/wRQCl1SClV64k0GeK0IpwIgY3OQvoTA8JidTul83kDuxkrTTOcLGenNaRNr8ginQleE4dWMG/y+VSUlx7hr/a6l56O66MfsA34k4gsEpE/iEg7zyTKAKeZbrUHGpJel2p0NlcFmo6rIExdUad0nlO9zcxU0wwnI6FAxPZ8v8L9dDdKMlkGOAhXYDquj0JgGPAdpdR8EfkdMBn4afxJInIjcCNAr1693JbTEbsBl5buiR0VKbrnbrgt0omBDVtX1C6db31use25xkrLnlxdSU5p36QUpdGCQML9wuAKzCT0MIgY93Qs6g3ABqXU/Nj3F7AU92EopaYqpaqUUlXdunVzU8aMcWodf/ulIcybfH7SwuHGegbptM75MEDnx5ojbQk3ellOad/SCwpinZCwrBHS4tpYM+WSpHoiiE04UlrUSqnNIrJeRAYopVYAFwAfeyaRC+QSmO+GAk3n+TrOPMsUHSdlhBk3elnJ8iSomab5YJTEE8TEn3SjPr4DPB2L+PgUuM4ziVwim0I5fVENERGabCYBZapAUz0/H5Rcy/v9YsYyauuscYGSqFmQMVv8MhL8JkijxKuoJL8bvbQUtVJqMWA7YyZfaOl22ilpLxSojhUqW+obm1v/33WgQTv/Y1hwS6HptkZLUEZJrr5xnUJPQz0z0U2c5vIXiGTty0uV0bpVqGwI26CozuRDL8uOoIySXMqmbgOgRlHHcOpeNiuVtZLWKaO9It/8j0HipNAARkx5UwvLLluCMEpyKZu6GSB5p6jT6a7YneO2H023jPYKp3SLiNB38szQKpagSLXDd742+F6QS53WzQDJq5GfdMKbnM45b2A3V0NudMtor7ALVQIrblf3iTxhICyhbTqSSxidbqGneaWo0ynUTuf8Zf56V/ej0y2jvSJx1qLdDDijWLKnrTT4XpDLHp1BxEonI7SuDzv3RTqFOtnMrWkLa1ybBJCvA0N2xHfX+06eaXtOTW1dm3WF5BI9oHu8vU6REXZk6xvXLSorlIrayW9XXhZll806H/GF2qngg7s+ZDcyWvdKYEey9I13hUDb8LHm6mPOtcH3sgzlu/9cp6isULo+nNwXSpGyu+LkU23BzS5lulNS7QjTok3xpEpfaFuukFx9zLl0370uQ8Z/7h+htKidlOnuugb+50tDUsYuA3zvrx+6MgPRK8IaNZLYk/B7t2bdcGu2YTZ57nUZMv5z/wilok7mt0unULf8rrMPOcyVID4PRkx5U2sfq9ek447zCq/LkO7+83wilK4PN0Zkc+lS+oGbUSNBrgWs2+i5n0xfVMO+g41HHI8WiC/v73XkUVvOW78JpUXt1oisToMFibgVNRL0gI9uo+e5ksng3L2zVtDQfKTzp11RoS/v73XkUb7lrc6IsvHT5kpVVZVasGCB6/fVEa9H1XO9t5Prwa/dk/OJxEYPLMXn1BPrO3mmrY9egDVTLvFO0DjCGDnUVhGRhUop28XvQmlR64LX1qobFn8YfN1hUSaZDs7p4MPVudfoNWEpV+kQSh+1LmQbnuSnz1j3GZJhCkPMtNEzPtzgCFO5SgftLWqdW8VsrFW/fcY6zZC0y0svQsi8KjOZWsjGh+sO2eRnWMNbndBaUaer1IJS5tl0bf0uQLooC6e8tFsDHLJ3zXjZEGbT6LVl14MbZJufYXD5ZYLWijodpRZkVEM2FTeIAqSDsnDKywKXtj5L9Rw3GkJdGr18IF3jKtv81GF8wE20VtTpKLUguzjZVNx8K0DpkmwxrNJogWuuGS8aQp3db7kSxLtlYlylm5+J73HewG5MW1ijhcvPDbRW1OkotaC7OJlaqzr5jP3EKS8r4nzVbigLtxvCoOPQnWRyY7Gvmto6BFpDCOPfDQ7ftLhTWZSfjz/ZlXfOxLhKJz/t8mjawhouP62COdXb8qKB1VpRp6PUwmahttXuc7K8dNM143ZDqNuglBsbtsZfn+h0qmto4hczlrG/vvGwyTq7DjRw+wsfpv2cZGRiXKWTn055NKd6W97MFdBaUaej1MJooergM/YbLxqoZJalW89JpVT8dh3k2nA4beIcT4sVnUhDk3KlgcrEuEonP4PuVfuB1ooaUiu1tmqhhhE3G6hUlmX8YPO9s1Zw63OLsyobyZRKEG6RXJVSrsrLDeWXqXGVqtzo0Kv2usHWXlGnQ1u0UNs6bkYEJatkyZRKEG6RXJVSso0dwHq3kmjEdsW/TJ6TDLeNq6B71X402PkxM7F+Lyx6Clb8PWhJDD6Ra0RQC6lmsCVbZTGILneusx3trm/Z5bLl3X4+/mSikSP3vnRz1b9cN9WIn9kLpFwJ08vZwH5soKCtRZ12V2LnGph6LhyshbKucPsqsNlgNefnGLTCrYigdKxipx5bEF3uXK3RTK73KuojF5ys17svG+w4cOi1xetHg62lorZL2F/8bSG9173Aul5f5N7XV7UWssd7vkz/+r0w7Fr44M+wdxN0ODbr5wQdemVID7cignKpZEF1uXN19aVzva7uxGzcTV67qPxosLV0fdgl7Iim9xm6+Ge89uLTrd3UrbV76bRyGpu6j4Ihk6wTN1khRChluUP2bcvoOWbPt3CQzCUBViO8v/7IRfvtlLkd6VSyVDIESZCbRXhJNg2r1xavH4tvaWlR2yXgMbIDgNPVh8yiEoALIh/QVXbz3zvO4r7uJwNiKeoBF8HSafDSt+DMb8G4u9J+TrLjBr1wsvrs1o0G++57rlZxtpanncsN3Blgy+eeYjbWq9cWrx+RZ1oqaruE7S61AIyMtMycUlxX+Hc2qs68uPck7ituD137w6YlcGg/vP5T67RlL8KYO2HRk9DvXOjUO+lzWo4bwotTrHCZzc4qQYR32inS25//EMSKVW45lq1y1W2Sjptk07Amu8atMSqvXUVaKmq7hK2I7AJgQGQDR7OLoZFVnBGp5icN13FMeTvrpB6nwrp3YO49sHcjDP86vP+/8PpP4N2HYMjVMPGhpM8pjRZw+5gT4YMnoXYdnHgRVAzLaIDSECyZ9pTctIqznXRit2VXtso1n3uK2TSsTtcAoel5pKWoRWQtsBdoAhqdtotxC7uEPaPdIer3llN8qJYvFrzNlwve5JPmCqZHRnNnS2va41T46Hn49/0w9Ktw4c8tS/rdmHKufgWafgsF0SOeU1NbR4EIdQ2N7J/5Y2iabl3z1r3QdxSM+C70GArtuuT0bibKxHv86Cnl4l7IRGFmo1zzoadot8hSLut22DXGI6a8GZqeRyaDiecppYZ4raRbSIyz7Nq8g+L+51Ff1JkfRv/CcbKV3xffwJ2XDfk8UXucav3teByMuxuKj4ITx1rHRt5qhfB9OveI57QMBjSrJu4sfIxJTdP5S/NoZl70bxg3BTZ/BE9dDvceDzO/B/X7snonL3adyNdBo1zwY3Anl4HoTBRmNso17DvL2NWTp979zPXdWsLU89DS9XEESsHezdChguLx98KOVRSccjm/7dr/8PMqqmDAJXD29ywlDTD6V1D5JTjhQnjvD/Dxi9D/ws+vWfYiha88weOykQ5FBzgp8hmPNI5nSuNVVMzZwiWTvwFDr4YNC6B6Jrz/B1g7D67/O5SWZ/QabvsO83nQKBf88Dsnq+Spek23jx3Arc8ttt34Np5slWsu759Lj8+t3mI665G4YfmGqeeR1i7kIrIG2IW12NajSqmpyc53fRfyulr4TW8Y82v4j29nf59pX4eVr8Oty6CoHcy5C966h+2qA6vVsSiEWU1V/KlpHCD2u0Wvmg3PXAW9/wMmvQCFRWk/3nlXasWaX4+Df/wUGuvh4nshUmBzpkX8MpV2eLbDeGM9RKIQ0TKq01ecdncvL41S39iccqfyPpNnOt5bIBC3WKa7rGd6bbqK3KmeJJLrbu65vK8XuLEL+Qil1EYRORr4h4hUK6XeSnjIjcCNAL169cpJ4CPYu8n6e9Qxud3n9Bvho7/Cv+6zFM5b98DQq7ns44l8tvvQEafbtqwnXAgTHoDp/wV/Hg/DvmqFAp75Teg/+sjzleLV+ct4YM4qoIB4b9M5kQ+5o/BPlEaa4ImTYN3b1g/NDTDiFsvl8vFLoJqg6wA447+Y/snBpFtYQeo9G1NWlk0fQmknOLgH3nkICgotJb3sReh8PIy5w/5d02XHatiyFNp1s9xVRe2yv1cO5GIBOg1Ei5BWr6kiyfrcQS3NmUuPL9W1mfT+Uq1HEn9eLoRpQbe0FLVSamPs71YReRE4HXgr4ZypwFSwLGpXpWxR1GnOOHTkuOFw6ldg3v2W8ht6NUx4kNt6bcws5GfIl60ByZdvtmK1C4qsaJOvvQI94xrEXWvZ8pdvcfHWt7kY2FDcldlNwyigmdMj1QyIbGC1Oha69oN186weQ90uqyH54AnrHu27Q3EH+HgGvPswK5qvpq5hJABfKXiTKwrmUkwDbzefwkONl7Kb9o4FOLGybKrdT/TFG9g3u4b2o74D3U+GhX+GJc9+flFxB4gUQtMhOPUqy+3z9P+h+sRv8J/rRrN3906+1H4xX+/0AUfvWQoDL7EU8KYPoeI0KySyy/GweSlsWgw1C2HlP2hdCbmwFDr3g11rLBdXWRdLjh6VLNjfjReXbOf9fV3Z3+EEbh830LfV91LhVMlvfW6x7fmJjaedogc4cKiR6YtqAlEWGftslbKioQ7tp/PuZYwtWEFX2c0hCpnffBLVzcexp/YAfSfPJCJCkaqjA43UU2R9Ghp44bXXmXioAZrqobAECov53aADPPH+FvY0FVBPEc0qQqnUU8ZBIig+UT05WNiBH599NDQ1tAYHANBwEDYvgaL2UFgM+7fDge3W37pdlgFy1DGWcdDlBGuQccixlvGw4X3Y8Ta8AbQ/Bo4+ySr7+zaDFFj32bXWqhON9bD9E9i+EhoPWmVYBKJlcOWfXc+blK4PEWkHRJRSe2P//wP4lVLKcQUk110fi5+B6d+A7y6yEiQX9m6B358JPYfDVU+3ZnIq68r2975NsGudFb/92FirIEx8BAZeDCteg2k3cOBQI1MbL2avKuPsyEecEVlOPVFWqON4u+hs+o/7JhOq+lmx30XtQCnenj2dN95bxOr9pSwtGoKKFNC1bg13lT7N8OYPmd88kCIaGRpZxUfNfdiljmJkZCkHKGaWOpNep49neNUZVoHqUGFZxMA37nqQK+v+SqkcYpPqTJRGvlAwn7VSQR8VG5iRCIy8zVK2jXVw2tegpNyqlJEINB5i3RM30fuzv3FAFROlkag08ZnqjvSs4rhtc61C3PVE2FZtNYitiBXHPvgKGHAx7N9muZJ2roEuJ1hy7t0Mmz9CbfsEobn1yhrVhZfVOfQb8w3GjDwjtzKAs+siV4s2k/tOX1Rz2HoaLXjS/a7fB5/OsRrK2vVWeevcz1Jch/ZBv1Fc89xaTt4/n3LZSxn1lEk9HThANBrl3Asuhk59rQZ73b8tw2JbtWWkNH3eG21QBRTQTEQ+1yvbVQeaiLTOhQCoU0UU03DYeVkhEejYEzr0tMrrthXQcCD964s7QHMTNOyP3S/mclQOPdZIITQ3Ws/t1Mcq54XFlgKXiGVYfeW57F4liesjHUXdD3gx9rUQeEYp9etk17iuqN+6D968A368GaIuOPoP7rFa3DR9rWn5snathb9eA5s+ZDud6Kxq+STSj+vrbmYjXY+4p5N/zWlWnXVNM98oeJnxBe/QSIRpTefw56YxKCKcKOv5btksxsp8oo37P7+oU1/Lr1+ziOZFT7OFTnymjuZ42UhX2cPvGydwb+OXWHNLH0txdqiAbicmTY8Rd7/Bf+z7OyfKBg5SxOymYXyojqeivIx5t50JKKvRqdtlKYYdn1rWScVpUFSWRorDeXe/RmTPeopp4JTIGsZF3mdU5EMKREGfs63Y9kjUklnEWpCr+yDoPdJqZPZthn1bYP8OK9qncz84prK10bLzg3ZiD0Miq/nTD2/I2s2WsqzU1VqWY7TESssjFLuil2zl2KOiPDv5K5aLaMMCy+Jr381SEp/OhT010FBnleNOvS3XWLcTrfyOllo9mn1b4eBu2PkprH3bUmSRQijvZfVkdq62rEEpOEwxHVDF1FFEHcXso4ye7RXt96//XMSi9tDrTMtt1dwIJR15b3dHbp9fxrqGjrSjjhGRpfSQnRRziH6yiQKa+VQdy0GilHKIctlHHcXsL+nBj779DSgt59VFa/nDnI+p3bOXnkdFuP7MHpzbr4P1jGhZzJBpgq3LrYZGBHbXWL2xPRut9+7UF/qNsiztpkNWuWjXxfpbWg4Hdlplpn6vdZ9da6w06Xoi9DrL6v1JAezfaqU9wFE9QDVDSUerQWiqt5RyYXFWZcSJnBR1NriuqGf+txUfPXmde/fMgHStpBkL1vDJjHs4tnkTuziKBxq/SD3FtgMjTpab07Piid/nDhIUQWM9bFvBewvf55+LP+HSQ68wILKBxoJSXmIUP9t/BfspRWiml2xlnepORXlZRvsWOg+K5ja4k+oZPdjBFQVzua1iOexYCc2NHCzqxP76JjqqPRRKzAIvLLWUkh1SAN1P5oXN3dnXaE1iKqKB7lLLyMhSiiVm3XbpD8edYVXuQ/vhwA6rcWh3tDURqiBqVda6WqshkIhV2dt35931dcxY9Blr90cpb1/G109qYGjJZqvR2rDAUigVpwHwwaebKKGBYg5RLA10YD8dJCZ7i/UGlnJvPBhLiFMt+aKlUL/H6pFsX2n/zsUdLEXe83Q45TIrMirWSNDcZN1fKVg5Cw7u5rX6U7lz7o4jy8H+7ZYLUik4elBrgxdPfK8zHa0SX251G9gLAjcGE4Nl7yarVfOZVNEVib6738xeQ82hLxxxnp1idfJ/pxPDqbAUva1SLSxm+uYu/HB+BXUNx/AoIxgo69lQeBzjT+tL88IaaGhCEWGdOobSaAHnDeyWkb/Wi7CmRNdSeVn0iMXrN9GFaUdN4rZvng9Njby0uIbJ05dT19BEEQ2cKOs5P7qUCb2LOWFgpVVmyrpaoZpbP2b50g9YuGojx9es4PzIO6293EYi7FAdeb75XPqOvJIR7TdbFujKWZaSjpZaCm/5jJgkcTlaWGK5hpoOQd1OAM6MfSgCDgEfAtF2lsU/6vuxnsYHUBClsbAdNQ2WL7ZeRdmvSliuenFUaQk/Gi5wzGArwqjjcZZ13NxkP+mquRl2fwa1n1kydz/Fcgkkm1EbKfg8umjQpQBcBFx0ls257bpanyTETypxMjgKRGhWynbrtLBMPgmC8CjqDv4q6mQuiBYSFZOTkk2qWG3umcqiTuVHjS/0TRSwTPWBBphTvY27Lxt8hOWcaSVxe3lPu4G9aESIFkjr2hdHPKOgkHv+sbr1mkNEWar6sfRQP6bVlDLvqwn+4I0d+eHH7ZPmZ3lplMVjx1hfRnz3yBP2brH8ueW9rS64aj7cFdfSrW44YFnDdbUxf31/S9HauNo2OlmSlwyGxLRPFrcfiflMO/VxPsdHnMqIk4UcpsknQaC/om5qgK3VMHSSr49NFXRvp5iclGwmA1RO0QDJnptIskJvN5U23UiFFtwOa3Ja+6K8NEq74kJXNjVNZxLFbodNXVs5qjvQPfbFpuqUdbY+GRCmELFMyPS9dJ58osOyD/or6o2LrBHZPiP9fWySlrzCIbPcsDQTC3jH0igiUHugIe1Ckmmhz6aSuLlamFNa765rYPHPxzhel4nc6Vhmidf5VUF1XaQ/VzJ5r6D3PXRCl9m/+ivqtf+y/vb2V1FnYx27ZR3lWnEzLfRBV5JsralM5E5nU9f463SpoDrjZkNmV3fOG9gtpx3k3UAX37k2itox09f8C44+OedV6zIlW+Wlg3WUaYMR1JrM8b2GpP5oBzKR2y4/W4YE7XpIulRQXfGiIYuvO7o0lLr4zrVQ1E6ZIs2HuHT9fGvJUp8Ju+8w0wbDzwYmMb9r6xqIRoROZdGMXDyQvtyZ5qcuFVRX3GzI7Iw0LxvKTHoCuvjOtVDU985aQXNDHd8v/BsLm/vzRvNp1DU08eqs17i04QD0PTsQuXSwjvMRp8HDsqJCFv3M2SedK5nkZ7IKqsPgEgQ7yOVWQ+ZkpDkN/ObaUGZqqQftFmxBi2XQNtbWUU+Uywve4tKCf7ce77s/tlFtr/8ISDKDF4TBWnVa07kl5tzttZEzxYu1zTMhl02B43GynAsc4r9ztWQzXUd84lA9NjDWQlFbiS/Maz6FsyLLaJlMMKJ4tTUDy2f/tMFb3KrkXuJUQedUb9Ni5/pcNi5wA7c2J3BqnJuU8mTzg2yMhMRNTILoPWmhqFsy/d/NJ9NN9jBA1lMajXB64SprGq8hbcKw40tYdiCxq6CpNgzwK+2D7pW4ZWk6Nc4t93Pbkg2DkWCHFj7qlsT/89/3QP2jXNJuBSeNqqT4jV1w3OmByqaLPzIdvBgp9+L9wzxQ6+S77lga9TVKQYdBLjfGcJL5gL0YI9LF55wp+i3KdP8wa8nLkydaS5t+cz4cPdBV+dIlbAvFuL10Z9je3w+c0qQkGjlibRLwbiMAv/PGS4PFb2NIV+MrXIsy9TsXljzXunwiXZMvueklYYuldbs7HLb394NcNwzwWg6vlLSXvYXE+GmvJ7mEMZpLP0V92rXUL/kbxavfYE7Tqfzknn8G1uIF7QfMFLe7w2F7f7+wq+hOqyx66YrwS+H41WDrMslFR7QYTIxn+uaujKq7j981fpGHGycEFv4E7g48+DHQ5PYgXVADL2EYEE3ELu0FS9mE5R2c8KvBDjqSRWe0U9T3zlrB5oZS/qfxCt5TJwHBZZZbis+vmFe3Yz6DiM4IOj44W+LTHg5fgzws7+CE2w22U0NsenDOaOX6mL6oJu1F+v3ALT+gn75eN7vDbrx/pgM3YfaLt6S93aBuWN7BDjcjJZK5N3SIZNEVbRR1SwY6EVRmuTHQEWZLIRfFn43PMZ200nXUvoUw57cdbg5cJmuIwxo65wfaKOpkC7vrkFm5DHSEzVJwSxFmYx2nSqswDDiFLb/Twa2eWqpNLSDY+HpdjQBtfNTJrA0d4nZzGegIy0w8cNdHnI1lmSqtwjDgFKb89ptU/u4gp2unU/aDGujWRlEnm0oatJKG3LqzqQb5dIpycFMRZjMIlSqtnNK7JbpChzQMaiEfncqREzo3YqnKfpAD3dq4PnT3T+XanXXqOurWlXfTv+rF5gtO+dASCgfBp2HLc4OwBnUpR0744d7I1n2RquwHOdCtjUWty3KCTnhlCQTVlXeyvtwMxfIiT53ilRMXQtDNHeI1YXAJteCleyMXqzdV2Q9ykFgbixpSWyFBOvq9sgSCyPxk1pfbPZtMLctUeWyXDzqFdAaFSQOLXKzeVGU/yEFirRR1MnTo2nnRnXUz89NtyJIV5pYFhIJoENPN48R8cFqMKsxRFpkwfVGNba8C2kYaxJd7pyXm0h1LAueyH6R7NjSKOswTIZLhVuZn0pClsuKDWrQm2zzWYXwjyN7evbNW2CooAW3GeLzCbhVBO3IdS2r5DYIxYkKjqHWbROBWxQxi9qOucb7Z5nHQ8bdB9/ac0kf59PwgSTb/ogU3G+2gjJjQKGqdlIvbFdONzM9EyelggdqRSx4HuXRl0L09p3SraANuj2SNuMBhjbauk1nSQZuoj1ToFH+p4wh7JtEafkTYZBPTq1MeZ4Ifvb1k6RnWdHODZPMv4qNKwrrYVwuhsaiD7t7Go5sbBjK3kr20QLPtceiUx5ngdW8vVXqGNd3cIN1yH3SvJ1dCo6hBn50ZdHLDtKBTZc2lUuiSx5ngtSspnfQMY7q5QbrlXkfjKhPSVtQiUgAsAGqUUl/wTiT9caqY5w3sxogpbwa2r5wuldXNShEGv6LXjWTYlYzXpFPudTSuMiETi/pmYDnQwSNZQoNdxTxvYDemLazxZOQ/6KiCTEm3UqRSwm69tx/K3stGMuxKRgd0HUBPl7QGE0WkJ3AJ8AdvxQkPidNg51Rv82yAUcfBy2SkM7iVzuCOG+8d9kEkaNuDhW6h+xIVqUjXov4t8H3gKKcTRORG4EaAXr165SxY2PCyexq2rm86roB0/K5uvHfYB5EgeXpm21sIg0vJbXRxDWZDSkUtIl8AtiqlForIuU7nKaWmAlMBqqqqnGZy5i1edk/D2PVNVSnSUcJuvHc2yl5HJWaXntm6hvx0pemYlmEkHdfHCGCCiKwFngXOF5GnPJUqhHjZPc3Hrm86cd9uvHemqwGGyVWSrWvIL1ea12kZhvW33SKlolZK/VAp1VMp1Qe4CnhTKXW155J5hFeZ66UPLNW9w1hg01HCmaapXTo4LYvastFAYlrpMh6QTp5m6xryy5XmZVqGqUF1g1DFUeeK110+L31gTvd2452C6J6mG9KWbpo6pcPdlw3m7ssGc++sFdTU1h22ypxdWumwXGi6eZqta8gvV5qXDUI+jD1kQkZTyJVS/wxzDLUu1pKb5PpOQVomiZEzuVSwVBV33uTzqSgvTbrBQMtyoXb4OR6Qbp5m6xryy5Xm5iYUieg0wO5HjzY0a324gU6Z6xa5vlO+NF7ppEM6Wy3psFxounmarbvNr1A1LxsELxuBTPDL0GlTrg+3unw6jWTn+k75MoswnXRIdY4uy4VmkqfZutv8CFXzcsamLhNY/HLBtCmL2o0WXrdBjFzfyS3LJOh0SScdUp2TbCU2Pwds8ynKx033VuJ9dZjA4lcvvU1Z1G608LoNYuT6Tm5ZJkGnSzrpkOqcZGu4+DmFX6cFtnRGhwksfg3MilLuz02pqqpSCxYscP2+OtB38kxHP+aaKZf4LU5S0nVFuOGyCFO6JMMuLVoiRhKpKC9t3WPS4A86uR1b5LFr3LOx7kVkoVKqyu63NmVRu0FYZglmErbnhmXiR7oEtbjSrc8ttj03zIPQYUTHxcn86v20KR+1G4TFf+h3NIfX6RKkD1yXCIO2jq4RSl754eMxijpD7AYxLj+tgntnrdBqZqDfoYheD+4EWUnD0jjnO/kYXpsuxvWRBfHdYx27YxCMi8bLwZ0gK6kZ3NMDP8u0br5wo6hzJOhoByd0iTN1i6DHBtxqhHRTAGHCrzKto/FlFHWO+GnpZVLJvbACg1Qy+dDw6KgA0kGXxsWvno2OxpdR1Dnil6WXTSV30xURtJLJB/eDjgogFUHneyJ+xE7r6AvP68FEP2aT+TXQFPSId9DP18WqywUdFUAqgs73INAxyidvFbVf4Vx+TWUNupIH+fygp6e7hY4KIBVBlzuvSGbE6Rjlk7euDz+7mX50x4IeTAvy+WF0GdgRRj970OXOC1K5c3R0s+WtRZ1vlkDQrXyQz8+XvNRlIaFMCLrceUE67hw/JrFkQt5a1PlmCQTdygf5/HzKSx0WEsqEoMudF4Sx4c9bRR3GbmYqgq7kQT0/H/MyTARd7twmjA1/3ro+wtjNNNhj8tLgJmF055hlTg0GQ5tDx3BPs8ypwWAwxBE2d07euj4MBoMhXzCK2mAwGDTHKGqDwWDQHOOj1gwdBzkMBkOwtBlFHQYFqNtKZQZDNoShroWNNuH6CMuiPm1xpTJDfhGWuhY22oSiDosCDOPUVoMhnrDUtbDRJhR1WBRgGJfBNBjiCUtdCxttQlGHRQGGcWqrwRBPWOpa2GgTijosCtCsaWEIO2Gpa2EjZdSHiJQAbwHFsfNfUEr93GvB3CRMSzWGbWqrwRBPmOpamEi5KJOICNBOKbVPRKLA28DNSql3na4xizIZDAZDZuS0KJOyNPm+2Ndo7OP+knsGg8FgsCUtH7WIFIjIYmAr8A+l1HxPpTIYDAZDK2kpaqVUk1JqCNATOF1ETkk8R0RuFJEFIrJg27ZtLotpMBgMbZeMppArpWpF5J/AOGBpwm9Tgalg+ajdEjATzNRVg8GQj6QT9dENaIgp6VLgQuA3nkuWIWadDIMhO4yBoz/puD56AHNEZAnwPpaP+hVvxcocM3XVYMgcszZHOEgn6mMJMNQHWXLCTF01GDInmYFjrGp9yJuZiWbqqsGQOcbACQd5o6jN1FWDIXOMgRMO8kZRm3UyDIbMMQZOOMirHV7MOhkGQ2aYtTnCQV4paoPBkDnGwNGfvHF9GAwGQ75iFLXBYDBojnF9GDLGzGQzGPzFKGpDRpip+gaD/xjXhyEjzFR9g8F/jKI2ZISZyWYw+I9R1IaMMDPZDAb/MYrakBFmJpvB4D9mMNGQEWYmm8HgP0ZRGzLGzGQzGPzFuD4MBoNBc4yiNhgMBs0xitpgMBg0xyhqg8Fg0ByjqA0Gg0FzRCnl/k1FtgHrMrysK7DddWHcQVfZjFyZYeTKHF1ly0e5eiulutn94ImizgYRWaCUqgpaDjt0lc3IlRlGrszRVba2JpdxfRgMBoPmGEVtMBgMmqOTop4atABJ0FU2I1dmGLkyR1fZ2pRc2vioDQaDwWCPTha1wWAwGGzQQlGLyDgRWSEiq0RkcoByHCcic0RkuYgsE5GbY8d/ISI1IrI49rk4ANnWishHsecviB3rLCL/EJGVsb+dfJZpQFyaLBaRPSJyS1DpJSKPichWEVkad8wxjUTkh7Eyt0JExvos170iUi0iS0TkRREpjx3vIyJ1cWn3iM9yOeZdwOn1XJxMa0Vkcey4n+nlpB+8L2NKqUA/QAGwGugHFAEfAoMCkqUHMCz2/1HAJ8Ag4BfAfwecTmuBrgnH7gEmx/6fDPwm4HzcDPQOKr2Ac4BhwNJUaRTL1w+BYqBvrAwW+CjXGKAw9v9v4uTqE39eAOllm3dBp1fC7/8X+FkA6eWkHzwvYzpY1KcDq5RSnyqlDgHPApcGIYhSapNS6oPY/3uB5YDO63leCvw59v+fgYnBicIFwGqlVKYTnVxDKfUWsDPhsFMaXQo8q5SqV0qtAVZhlUVf5FJKva6Uaox9fRfo6cWzM5UrCYGmVwsiIsCVwF+8eHYykugHz8uYDoq6Algf930DGihHEekDDAXmxw59O9ZNfcxvF0MMBbwuIgtF5MbYse5KqU1gFSLg6ADkauEqDq88QadXC05ppFO5ux54Le57XxFZJCJzReTsAOSxyztd0utsYItSamXcMd/TK0E/eF7GdFDUYnMs0FAUEWkPTANuUUrtAR4GjgeGAJuwul5+M0IpNQy4CPiWiJwTgAy2iEgRMAF4PnZIh/RKhRblTkR+DDQCT8cObQJ6KaWGArcBz4hIBx9Fcso7LdIL+DKHGwS+p5eNfnA81eZYVmmmg6LeABwX970nsDEgWRCRKFYmPK2U+huAUmqLUqpJKdUM/C8edfmSoZTaGPu7FXgxJsMWEekRk7sHsNVvuWJcBHyglNoSkzHw9IrDKY0CL3cici3wBWCSijk1Y93kHbH/F2L5NU/0S6YkeadDehUClwHPtRzzO73s9AM+lDEdFPX7QH8R6RuzzK4CZgQhSMz/9UdguVLq/8Ud7xF32heBpYnXeixXOxE5quV/rIGopVjpdG3stGuBl/yUK47DrJyg0ysBpzSaAVwlIsUi0hfoD7znl1AiMg74ATBBKXUg7ng3ESmI/d8vJtenPsrllHeBpleMC4FqpdSGlgN+ppeTfsCPMubHaGkao6kXY42grgZ+HKAcI7G6JkuAxbHPxcCTwEex4zOAHj7L1Q9r9PhDYFlLGgFdgDeAlbG/nQNIszJgB9Ax7lgg6YXVWGwCGrCsmf9MlkbAj2NlbgVwkc9yrcLyX7aUs0di514ey+MPgQ+A8T7L5Zh3QaZX7PjjwH8lnOtnejnpB8/LmJmZaDAYDJqjg+vDYDAYDEkwitpgMBg0xyhqg8Fg0ByjqA0Gg0FzjKI2GAwGzTGK2mAwGDTHKGqDwWDQHKOoDQaDQXP+P4zewyNpLdWUAAAAAElFTkSuQmCC",
"text/plain": [
"