{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Vector data cubes\n", "\n", "Exploration of possible implementation of vector data cubes in Python based on xarray and geopandas. The goal is to mimic what stars package in R can do - https://r-spatial.github.io/stars/articles/stars1.html#vector-data-cube-example\n", "\n", "![vector data cube illustration](https://raw.githubusercontent.com/r-spatial/stars/master/images/cube3.png)\n", "\n", "Because we need to use geometries as an index, we need to ensure that geopandas is using shapely 2.0 (beta) as a geometry engine. Shapely 1.8 geometries are not hashable, while shapely 2.0 are." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import geopandas\n", "import pandas\n", "import xarray\n", "import numpy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We want to recreate [this example from stars documentation](https://r-spatial.github.io/stars/articles/stars1.html#vector-data-cube-example).\n", "\n", "Load geometries using geopandas:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "nc = geopandas.read_file(\"https://github.com/r-spatial/sf/raw/main/inst/gpkg/nc.gpkg\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get the GeometryArray. Note that this also contains CRS information at this point." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "origin = destination = nc.geometry.array" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create dimensions and dummy data." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "mode = [\"car\", \"bike\", \"foot\"]\n", "day = pandas.date_range(\"2015-01-01\", periods=100)\n", "hours = range(24)\n", "data = numpy.random.randint(1, 100, size=(3, 100, 24, 100, 100))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Create `xarray.DataArray`." ] }, { "cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray (mode: 3, day: 100, time: 24, origin: 100, destination: 100)>\n",
       "array([[[[[ 3, 22, 90, ...,  1, 39, 76],\n",
       "          [31, 59,  8, ..., 62, 20, 45],\n",
       "          [23,  1, 21, ..., 88, 37, 35],\n",
       "          ...,\n",
       "          [ 4, 43, 45, ..., 55, 44, 46],\n",
       "          [22, 78, 26, ..., 41,  3, 65],\n",
       "          [98, 35, 77, ..., 28, 72, 47]],\n",
       "\n",
       "         [[ 5,  4, 87, ..., 45, 28, 85],\n",
       "          [58, 81,  9, ..., 29, 35, 46],\n",
       "          [17, 11, 99, ..., 77, 29, 88],\n",
       "          ...,\n",
       "          [30, 37, 17, ...,  1, 69, 82],\n",
       "          [94, 21,  4, ..., 73, 10, 63],\n",
       "          [98,  3, 24, ..., 19, 88, 88]],\n",
       "\n",
       "         [[10, 29, 63, ..., 75,  6, 70],\n",
       "          [23, 20, 13, ...,  2, 44, 38],\n",
       "          [60, 19, 98, ..., 40, 11, 17],\n",
       "          ...,\n",
       "...\n",
       "          ...,\n",
       "          [19, 97, 35, ..., 14, 41, 54],\n",
       "          [ 7, 16, 64, ..., 29, 75, 74],\n",
       "          [66, 62, 62, ..., 67, 39, 36]],\n",
       "\n",
       "         [[61, 63, 53, ..., 63, 40, 85],\n",
       "          [78, 61, 62, ..., 32, 59, 43],\n",
       "          [41, 68, 64, ...,  6, 51,  6],\n",
       "          ...,\n",
       "          [51, 89, 55, ..., 82, 61, 63],\n",
       "          [46, 20, 85, ..., 76, 41, 52],\n",
       "          [77, 68, 81, ..., 40, 75, 99]],\n",
       "\n",
       "         [[53, 80, 42, ..., 77,  1,  1],\n",
       "          [78, 76,  2, ..., 51, 64, 11],\n",
       "          [ 2, 76, 35, ..., 87, 55, 91],\n",
       "          ...,\n",
       "          [31, 23,  7, ..., 69,  3, 91],\n",
       "          [24, 51, 44, ..., 48, 85, 20],\n",
       "          [81, 49, 65, ..., 46, 77, 31]]]]])\n",
       "Coordinates:\n",
       "  * mode         (mode) <U4 'car' 'bike' 'foot'\n",
       "  * day          (day) datetime64[ns] 2015-01-01 2015-01-02 ... 2015-04-10\n",
       "  * time         (time) int64 0 1 2 3 4 5 6 7 8 9 ... 15 16 17 18 19 20 21 22 23\n",
       "  * origin       (origin) object MULTIPOLYGON (((-81.4727554321289 36.2343559...\n",
       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...
" ], "text/plain": [ "\n", "array([[[[[ 3, 22, 90, ..., 1, 39, 76],\n", " [31, 59, 8, ..., 62, 20, 45],\n", " [23, 1, 21, ..., 88, 37, 35],\n", " ...,\n", " [ 4, 43, 45, ..., 55, 44, 46],\n", " [22, 78, 26, ..., 41, 3, 65],\n", " [98, 35, 77, ..., 28, 72, 47]],\n", "\n", " [[ 5, 4, 87, ..., 45, 28, 85],\n", " [58, 81, 9, ..., 29, 35, 46],\n", " [17, 11, 99, ..., 77, 29, 88],\n", " ...,\n", " [30, 37, 17, ..., 1, 69, 82],\n", " [94, 21, 4, ..., 73, 10, 63],\n", " [98, 3, 24, ..., 19, 88, 88]],\n", "\n", " [[10, 29, 63, ..., 75, 6, 70],\n", " [23, 20, 13, ..., 2, 44, 38],\n", " [60, 19, 98, ..., 40, 11, 17],\n", " ...,\n", "...\n", " ...,\n", " [19, 97, 35, ..., 14, 41, 54],\n", " [ 7, 16, 64, ..., 29, 75, 74],\n", " [66, 62, 62, ..., 67, 39, 36]],\n", "\n", " [[61, 63, 53, ..., 63, 40, 85],\n", " [78, 61, 62, ..., 32, 59, 43],\n", " [41, 68, 64, ..., 6, 51, 6],\n", " ...,\n", " [51, 89, 55, ..., 82, 61, 63],\n", " [46, 20, 85, ..., 76, 41, 52],\n", " [77, 68, 81, ..., 40, 75, 99]],\n", "\n", " [[53, 80, 42, ..., 77, 1, 1],\n", " [78, 76, 2, ..., 51, 64, 11],\n", " [ 2, 76, 35, ..., 87, 55, 91],\n", " ...,\n", " [31, 23, 7, ..., 69, 3, 91],\n", " [24, 51, 44, ..., 48, 85, 20],\n", " [81, 49, 65, ..., 46, 77, 31]]]]])\n", "Coordinates:\n", " * mode (mode) \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray (destination: 100)>\n",
       "array([61, 37, 67,  7,  8, 68, 28, 23, 79, 36, 56, 37, 90, 33, 94, 91, 20,\n",
       "       38, 29, 21, 93, 75, 43, 67, 57,  3, 73, 69, 84, 79, 49, 72, 74, 87,\n",
       "       94, 31, 11, 86, 42, 72, 53, 93, 19, 42, 81, 36,  7, 80, 24, 21, 89,\n",
       "        1, 82, 49, 21,  9, 96, 13, 52, 79, 80, 41, 71, 53, 30, 96,  4, 83,\n",
       "       95, 84,  4, 45, 98, 11, 21,  1, 95, 21, 69, 79, 35, 14, 16, 73, 67,\n",
       "        9,  1, 47, 28, 17, 41, 59, 98, 77, 90,  6, 99,  5, 53, 24])\n",
       "Coordinates:\n",
       "    mode         <U4 'car'\n",
       "    day          datetime64[ns] 2015-01-01\n",
       "    time         int64 12\n",
       "    origin       object MULTIPOLYGON (((-81.4727554321289 36.23435592651367, ...\n",
       "  * destination  (destination) object MULTIPOLYGON (((-81.4727554321289 36.23...
" ], "text/plain": [ "\n", "array([61, 37, 67, 7, 8, 68, 28, 23, 79, 36, 56, 37, 90, 33, 94, 91, 20,\n", " 38, 29, 21, 93, 75, 43, 67, 57, 3, 73, 69, 84, 79, 49, 72, 74, 87,\n", " 94, 31, 11, 86, 42, 72, 53, 93, 19, 42, 81, 36, 7, 80, 24, 21, 89,\n", " 1, 82, 49, 21, 9, 96, 13, 52, 79, 80, 41, 71, 53, 30, 96, 4, 83,\n", " 95, 84, 4, 45, 98, 11, 21, 1, 95, 21, 69, 79, 35, 14, 16, 73, 67,\n", " 9, 1, 47, 28, 17, 41, 59, 98, 77, 90, 6, 99, 5, 53, 24])\n", "Coordinates:\n", " mode " ] }, "execution_count": 34, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAh8AAADDCAYAAAAiPnOsAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/av/WaAAAACXBIWXMAAA9hAAAPYQGoP6dpAACFJUlEQVR4nO2dd3wc1dW/n5nZvpJWvVmymm3ZcjfuVAPGdENIgFBCAqmEN5UkBN4QICGk/9IJSXgpIQFCIJQQDAFswLhg3HuRrN7LSqvtO3N/f8iWLattlYw9jz/r1c7cuXN2Z3fmO+eee44khBDo6Ojo6Ojo6IwR8ngboKOjo6Ojo3N6oYsPHR0dHR0dnTFFFx86Ojo6Ojo6Y4ouPnR0dHR0dHTGFF186Ojo6Ojo6IwpuvjQ0dHR0dHRGVN08aGjo6Ojo6MzpujiQ0dHR0dHR2dMMYy3ASeiaRqNjY0kJycjSdJ4m6Ojo6Ojo6MTBkIIXC4X+fn5yPLIvo2TTnw0NjZSWFg43mbo6Ojo6OjoREFdXR0FBQUjtjnpxEdycjLQZ3xKSso4W6Ojo6Ojo6MTDj09PRQWFvZfx0fipBMfR4daUlJSdPGho6Ojo6PzESOckAk94FRHR0dHR0dnTNHFh46Ojo6Ojs6YoosPHR0dHR0dnTHlpIv50NEZD/xqiA/aarAbTCQZzdgNJuxHno2yMt7m6ejo6JxS6OJDRwf4V80O7t3y6pDrzLIBu9GE3WAeKE6OCJSkI8tzgg5CPQJFkTHIMgZFRpElDIqCIksYjzwr/etkjEeejy4zHLfu2DIJo1nGJBtQJFnPf6Ojo/ORRxcfOqc9qtB49MD6Ydf7tRB+f4hOv2fEfhb3TGLbruZ4m8esealUZu8HQELCKCuYZANGScEoGzDKBkyycuTZcGz9kYdZNmCWjZiVI4/jXlsUI1bVitZhQpYlZEnqf1Zkue/1ccvkI+JJkiQUWep7HrC+TzA5bBYUWUYIgSYEmnbkWWhHXtP/txAC9cj6o38LMfAzOFFvHX0tSRAyHe1bIOjbVuO4ZUIMeA30/201GLEbTSQZTdgkIwpH33NfxL4sSUhH35skgUT/3yOJQE0E0TQvGkEkDMiSAQkDkqQASlQCUggNgZ+AGiIoVFQRJCRCqCJESAv2PYtjzyERQtWOLguiCQ0QiP7/6f9L6jHT+pZyZD99n6EQAo4+99twwoERx/95fLvj7RZDtoe+pJJqSENT+55VVSUU7Hu2z06ic44Pmb7PeuhnGUWSkOg7JvZqB+7G0ID9iuPegziy4KhJQoj+Zf1tERgMCkEjBEPqgO/NgPd73N9Wgwz7WtECat/3TRP9+xGahoD+Zf1IEpIEd/zket7+5yYO76mnoCyHiVPymDgll8z8NFa/sAlJklj2sQWcaujiQ+e0542GfdT0dsbcj7BqcbBmMNLEXvAd2QeCgBYioIXi1v9kXynb33LHrT+AmcW57KiLvxA7kckFmWxOb4xLX+e6iti5qymibeTjBIp05PUD//NPzOaGUbZUBgkSWTIioSBJBvJMKahqHUL4jzwCQBCAQ9oVvO2M/ft6PLm+Kbx+f29c+4wV8RRMfKyUfdJon+UxZtSXs/sf8fls8vNS8VugqdMVVvvS3DRC71TS09AV0X7OvGwuK25Yyl/uf4E/ff/5/uVLLp7Fhdct5tEH/sV5V88/5TyeesCpzmnPGw374tLPRmMVc6bnxqWvo5SWJXPIF/7JNxoM6kf3HsRkil88jhYSozc6cRshCGkawZBKIKTiC4aQ5e4wtlTRhB9VuAlpPYS0LgJqK2bZSI7RSiDwAarahKZ1IoSbo8IDIF0JRGznaDhkSMk4uS5ukoDaT1cx6YlU5jUWM9GQOeo21Tm1yHJ83kdjkxNPs4fpE3PCal/V3EXP/AkULC4Lq31ymp05Z5cjNIHRZOBLD17Lys+e17++cHIu88+fTk+Xm+p98RHYJxMf3bOOjk6cSDVZ49ZXXWpH3PoymWSSZ3pp8cWtyyGRQwkIqB2juzTFGL/7JzUUH8+Voow8PDccFiWDZLGXUHDkdsmyC4jf51tsSWFG6jPMeAWU0GSc9RM5uNHBhpdVmg+rcdtPNEhA09stNL3dgmSQKHm8kMNqy7Dt83pyqCK6z38ovN4AjYfamT8xh/3dXbi8Iws/p8uLK1lm3jXzOPD8lmGPktlq5PqvXsyVt52HwXjs92dLPnYuKp1egMls5Jwr57Hu1W2UTJsQj7d00qCLD53TnhSjJW59NYguFk+fxLbdsQ05mE0S0y6VOZhgrweAlAjxcWJsQIKQDfG7CAcDsV9oDYqMLEcnYjJMhRCsHbWdiXYgK6p9DEW60Qj02awa9pNcvJ95xTDvOlC0InpbSmivdmBP95OU7kENGXA22Vn7DzM7342/F+ZEFJvChBV5iJCGwWcG49DtcjyZND0X7IutiCNJdjN7ttWRlGxh0tRcttY2jSiuVVVjU1c3cz61mPrnNhPyDlaTfm+QP9/3PKv+9j5f/OEnmHfuNAAWLZ+JLclCT5ebigWlACy7ZiG/v+sZbrzzsri+r/FGFx86pzUBNcTLtTvj2ueGlEMsmTOJrduiFyAGk0KVP7L4g2gRofh7KcZsfDqOuikYjN3zYTZGb5ABL2FF8mj1QCbx8n4kjWCyKtdgzauhMG/gcmseXDt9FjuXFcfFhpEoe3ASW9Orj7wa3rNoezOTVk/8PI9HMRn7LpO9Lh/7N9VSUZqFywp1HSMPr22rb6P4qlmI9w/TXTt0HErdwWbuue63LLlkNp+77xrK5xVTPq94QJuKBaV4e/00VLUyoTQ7Lu/pZECP+dA5rXm6ajMNnnDG6CNAgvW2Q6PGf5gNChMzU0+cAACAuzdEqTlviDXxRwueXGP9kaAq8QvyDQRjD+I1GRJ/ShXCQ47RFrf+rFF6aoR5B9/6u5ei6Ym9h/WljX5cjKqB2qr4BuEepa6xk8KJGf2vq6va6NzTxpz0TKbkpo/o5atucdI1O4/CMyeNuI/1r23nC+c8wOt/XzdonSzLLLtmAev+sy3q93AyoosPndOW3qCfP+xdm7D+1aSBbvw0u4U5JfmU5WWQk5aM0WSgpsPJxKxUZhUNFioWrz1hth1PKDA2QySJICTFT3z4/XEQHzF4PiLxZBSY4zdUaJFGCTIZgeTi/3Lbn1/me690MX+FKW42HU+vcfSgJ4tqRpETdDmTJAz2gWM9QsCBPY3Ub2um3JLC7IKcYUVIj9vH/iSZtOKRA2YrFpShqn3f57qDzQOmKC/72ALef3VbbO/jJEMXHzqnLU/v20xXIH7BaSdyyNBCdmqfgCjLzQBZYmt1I4eaO2h2unB5/QDUtjvZXtvM5LxjJ6c5JfmozVYQifdKJEJ8jNWoi1+OX0Ckzx/9RfgoJuPYvPEcU/yOmVHyx9hDCGPGe6z8/ktxsed4JizPpS7UPmo7l8lN6fU2JpVnYTTEP4bpYHUbkyYPPeulrqaDgx/WMT83B5t56ICUQFDFkjHyzcTO9QdJz+6r5P7vx98leJwYnjglD03TaGtIjHdnPNDFh85pS9P6HpYcKGFaML7TY4/SJTzYp8nMLc2nuq2Lzl7viO29wSAFGQ6yU5LYWdfMBzuamVw3A5M2TIRdnPD7x3dGQyx4ReyCAQAB3sB4D7uEv22qEr8pUArxyfEiDRcJGgOWT6SE3XZfdhUHL9tD2V0mPnXXfD57y1ksWVhGSnIcvERCYEm2II0wjXfP9nqMjT4WnBggc4SGfDtll80acpgVQFM1fvT5v7DxjZ28/c+NeHoHHuNl1yxk3Wvbo30HJx16wKnOacvimcW89vBetL2C+aVFhMpVthnro4rjyyaZMncuvc4g9iQDwibwGQJIzQa2VIc3R3+oALate1tJrc2nYraNPcl7IzcsDPxxmOUxmLHxAMRLfNjk+Fw4TYnVif3YpG7AHJ/ORHhJtEbD3bgoLv0cJbkkiT2W+kFZUUdjr6eWvdSSnpfElfMX8428T+FpU9mxu44nn15HW0eYydSEoGRiJo4kK02NTnbsqGP2rEJ2bht+RlIwqLLnwxqml2YjUozsaWrrX9fj9rMZPxWfWoR/fTWSLGHNsNOyo4HAEaER9Ie471MPA+BxeUnNTO7f/tyrzuBnX36clZ9dFtkHcpKiez50Tguq6ttZu7WKLXvr2F/dSl1LF2dUFHLbVYsBqKxqp+a1LmZvKSBfc0Tcf6kvhy07mjhQ287WPc1s+7CFfRu62FvdNvrGo+B0+9i0vpPCUH7MfQ2Fxxcn78E44NJiHTLowxqnu3ZTTFN/w99W0eKTPTbNYCGkVselr9cezhi9UQTkf34CIRG9MO4M9PL44Te5bv1DPN79GsVLHPzlt5/mrMUjB3/abRqf+qREZrKdmqp2duyoo629T6Dt3NNASenI05yFgMOVrTgPd5FkHhwHs6e+ncrCJA5NsLPTAtkfm0PhWYNtWv3ChwNep2c7MBgNONvjIxbHG93zoXNK09Lh4pF/vs9/3t8TVuqJ2vou0lw2is5Oo0YJP01ytzW+6clPRAiwt2ZAfvwzHXp88c/VMFbzZ3riJD5skhF3eBNdR8Q4Rp4PoTVilEoIitgCbmclJRGxa2EYCqfCtrf6/jaYYM4FJqafpeJ1Sfzjx5F9trJZ5vDEduJwSFCFxtut23m7dTuTkvK46XPns2BeCb//y2oCgRBWi+CqKz2YTCo2q8q8M96mpXEuT3WVDOpL0wTOXj82mwmPZ5SEY04PJo+Bs+dNZEd7O929Qw+V7apvIzsniaWXzuaD/+4iFOwTXM8//CZXfOYcUtKT+tued/V8Nry+g4tvPDOGT+TkQBcfOqckLo+PJ17+gGff2EogGNndU1e3h5Q1Fiafl8VBJTzPRbPoIdGXXCkRiUg1CU8cYh0G9TsG6sNqMtBDfISTRYpP5IPREMuFPJIPTWOCyUa1P7Z6LDmGdtQ4Hf4FN77Awo9PR5JDCGMl4rhjYzB+jL//IPwdlX22hG2h+vgYdhyHepv4wZ6n+f1ZX+aR6TezdfcjzJz1AoJjniQBtLaeN2wfHZ29TCvPY9/uhlGP2GWXz+Hzn19Gt9vHA394jQ93Dz1k09rVy/nfupqv/OxG1vxrE0/97N/0dnt57ndvcNu9H+tvt/TS2fzyq389JcSHPuyic8oR0kI8XPNb1lXtjFh4HKWn14fnzSDzfRPJ04YPelvcM5mFXZMoaIhfxsnhkCMVH0IiS2QgC5lMkU65fyqpDHwvNuKXWn6ssVviFPMAmON0H2aKQXxEsmXIsIhqf2zu93SjFU3dFVMfA9EQ5p1oxr0DhAfA1EteYsbZ4buFupckbihQFRoP7nmGvIIkZs15eoDwOMq+vUlDbHmMvfubmDGrcNR9rV93kN5eP1lpSfz6ux/ny588G0UZ+rL7wn+34chIYumlc/C6+zx6L//fO3Q0O/vbWO0WrHYzvd2Jm6U3VuieD51TjhcaXmBv7x5K551BZQwxmh5fkMrX+zImnjW7jAMTW2iV++40M0misCWbbYebxiqTOMiDd1SmlmJw2/F7wJIsQNaotVSTKpLxH07jQEsXEll0AAfpYm5pEdkZfjRJJSgFcCtuKubnHCszrvW5loV2pNz4kb81TRxZ3vesqX3PqqohNIGq9i1TNQ1NFX0VXgVoCfSA2CzxG+Mwivi4lYyGWIZBwvuwaoNLWNN5NhbJDTgJiBY0Rp5JNRQXOPxooQQXDupH5fxPedm9VkGMMn18woo8DoxQvyUemGQDPb4taNrgoVW/50zefHP0nCUHDrWQX5BGY/3ww7PNzd387z3/4Be/vBGr1cRNVyxkXsVEvv+7V6lvcQ5ou377YRpbu2mraiW/JIu6gy0EfEGe+dUqvvzj6/vbnX3lGXzw5i7Ov2Zh+G/4JEQSYsxOnWHR09ODw+Ggu7ublJTwp1np6ABU9Vbxw70/RCAIrVnAgR3xO7majArlC7JxpwdpP+inuWtsS5CfuTAfn80LPiPuXkF7t5cm52AbLAYDQVVFDeOnXVGcyTalES3Op4EVO61U/+5DZFnCYDZiNBswmI49FJOCYjQgGxVko4xiNCAZ5OMeChhkhCL1PcsSQpHQZAmhyKgSpEw1EjyzARkDiqT0lacXCiCDkNGEjKZJaEJC1SRCmkRQFSgihFnyERISIQ0CqiBXddB56Li7Sel4KSAGpIv3hRQ0TUJVIaRCSBWoKqhBWDhV4rzFv4rqM5tonYka3DBimwOBC/np4WK86kDPQIrRSrbZTqrJhN0gsCgevNrmYYe/FiRnMFF5PSo7Y0HRivj1DfNorR3eIzn70Rl4M0OoCIJaCJ8aQBN9x0BGQpZkJCRkqc9b1ODpIKSpKJKCIsnIkoQi9z3LyP3LZPqei+zZ3F1xPV7PI7R3/2TQ/h/945fZtCm8+3JJgsllOSiSxL69jcPKx6VLJ3Pf/R/r93q4vQF++cTb/Ofd3QPa3XzlQm6//myEEGx7bz93X/sbFIPMn9+/j7yivjxAakjlD3f/g//56SfDsnEsieT6rXs+dE4ZglqQ/6v+PwQCCYnaQ/EJRjxKIKiyc10TFTPzae6Kfw2J0eh1SXy4d/QgWF8o/LH1PdXtnDGvgE3OulhMG4Ti7xMzmiYIeAMERqkGGg1Tbi6mdX4Es4nkvsf5mRmgDc5sW1IcXjft6mRcIeeQ6yYnh9nJkAx96QpqNv7Y+kXcIcE+VzfBIaYX9wS99AQHej/mp59DmrkSgR9VeNCOVHu1yQZKjVuJUxHfyJBdtNWFOP69WnMsWHMsmBwmUpeksNG4H6KoeBASHBu7GmG09fbJl2MzWGj3vTPk+qTk8D8YIfo8IADTpk+gprIV3xCzx9atO8jDD7/FHXcsB8BuNfG9L17MollF/PTRN3Ef+X28snonn71mCSajgbnnTKViQSl7NlXxt5+/yp2/vQUAxaBgddhoqe8kpyA9bFtPNvSYD51Thn83/ZsGb18VWIFgyaXxiwk4no5G1zj9cBLjpFQSECgr+RJ/ZZPM0R2FWN6tIiUNKzwAZCn+x+it3pvY0NnBzp5OghFMPf2ws5H/Nll5symV1c351LgWIWlnMd+aR0hrjbud4WHimjtNmI8Lqch7aAKNd7qp/lwX22bUJNyCEnsumtaL1//hkOuLi6Pzlu7d38SX/2c5hYVDC4J/vfAhz//zgwHLLlo6jb/++FPMmNyXmMzp8vLe5sr+9Zd/5lwA3v7nB9TsOzbTrbRiAk//5o2o7DxZ0MWHzilBraeWV5teHbCsJX8T5329nYu/KMjIip+TLzMjKe7DFOGQqD1qCehZ8iY+a6oyRqnMj8dsmDDiepnoRZeKjJ8iQuLYd7VbnchzjfERctXuDmraLNzzehlP7/42P1r7v/xo7f+yr/3KuPQfDirNzLr6eb73xhq+9lgQgwmM7rG7DNkUMzmWNDz+9xluHm9eXvSBvG0dvfzu97ewdOnkIdc//PBbbNxYOWBZXpaDh++9ns9cvZg5UyewePaxKb5nXT6XtKwUhBD89Wf/7l8+c1EZ/33uAxoOx55HaLzQh110PvKoQuWxw4+hDnFX2Ko2gamJ8pvSqPxHMS0Nsc0rnFSUybae1rErXnIciQrPSoT4EN74T989EdkUvecj6ncspQLD19eQY/gs6z3beKLlC+x3dZFnSaLAasQZFPi0+NTzSMLOgUpBtz/Ai/sA+u7wNWZy5+J3sRiccdlPOGj0kF7+Gvf9dxq1HiOqP5NdPaPXcImVEnsukiTh8Q0fW5Oe4QRGFpnDsXtPA7fcdCb3P3ANb725i/37m6iubqe6up2uLjfTpxdQUTE4WaBBkfn8J85EiIHxRUaTgYtvOpOn/99rvP/qNvZvraZ8bjFZ+Wlk5jl46v+t4ju/uTkqW8cbXXzonJTs7j7Iy41vo0gyuZYsFqTPZFpK2ZBtVzWvotpTPWJ/TrWLkorJMYmPovw0qlQXQW08BssT5/kIxXC3PhxjIT4kU5QCMIahkYA28qwYWdKiPlCNoWXs6emLJar39lAf+QSWYRECkruLqfQOvlPe2eLiS/+5laWFFnLsXjKsLtKt7TT15rG/Ix1fCEwKLCveTciQRotBZaq1g3z5dRQptuOsGfdS4NjLjcLI74I30+CNItgjAkxK3wwps6li2DY2ewswPar+9+xtRNMEsiyx/KKZLL9oZv+67m4PFosR8zDF54ABwuMol37qbJ79zetoqsaTP3mFB5/5HwBmLizj7X9t5rovX0hx+dD1ZE5mdPGhM64IITjYW82Wrj3IkoxRMvBh10729Ax0Tb7W9C4/nnUnRfZjdw3OgJN/N/2bNW1rwtqXLXXwRdZiFuTlQE42pDpU/H4JV69Mj0uiuwc6nYLUpGSyClPY3NyIGLPcnYOJt+fDbFSYXJJOgz/+J3zVm/iU7VIMno9ocakjK4Joh12EgDfapwCJCWQuF9N4r2F4F70nGOLNqqMzp2zAxCN/HxuCeLdmIteeG6LOe5DdLsg0Xc3l6VWkSptjtk+RgsxOkWmIo+Aaiq1dh6h1tzLBugIJ06B8JACKIfoswm63n5radkqKB+f9cThsUfWZmZfKmZfO4b1XtrBlzV52rDvArKVTmLGojLde+JC//uI1vvenW6O2ebzQxYfOuOBXA7zW/C5vtayn3jt6nQqf5udHe//I7ZNuwBV00+yr5Y2W14aM/B8Ob8Ehbv9mFpn2XaC4UOkiJEbPa2mWCtm880z2dBnxxHcCTUTEU3rMnpHDdn8jH4bq45LC+kRC7sSLj2hjPqIVH0JAZ3DkMXa/iC7IuSF4AXtdiREe2VI2HxxwxtxPSZaZ+uC+/tftASdPtWRzSebHyFO2YeFwVAG3QkCvVkqOqQuDpMRUzyUcnq9/j6+XX4Pdch69vsFBm4IObvtcN8/8PQW3e/hvi91morAwg4IJaRQWpFNQkE5hQToT8tPibvOnvnMFrfWd7N9azc51B/vFB8C6VTs4uKOOyWEkPTuZiEh8PPzwwzz88MNUV1cDMH36dO69914uueSS/jZ79+7lO9/5Du+88w6apjF9+nT+8Y9/MHHixGF61Tnd6PDuZU3zc2x0dlLvDRDu5aDV38F9u38LwGxHSUTCA6BbdbJGdXKGyKVUdhJSw0uo7Rd1zJjxDDNnKJikQvzeAl5ZW8KHB8bWCxLPINda0YUvXnm1hyDojv/U2kFEO+wSJWZDLoHAyOqzJaiSY4gsJEgIWNVRxkixJNFiEAq9Tan4VWfMfS0oV6nTBn4HQyLEK20uoAyDVE6+JY1LHW9jkUa/oVCFwruuW1nb6cUVOurySHyg8htNH3Jb6SUk264cUnwALFj8OAsXJxHwVdDTnUVraxo9zjlkpE1n4hHBkZpqG3KYJBEUTMrhl6/eyVv/2Ijb1fdZTSjJIi0rma42F0/+/FV+8OQXx8SWeBGR+CgoKODHP/4xkyb1VeB74oknWLlyJVu3bmX69OlUVlZy1llncdttt3H//ffjcDjYu3cvFoslIcbrnBwIoaGhokijZ5wMBPfT03klc429zM0Cn7yEew5HPl7pUaM/SW3u7mKnlMcV2cUEguvD3k6g4hfVYKlm5YUbuXjxOWzYmcubW+NcdEUCRZb6MomKY9Isnp6PVKOVJm/iqmP63Yl3EUlRVpCN9nJhVHKAkV3y9d56kpKuwi5peLem0bDejM8tYUkSTD13HxlT30U6YbTIpZVzwBV/4QFQ5K9gXWfs2ULTbAot7B+xTUiEqPW24XPkYRkiZfnxqCETr2/4FmuzDxKSEy84jsenBfh3w3qun3gRkmRBCB+KnImqDQx4FfRitHxAhgUycgCewmycSXrKl0ixXTWmNgPIsszy65f0v5YkiRkLy3jv1W18uGYfuz6oYsbC0jG3K1piznCanp7Oz372M2677Tauv/56jEYjf/3rX6PuT89w+tHBE2zmYPdzVPW8BMDcrG9SnHzJsO1VzUlD66WEQof7lwmlnDurIg/usitWjHJn1MGDEnBtjg/PcW7k6JD5z9s3sHFf7NMFZRlmlOZxuNtJl2fg4Hdf1kYJs8GAyaBgNhgwKjIGRcGoyCiyjEGWcTgkkotagSPpH5FAHP9337PBaKBXqCiS4OX98Y/5mP3NQ/h6EytAzvh/M6icfXj0hidwSXYa/tC6iLezmpawx3Uo7PbGp8/mw+cHipWMnCTOWK4w+dwtpEzchiRBe2g+P60qjtie0SiRitmwKz41QK5eaKHZsDWstmel5TLP/AyyNLRnUg2ZePrlb7Fuh4spZRlIczposDbExc5wyTQ7eHrJPYRCu/D4NxIM1SEh09X752G3MSi5pCffgSPpBmTp5LihdjnduJxeNr+zl/de3cZPnr1jzLwxQzEmGU5VVeW5557D7XazZMkSNE3j1Vdf5dvf/jYrVqxg69atlJSU8N3vfperrroq2t3onGQIIejw7eRA99PU965GHOcm3d7+KyYmXYQ8RPlVIUK0dnxxgPAAkLUGSq2zKbPKuDVY5wzPXe9WvcxOnkaNZ09U7+PSLAee4Laoth2IxmXnP8NFC8/lL6/k09R17IcvSTBlajamJAOSoL9OimSUaNJ6STGYsIdMCEXQJXwcdnXhkgKDhAf0DbloqiCoBmCEa/oZDjtt6gh3qEfNO27E5app0/nPfhcBLTwhd/tUN9mWQP/QghBQ5bLxbouR6l4ZNJFw4QEQhqNtSAySOtJHOCxeNUKhK8GZX8mk4S2o3t13V93R0ssbT8EbT5VRUHYG8y4KkLQo/h4oG1YOVMbHV2ZUwGOpDDs2aG1XM+s4n0yzg0yThWnWdibIryNL6gDhAXCgsgMqYfq02ZhKfRx2VBKSEz9Tqt3fzZrW7VyYOw+LaTZwNKBb0NX7lwFtT0bRcZTkVDvJqXbyis7i3Ve2sm3tAeaeXT7eZoVFxJ6PnTt3smTJEnw+H0lJSfz973/n0ksvpbm5mby8PGw2Gz/84Q9ZtmwZq1at4u6772b16tWce+65Q/bn9/vx+4+dCnp6eigsLNQ9HychDe532NP5GJ3+3UOun572OWZkfH7IdR3OB+jufXiE3s1IkolnOy9jY0+AG3IUHMYQjzdKeIeZQGBTrBjljmHH14XGIBc3QEWSg0mmtQOEUzwwS/k88vwFdPvtTJjk4KC/kzZvZHeeczJy2b8z+sRB5yxMojVtW8TbFVvL+Neu8ITf6xe8AtrQHgdJysIdKueuaYnJLns88/84k0NTqiLe7qpcKz2BLRFv1y0q6AxEcGxUCRRB8ua5rH2oZ8SmRbNyMS2xsH9KKz3m2Kd8lHhmsuFwfIqzrZhloycpthkt56TnMt2wiude+ny/8BgKu8XE1HNt7M8a+hwTT8qTC3l4/lcHeAqEEPiDexDCgyY8CBHCZjnzpBMdQ3F4XyO/u/s5fv78V8bN+5FQz0d5eTnbtm3D6XTy/PPPc8stt/DOO++QmpoKwMqVK/n6178OwJw5c1i3bh1//OMfhxUfDz30EPfff3+kZuiMMb3BetY2fYuRIg9CYvCFVtW66HB+n17Pc6PswY8Qfq5Lf43z05aQyWrAzwTLtRzyDH0nlGyw4RdDzBBQJdJfX8r+DV1MvzCVlrPW9YsQu2xguvUAvhhiRoZ9B6KRvJkp7Kx2UuOM7m52e0czhWkptHWFFwx7ImZLdHe7Tf465udP5sPGHhQ0Li1U+Z/Jz9LnKjGAZERgQKCANnwdGCHasCsCmB2VHZEgRem3lYg8fbaEia5AhEmwlCPHwjj69NuaHc2wA4wGmfkLJ6AulNlT1ITfGLkXIF/OY2OchAcILOkN9MQYP/xuZzPvanOR3CO7q9y+ALUfGJk6dyb78nfGttNRKLBlDlomSRIW08jDwG5/AH8oRLp95KmzQgjcIScWJQmDHL8KzMNRMjWfqXOL2PjmbhYvn5Hw/cVKxD9fk8nUH3A6f/58Nm3axK9//Wt++9vfYjAYqKgYmLxl2rRprF07uIjTUb773e/yjW98o//1Uc+HzsmFO9jIaCGPLZ4PqOt9ixRjMT61i9reVbgDlRSINchhCnEhuslkVf/rVGX4WAqBQMGAeoI/OGfjIjat7qtd8eErbRTuXkhGiRFTXoDZc+rxqUPXdIgVIWBjkzemwFCjouAZojBV2Nubo8sz4dcCYN3NV8/I5TLH0wjhZIBP9MjfYR1GzRqVDRETpfjQROSeBYuhEOEf2XsxFEKDxpdsQHgVkNWQxqF19bAOctJtqN9MosHhjGCHEq6mFASjFyAMhyWT7LQGYo2LOoIssJ69l5TmifSMEJDc1uVmotsen30OQ0VKEd+eet2QHgK3P0CDs4cGZw+Nzh4aunpocHb3L3N6fKTaLDz7+U8yMT11wLZ+1UO1eweVvR9S1buF7mDfeciupJJszCTlyCPVmMOctBWYlehyfwzHDV+7mB/f8QQLL6hAlk/u6ikx5/kQQuD3+zGZTCxYsID9+weONx84cICioqJhtzebzZjNiXfR6oTH3q4n6PTtxSjbmZJ6PanmyQQ1N9s7fjfqts7AQdY13zVoud12LhladOW75yQH+XAYJ0Krv5OpyZNp8u3tH3px1J2Js9fO8cma6qq6qDvind/8TApXfuYa8hY/H5U9IxFQc3H6Y4t1mJmWw5666O9aJUOQWBKWvutsZnlaOUZ1Y9R9CHVsxEe0s13UITx0oyFL6UDk4iPl/UXs2B5dEbeeTg8ZvzOQ/nU7nbbwPGFTmMp77fHLF1JQ4KQ+jrOmPeZOZq2YwNoXRm4nGYaW8EKALCSEHL3ETzHYuKXkIj5sreJwrXdIcTEaTo+P2//2En+77Vp8chOVri1U9W6mzrMHbYjhXLfqxK06afYdC1j+oPNlLsv/CqVJc6N+LydiT7ZwwTULWPvqds65In79JoKIxMfdd9/NJZdcQmFhIS6Xi2eeeYY1a9awalXfneq3vvUtrrvuOs4555z+mI9XXnmFNWvWJMJ2nTjT4dvNjuNERrXr3xSnXIEn2EyXf28MPUd/ophmfJ6flM5lTU85r7UPPgvuc9Uw2zGdGm/fGLGc1M32bh9ZaVbcXYPvcEMhje1r0rEkX4RjyrvIxugqWA6FPzjYjRspoe7YhoM0JTbxAdCtphPLOwn6rZjt5oRPtxVKdN+roBb5kFiQ8NzmZmcmkibhT+nCtnohG/4vtnolHU09pP48ifLlJeyd04zTOrzXxi7Z2HEoPA9LOEzOMVMfiOV3PzStWTtZtOhsNm4c/rPZ/aGTvPoZmM0KRrOMwQyyCTrqgtQ1djP37HT2FmyPav89IQ93bf0Lha1z+GBf9PFVlW2dPF1zPx3arujsCLbxdM33OCPtMi7O/9KIbd09Xg7va8LZ7mLpxTNH9Gqce8Vc/vC95znzklkohjinAYgjEYmPlpYWbr75ZpqamnA4HMyaNYtVq1axfPlyAK6++mr++Mc/8tBDD/GVr3yF8vJynn/+ec4666yEGK8TX0xy8oDXAo3DR6bRxoJfSGjCgixFd6E3qFtZnlzHqvazh0xvvr37MNMss+nckU5zj8Dj7yR7cR6HXxs6GLHyoJPKXzo4Y8mnWHjrn5FjuIs6Hpv5ELI0J+pkYDk2Owcq24glAXgwys/4eNQYi13X7MolaDChpJox20xYbSZMFhMmswGjyYDBqKAoMrIsHXF7C4SqoWkaalBFDYYI+IIEfUEC3gA+tx+Py4sWOkGYRXVeFUhRvD+3Gt7tf9cfp3DwwxYgBYhP3IWztRfn33rJejcV5Q6FDvvQAiPbXUaVP16xHjBvSpDaBKXgcM38AOmDUob7qXj9Qapqhh86Crmi/45KmoR5TwUftMdWETbLIdOhxR4Y2+av6f9bVTWaato5vLfx2GNfIy11x/LALLloBt/85Y3YU4b2MEqSxBW3nMW7/97GsqvOiNm+RBGR+Hj00UdHbXPrrbdy660fvTzzOmAz5hJjzc8hafDupl2ZRLElixR1dVR9CM01olXVvc3s3nasxa72NkYbNd68vovk1FupuObRuBSpNcg+CpLt1PZEd/dZak1nm4i+rgSAT0QvPq7OMrIwaS+SGtvdrt/Td1IUQuBz+/HF7AGRwG7DZFKw2MxYrCZMFiOiWSFyF42ER8zBwPBxaCciBHQERh/KMLocVG6JboglHNrqnGT+2sGE9GzExUZ2Fh/LjTFByo9jkClkJRto1OLv9TiK1Z82rPAIBzm773suaRLFvSUYmxw01/jIWhykMuXY0EZKKIXc7nwaHXX0GvqGrkp8Faxrjz23TVu3Spo0lS4R2+c0OXkhAK0NXXz+/Ifwe0cWuuvf2MVXr/gl//unW4ctKDdxci47N1QSDIQwmk7OKiond0SKzpjiCtRGdVcYDn7VyX73QYLylKi2l+RkhvIISJpEvnMGvo0D57aHmbKCNa/1sPMfnyXoik9VSCUGEdPYEPsJ0R1myvih2OoSyFoNsRZ78bsTE8MVCqj0Oj20NzlpPNxGUIvOznc72kk2n0eyeREp5vkkm2ZiN5ZikjORxEB3iiQUrIZ8DCg4DOmkG7PINOWQY84n3VNEUs0kkmrKsNeUYXhnBlq4X7woaW/spnZXC82/aWZqU07fQiHhak6O6y3D2dMVQiJx+TZsvYMLr0WC1+Km0FNI2b457HzZz5ZNrTS29rDzFR9F7mMxhnnVZWx5zYX2xgRKvWWU9cxm46bIY3eGRuLVdwpjPmdOSV4EQF1ly6jC4ygNh9v42pX/j+7O4W90zv/YfLa+N3JW2vHk5JREOmOOKgJsaLk37rkvTqRJK2QiByLeTkh9c8aFgIL2uQS7zbicgsZWN5W9Pjhh+mRpRipthDe2/94b3XzwznmcfbGZsvNfxZAUnTtWE9DQG12Ohoq0LA7vji3FdpbRSurL08nLliDDjyuznTZHU1jbTrObyDcDUhqI2E7Ofq+JETOhxQtbCOMJCe2Oy9jAiWL1eO9WZ1Ci3T/UiTkZi2xHkiQCmg9VhCDkATyDHIKe5xex86XjBePo9UziRTAQovvXXaR+30qWpZi1sQaZCsFUYzp2YaBT8tPyDOSUnIlU3EVb1iFUQ3xr9WjOJMAZ9fZ7/nX0PDXQ06QJgc2VAnawqlb2bO37TbV0ugkdSKW2K77F+xZXBBExBFllmApIN08AoL4yMq+Z3xugvcmJIz1pyPVWu5nkVBt+XwCzxRS1jYlCFx86AOzqeITuwMGE76fFt58c2wzMEQZpCczkmBTkg2fw3sZOGEJYpFjNFKU6SArIHPpvdUT9+/0h3nwpxNrXl3PuZUZKlr2EYnVG1IeqpRAYIX+IpAkm+E1keWS2Zhy7OCebzATbYhd9FZKDvZtrj1tiZPa5C3CfUUdbUvOgoaXLM03kmQMUmpqwiMNoWgsixmBVAF9v4nMaABTMbsYm10e1bbt/+O182ujeozQK2PiSk1jic2LF2+unYu9k3kqJPsg0XTYzJ5RBb6WLxob2/gnBNUBNdd/fOTkzsHxmO1qMNViKrVbmp0js6pWo60xcIKRwG7CqVtI3TabN3yc2plbksrkrfsNSADOKFYxZ7xJp4tvjOer1gD5vRqR0tvZQNkJakqnziqk90EzRMMMz44kuPnRo9W5hnzP6ejyjocjpmJQJaFIyfk1jnacHj5qM3ZDMdEuILPH2qH00tttQN85i/cHhvQP5TdC0oWbY9eHg8wV5/fkgSa9fzrXf2YU1f1vY2xqVHqakp9Dc0stEt4Esp8DaFkQ09uKu7aK9qo2gP0gvMO2hs9nrCGA1GCgMJnO4Nfa8DO7qwReh7e+0wDsmUlPLKVlp5HDxQNFXpvwL1JgnyAzA5x6bCHtNjj4fSqx41xUQy517rAhgwpWTaJqQRE97ZJ4qWcB8QzZJLVC5t4UD6vBJ4wDa2z3kSWpUMssiSWjAymwLdvFfELDYDmev3Mf5C+awak0muw/Fd6h35+ZOQhszaNOOeTlCdimamdLDYjbCuUv30K3G9h2cnHKc+KiKPF6os3XkNyVJEhm5DkJBFYPx5Jr5oouP0wivN8Dh6jYyMpLJzEhCUWT27G1k/f53SJoR37Fqg5yBYiinLdBJt78VGFw4yh9wU6uUkhXGt/BXr67kUNMowyhTk0kLanQ1xn6W6e3186//N4uP392FKa0GIUCETIhgMlrIjCG5cZAnwVdVStZ3tyAOd+AHRronL3qjjepPpDFFyuBAQ2zTMQGKzSnU1w7fj9PpZesTXmYsmY9vST0dtlZ29Wqcbbcioki6NRJ+z9iEkqlSMN6x0WFhkZNY/3Acr2QRkn1uIW2Lk1jd0wntLtIsFrp8owcalyjJTPIk07Snncbu8IbjoG8GRpprAs6UyIq/LU5NpsK8Ggkjqhj42w2KGjLyarjxk2BkCrU1M/ndk2lDzmaLFH9gYKyKADyh+ArVa87R6Faj87odxaakMMF6LFYtGs9HV9voQ8tJDhuqGs/bi/igi49TmEAgxI6ddbS1u3h/3SE2bT5M4MgPU5YlMjOS6Oh0I0nw6R9OwBDhyWUQQsFqmkWvZqTeV4Xwj54ZsVf1IhRGnG3S483nUNPoruV9ze3YJhiZNaWYg2uqYz6NdXV6efJ/l6IoZ+HzBgcEE85ZmMaCa9/AlFYNgGd3Bc9/xkBvmGPKlav3cW5oCtvOUeLiuS8J2ggn5n7X+lZYb8JgmIhyeTqrzpNZkfxs7AYchzd+qSZGJDQWcSVDkOzLJ+Ae+5N5+pwcfBdl8Z6rDXqOeQBL09LZ3DT0LKkkycgZIpNQjY+aw+3sjdJbY95VDEvDPT8Izk9Podi4Gk0Mjsc6kSAHyCs6wMwpt7HjAMR7KKtiZh6b2uIbj2MynphXOXImJS/sL8IZ8AVprY/c+zma5+MoiiKjhtSTKu+HLj5OUdraXdx7/7/Yt3/oOxxNE7QeUc0OhwU5KfopniZlIkIupMnfgNcT2bBHR6CJGvOVFPPykOt7fbn89KUvMVSMx1B4AkE2BFqoWDER/7ZOultiuxL6fUOfYrZ90MXurYspKLqQ0gKV9b/4AJ87ssyZle8doPhwO+4bZ3LIGMOdmRC0H4rsxBUKaWx+sR2bZzreKzKx2mL3vhzFF/2Em4gIjpP4oNfGaBfUeJJc6sD4sSI2uFsQrsF3x55ggGy7nVb3kQ9eCOYYs8jqMFC1u4XKQGx36AC7NzjJW2AhOEpSvsk2G8tSmwmE3g97xtlRrr/uUT6JCUVKZf/eZTz6XHwKi3bI8f+exMNDE9KO2dVU006ENV4B6ApTfAC4ur2kZgwdnDoe6OLjFKKzy01lZSuVVa0898ImOjvDuwp0d/t4+cFbmDRNYdqlL6LYRr97l7FjMk2nKxSgzlcLYd13D80WVzVNlsuZZmonVWzoX97ZO4kfPH8z1a2RZ6Tc09yGvcjEjClFHHqvJiFhgcGgyuFD3Rh390ady6KjvpMJrxyi8uqJiCiTjUyzptPWGl0w3XtvdLDpnVs4b5mBs5e9SkZmdNkaj8fbm3ivgKRASMR3Bka4VL9kYizEhyXLhuOGSawPtRF0D39897a3Y5Rl5uXlke5U6Piwk7a2FuInJ8HrC2J7+wy6l6/rCxo5gRSDgauyQFHfInBiMrgIEAQIiVbKpj7L0/edRZJczarNZ/LHV5JQtciH86xWI9WdzqjtSSR7e97n/EArDlM29VHEe0D4ng/o83Y7212kZiaP3ngM0MXHKcI77+3nvh+8GPX2zU29NDeB3XEpJcuHDj4VAqzGafhJpcF3mFDw0JDtoqHJV0tnwMbZ8gR21V3E5qpiPjgYoDeG22i3P8BGfyszyzNp2h/PU/ExJmfaOPTC1pj6aNjTwLlZyfhnZbM9WcUTYcbVXI8ppjJiPn+IVatCvP76BcydczU2u8Djlllx6fsUlb4VeX/uxE7XBrCmjs2MmhNJkXLY+Eo3iZzlYkwyknVjOZvM3fQGwhOVQU1jS1MTFoOB2fkpEEYsQKTs2tzGGclLaVnyfv8wqQxcmpVErrKBkOqMWwjOBMsU0qVNCK2dKxbs5bIFRazeuYzf/CsNfyC8z16SBNdeEMCQ0ktIk1A1jjxLx54FqKpESIPK3Xb8QQX1SMIBVRKoEoQk0ddOEoQQhAARh6+4QOODzpdYnvs5Gqqim94fTszHUVLS7OzeVKWLD5340dLYxVvPf8ickmy6vAFqmp1R9/Xe6yoTz0rFnpxHSErHIKlAiKAw0RnooN7bQrzSRp+IiXRu+d05R17Fz3dvK0mBBIiPrGQzXRv2xKWvQ+/sg3f2kWM2kvrlpWy1hTcMIyFoiKE+xfEIAVu2HruTWrw0N6p+PL2JS051FFva+IgP+UARcZ02cQLZZxWwf6GRXd42iGIkzhcKYUpNXE6HzWvaWJS8iMaZGzkjJZn5SQfwhyoJxdnZ1R3sIPe4tyFTwwUzH+f8mblsOHAhv3wuE5d3aE/IlAKVz1/WREXh+8hDBLoPx3/2fYwnfzBwmFo58jgxbd7bD4I1OY2kNBu2NAvWNBPWNAMmh0LKx9oIlgzcr/HtGVhqJ3LeV2ajJblxBTtwhTrwhLrRhEr94eg9H0KIIavzDoXJbKStsYus/LSo9hdPdPHxEefQ3ka+d/sTdLX3xTZIksS8xaVsqY7uguTxBKj94GMYlr5IokTGcIR6ihPSb5W7Z4iUU7GRmWSCrQfobIxP6fKjBP1Bev64gYovL2GPefRhhb5h4sTchb/wfBFlU/KxRxgP5O1N/HCIJW08Aucktj+emPcmgIJPV/CuvYuQN7aZGb70xOYead7RxW0XqKi8hj9BOrNXbUdI+XDC4JFEM0umPMVz96SzrXo5v3guj7ZuGZNB49Mr3Fw0bydJ5g+IZvL4Rbe8zLvPXUb17vCEgNflw+vywXGpdSRJ4vzPDBTGQoPdv3fReGA1qx/9gOu+vZJpS2Yhd7lxdbl5oes/uOpbqajoy1gr6CtNoGmCQ4fa0UZIJOL3Bgj4g2EnEZs0s4BXnljLlZ8+O6z2iUQXHx8xVFXjlWc28MG7+3F1e6mtbMXvO3ayEkLEfPLa8rbKgsUSUpwKroWL0eIG4h8Q1dbjZvb0bBrCPKmMRkaSCXnHIdobY8tIOhx+jx/Dnz+g7AsLqBwtEFWSyC1KpWdnfKfLAjQ3u/nzH27jy1/7BUZTeMG0QvQlv0o0FsfYi4/M4CQ27o1/rIcxyYjtS1N5u7c1LglXtna3MCXHTk9LfCN/JUWiZMUEli99B5X1ce17iL2BNJJ3q5M5xc/y128l0+ZaSGbydqQYo1wMxhAls6xUx1ArrnTBBNT0gckaTRum0nigLxC/1+nm0bv/HnZ/My+ax949w98ELv/Ewoiyl0qSxP6tNbRcOJ2cgvSwt0sEem2XjxD7d9bxjZsf4Y8/fpUt6w5xcHfDAOFxFE+Pl2Rb9K7XtjY3gQ2XAiC0scvg6DXup2RCYr6S5onxETXpdiPKzkO018c3TfOJuJ0eLI9tpVAd/f5ASmD8w779Pfz8oW/w3uov4/OO7qoN+u0Jr28CYEmJv/iQMWCSbMOu716XEfd9plZk0PulEj7ojV9BupCmkb4wvrbmTc/AcHUaa9QGfryxgnbP/Lj2PxSHAyEkpWiUVi6ykt+KWXgcZebZzpi2L7xwYDlLIaDq0egDRCwGmLm4bMh1VruZT3/7soj7zClM55nf/jdqm+KFLj5Ocnp7vLzyzAa+fO3v+OoNf2T/ztGnzdVVtZGTNlpN15F581kv+373Md793grk5hHy98aRoPCwdGFi7poP9TgHjk5IImLPTorVgHF3Fe11iRUeR+lu7aEwjHiy9jjX3TiRmloXf3vKwBN/uQNVHfmiH/Q4EmrLUUwp8T112ZQ0kowFeDTItUwn01wyoOqqUbKw+fH4BnHmryxj14VWqntjLyh4Ijv87RgssQs0c5KRwmsK+KCgk0pX3xCjKxDgJxsW0uGdHXP/I9EVbKQmNLZTQ5deuYrbHsomuyi6mAhvExj3FiNCfd9P85bJVH4Q3VTn5PQkvvbw5/jpP/6HX774NZasmDlg/XVfvpD0nMh/bxNKstixPvGlNEZDEtFMLk4gPT09OBwOuru7SUmJzzzvjyJCCP797Eb+8otVQ3o3RiNnTgG1LfEJjHM4LJzz9V1oabWjN45lP9pMfvenfBIVwzAPB7Xbm0nJULjgG346Hduwy+n4t07jrT+PLnqmSwH2vhODTzZC7Gl2qm6fQ1Aa+SdqkmSyNwUIxTvqbwguvdTCldf8ZNj1nQ1TeeCqqQm3Y/Ed2YirYptldJRMcxFdARdudaC4KLVPpdG7s/+1Q8rB2FRA47sGdr3UQ9Ab3alTUiRybp/BO6H4eTuGojQljXxhQ23wYbObERo07G7D7Rx66EgAJpOC0WIg4A9RsDCbPek9tPuGHnJLt1j49qL3ybBtS9ybQKLMNp1UUYPQos9FFClCwOGdc/m/e4o5tDXyBGXWZAu5kzM4vCX6xI0/ePkuFl9+xoBltQeb+ecf32bvlmp+/9q3MFki93ru21rNP/7wFvf++baobRuOSK7fuvgYQ9SQSntrD21N3bQ1d9Pa5MRqN7Nk2TSyco8pWL8vyO9++DL/fWlL1PuqWFDC1jgODRQUpjD7zhfj1t9wGLrO5flXrfTEOUtmZprMx67swWepwaUNDMZNo5CNv51Ey/5jJ9mjP4qjMqgs08LhFxI9zj2Q8k/M563y8IbPFrWmUFOVmOnEJ/K1r9cydca/hlzXfGA+P76xIOE2nH1XLoELN0e9vRBgkIxkmIuo99ajDVHNWUYmy5xGd7B50LayOIeuToGxykz7W700b+oiHNFszU0idEsRu3rG5lgBpBhMGLsk3P4AsiRRnJ6KhERAVQmoKn41hD8YwhcK9Xt7Fk4rYE3P6AkDbQYD/zO/iUlpQycJjBcTLOXkiA8Tuo+hCPgsbH3rXt5+uorta8buxuPj37iCL/z8U8Ouj6VWi8vp4cVH3+Hmb14SrXnDEsn1Ww84TRBCCN5/czfvv7WH1iYnbU3dtLf2oA2RY/8PP3qFSdPysdpMdDs9dLb20OuKLbBtz6bDFMwuoD6CJDQj0dI8NmkrQ2nv8PGbrZg8c3h/fRo7D4x+N282gcUk0d07vI5efn6QNtOHQwb0dVFHxwWzCZ5lJ1ey0S0FaPT1IklQaE6mS/VzICBhnXAhPl8Qo0HGqCgYFRmjIqNIEgZZQkFCkfrGMhUBkhDIGkiaxpEkA33PQRURUhFBFS0Q6nv4g2j+ECFfENUXIOQNUllqgzCTONuzbVAVVtOYeWf1zGHFR8Ab23BfuBjsEOlgk1lOIqj5MckWkoz5CAQdgfYhhQeAhoYk2bHIKfi0Y78js7KY99tr+g70pL5HxpdTyHVloe2SqHu1k976wQHAWWdOYP8iE21jKDwAZiRls6W57w5cE4KqjtFnaPnCTGThCYX4yYYsbpz+Rc4p/CNynJ2WRsnCBEspBgwYjOcRCqyHMcxsa7L4WHLln7jktlW01Hp45M4nWfvCxoTuc+qiydz20A0jtomlSJyn10fx1PGvcquLjwTQVNfJH370CpvWHgh7m0N74+tSLJuez56W7pGLpkSA0SAjRNy6G5Gg8BK0rmfO+bDs/FlUVxZy8LDEwcNBQqqEQREU5BrIy5UoLAhiyt0BSBzePI9tu1WyMhQy0sCRAsnJKo40H6GU3SPOJHAHVerdvdQycPx9d6DvRGdGIcej4AvEcW6hETAawG4ALINWT8vJorozvAJgvdaxqzWybXs3zs7JpKYPHjf2u61jYoNh+LjQQQgBudYp1HgOk2cpoifYSb03PKXW7KvFKJkotM2k01+NJJXy3hB1QrpCPXRZe2ABsAAKlWxSWx14PghSvaqd3E+U8569i6A3shT8sZJptrK7OvJhAzHKUN+J/G13kCrnN7iwaAdBzUJZ2r8j3ueJSMjMsOUihfou9mowiaT0P4PWid/7EkH/e4QrzsPDhGIoQjaUoiglKIYSFEMpiqEUJDu5xUkYzYm9ZCal2rnn6a9hMCZuP7UHmymZlp+w/sNFFx9xYsOavWx67wC1la3s31VPIFET4Idh0lSJa2/tISXVi6oa2LwrmT0vxE8pzFiQPCbC4yh2JQu/sFIXPIBSfICpxTD9fANmKQmv1oM4oiSCQPDIdTd77ptcNHdwXy4YdQqjKzByXI2GID3NTmNL/IMDh6Os0E3hNA1NA+1IhkZNk1CPZGVUNYGVECnGICmGOuba/fj8Cl6vhM8PXq+M1wter4THCx4P+P0QU0yNEJgMCts3XcqSi1850lffQyBBioOpFw4cKhqQAEkaIiHSkZdaEDzVFmRFRjH0PWRJQpIl5CPbSBIgwPkmyB/ORzIAskCS+1KuS0fcTpJybHnqWRI1oi92o8kXWe0hgKAIUOXeR6phKhs6ehFh5PFsUFtpyGiFSyB7ZQbdzX6Ch8e+GN1kczpbQ5Hf2KhR5Cpd39DN+oYiFEnimwuvZHJ6bEMxJbYKpNCxoU4heunt/BzJ6X8hJeOvaGonAd9/8HtfJBTYSCwljiU5nbScrUjS8B4FVVXZvjqxQy93/t/t5BZnJ3QfLfWdzDsn8XFZo6GLjzjw4fsHeeBrfx9ySCXxCL5wZ4hzlr/E8VlBi0sl5s5Zyd0PFBAKxa4ass+oiXv1ciFAdc+ktXUKGWlNmFLX9QscqyGTFs/hAe1VEcIjnHG2os8OV2BkJ34QjcAUGakVxipKKjOrCS377RHbrEjNxa6u6XsxhPAajASY6MvZaAJMCNH30IQRTTv6MKCqRtb9eyabXxP4XV7cTg+uzl4C3gD/2WJGlOdD9vE5CASUbmHq96J4s4CEwmsXleDpiJ9wlySJM/+aCjGGjyUrZWzq9BHQIretNdAB6R1c5JjJG9vin/Z8OPJtyeyMwusBfcMz0aIKwYsHS/nWoqi7ACBJGup86sfVeStJab/FbL0ci/0mLPabUNVGAt5/4/e+hBrcHvG+hNZJnxdlaPHRWtfOTz71WzpjyB49Glf9zyWcedXChPV/FAkJRRn/ia7jb8FHGJ83wIt/W8eP7nx6XIRHbgH84R+VnLP87wxORy4onfgiX/li7HfqqalWtMLog1+PR6gWQq7ZaJpCW+1V/N/rE/j3VjdPvJ3CmrXX4us8F7uSi1uNLVFaJEjCFtbJdk+og+lzJ4yBRX0cPjR6/ESVXzA4+fNICPrGzHvoyx7ZiCRVI8sHMCi7MRm3YTFvwmZdj6+riRd/Ws3BDQep3V1PR0MnAW+fSPP2+tnynRI8/z4Tcbg0LoJMoDLn4/FNfHTGTWW4S2Or6ppkKGVLl4pPi21Kc52yk0vnJxHLHXokFMjJhLTozktajDZWOd30+kti6sOoDXfcgvR23Y7P84/+JYqSjzXp86RmvUpq9rtYk+9EMUyOaH+aOrRQe++FjXxxzp3seCc+pRSGYvIZpXzupzcnrP/jUQwnx2Vf93xEScD1PDs3fkiKJcTvH2vEkVJNZeVctm9JobVFoaVJoqleo7VZI95TR7/3yx6mTNuKLNcy2njC0vmvkvn989jwYSYvv2qKyhaXy4/iLEbLqOxfJtQkejuW4PXayczZjmw5PGg7TUhImhVQUD2TqKmbwoZKD+5AkBRLESbDQFsOtvYwo8RB0O9EkpwR2xktkgj/tnidpZHJGcm0dCT+Drb6sETR0pHbHPS2MCVtPubQ+3Hddyhg5JFvLsTnHn46aPWOZqp39P198QNLsV22Lub9TjxPYt0jMXfTz5a/V3FmyUx6L945euMhSFKK2d4lYhYeRzksdnH5ogr+vTGxsR8lSansqA4vXmgoVBHbzVRI03i3/nIuLfttVNvblFSEVjdCCw238xsIzYU1aeCUUcVQii35a1iTvooa2ovf+xIB78to6kj9gaY2oBgGJjX72w+f5/F7n4nqPYSLLcXK9579BiZz4usVCSF08fFRRYgQas+PEJ7HmXVC7q0pk5qYMunELQxoZNDdXUJVZS6bNth5902N+Utkrr2pmeqqVF5/xcKOLeH92GVFMHX6asAZpsWdTCl5gSklcP7Zl/L175YiRGQCRFU1NvyxgkU3ZKAVbaGpeiX/3enDH1LpKy1eTlnWAlKsMpoGKVaYUraLhoZprNp5/EX6mBemxzf0ydyWsZsIzYsZoYafyMgnVOQKM7yXePGx4sqOsAq3V/lNTItzws+Xf/dJKreHf/F668FWLq8oRiqpjmm/5pL4zgTRVI21PzjMGedaUK2RzSCzKxPZ4ZTxqPFNqV4Z2sP8khl8eDjO88mPI0Oz0ET0Xs9QHFxZr1b6OXdiAXZj5J6nLGM+iJHFAoCn5/uEgtsxWy7DaD4bST4WiSxJEgZjBQZjBbbkuwgFtxDwvojf+2/EcdPtZTkPk/UKZMNgT81Fnz6PN596h/oD0Qu50fjmX75EXmlOwvo/ns6WHjLzxr+oHOh5PiJCCC/Bzi8gAmtj7MlIn+u1b/y407mEz1xbHNaWxZMkHvz9k0RbBKKt60Ju/8bUiAUI9CUgmnVuIS/7GhAJij5dMctBRulLCel7OIRvFk98GFmmwCvas2BfJ0KREZKEKkmoMn3PUt+RDQJBIQgiCApBQBP4Na3vxB7G5/flb9XiMuwYtd2cpAmUSm9FZP9I7F23jJ9+JvL7kvzJmSx+fBdYor9YSyi8dmEJnq74Bmwv/PQklKsa8GWHJ25sSiG7u030huJfMweg0DKBN943kKiEetluCz3e6KekTsnPpNrUjTMY+bG8rDYdvw1UBeZY26i4+vGIthcCXnr/TuZPaOTSKX/BIIebntyM0bwUk2U5RvMFKIahh0iFCBEKbCDo34DRfC4G0xlI0vDegM7mLv79x/+y6/197F1/AJ8nPlN9jWYjn7rvWq7/zlVx6S8ctq49QMm0fFIzEpM5Vs/zkSiEQATi4d4eGM+gauEfhupDgsqDl1M2ObpI8qy0N7n15ok8+mTk+RhKirP48L/VnDc7n7WODoIJGLveWe3i/81PxyhpKJKKARUZFQ2F+qCN97o7Yh6PPhFNi3x66Fup7czdUEVrVfhZKg1HHpIkYbQYMVmMGM1GDGYDRrMRxWTof9gyrfhFeLMUUhXiUpAMoKc9j0e+5WBwDNHoNB5sp+E3i5nw7TVR71+gMvsT6az/U3yzf37w+CFsL1iY+/AEPCUjZ520Kfns7bbQG0pcbps6XwOXnDGD1zbH3/uRbDTFJDwADjS2M29yPmuDo3sfjudMVxq1rxwiFOwTDO8lmym/OBnFGp6nMKTa+GDvLaze18PqfUm8sON/+crZHzAv77UwtvYT9K8m6F8NgGKowGS5EKPlQgzGOf0CQ5IMGM1nYTSfFZZN6blpfOq+a/vsC4ao3FbNrrX72Ll2L7vX7sPZFnkupfOuW8pnf3wTOUVZEW8bC73dnoQJj0jRPR8RIIQg0LoQtPjW9vjJgzey7p3wrx5GE/z5Xx9gNOyNco8Sja0Xc/+Py+joDO/Oa+aMAnbuOuY+LZuaw6acHrxhJiMKl6ev20mebdWw653yBfyjIz6J047i7z2Hp7dF/jOY4TQR/N77CQk2vuRnxZjPrERLcTLCTVk/c5ImUCqvBxFbLMGfv/Vp1r0cfUpogKv+akGu2Bf19qYDi3jmk8NX8oyFWZeVotwzfPyHTc5jnyuZ7mDihkSOZ0JwJm/uiO8Q3qTkdOoPj55IbDTmTclnrfuI+BCCTGEmL2hhl7EbcSSbmFGTmOVLIadDgVoPtVsbUE9I83/1182UXPzckPsQAhrblpOX9V9kCf702p3sbhj8+z5vShJfWvwc+cnRBX1KciYm8wUYLRdgMp+DJMfnAiyEoOFgEzvf28uu9/exa+0+Gg8NP8OofEEZX/zlp5lx5vhMdV33+g6WrpiVsP719OoJQvO9Q7DrM3Hs0cDa91fys/sjr0B714/dzJz7z5j2Lsjip7/5BJu2jBwwMGtGATt2DR63nViSyd5iH041PsF4n13Qy00zHh613S71EtbFMUtkT9cKXtgd3UX7zAYTgUd34GxOTP6PZd+cTMptG8JqW2jJYKG5CqFFNz59YNPZPHTT4GRnkZCcbmf5izVgj/7inRws44nF8RW1R5EVmQUvGghkDL44W+UcDrpS6QqO3XRYGYl09wze2xe/fS5On8C2A7ElLSxRbKSuacfr8mOxm3D3+AgcqTE14aZy3szsZHFvKsbX2mhvcI7YV2ZuEiu+3EHG9I0Y7MfahlQzqzbdzn93uTlrSjKzi+v4y9sZR2LJBmNQZG6YL3HD7D9hN468z5ExYTQvwWi+AIv9k0hSfBPjdTR1sfuIENm1di+V26pJz0vjtodu5PwbzkKWxy/gc+Obu1l0YeIKheriIwEIIQh2XIsIRl9P4niCaimP/WUxrz4f3cd/36+dTJ4ae2xEMFTBzV88f9hcIFOn5rF/f9OwUylz8x3UT4WWGMfGC1JUnvr4X0G0jd5YslMvzuZNZxeBGKPyhYB1ey/nYGf04sEekjhvk8qBv26KyZajSJLE+XeVkXxWM1pxTUTJ3TKMyZxr2UqkmR9DQQMPXPNx6vbH5tW75P4JWC+PbWgyUXEfR8ksTGPabWm4l+1FGPv2YZGzqOrNoCMQX69aOBgkhY7KSTR2RT9UMictB8kDwZDKoZbYPbPn7JGo2T28iM2emEZrbWTelbzFBcy7zEfR7P8A8NTqm9hZH/nnnW638KWzGlle+gSKHNvv32BaSkr6Y0hy4soCeFxeZEXGYotkSnxi2LR6DwuWVSSsf118JAAh/ASap8Xcjy8wjX/9cy7PPB7eFFxJkZmysATZYkQg0ERfAqBFC1q5dPmvY7YHYM/BT3Dvg7kDlmVmJpGXm8q+/U0EgyPfhWZkJNEzx0x11HeMgn/duIk00zsRbSXJE3i9dyrV/hguGIGpPP5B7OOuS5pMtP/gvZj7AZhxUQllv4q+iNaVaSkYQh9EtM3qv1/Lk/fH5qovmJrFgic2IxliH4ZqemQB6/8UhhCNgfQ8BxW3ZqJd1EpVKI12/9hlrz2eJIMdpXMSGw5Gtv/ylAwyZCsI2FrVENfkd0u2qjQeis/nLxlk8j89g3Xdx+J4DIpMKMbhyml5yXz1rHeZPkoivtEwmBaQnP4Estx3vTns3oPDmEmqMRM5nDHPjwjdnb20NTqZNCNxhR/1gNME4PK4kUJmTIbYArmeeWou/3pGMJzwMJgMODKTsDms2DOSaHC62VE9+CTQ2mnj4gsgHh68isnP8ddHzuCxv19ER2carW09NDV3094enuu8o6OXlE0h8ufbaAxFNnyhSBpPXbs7YuEBILQGVthdbDAsYbu7M+LtAerbpwCxj4/7zPGbtVByfWxTOyv9SZRHMPW2pyOX534Z+6yOeV8NxkV4ABSdJ7P+T3Hpalg6m7rZ+vsgjilZtGc5E7uzYSizF1Pj7kJy1JKTkkVLT3hDmPPScqk83E5NKLrv/WiE/PEb9ipcOYX3ugcGEMcqPAD2Nrn44nNzuajiLL6w8G9k2ytH32gIQoFN9HR8kpSMpzjsaeDPVfcBYFOS+OqUX+AwZsRs68lA4+E2yqYnvuJ0uJw6si6BeANBfvr8B/xg1TfQYry7WHp+N1OXljHt3ClMWzaVSeeWU7C0jPQ5hRjKMnFn2GgUGoecbrZXttDeMbQAaGv3sL9qeWzGHIfVvBkBbN9ZR1MU8Qs9PT5scuRJci4uD44YYDoaQvSwyLSa8xyZkW+r2VhXG5+x9q2pPsrOjz2ILLsoHWlBbPUjdnsakOTwC0f961cX442xijKAIdcZcx9HMU6pRzElNuGLwWQg98EJtI6j8NjZ3UBPyEN3sJfJ0zyDEu8dT7LFwJJyBysW2nH2dA8bGxEPQnEsoNj+fj3ptsQVHHxjj5sbnrqGJ7d9A28oukBSNbidrvabeLnh0f5lHrWXfzc+Fi8zxwVN06jcXc9zD7/Fq0+tw2RJfCKzcDntPR8vrt9NVXMHN5w3l9y05AHrunq9PPvuNp59dztdvX13hvOLb2HljCeG7U8I8AWT8KrJJBk7B3hKml1F/P6pCvYfjo8785VVC5n65f/GreBba2tsWnRi0MqhCBMbTUiJR7BqgCnya6Smr+DFzvDvBD2uJXhC8QmWRZZY9XEHC849B+f33426m5k3pCMpsfrPJdZ5izjDXoQptJGR5uEe3rGINf+IbXbLUTSXNW53M26pgWv+vJh/3BJdbZJwkCSw7zOjFMioxujvxIWA7J5Msqoz0WoVhCpAA00FNIFsljDkSgSz/HSnd9Oe3InDlkJl78Cg6TpfC4vPyMOiOgj4ZbpdGh6/RkGOAdXSQ62vjkah0eiH3Hnp9K5NoceVmLiYgD9+5Q0CPT4KrEl0ehKTMwXAH1L58/sSL+38Gt9evp+FOc9FfF5c6y2nxV87YNnO7g0cdG1ncvLsOFqbWDpbutny3n62vLefre/tx3nEg33O5XPG17ATOG3FRyAY4if/XMPz7/dNu/vb6i3MLsknEFLp9flxefx0u32DaiP8c9MElpcnYzMeu2MOaQpbGy5k9f5ZrN4Nnb19d5FWk4HFU6xYjBrbqgU9nhCmg/ELaPvgwx5eW/0lLj3/YTQNenoLSU2JbF7+UYSAltbYTmQtO9phZkR75foZr8Wt1EUOG4Ap4TVWs1i1Pz77PZ4DqUFiqkmpxUdJNgW6+HcAck0LOcOuYA5tAAbeKWuqxFM/mATEJ6fGocezmPKzaKd/D0adsYnpl81m96uJGVoI+kPs+O0ucldlYbjLTnNe5MNvhS0T6PmzRHOLi2bCm2VkMDowfs+OJ3mwsKrxNQFNfT7pI3nvDgnghOt2c6CTaWdb2PqmQiAQ/7C9gC/6c0HB5BxSchy4WnvwmaFtsYMdHYmN3znKzDIT78rNdPmvY5n5P5ilwZ7NkDDjlXLwiEzcpOLR7Dg1K+u7BwoPk2xhWfbHKLbHHuuXaNoau3jx0XfY8u5+qvcP/T0sqRi7ulThcFoGnLb3uPnaIy+zqya6u6qMZAsfX+QlpMo0ddtYtz+E0z16LEhRZio96+KbPAmgsCCJXneI7m4/v/vZ62RnRHZVDYZMPPmPr/PaG7HXm2g7x0IggoxXb33mBRSiG6sdhDKJP7WHlzq4selK3qiMPdZjEEIw41vb8UeZBXH6BcVM+m18ZlQdT7bJwXy7CUtoPUdnwqz95zU8ek98Z3dc+YgDw/zIq4oOR5JWyDPnGgl4Elu4MTkjiawf5FJbPPRFUgjIcKeS3p6OqdlCqAHaD7poGyInxWhM/1ghW86Oj/ItMudhdKZxYJ+Pts44efGEoOif0edZyS3OpMHZp5ZSJ6Ygzsxgd3t80+afiFGRuWJ5DyHHse9eujGDKTY7blWiV1PpDflxhdz4tJHPcxIS89LOY0XuDaQYT45U5COx7f0DPPTlJ+jp7EuKZ7GZKJqSS1F5HsXleUyckkvxlDzSc1KQEpSZ+ih6wOkovLXtUNTCA6DD5eORNyX6btnDy4I4MTMV75bE3MHV1R+LC/F6I0sT3u2awC//8Gl2743PRSjPaKNmlARNl1WUc/XsCj7/zIu4QrNJNVQBBkymeShKIQI/CD9CBBCEQAQRIgQEjjwHj3sO9r8OSFmEO8W0KVGpHCSJ1PxUWg5Fd/Ku297OoPJAcaA10M0boXJqfZ9haXIXKXubePKB8MRm5rx0rN8twdoIB746WBgJIKXMgauym7pXUimZHz+7e+U6PvboQp75ZPxEe2ZBBjnlBSiyRHdzJzU7a3F19OK6/RClS4owXmMllG7A3mJHNCr01gVoreqmze2njdimsRpNCj2T45fArMbfBNYm5Lkw35yD3ZfO1k3emIZjzFJshYKaazqw5Djw+YI4a3sQtT2cdUkJe0zuhA29LF9gJeQYONusM9jBhu7IjleRbSpX5H+GAltZPM2LOz6PH4NRYcf6Sla/+CHXfG7ZEbGRS9aEtHHNJRIup6X4qGlNwB3vKPiCIYL+xIzPHk8wGFlg168fiZ/wACjRkqhh4MnVqCjMzMum3e3hk/NmceviM5AkiX999kbKszNBewAkC/JxRaGi5d5slYDmw695Cahe/JoXTaiII/9AIITG8rQkevwh3MHAkUcQd+DIczCAJxigNxjAEwz2Px9t6w0NfxzTAwpdDdF/v9InJI/eKEKEAOSzebezHVW0UOkGW0YBc/6ax84vHcLX4e+bfCVg6udLMTkMdO9ykbbIgbYkhMcs2NLZiDJRpvisPNrW9rl1hQRlt1dgmGqjsshJ+r1GOqq9hFtIXXMVEuwuxGBxIVs7wdwGsgrebFRvBpISRE7dT2DKB9z44nyeu6EzJg+IPc1O6fzJ7N9RT+fWYy72rGlF5BWk4e504Wl143g7g23724H45/yQZZnazMRkb63zt4DUwsypC/G4JLSQwOX2I0mQkmKhxeumy+fDajAwwZLMjsq+43jivXCa3YJkkBGhKD9rIcgvSKPqUGt//9WvHSbZolBxUTH75F7a3fGt6muasD6mCgNGyczHC29nlmNpwr0D0SCEoGpHDR++vp2Oxk7Ou24pFUvKmXdO3+OjiC4+xgBFlslxKdSPwQDXoeqpTCoJb9qqs6eQHbvie4Jt2dxKwVnJpCbbmJjm4PzJpSybXEqyZXCCnWk5R/JrKOlx278sKVgUOxbF3le/bxhKYsiurGoantARsXLk2RMKokgSZq+GmHURPrcPb68Pn9uPz+0/8vexZV5332uf24+v14fb5aanvZeKz0Se7XYkJKx0aUvZ5RwYC+TRvBxMq2LCPzKQkekUTpJkG83qkfHii6D7SOXkCfIUwIMqNPzfTifzW1kYNBnbm362XdAFdEEILA9kseDAwDtbTZMQPWVoPgchXzIBdxKtdWkcPBSktePoBSgTjkgWSWJAvoqczHnMmCmRnt/G/PunseXBzfickV24DGYjFWdXUHWwld2bawat72jpoaPl2O8gLYH3CH5fkNLGCeybONiOeFHX5eVQ2wnnuOMctN1Ac28vWfk2egMBStPSUV0h7DYTPkVjX1sbBf8ziYnb/DSvjtxO2SBjNA72nqg+leqXK7EbZaZeVMIhs4dmV+xeoHNmpKEZYnNlBoWfZl8tsxxLY7YnnvR0unj9/1bzxhNrKJyaz7nXnsnHv3E5iiHOZazHgdNSfNQmUHwYJImyJAdZmKA7QGtNF0KESJ8okzUlD03VUIMaitXI9v2xpUAeir88YUJVv8RlF46eptyRXEd6mpXOrthdoZmZSVxwXgUXXlBBWWn2SXn3EC8UWSbZZCbZNEzGwtLI+9za+Rzvtz4MVMdi2gAUuZh9vXk0+4cPQu7Qjv0WutXhTuDHBFH98XklLhjYqinkpG1qEcFdKxEa+NwWtm6WjhMZ/Xsa1p4TI9Ba2t20rAawAt1UfHIRdQ+vDrsWbPHMifQGJXZtqR298REMCZ6O6LXFpyrqcPR4w4v9aDvifdjd2opRkQl5j5VsrHH1UFMG8+fOQv1nDb214c1iS06zY8lNZf/e4YNv1aBG9auVGA0SZ11YwuEkHw3d0YuHkNr3vYn1lLO69XkyTDnMTz8/to7iwMEtVbz8h9fxuLwsvXIBv173ILbkxE1XHg9OO/ERDKk0diYmhfKZydlUr6uhK+AalLaqu2vgCTjFYaVoQho1Mbjoh+OpZ83MmzWDvOxdI7bbue9y3J7op9R9/Or55OSkUFaazayZhSjKyT/OeLIyzXExO7tepCcYXV2WE5HlRazrchOIQxHEYAQzcGqdKmv6SwYcE7WFBanYkozs3xfbrIc9DR2ccf0iDj6zcdS2JbOLaWn34gvzYnyUUKzJfEbAYjVSl5aYYRfouwh3+yLP2RIcJunXhz2tTP54PvxydPFhtplInphBXU143zkRElSvqkJWJM48v4i61BC1zshzDK3b28WljoswFL0R8bYnkmKMnxc2Gj58Yzt/vf8fdLe7+PojX2D2eYmrwzLenHZXi4aObtQEnFyWpGRxcE0VwUB4iX96ur20H+5gamlMEzOHJBhUueeHF/PhjmuGbbPuw0/y4M/L8EcRhyLLEpkZSay8ci4f/9gC5s4p0oVHjBhlK6mmgdkHFcnIjNQrkCO5RxAGPOIC3unoIqDFZ/aDe4QYlxOpsTQx69J0ZszNISPNTklxOrMuyaDu/EPsX7yHqRWxp7KXkkYvflc6t5Tmdk/EwgPAkMBgvbKVuWgx53EZHgU57ue3gz1dpFWMnMRPNshMmFEYtvA4HqEKav5bjfpcPWeq6ZSkpQ7b1mxQyElJwm4eODz5nw1ujL2xJfmTUSi2j0+1WYBtq3dx/zU/Y/Z50/nT9p+f0sIDIvR8PPzwwzz88MNUV1cDMH36dO69914uueQSAD796U/zxBMDE3AtWrSIDRvCq8o5FiQi3mNRShaVqw9HvF3AH+LwzgZmz53I9gPxueM9Sk+Pnxdfncr8Iaonv/HO5/jTYxZGSj51lOKiTObOKWLu7IkUFWXgcNhITrIgy6fusMpYowmVNxsfotZ9rDDd5ORlLMn+HCnGXDItk1nT/MtR+1GkHKr8U6j1RJfrZTi6AuHHWKiyyo7s/ZANzGbQ3BBTsQrRVUU/tg/DyOKg7IwyGppc/VVYI+4/EL8EW8eTmZ/CloUJSC5zHJqksaAsg/UH45tXw5A0fCySACYvmsSe3bElq5OAmrdrEMC8Gyaxpe2Yh8igyLx71+dJPZIpVdU0Drd3sbO+mZ11zeysb0YapQbVaBTaJmOSY6vqHC2Hd9XyzE9e5NfvP0jprKJxsWGsiUh8FBQU8OMf/5hJk/omAz7xxBOsXLmSrVu3Mn16n0q7+OKLeeyxYylpTab4BtDFSk2rM259TU1OJa1di0p4HEXTBHs21zB3fhFbh0kOEy0rLxl8oqtvOuOI8Bia/LxU5syeyLy5fYIjPT2GyEydsPCrLpq8fUNkWZYpnJ39ZfJtx7K1zUi9nHbfIXY5Xx62D0WewyYnuNU4u/SFRJs/fvOSd2TtZ9r0Cvbujv7iuLWxg4mfPZskowEZCaMQdLy2E2dNO5MXTKK2rptgDOnB3REGtIZLwBtCFjKalNicJSQ7495lyD1YkAlg0pyJeFQRs/A4HglwvlJLxeUT8GsqeT4Tvsoevv70q+Q5kslLTSY/NYX8tBTOLS/h6nl9155na3/DVmd0+7TINs7PGd5TnGiKKgp46LV7TulYuROJSHxcccUVA14/+OCDPPzww2zYsKFffJjNZnJzc4fa/KQgXp6PKckOOtfUE6/MHXu21JJfmkFjS/ziUf78RCkTCyrIyTx2q5mfs5lzz1rBe+t60TRBWpqduXMmMm9OEWfMLSY3N7I8ITqxYzWkcnnBQ7T5DjDVcRHSEJU0z8n5H/xaLwd7jlXwFAIMSiFebTIfdDQcmUocX1KMqdSKOF4sJfBWdENs5WuobR0YG5C0tJiZ10/n4BN7Yq5L0tbQhTHJGvYQarj0dHkor5vI3qLquPZ7IgHJB8R3NsSJ4qOoIh9hNXGgMjGZSwOuAIGn+27qqukTOt3FKhurBnr1SjLT+MftN2A3mzDK0ZWsL7BO4oair5NuyonN6Bj4KOTliDdRB5yqqspzzz2H2+1myZIl/cvXrFlDdnY2qampnHvuuTz44INkZ8c/riFa4iU+ckOmIxMR44OmCbprncyaOYEdcRqC6ezy8oOfXckD3/WQnloNgGLI5pt37OS2T19GSFtEcVHmaaW2T1YyLaVkWoafJiNLChfm3UVQ83LYtR5ZWcKOngDdQRdQnzC7bEoq8aj6ezxp7WlUEx8PTWqulazFDvYa2njR38ycn0zD+IsqvPXhJf8bCq/bz/SZhezaH//AUONmIyTYq+6TvECcPZZHro15pVnYsh0c2BdfL+1oSECmzUavf2AMz+H2Lu578U1+eu0lnJl5KV7Vxa7ujaMKcYcxk1J7BaVJ05mbeg6GMIti+kONmA3hF23UGZ6I06vv3LmTJUuW4PP5SEpK4u9//zuXXnopAM8++yxJSUkUFRVx+PBhvve97xEKhdi8eTNm89Cq1O/34/cfm3rW09NDYWFhwtKrL7/nT7R1R39igr4Yj8MxDLWMhM1uxpRlJxBSMSgKBoOMYpAxGhQUg9y/zKDIGIxK3/ORdcfa9C0zGPqeS4o0lp3dhtG8CFkp0cXGR5iQFuD/Dj/BW63xT8E+FBPM01jfHt8LzZy66Wx7K7aMpdmlySgLTezwNKGdcArLMNmY/qag5dXYcmlMWTSJ/QfiK0AMRgV+EMJtTlyRNYMwcHhzfL3P0xzpTOo0s+fFA4OmQ48V0icKqBtmNsz3V17AdQv7AtxafHW83fo8O5zv94uQDFMuJfYKSpIqKLFXkG6K7Ia42/chtc7f4g7sZVHhhiG9kzqRpVePWHwEAgFqa2txOp08//zz/OUvf+Gdd96hoqJiUNumpiaKiop45pln+NjHPjZkf/fddx/333//oOWJEh9t3b0IcaSWmTiS8/K4ZwR9uTAFyJKEJEtIHPlbkpAkkEMCRZb7X4fH6A0VRcJqM+vBnDoj4g55+P6un1PnjX+emBPJN89kQ3t8A1glTaLonclU10TnUZkwI40DkzvpHaEicbLBTPmfXLj2O6O0EuwpVoyZKXR1xTcGZPrKQracuz+cU0LUtG4rwh+K77ARwCJHLj3/rCEvPcCtn+0mNdVDfX0KBw9Y2bFdpvKQIFFvTL0mj6aeoZOSGRWFp79wHRUTjg2dtPkaaPQdptg+DYcxY9T+hVDRhA9N+I88B/CF6qnv/hNO37r+drPz/kmKec5A2zQ3imyP7o2dQiRUfJzIhRdeSFlZGY888siQ6ydPnsxnP/tZvvOd7wy5fqw9Hzo6pwId/i6+t+undAQSm6033TCDbV3xH9Yp7ylm/wvReSBLP57Dev/oScOmWjLgzv2IYPQxKyUz8jnc4kLEOevp9E8UsmVJ4gSId++UuKcwP8rsbCN/W/Ew0omldgGBjUCokO6eXHbtzuXJx610xKmmnPeqHDpGeE+F6Q7++eUbh8ymPBKVHT+gyfW3vjpSYVDouJ3itG/0v+7yvs/BjnuYl/9vDPLpHaAfifiI2XckhBggHo6no6ODuro68vLyht3ebDaTkpIy4KGjozMyGeY0vlH+hYTvpycYp0qpJ9Bgi76wo18J745+n6+D/C9Gn7fB4lAousfLokc9TLssvoHYu5+rY97qcgxqYtJkJ5kTl6V1e2uQ3T0XDrlOwoPZsJ/s9Hc4/+xneezRp3n08So+fp2KogwUgZIU2X1vYJR8M3Wd3dz9z9c51NqB2x/+9zbDflHYwgOgyfU0u1u+wP62b3Kw/W52t3wWCQlNRJ7c7XQmooDTu+++m0suuYTCwkJcLhfPPPMMa9asYdWqVfT29nLfffdxzTXXkJeXR3V1NXfffTeZmZlcffXVibJfR+e0xSwnfhp7hz9+FViPp9fgIy87n6bWyKfxmuskGN2LDoC3MLqLcFaFnUk/8lPPgb50OB9vZtmlFWz/maCzOj7p0Xe/UkfhnjxCNwaoS4tf1V4Auzmxyav/542JvP1xGxKjlaf3k532LrfcADffkM3BqgW8viqTyy5vpXTiagTJ+IPZvL16Nn/8vZGRXEEBdXTR+dbeSt7aWwlAssVMriOJXEcyuY5k8hzJA17nOpKxGA04zAswK/n41fCGMUNaF53et/pfm5QcZuY8iUkZORGbzkAi+oa2tLRw880309TUhMPhYNasWaxatYrly5fj9XrZuXMnTz75JE6nk7y8PJYtW8azzz5LcnL8K3Xq6JzudPgTO+RilMx0BxPjule0vuDoaIhkoLjF7BupvuCQlK/MwHJbA23aQGFUb9tD/v0WytZWsOnPzgh7HZrGyk4MP5Ip/+pE9heEX39mNBQlsXFj7pCERGQiTKaV8tJXKb/92DIJH1ZjG5ddtJsli5fwv/9bTt0QccICIo5hcfn8uHx+DrYMn3U1zWYl15HMhTOuYFbJk2giskBgg5zGzJwnsBgLI9pOJ0Lx8eijjw67zmq18vrrr8dskI6OTni0B+KVZWZokg1pQAJcyZrE1O1T2d0Y+d1+0fwMtmQ3hZOcF4A6fzeTzDKaP7wNFn07F+eZe3BrQ1/o/JqP4NLtTG2Yyb7/xCcnTyik0fqIi5xvpdGSErugVDBwqCXyGimRcHZhCIhvQGt6ynp+/5u93PO9q9i5Y+A6oykxw1NOj4eFZT0UZtYjYYEh4liGQ5GSmZnzODbTpP5lQmj6TJgwOe0Ky+nonCok2vNhkZOJu/gQMHv/dLZvj3wK68S5GWyf0EpgGGEwFEGhkXZGDh3rRp4uLBvgnN/nUJ+zk9FytWmoFFzmZ99/wjZjVDy9flIeS8PyJTc+U2xxNhPUUqoT5LE6ysL8xEwVlnBy/Q1d7NyRNmC5YonnpUqwbFqIi2Y1kOXYjCr64o/UCDxqsmRjRs6jJJmP1V9x+XdgNZZgkHRPfzjo4kNH5yNKoj0fEtFljByJclcx2zdGlzvDXR7C74582knKolRCdV0EfSpBr0bQHQJxbFgiKdfEvN+ZqTeMXAX6eNqTapBNeYiAIDc3hbQUC4rQ+oSLxLHn4zm6/OiL415LEgQ9fiz/z0rtVwVeS/T1ZQJOO4wSixErFemJ++5NK9+OPelc3L3HPkDjMOLDbBBkJgsykiDVppJq10ixqCRZgySZQ9jMQawmP2aDH5PBj0HxYjQ0EhJ96eAjERxHkTBRkf1HUizz+vrQvDT0PEZ20koMsi48wkUXHzo6H1GCWog0Y98sDNGfTkkcyWNz9Kx6XA6bAUsYlAXyaJtjSw0YJIWQiJ97Xd5vA8KfYivJkDvZAZKENRSdGMpf3o5y7sBcJTIyimRAkRRUodIoIvM2+DQPi//kpfOBNPavrST6uTuDKfxFHg3flPBYIveAKBjYXe+MozVDU5AU/+yvRzHK1Tzzt3q+9s2bqDzUt8zd6WFWVSoGo2Dh2V7mX/Q8mnAjIow7AQjFlFxCYVr2b0mzLu2zK3CAw10/ZVL697EYJsTS8WmHLj50dD6ifG3KZ8dsX5rQCAkNVdMICRVVaKjiyN/acX8fXa6pfe3F0WeNoBpCpJvQlgtUTaBpGpp27G8hjv7d9zqkCmoynDx6+ANUIYjiOnOEwcGXGhqaCBCM4ULUJhpxfEqGtdH3MRR1O5oo+GkOzd+S6LVG9qbHYsgFIMVQO+rwVCwE1YlUH9Y4mg3CbtNYuqSdBVe+RYjmqDwWsSNRnvlzMmwXANDpeYeqrh8xPfvPesBpFOjiQ0dHZ1RkScYkybFnBoqidtf5JZO54/1/0e6LLinZthYbeWmjt4uG7qJ6Zq2cyY6XquPab/3uFvIfykb+jkyPPfz4irEYcsm2CxQRT1/PYIxKFX94pJhKn4I9Zxch0edpiXOut7CxGoopSvs6WfbLAGjp/Rc1zl8zM+dxrMaJ42TVRxs9LFdHR+ekZmH2RF5ecStzM6Jza6ebJaQE5jIXN9RhtsY/50rj/lbSfwTJHktY7RWhsKfBGXc7TuS8idHHo0RCftbbzCpuQRm3iq8S6dYLmJHzGGdMeKNfeDR0/x/VXT8/IjyKx8m2jz66+NDR0TnpybUl8/QFN3HT5HkRbbc0z47dvmXUKqex4E11MvvTBQnpu/lQO1PWZ4XVNlXKwGZMXGbTo8zPTUziuaFwqDuYm/MXcpOuH7N9GuQ0ClI+z4IJq5me8whp1rP7p8+2uV+lvucvzMz9qy48YkQfdtHR0flIYFIUHph/MbPS8/nfTa+NOuV2TqaN3NQd+LXEO+udl+1m4nvF1O6KfyBmy2stKOdJqMrIAqqDVixT2zhLm8TOvRqZNsHHp3SRaw+QYfGRZvZiM/p4p34iz+6z4TCDSRasrR86h4aEIM0CDjOkWgQOMySbNWZnDp+0KxHIoV1MzvwhmfaLOdh+d9iZSMPoGUWyo8h9D5OSSU7SNWTZLkOWBwc3u/y7jng8nsRmLI2TDacvMReWizeRFKbR0dE5PdnV2cyX1j5Pg3voZFqTUy3MyT+AR408fXu02DsyqLxdodcZ/xwYFZ+cwparGsJubxN2zsoI8p3c50dt66MvSZaCH1nyI+MD4UNi6Nk2wjCDUHBb2LbEimz9OMbUnwIQ0lwc7voZPb7NKLLtmHjoFxFHlyWhSLYj62x9r09oL0tmpDDLkgdCrexu/RxTMn+O3TQ5kW/3I82YVrWNN7r40NHRCYcuv4evrnuJtc2HBywvTDJxVlE1PaHEJmEbisw1M9jw0BD5weNA+YpS2q8I0JA18vtKN8rcXRRgjvk9JJGYIZKQnI9Q45cOfiQk01JMGU+Nyb6GQgjB3rbbKXR8kWTz7HGz46PAmFa11dHR0RkP0sw2Hjv3Or5UsbR/WabFwLnFDeMiPAA8OxMX2Lr/9SocT43cv02GJ6fsZ67ptYQJDwBFyU1Y38cjmRZjTPvdmOxrONo9/yHJNFMXHnFGFx86OjofWRRZ5luzz+Phs64h12blssmddAXjWyE2XMyuZPa8Xjd6wxio2lDHlKrhLvyCP5Y7sYp9CbUBQAruBCk9cf0rJRhSf4cx/W9IcmrC9jMaQbWbdvcqCh1fGDcbTlV08aGjo/ORZ0VhOY8vu5T2QGIv/iOR9E4JoWB8i60NuZ8dQweIfr9UYoL0XsL330cQg7E8/t3KuRgcD2HMeh3FemnYMRmJoq77D5Sk34UkJaaw3emMLj50dHROCaaklHNryVcwSvHPuTEqqsShFxJba+cota/Vk+UcWEPkE9kmzrbEsdJdGEih/SCFl4Nk9M4cKMl3Ycp+G8V2HZI0/hMxD3XcS7J5lp42PUHo4kNHR+eUYW7aIr465X9JMaSO6X4zdk2moyGxZeyP0uv0oH2nm/LKvuGXmUlGvpD1JlKcS9yPiuhBMcYaB2FBsX8JU/Y7GJI+jxQvMRMH8lM+059YTCf+6OJDR0fnlKLIXsadUx+gwFo0ZvvsfWtsvS1up5eWH9ZR5knmJ0XbkMTYCJ8TkdUGhqqdMzoKsu0GTNmrMaR8C0k++WY22owl423CKY0uPnR0dE450kwZfG3Kvcx0RJYRNVoU49jHJvg8Ab5i6cAsqsZ83/1ozcjGORFtIlsux5j1BkbHD5GUKIr96JwS6OJDR0fnlMSsWPhs6de5IDvxrnNzzvicSnPzx2dK8fEoIrykapLpbIyZL2NM+w2yQfcqnO7o4kNHR+eURZZkriq4gRsmfg6ZxM1YkGzjk6vx7ZezEYxznIR6CMkw/MwXyTgbY/pTmDKeQDbOGEPDdE5mdPGho6NzyrMk8zzumHwXNsUe975NHhv7n2mLe7/h8MQvPHz5mnOpqlqGGMdSXYo0+HOVlFIMqX/AmPECsnnpEFvpnM7o6dV1dHROG1p9zfyx8me0+Ztj7ksOGEnbMJn2V4MYJDP7tlQPyEshhMBqN2NLtmC1m1FVjeaaxBVlmzbPzJ0/7SY3Z13C9jE8EpphNmpoL0jJGJK/imy95qSYMqszdkRy/da/GTo6OqcN2ZZcvll+P49W/ZqDvXui6iPj9el49kPVu83s7zpWVyY1K4WUdDsBXxCPy4e7x4P3yAPAbDUxY1EpmoDDuxvwuv1xeU9HCYUEKQ5nXPsMH4Ec2oZs+xQk3YksnzxTZnVOTnTPh46OzmmHKkI8W/sY6zvWRLRdalUhmz7XHvP+TRYjk2YX4XH7qd7bGHMmz5u+buf6W99Boidm2yLGOKevCm7S7UiWi8d+/zonDbrnQ0dHR2cEFMnAJyd+lhxLPi81PI0gvHsw19/tQOziI+ALsmfjIQDySrLIzE/n8J5G3K7wZo4cxWKT+NnfNUrLXonZpqgxzkNOuWv89q/zkUQXHzo6OqclkiRxQc5lZJtzebz69wS04YdBkhtz8D6Wxt7VNXG3o+lwG02H2zCaDUw7oxivJ8DhPQ2jekOmzC2ifKZGft5zcbXHFXLgDGZRaD00ckPJjpTyfSTrVXHdv87pgT7soqOjc9pT76nhkcqf4wwOXZ8l5ckZbH6icszsyZmYSXZhBgLQVI1QUCXgD+L3BsjIcdBS005bQ5+tucVpfOUHLcyZvwGAkCrT480h2dKB0RAY1LcqZH5b9XmSjQGuyXuZfzZeSY7FSa6lnRfqzuCVaj9ZNjOvnvNnbIbeoQ00zEBK/SWSoThBn4DOR5FIrt+6+NDR0dEBuoNd/Knyl9R6BmYMTf/PdDb+YhyziIbJ8usmMKHYz8aucjbvayc5yczZ842cN2838ya9jUEJARDUDCx943Y6fD5STCZ6AoMFCsD5BSn8bu7vsSgnDAXZbkNK/jrSeBTw0zmp0cWHjo6OThQEND9PVT/CVudGAMw9yey7VhAKjnHRtiiZcuU8tjndg5Y7ks2cs9BA2ZxO3nBns6beg08NjdrfwpwU/rTgL6QYnCBnIDl+imQ+OwGW65wK6OJDR0dHJ0o0ofFq4z/Z/O99ND0doGZ77DlBxoLMKbm0FmXi8weHbZPqsLJvaYAAWtj9Ts9I4a/n7CUt824kJSsepuqcokRy/dYznOro6OgchyzJXDHhWs5yLKdxb+wzW8YC2aggzS0eUXgAOLu9TJIcEfVdlpZPStYvdOGhE1d08aGjo6MzBBdev4SH/vUNUtKTxtuUUZl83WJqmsIrMpd1QEEK0+F99eQKfnn+pRhk/VKhE1/0b5SOjo7OMMxYOoVf/fe7TJh08pZ+LzlvGltrw/fQHD7UxvLW/FHbXTNlOj9fdrEuPHQSgv6t0tHR0RmB/NIcfvXfu5lzztTxNmUQSdkpNCSZI95u35ZGzgkNL0CunTqDn563AkUXHjoJQk8ypqOjozMKyal2fvjPr/Hbbz7F639dOy42SLJEamEGSbkOLJnJKEkWumRobQhvuOVEzAF5yCvAJ6fN4sFzliPHmPJdR2ckdPGho6OjEwYGo4Gv/foWCibl8n/3Pc9YThR0FKTRM6OABs+RnBwub98jBhS/ANvAZZeXlevCQ2dM0H1qOjo6OmEiSRKf+MrF/O+TX8JsHZskW5JBJve8CtyeoZOBRYvDMHi45uPlM3ThoTMm6J4PHR0dnQg58/J5/OzVb3PfDb+ls7k75v7y5xaRNL0ATZaQuz00bqikt6Wb1IkZWM6cElFAabg4TBYcZgvdfh8AsiRxRu7ogag6OvFA93zo6OjoRMGUucX8+r/3UDqjMKZ+LKk22ieksa22nR3VbWzrctM2NZecTy7GNaOAqoah683EitAESyccs316ZjbJpsiDV3V0okEXHzo6OjpRklWQzi9e+w6LVsyKuo+ChWU4ewbGbwgBhxu6cPX6YjVxWFRNcFZBcf/rRXkFCduXjs6J6OJDR0dHJwasSRbu/dsdXP2lC6PaXs6JLONovNA0wVkTivpfL8yLzYOjoxMJEYmPhx9+mFmzZpGSkkJKSgpLlizhtddeG7LtF77wBSRJ4le/+lU87NTR0dE5aVEUmS/86Hru+PmNyEpk93TtvvgGkoaLpmkUOVIpTHZgMRhYkDdhXOzQOT2JKOC0oKCAH//4x0yaNAmAJ554gpUrV7J161amT5/e3+7FF19k48aN5OfrwUs6OjqnD5fftozc4ix+9JlH8IQxFTZv9kQqW2IPWI0G7chU4aevvJYcWxJGRRkXO3ROTyKS6FdccQWXXnopU6ZMYcqUKTz44IMkJSWxYcOG/jYNDQ3ccccd/O1vf8NoNMbdYB0dHZ2TmfkXzOCXr99FzsTMEdulFWfSNTFjjKwajKb1iY+CZIcuPHTGnKhjPlRV5ZlnnsHtdrNkyRKgz4138803861vfWuAJ2Qk/H4/PT09Ax46Ojo6H2WKp03g12/ezdQFpYPWOQrSmHLzWbSWZdPZ7RkH6/o4Kj50dMaDiMXHzp07SUpKwmw288UvfpF//etfVFRUAPCTn/wEg8HAV77ylbD7e+ihh3A4HP2PwkI96ElHR+ejT2pWCj956U7OvXrBgOW5505jW00bwaA6Tpb1Ict6MjGd8SPiJGPl5eVs27YNp9PJ888/zy233MI777yD1+vl17/+NVu2bEGKIEPed7/7Xb7xjW/0v+7p6dEFiI6OzimB2WriO3/5HBMm5fD3n/0bAfSOt1FATlYK/3Pb+eNths5pjCRiLFBw4YUXUlZWxrRp0/jGN76BfFwVRFVVkWWZwsJCqqurw+qvp6cHh8NBd3c3KSkpsZimo6Ojc9Lw5jPrefW/O9jcPr5Dy5npSfz2R9czIS9tXO3QOfWI5Podc3p1IQR+v5+bb76ZCy8cOM99xYoV3HzzzXzmM5+JdTc6Ojo6H2kuvH4JGdMnsP+hF+l1+8fFhvRUG7/6wXW68NAZdyISH3fffTeXXHIJhYWFuFwunnnmGdasWcOqVavIyMggI2Ng5LbRaCQ3N5fy8vK4Gq2jo6PzUWTuzIn84Sc38u0H/klza+I9IFaLkSS7hSS7meQkC3d+aTkTC9ITvl8dndGISHy0tLRw880309TUhMPhYNasWaxatYrly5cnyj4dHR2dU4riwgz+8JMbueuHL3CgsmXYdhazEbvNdEQ8mLDbzNhtZpLsfQ+7zYzdbiZpiOVJdjM2qwklwoRnOjpjRcwxH/FGj/nQ0dE5HfB4A7y7/kC/qLDbzSQfFRU2EwaDnntD56PFmMZ86Ojo6OhEjs1q4uLzZ4y3GTo644Luk9PR0dHR0dEZU3TxoaOjo6OjozOm6OJDR0dHR0dHZ0zRxYeOjo6Ojo7OmHLSBZwenXyjF5jT0dHR0dH56HD0uh3OJNqTTny4XC4Avb6Ljo6Ojo7ORxCXy4XD4RixzUmX50PTNBobG0lOTo6oQN1YcLToXV1dnZ6D5CRCPy4nJ/pxOTnRj8vJx6lyTIQQuFwu8vPzB9R5G4qTzvMhyzIFBQXjbcaIpKSkfKS/IKcq+nE5OdGPy8mJflxOPk6FYzKax+MoesCpjo6Ojo6Ozpiiiw8dHR0dHR2dMUUXHxFgNpv5/ve/j9lsHm9TdI5DPy4nJ/pxOTnRj8vJx+l4TE66gFMdHR0dHR2dUxvd86Gjo6Ojo6MzpujiQ0dHR0dHR2dM0cWHjo6Ojo6Ozpiiiw8dHR0dHR2dMUUXH2Fy4MABVq5cSWZmJikpKZx55pmsXr16yLYdHR0UFBQgSRJOp3NsDT3NGO24bN++nU9+8pMUFhZitVqZNm0av/71r8fR4tODcH4vtbW1XHHFFdjtdjIzM/nKV75CIBAYJ4tPbdasWYMkSUM+Nm3a1N9u06ZNXHDBBaSmppKWlsZFF13Etm3bxs/wU5xwjwvA448/zqxZs7BYLOTm5nLHHXeMk9XxQRcfYXLZZZcRCoV4++232bx5M3PmzOHyyy+nubl5UNvbbruNWbNmjYOVpx+jHZfNmzeTlZXFU089xe7du7nnnnv47ne/y+9+97txtvzUZrTjoqoql112GW63m7Vr1/LMM8/w/PPP881vfnOcLT81Wbp0KU1NTQMen/3sZykuLmb+/PlAXz2OFStWMHHiRDZu3MjatWtJSUlhxYoVBIPBcX4HpybhHBeAX/7yl9xzzz3cdddd7N69m7feeosVK1aMo+VxQOiMSltbmwDEu+++27+sp6dHAOLNN98c0PYPf/iDOPfcc8Vbb70lANHV1TXG1p4+RHJcjuf2228Xy5YtGwsTT0vCOS7/+c9/hCzLoqGhob/N008/Lcxms+ju7h5zm083AoGAyM7OFg888ED/sk2bNglA1NbW9i/bsWOHAMShQ4fGw8zTjqGOS2dnp7BarSOe0z6K6J6PMMjIyGDatGk8+eSTuN1uQqEQjzzyCDk5OZxxxhn97fbs2cMDDzzAk08+OWpRHZ3YCfe4nEh3dzfp6eljaOnpRTjHZf369cyYMYP8/Pz+7VasWIHf72fz5s3jZfppw8svv0x7ezuf/vSn+5eVl5eTmZnJo48+SiAQwOv18uijjzJ9+nSKiorGz9jTiKGOy3//+180TaOhoYFp06ZRUFDAtddeS11d3fgZGg/GW/18VKivrxdnnHGGkCRJKIoi8vPzxdatW/vX+3w+MWvWLPHXv/5VCCHE6tWrdc/HGDDacTmRdevWCaPRKN54442xM/I0ZLTj8rnPfU4sX7580HYmk0n8/e9/H0NLT08uueQScckllwxavmvXLlFWViZkWRayLIupU6eKmpqacbDw9GSo4/LQQw8Jo9EoysvLxapVq8T69evFBRdcIMrLy4Xf7x8nS2PntL49v++++4YN9jn6+PDDDxFCcPvtt5Odnc17773HBx98wMqVK7n88stpamoC4Lvf/S7Tpk3jpptuGud39dEnnsfleHbv3s3KlSu59957Wb58+Ti8s4828T4ukiQN2ocQYsjlOkMT7jE5nvr6el5//XVuu+22Acu9Xi+33norZ555Jhs2bOD9999n+vTpXHrppXi93rF8Wx954nlcNE0jGAzym9/8hhUrVrB48WKefvppDh48OOykh48Cp3V69fb2dtrb20dsU1xczPvvv89FF11EV1fXgHLHkydP5rbbbuOuu+5izpw57Ny5s//EKYRA0zQUReGee+7h/vvvT+h7OZWI53E5yp49e1i2bBmf/exnefDBBxNm+6lMPI/Lvffey0svvcT27dv713d1dZGens7bb7/NsmXLEvY+TiXCPSYWi6X/9Q9+8AN++9vf0tDQgNFo7F/+6KOPcvfdd9PU1NQ/bBwIBEhLS+PRRx/l+uuvT8ybOAWJ53F57LHHuPXWW6mrq6OgoKB/eU5ODj/84Q/53Oc+F/83MAYYxtuA8SQzM5PMzMxR23k8HoBBcRyyLKNpGgDPP//8gLuDTZs2ceutt/Lee+9RVlYWR6tPfeJ5XKDP43H++edzyy236MIjBuJ5XJYsWcKDDz5IU1MTeXl5ALzxxhuYzeYR43V0BhLuMTmKEILHHnuMT33qUwMucNB33GRZHuB5Ovr6+N+TzujE87iceeaZAOzfv79ffHR2dtLe3v7RjsUZtwGfjxBtbW0iIyNDfOxjHxPbtm0T+/fvF3feeacwGo1i27ZtQ26jx3wknnCOy65du0RWVpa48cYbRVNTU/+jtbV1nK0/dQnnuIRCITFjxgxxwQUXiC1btog333xTFBQUiDvuuGOcrT+1efPNNwUg9uzZM2jd3r17hdlsFl/60pfEnj17xK5du8RNN90kHA6HaGxsHAdrTx9GOi5CCLFy5Uoxffp08f7774udO3eKyy+/XFRUVIhAIDDGlsYPXXyEyaZNm8RFF10k0tPTRXJysli8eLH4z3/+M2x7XXyMDaMdl+9///sCGPQoKioaP6NPA8L5vdTU1IjLLrtMWK1WkZ6eLu644w7h8/nGyeLTg09+8pNi6dKlw65/4403xJlnnikcDodIS0sT559/vli/fv0YWnh6Mtpx6e7uFrfeeqtITU0V6enp4uqrrx4wJfqjyGkd86Gjo6Ojo6Mz9pzWs110dHR0dHR0xh5dfOjo6Ojo6OiMKbr40NHR0dHR0RlTdPGho6Ojo6OjM6bo4kNHR0dHR0dnTNHFh46Ojo6Ojs6YoosPHR0dHR0dnTFFFx86Ojo6Ojo6Y4ouPnR0dHR0dHTGFF186Ojo6Ojo6IwpuvjQ0dHR0dHRGVN08aGjo6Ojo6Mzpvx/odxoeG7hXlYAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "gdf = s.to_pandas().reset_index().set_geometry('destination')\n", "gdf.plot(0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The issue is that we have lost the CRS along the way. CRS is stored on a GeometryArray level but when we create DataArray, xarray uses only the underlying numpy array and ignores our custom attribute. " ] }, { "cell_type": "code", "execution_count": 36, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 36, "metadata": {}, "output_type": "execute_result" } ], "source": [ "gdf.crs is None" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## What needs to be done\n", "\n", "It seems that the it will not be very complex to implement convenient vector data cubes. We need to figure out:\n", "\n", "- CRS handling\n", " - pyproj.CRS class storing the info on the GeometryArray level should be stored on a DataArray level as well (per dimension) and retained all the way back to a geopandas object.\n", " - It may require some custom logic during I/O and figuring out in which case we lose the info. Hopefully, most of it should be doable using xarray attrs. But some custom `to_geopandas` method may be needed as `to_pandas` still uses geometry as an index and that cannot contain CRS. So we need to convert to DataFrame/Series, reset index and create GeoDataFrame to which we pass a CRS stored on a DataArray level.\n", "- Figure out which geometric ops need to be available on an xarray level\n", " - I have to check stars docs to see how the actual data cube interacts with the vector data apart from using them as an index.\n", "- ... and probably more stuff I missed now." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3.10.6 ('shapely2')", "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.10.6" }, "orig_nbformat": 4, "vscode": { "interpreter": { "hash": "716119b7942b628d25d105f893f626cc168d483792084e8427e38a31c13f35ec" } } }, "nbformat": 4, "nbformat_minor": 2 }