{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Averaging over a region\n", "\n", "Although this may not sound like _real_ regridding, averaging a gridded field\n", "over a region is supported by `ESMF`. This works because the `conservative`\n", "regridding method preserves the areal average of the input field. That is, _the\n", "value at each output grid cell is the average input value over the output grid\n", "area_. Instead of mapping the input field unto rectangular outputs cells, it's\n", "mapped unto an irregular mesh defined by an outer polygon. In other words,\n", "applying the regridding weights computes the exact areal-average of the input\n", "grid over each polygon.\n", "\n", "This process relies on converting `shapely.Polygon` and `shapely.MultiPolygon`\n", "objects into `ESMF.Mesh` objects. However, ESMF meshes do not support all\n", "features that come with shapely's (Multi)Polyons. Indeed, mesh elements do not\n", "support interior holes, or multiple non-touching parts, as do `shapely` objects.\n", "The `xesmf.SpatialAverager` class works around these issues by computing\n", "independent weights for interior holes and multi-part geometries, before\n", "combining the weights.\n", "\n", "Transforming polygons into a `ESMF.Mesh` is a slow process. Users looking for\n", "faster (but approximate) methods may want to explore\n", "[regionmask](https://regionmask.readthedocs.io/) or\n", "[clisops](https://clisops.readthedocs.io).\n", "\n", "Also, note that low resolution polygons such as large boxes might get distorted\n", "when mapped on a sphere. Make sure polygon segments are at sufficiently high\n", "resolution, on the order of 1°. The `shapely` package (v2) has a\n", "[`segmentize` function](https://shapely.readthedocs.io/en/latest/reference/shapely.segmentize.html)\n", "to add vertices to polygon segments.\n", "\n", "The following example shows just how simple it is to compute the average over\n", "different countries. The notebook used `geopandas`, a simple and efficient\n", "container for geometries, and `descartes` for plotting maps. Make sure both\n", "packages are installed, as they are not `xesmf` dependencies.\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import geopandas as gpd\n", "import pandas as pd\n", "from shapely.geometry import Polygon, MultiPolygon\n", "import numpy as np\n", "import xarray as xr\n", "import xesmf as xe\n", "import warnings\n", "\n", "warnings.filterwarnings(\"ignore\")\n", "xr.set_options(display_style='text')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simple example\n", "\n", "In this example we'll create a synthetic global field, then compute its average\n", "over six countries.\n", "\n", "### Download country outlines\n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
idnamegeometry
5032ArgentinaMULTIPOLYGON (((-67.19287 -22.82225, -67.02727...
9156ChinaMULTIPOLYGON (((77.79858 35.49614, 77.66178 35...
37710South AfricaMULTIPOLYGON (((19.98200 -24.75230, 20.10800 -...
67724SpainMULTIPOLYGON (((-5.34065 35.84736, -5.37665 35...
98466MaliPOLYGON ((-12.26352 14.77561, -12.13752 14.784...
155484MexicoMULTIPOLYGON (((-97.13797 25.96581, -97.16677 ...
\n", "
" ], "text/plain": [ " id name geometry\n", "5 032 Argentina MULTIPOLYGON (((-67.19287 -22.82225, -67.02727...\n", "9 156 China MULTIPOLYGON (((77.79858 35.49614, 77.66178 35...\n", "37 710 South Africa MULTIPOLYGON (((19.98200 -24.75230, 20.10800 -...\n", "67 724 Spain MULTIPOLYGON (((-5.34065 35.84736, -5.37665 35...\n", "98 466 Mali POLYGON ((-12.26352 14.77561, -12.13752 14.784...\n", "155 484 Mexico MULTIPOLYGON (((-97.13797 25.96581, -97.16677 ..." ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Load some polygons from the internet\n", "regs = gpd.read_file(\n", " \"https://cdn.jsdelivr.net/npm/world-atlas@2/countries-10m.json\"\n", ")\n", "\n", "# Select a few countries for the sake of the example\n", "regs = regs.iloc[[5, 9, 37, 67, 98, 155]]\n", "\n", "# Simplify the geometries to a 0.02 deg tolerance, which is 1/100 of our grid.\n", "# The simpler the polygons, the faster the averaging, but we lose some precision.\n", "regs[\"geometry\"] = regs.simplify(tolerance=0.02, preserve_topology=True)\n", "regs" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
<xarray.Dataset>\n",
       "Dimensions:  (x: 180, x_b: 181, y: 90, y_b: 91)\n",
       "Coordinates:\n",
       "    lon      (y, x) float64 -179.0 -177.0 -175.0 -173.0 ... 175.0 177.0 179.0\n",
       "    lat      (y, x) float64 -89.0 -89.0 -89.0 -89.0 ... 89.0 89.0 89.0 89.0\n",
       "    lon_b    (y_b, x_b) int64 -180 -178 -176 -174 -172 ... 172 174 176 178 180\n",
       "    lat_b    (y_b, x_b) int64 -90 -90 -90 -90 -90 -90 -90 ... 90 90 90 90 90 90\n",
       "Dimensions without coordinates: x, x_b, y, y_b\n",
       "Data variables:\n",
       "    field    (y, x) float64 2.0 2.0 2.0 2.0 2.0 2.0 ... 2.0 2.0 2.0 2.0 2.0 2.0
" ], "text/plain": [ "\n", "Dimensions: (x: 180, x_b: 181, y: 90, y_b: 91)\n", "Coordinates:\n", " lon (y, x) float64 -179.0 -177.0 -175.0 -173.0 ... 175.0 177.0 179.0\n", " lat (y, x) float64 -89.0 -89.0 -89.0 -89.0 ... 89.0 89.0 89.0 89.0\n", " lon_b (y_b, x_b) int64 -180 -178 -176 -174 -172 ... 172 174 176 178 180\n", " lat_b (y_b, x_b) int64 -90 -90 -90 -90 -90 -90 -90 ... 90 90 90 90 90 90\n", "Dimensions without coordinates: x, x_b, y, y_b\n", "Data variables:\n", " field (y, x) float64 2.0 2.0 2.0 2.0 2.0 2.0 ... 2.0 2.0 2.0 2.0 2.0 2.0" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Create synthetic global data\n", "ds = xe.util.grid_global(2, 2)\n", "ds = ds.assign(field=xe.data.wave_smooth(ds.lon, ds.lat))\n", "ds" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAADnCAYAAADvoXoeAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAABTYklEQVR4nO2deZwcR3n3v0/3zGp1H5Ys2ZJtSb4wvjDYxtgcNpfNFQMvJBwhBxAnBnJDAknIATlI8kJeknD5Bb9AwhECdnCIwTgEYwwGfOBbPmTJRrJkyZIlWddqZ7qe94/urq7u7Z6d2Z3dmR3VV5/9qKanp7u6++mqp371VJWoKh6Px+PxtCLodQY8Ho/H0//4ysLj8Xg84+IrC4/H4/GMi68sPB6PxzMuvrLweDwez7jUep0Bj8fjmWlcfNFc3flk1Na+t9116DpVvWSKszTl+MrC4/F4OmTnkxE/ue7YtvYNj3po6RRnZ1rwlYXH4/F0iAIG0+tsTCu+svB4PJ4OUZSGtidDDQq+svB4PJ4J4FsWHo/H42mJokSH2VRJvrLweDyeCWDwlYXH4/F4WqBA5CsLj8fj8YyHb1l4PB6PpyUKNHyfhcfj8XhaoaiXoTwej8czDgrR4VVX+MrC4/F4OiUewX144SsLj8fj6RghQnqdiWnFVxYej8fTIXEHt68sPB6Px9OCeJyFryw8Ho/HMw7Gtyw8Ho/H0wrfsvB4PB7PuChCdJitSu0rC4/H45kAXobyeDweT0sUYVTDXmdjWvGVhcfj8XRIPCjPy1Aej8fjGQffwe3xeDyelqgKkfqWhcfj8XjGwXSpZSEixwCfB1YQTzl1hap+tLDPe4A3Jx9rwCnAMlV9UkQeAfYCEdBU1bO7krECvrLweDyeDok7uLtWfDaB31fV20VkPnCbiFyvqvfZ86n+PfD3ACLyKuB3VfVJ5xgXqeqObmWoDF9ZeDweT4d0s4NbVbcCW5P0XhFZB6wE7qv4yRuBL3Xl5B1weIluHo/H0yUilbb+gKUicqvzd1nVMUVkNXAW8OOK7+cAlwBfczYr8G0Rua3VsSeLb1l4PB5Ph3Q4gntHO/0IIjKPuBL4HVV9qmK3VwE/KEhQF6jqFhE5ErheRO5X1RvbzVy7+JaFx+PxTACjQVt/7SAideKK4guqelWLXd9AQYJS1S3J/9uBq4FzJ3RB4+ArC4/H4+mQeCLBoK2/8RARAT4DrFPVj7TYbyHwAuDrzra5Sac4IjIXeClwz+SurhwvQ3k8Hk+HKEKje9N9XAC8BbhbRO5Itv0RcCyAqn4y2fYa4Nuqut/57XLg6ri+oQZ8UVW/1a2MufjKwuPxeDpEla4NylPVm2D8QRuq+lngs4VtG4Azu5KRcfCVhcfj8XSMdG1Q3kzBVxYej8fTIUr3WhYzBV9ZeDwezwTwix95PB6PpyWK+MWPPB6Px9MaBRrdmxtqRnB4Xa3H4/F0BfHrWXg8Ho+nNQptj84eFHxl4fF4PBPAtyw8Ho/H0xJV8S2LmUg4b67WlizpdTY8Hs8MYHTT5h2qumwyx4g7uLs23ceMYCAqi9qSJax89+/0OhtwmIXSVSLa6xz0F94u+somNv72ux+d/FH8Gtwej8fjGYe4g/vwcgJ8ZeHxeDwTwI/g9ng8Hk9L/AjuGcxUPLcpk1n7R75tnyl6L2bi++btwuEwtgvjWxYej8fjaYUqNIyvLGYmZZ7IJD21Mu9mQl5lq99MMI+d5mNCnpr7m7LzTeCYU+IxdnpMbxfZabptFxN8vn1hFx0Qy1C+spiRlBpbhbFUvlBtvGgqHb6QxX0rfjuZPLVD5fEr7pFKxbmLBUWbL2RHhUGrPHWDXttF2X6DYhda+DwOfWUXHeJHcHs8Ho+nJT50dqYiQDCOq+V8nXvG7ocq77CwLf1JS09Sy9NjftMNL7NDtEpGcG9FRashlwXXy2zx3lS+U4Xt5a3DqptRfb6O8HaRHb4bduG2MMZ5Rn1tF+PSPRlKRI4BPg+sAAxwhap+tLDPhcDXgY3JpqtU9QPJd5cAHwVC4NOq+qGuZKzAYFQWQGVgQlVhZrdreYHhWKi4+5Pfd9xCRLPTVH1nz5P7rqrNXb65LaoKs3SDVhQY+XLTuUe0LBzGKwxy3xczJRXpqn06ZbrtYrxKYlDswv2u4h73tV10QBfX4G4Cv6+qt4vIfOA2EbleVe8r7Pd9VX2lu0FEQuBjwEuAzcAtInJNyW8nzeBUFh6PxzNNxNFQ3ZkbSlW3AluT9F4RWQesBNop8M8F1qvqBgAR+TJwaZu/7YjBqCxEO5Kh4s9S/p2O3ahI3odox4uralEUvcqKfLQje0yYMZ6e4y2XeHWtAl7GyA9tnNNesnuRQrXXWHYzpkCGij9PvV1I2XEH0S7aeUb9ahfj0OGgvKUicqvz+QpVvaJsRxFZDZwF/Ljk6+eIyJ3AFuDdqnovcaWyydlnM/DsdjPWCYNRWQAalukBxc+O8auO3a+gU7sFhCZWOOEQSbfyKCkIxNlnzHm6UUAU7bryBUwOL9mHWKcuKSAm8GKqZMeqLAhEywsRm5kCnebD24VzwBafJ2IXEyys+8IuOqQDGWqHqp493k4iMg/4GvA7qvpU4evbgeNUdZ+IvBz4D+BEyq+yW+5DjsGoLAQkLHnJXRxBVTX+y3+XFBRlBYRz0Jw36Wi5uZfa9R61+FlKC4JioTCmZVJGOybh9imUHUKKXqPzf1kBgZb2fRYrEZWKz+4L7xYKzvZcweGcpMy7HZORVkyzXRRtwGa71E4Gzy7KHIu+tIsJoHQ3GkpE6sQVxRdU9aox53MqD1W9VkQ+LiJLiVsSxzi7riJueXSdwagsPB6PZ5rpYjSUAJ8B1qnqRyr2WQFsU1UVkXOJQ3p2AruBE0VkDfAY8AbgTV3JWIHBqSzcPosqz8pp5bqeomrmWZZ6kyavxrYlPRS9x/S4hvG9zMI1lEbMlH0uo8z5cTy6Yiikut5amTcZCGK96XHO4x6zTGIIKrxG16stSA85973N82cZab192uyiRHoaFLto1ybsMfvBLiaAqtDs3gjuC4C3AHeLyB3Jtj8Cjo3PpZ8EXgdcLiJN4CDwBo0NtCki7wKuIw6dvTLpy+g6g1FZiBKEZsxmzWnAmn0uvAiS/rRQQEhaGASaLxhSexzT1rbZyadVnHOQT7v7Otvdc5UVCpOa1qH4wrvbU/vXfAGRvhdi4oIh3l1zhcuYVnmVjBFkL78GFQWBfQ/zsoO4hUu6ucOb0Qu7GGMTyXEHzS6sI+HIcDPFLjqlWzKUqt7EONWaqv4z8M8V310LXNuVzLRgMCoLj8fjmUa63WcxExiIykIEglrWskg9RXHcIFVspIs6nYmq5Lwm6w0Zsd6NGMm8HtPKZXIoeo+p41MhN1Rub1OSqESy/9Tx4lJvTV0pwb20IHN1NM1fut2mBdsj3Oq9cSUDV2Jw0uKk3f2lMk3uZrTjRfaNXZRJTwNiF2mrM76x4+WhP+xiovjKYiYiEIQlb4vmJQe1L3yQlXGqqMkKDsoKFDSTHipeIve0pdFQyUsljlomhnhwP8l7YPK/sfuUShiFW+B8Lp8egfzLa19Sci+aOgWkLSyD/Pm0rBBVxuTJLYRcucEWtoGTcSctonEhkewvtuAwuciXMukht52C5NRLu6gq5AfNLpzr1Io89YNdTBa/+JHH4/F42qKL033MCAaishCUWhjZz3mvMfUCne1BhLFeo2QeiAFNvbgA690h4nTqKTk3pgobPVPwJgFJs1rwLMft7Cwep+LUJb50ehl2B3s9xvEaxflBkCVFQdOZDQr5ye51C6/NkQY0501qrlNTrGeZeZAiSpDzJjU75ARlKDfdC7uwju8g2oVtabo/apGhPrCLiaAKTb/40cxDBOq11pWFUSkvCFQxyduiEqCp8RrJQiHRrEmda+8qpa9miQRhCwjn5c/JDW0UClUDtFpRNrAq/2I6x3HTrsTkSBIEFfJCWV5KBlDhSglBsSBIdglMriAIZGwBEYh2pbKYNrsou1eDaBft5KPP7GKieBnK4/F4PC3xfRYzFBFlyJGhXC0x9SAjIzZSIzJBFsRiAsezMs4i7Gpj6wky6UDbVBvcjr6cd1hMJ9nOeZaOh1aMmBl3gBaFjkwp2R5k16CB4ykWPMicY+n047rH1nYcN1feCMg6SB1vUgLsWBkJ1HqNgdN5GQbGeophkPcgA1xv0rlmJ3/9YhdlrYmZahcGJdr3FI09T2I0Yvj44xGRvreLbqC+spg+ROQRYC8QAU1VPVtElgD/BqwGHgF+XlV3tTwOyqx6E8g/QKOSFQoitrAIjFrpIRJsREtEQJC0/w1BFlIZkW+e27QrW+RlgTHRKe5Lnr78kVMQRC0kibICokVzX8YkyM3148oHGjjadCG6xT18LqTS+aJMN9f03K4GXaKZuBJDEBpbQISO3DAmzdhCIRC3QCgvpXpiFyX2UCUZzVS72Pyvn+DAo+upL1lG48knAJh19CqO+c3fK61g+s0uJoPv4J5+LlLVHc7n9wLfUdUPich7k89/2JuseTyeKpr793Lg0fUsOf+FLHvJKzGB8sS3/5PdP7iBxs4d1I48otdZnDJUfZ9FP3ApcGGS/hxwA+NUFoEos8K4ZWEKnZhp7R+ZgCiJXghFaSaui4hiUjeG2IsEEDVI6loptlMTQ96bLMMZ3JV2SpZFt4gjNwSFTs2gNDJG2/IgLUWpwXqNkpumwbgRLUF2mbnjlHjQeXnC/VA4v/v7MZ2XmcQQJulaaAiSdBgotSQdSLaP603mZIcKD7LnduFOKVLRUpgpdjGy43EObnmEx2/8TwCOfMErwEAgwrJLfo7dP7iBxz79CZb+3GuYe+ppYzPSR3YxccTazeFCrysLBb4t8VP9VLIgyPJk5ShUdauIHFn2QxG5DLgMYPbyeQwnlQVkzUOjWSdUUwKiIH64kQms8TQkICot9INMm1ZQk2mo6WCtVn5FsXVdVhBIlBQGlEgPTkFQKj0U5IzxtOl8oaC2UNBArE6vYabxuxGiuYuVTJ7IhUpW3Yjk3DbSBc0VCulgyjAw1BJtOgyNffnrgckVBGFyM2qByfRr0Zw27RYMbiXRL3bhSndTZRcje3ew/sbPc2DXYxz3rFez/KTzIXRWdpugXdAY5dGvXkFQH2L+mlM48qJXEKigJuvnOO6d7+Hxq7/E45+/EoA1f/u3HHjoIVBl9mmn2Cz00i66ge+zmF4uUNUtSYVwvYjc3+4Pk4rlCoBFTzuy+4KkxzMD2bdzE3sff5htD9zErLlLWHn6S9l05zd59Lb/4Bn/688Ymr1gUscf2bmV5r7dAIzu2k44bz4rXnRpbp9Zy4/i2Hf+Ho3RAzzyF3/Cxj/MhIHj/vHvJ3X+fkHxMtS0oqpbkv+3i8jVxOvJbhORo5JWxVHA9vGOEwBzag0g9h7Th5jzIIPADqJpSkBDYi9LRGmW+MWq2RxARsV6PerO5eAOyoJ8Z3bReyyJegkK3qQbARNEatN5DzLZXpAzysh1WAq2J9PtvNRAMaEjkTjOZ9rINkI2KE2y87nTO+Su2Z4v/Y3adK7z0oluqYWGMEy9xijzJgNDPdFeaoGx0kNNCh5kmq5o4/TCLmxRUrCHsk7tbthFc/QgG77/BdQY5ixYzuozf47ZC47k2FNeyk+v/RA71v2QY067eFJ2EdZms+S083jynh8BcODR9VaCyw3kUwhnz2HNX/wNzYNPselDf8Os49dA4N6n3tvFhFHai/gaIHpWWYjIXCBIFiifC7wU+ABwDfDLwIeS/78+3rECMQyHDfs5lRuaJswKBQ2sHt0MAvvSpYVDipVhAsnebwXS/dxCsBUVMlTgRLdQKBRsQRBl2nRcKGTbs4JGq/XpoiZNLCvY0bIhNlw0LvCz7bnwUvsiZ/kxkt0jafeFKfR5pHp04GjQgfvyh4Z6EgpdDyIbFl2T8kKhFkQdyVC9tosyGWqydvHTaz/EyL44TuSsl7yHOfMT9bYR77NgyXFsvu96tm34MYuPPIkVJ1zAvMWr2rKLg/ueYN3X/w9D8xYxum8X8445kaf9yvupL14c7xdlduGGEgM88bV/Y99ddwCw/O1vzd+LHtpFN/DRUNPHcuDqeJEoasAXVfVbInIL8BUReRvwM+D1Pcyjp49pjjS46oVXln538utO4ahzj2bRSUcwd/k8pNtB9n3G8LylNEcPcM4r/5ygpBA78VlvYM2Zr6ZxaB8b7riau/7no5x83i8zvGAZQa1GWBsmnDMXt8ZTVUzUYMvt38I0Rlh1wWuZvXwVMmc4/n6cPB1Y/xD77rqDxRdfwsKLL0LCMN8PNoNR38E9fajqBuDMku07gRd1ciwRZW7tEBAvdWgjXYImzSSkI/Yg4/SoCVvKFXE+slh843TEaaRZwVM2rQMVnZglXmMQuWnN0s3Ma4zTyXGdTk1Uc+covzFOZ7fTMe1Os21qEKTXSdZcNwjpbNPGCe5yvWHXkxQtyYadjiHLpARq4/oDUSsr1MKsI7MeRtabHA6b1mscCiJqQeZNNg7tt6dasGo+z/uT5/DkgzvZdtcTrP/PB3ngq+vs96//90tZsHoR0Bu7KBtn0S272PvEBnY//gALjlhNGGHHgRQvZYhZDM2axennvZ1HH/wOj93/P0TNQ5holEbjILPnLWX+0jWMjuzh4N7tjOzdgaphwdEnc9rr38esBUsxocTjSyi3C2MMmz7/CQ4+usE2L+aedhqS61yfWrtI06HTwR1UabWTwMtQM5AAZXYiQ8V6dBLdgtBMjKShgU0HoowWZIb0t2HiLURBNt+Ma8gSqI2AiQdfpRbjvDmOEaUFgjsKtzRcMooLAEgKCDddVllUhUs6uAOuyOnRrrYsaC27Bts1IeTnJ3JXRHMKO3ciubwUlunRuRG5TgEbOPc4dMIf60Fko9uGwoih5OUfCpq2IKiLYf6RNc759TO45VN38dTmvRx7+gJOfPZi4ATUKB991lcAOOllxzCydRfL184FptMuxobLpvct3T4Ru1A1bL7n2+x7chO7dzyMSMDqp13MUceeR9Bozy5Wr30hq9e+EJJwWRM12bNrI9sfu4PZc5Zy9EkXMmfp0YiEaD3E1IAobxdEyuM3/Sc7br+BxWeez9wTT2HzVz9jz3P0W3+DOSechKkr6mZqiu3CVhaorSS6Hzrro6E8nhnF2Zedzi2fuos5RwxTn10jLSlv/9cHADju/OU8+M1NPPjNTbzrtilZx37a2frgjWy6/7+Zv/hYTnrmL7B0xamE6atsJlYoBmGNRUeexMIVJ8WHCcU6AmVH3HbnDWy5+Rr7ededP2TXnT9kyfkvZOkLLobZ9biCGVBUu1dZiMgxwOeBFcRV8RWq+tHCPm8mG2+2D7hcVe9MvnuEwkwYXclYgYF4nKEoc4JR+zkdQNXUkEbShq9pmPMgXdz4ew1TuQFMOmeQBBgb+B521MFdOiivpPMycDsymxA00+1OayLSrFPTOHJDfBEZjpSqNtJFSXUlCfNRK3a6CzRrikg262rs9Tn5d2SotoJMnI5MceLg3YgWV2IYCrPOy6EgssELs3IeZEQ9yC76qFMWMi+MpciR3SN8/x/uZO7SIYbnBLz0fWfwrJ9fiwni73tuF07LbCJ2MVSLW0hLlp3M8hVnJnZhsuNMg13s2XD3mEs7/vL3Uz9icXzaPrCLmkS4vSRd7+DuXsuiCfy+qt4uIvOB20TkelW9z9lnI/ACVd0lIi8jHjbwbOf74kwYXWcgKotADPNqI0DcZxElFt/QkFoiK9TVcKhKj0416NCZMygIsmmPA4M4I3ttPZALj8wMx7VJGw1VVkA4erS40lNTnQJCETdc0qY1H0ZbQixDpYWC4ykahTAtCJz7gGTrAoiTDqR8MrvidZacH5IywdGps0gXJUzusRvdMhRkEsNw2GBWosPMCpp2pH5NIup2AQh4+KZtfOnt3yNqGH52x24A3vHF81m8cnZ8bdogSsT2abWLkvuT68+agF2sWPUsHrr933j0wes58sgzmTt76bTbxYkv/TUe/NYnOPDEJrvf0LxFebnNPb+z7sl02UVdIjtYr5/7LJJByOlA5L0isg5YCdzn7PND5yc/AlZ15+ztc3h153sGkr+4+1UgsGfbCEuOncNbPvZM/vLuS2xFMYg89+K/YuHiNdxy0//m/nv/PRtVPk2E9Vkcccpz7OcT3/KegY84c1EEY4K2/oClInKr83dZ1XFFZDVwFvDjFqd/G/DNXHbimTBua3XsyTIYLQsyGcpoQJT4dQ0NaSTCaUNDO3Oo62W4A7SMCs1k6ofQmU5AJMyt/duO3FAVDVWUHtK4+aDpSg+KNJx0M/Mas9j6zIOskoLixesTDzKUbE3pIJtC2qjjMTj756b1cPJcFQ3VsqjKyQ1Zp7B7j8PAOB2WEUOO1zg7HLXp1JuM5YYs/ZF7X8LoSMTOR/dyYE+TebVs3E18nX1kF44tTNQuQg155lmXsW/f49zyk48yZ2gxa469aHrtItn3+J/7DYaPOArjBj/oODZhj5skp8guwqmMhmp/1x3t9COIyDzga8DvqOpTFftcRFxZPNfZPGYmDFW9sf3stceAVBaGuYkeHWlg+ywajjZ9yNRtk9TFkE0q19TAhuPF881kGmo2wRzkdBjHm8qFshakhnykS/zlU489xLZ7vsve7RtZcfz5rFx9PrNmL4oLhEYDEzUJtYYkhVbQVCsGizFZtJJqvk1sV3gTZxh2gCTatKmJPSa2qEz6JtKXN8RO0R1XUmm69XXauyP5G2KP60765qTdUbhDYSYrzA5HbUEwJxy1EsOsoJGXGzDMmQNHPH0WAUNAIkumYdS9touSe1ZlF/kBd44M2cjkyaBhwCjzZ8eD7zY++l3WHvWCabWLBUfHneEPX/NJTrv8b5BwVuV12lz0wC7iY5rujuLW7kZDiUiduKL4gqpeVbHPGcCngZclQwzirJTPhOEri0EiCGvs2RpH7Wx58Aa2PHgDcxccRb0+j907HwJg/ryjOe3UNzF7eHEvs+rpc1QjRg7tZrg2H5mmkcWN/Xtsev1XP8YJv/R703LevqFLdY/E+t1ngHWq+pGKfY4FrgLeoqoPOturZsLoOgNRWbjRUJGKXdWsoSGjGl9iiDKiYy+3oSGRnScotHLDqITWu3HnqsnHiVdkqNgJ7Mo4SVpVmb/omDE/Pf7US+HQKKtXPY9Fi1bz6Mbvcuddn2XV8nPYtWcjK484iyMXnQxteJA4HqQEaiN6Ag3QWvxFkAzFi29e1jOrTQjScRkhuTEX2XnJe5JF3PvkztfjLFqTevVu5+VQ0MzJDXOs3NAgHD1EYyRizmJhSDK5IcDxSp0ImLQ10XO7KGmBldkFJIPymlnaSlKRxi0K0o7vOH3e6e/gvof/g5tu/zAAa5c/n7Urnh+P5C6xi/2jT7J3ZBtLFxxPrT57wnYxf9lqe+iRJzbH11DW6iwyBXYxLE2bbscuukEXWxYXAG8B7haRO5JtfwQcG59HPwn8KXAE8PGkbygNkS2dCaNbGXMZiMpCUEeGEkduqBGYzACtrQTkIqZSbXo0qFEL4tCQmqObuovAB4Fi0rl0WuXJLRCcnVNteu/W9ex59F7OfsX7ufW/Pmh/d9fNn+R5z/9zhqiDwtqVFzKyfxcPPvJNZtXns2PXA1x46u8xFAxnBUHFOyABWcURaCY25/YPnKgnLcgNaZ4l32eR/tQtHFrdC6cgEGdyt9CdAC7IBlO5GvSsoMmtX3iQjbc8yfzFNXZs3IuI8L4vnM5wEPdL1KXpyA3lGeoru0jvq7Nzvs8i3zdl+6zcyLimQZKKY8HQMp5z4tvYf+hJdu19hHsf+y8efeJmjl54OkvnHU9kRnli73qWzj+erbvvYcf+DTY/Zxz7WlYsOz3+MAG7OPYFb+Jn3/tidj2uIzEOk7WLWcnzH5Ymc5L3vy5RR3YxURTsqoqTPpbqTVS7nuk+bwfeXrJ9AyUzYUwFA1FZzEQWHHUiCxevwew/MOa779/45zz3vPcxa2geIgGnrX0Np655NYEqP33oS9zxyL9z7tq3THuee8XeHaPs3HSAnY8YhmcH/Pz71vY6S33J3FlLmFtfxFELn86PN3yOTbtuZ9Ou2wEYrs1n61P3MlSby0VP+x2aOsotG/6Vu352FYsWHsfw0IKOz9c4+JStKE55y5909Vr6HgW62GcxExiIyiJAmSPJ3FAijsSQeYEjWs9+YMAEmQeZzhM0FDQZlfi37iIqoZhcR5yN4ChGwCTGM2YsgtuRqdkgu1BDhmrzOOfiP+GW6/4yd003/ehveNFzPpB4k3EcjzQNp698Bf+z7h+47u6/4sI1lzNcX1Ad8C0CYXydGmo2WKumY3YDkEiyaaKjbA4gd6CXOBMCjb3OLJYep1M0TqceJE7su+ZmCx2yUSzGepDDQYNX/97x8HvHMyc4ZL3GYTlE3ZEbQpN4kxXB4L22i7KxFVV2kZtd1u3gjoyVnqRhkCjeSZomliWT/WtRwAXH/iqooqocaOxi7tAS92YwK5zNc1b/Ct996B85dGAXs2Vux3Zxx1f+AoAjTjqHWfOWEBWvzaaFzFC6Zxep9JS3i0bOLlJJqutTlNO9cRYzhcGoLMQwLOl6FoGNhmhoyIgr1aQSLAGNZCTScNCgmUxydsjUchOSpcYr4oz+7NSZSF6a0sFXSdTL8KwFvODiD/G9696b++l3bv5TXnLW+21hQdMwFA3xkhN+n+vXf5h9I08wHMxN5olyNAT7ZotTiASk0nzuEgKQVLNuGjQpLMdEQDnyQk6XbpdUDRMnRFIMNUnlhoh6kp4VNh0pIcrkhqBhn3OcjjXrPdsP8abnxtLKq39xIZf/wSJmz0n7KuL/+8ou3G6OFnYBceEfpBWEUSRKn2dkKw6aBmmms/s58atGETXMlfnQaIyxix9s/DQAc2tL4mO4eRzHLn741XfbXddc8Iasoui0AJ1Cu8iFzmJKo94mxWFWWfhBeX3EqlVZ6PSSBScA8MDm68bsFySjjxcNHz3mu/2NXVOUu/5lztyAel246idr2bcn4pcufoxbbjrY62z1PUfNi5c4/eHDn6YZHWLD9pu4ff2X+OF9n+T+R/4LY5q5/b933Xv5/jV/kKsozv7VD09rnvsHsTMQj/c3KAxGywJlOGmeGhU7+Cp0q35nUemGhNYriad+SFfcihyPJj9XTToVAWRNZ23RaZb7yokgiqWHVNLJr6N80vEv4/jVL2bzozdx3NKz2b57Pfds/Cqr5p/B/OEjY+8ximhGseT21IEtLBlelbQslMhEbNz9EzYfuI+azOKFqy4jCJNHHIbWc1QcRzgIbAeqBlru3Royicm9luJ1VtwEcXYKCh3ErtyQee9RLm4+lRKGpeHIDaNWhpg/P6JWi73KD/6fI7j5f/bz13/wBOc9b5jf+uMlzF8YzGi7sGNrCh3cpC2OxC6A+H93ZF0uGioNt4q1oFMWP58DjV08cWAD33kgK/SXzFvNpm0/YdO2n1Rex9JjnsUJz30jTSeQo9/sIu7gzqKtur6axmHWshiQygLmJs3QBkF5oQBESTt3ljSsBl0vzCWTnxffjXqJj+Fq0y0pSA05ucHpC8Dq0XEETE1qrD36ecihiIMH43E3D2//Pmcd/eq4IGg2qSdrXD705M08+8jX2oLih098if1R3LJo6iH2HHycxcNHJfnJ7oWQBUbRNFZuECOZ5BU5MlRULkO1LTsUtOlsJbOi3OBMBuekZ0kW3ZJKDMPStNJDXSKOWllj15ZDrFpS48UvGuK8c5fx/vfs5uXnbOLj/+8InnXB7L6xCynYBoxvF/H2rILAONJTYhdx2pGh0gqkSBjajoOzlr6CdbtvZNNTd/K0ZS/kuOXnAbB9/wZ+uvHLuZ/NmXsMp53ziwwtWmIPX7SLjiTKabCL9FnXp0CC0i5FQ80UvAzVB6jT37Br90YOjDwJwPErnseSucexbd8D7Bl5PPebC1f8Krsaj3Hv7u+yt7GD72z7f+yPdvGMRS+1+ywcWjY9F9AHLD8qZNvWrHCcNz/g7/5pMb/y6/N425t28qE/PfzkuXYIJOTUZS/hkuPfzerF59jtyxaexJknvsF+Pvlpr+XZz34Hw7MX9SCX/Yq0+TcYDETLQlDqqTRAZD1IlwixTdWG1qz30QhCRkwcEVMPjPVA3MXe3Rjwjh+9On/kvfE0nn5k/04O7N/L4uFViBFCCZGkAnnmyv/Fw9tv4kc/+xcWD6/k5AXns3BoOcPMZk64kE3772HT/nvs6R7am0kHQaRgdeca8UzIgAjSTHsWDSRpCQIkdKWQJKOFaJ6WA/EqsJNIuGMTnPtaF2Onls51XkqDocSbHHamcoijXtJ0xMqjA3ZsaVCXIetNDg/B779nPsuOCPjLP3uK939w4Zh89cQuSu5lmV3YtB1nYbKpyJ0IKJrNuEWRpu3cUJqfKzxd+lDVCUmqtovQrmwFK1ecHXeYj2cXHUozU20XmQw1tkU5abwMNfMIRBhOjCsCsinkMk9zmIaNjBkOGowmL4I78VjNRHbCMVendg1ZHOlhjONQUaDGzfNM/89etvjznNlLCYPFyKEmS+YdgxyKYDTJk9Y4edEFrJn/TDbtuZObt/8bAGvmnMmZ817IzXuuBuCcuS/nlv3Xsj/azXMXvI654SK02cyWs2w2sY9bmllh0TR2SK4Yg5jA5jNfeKVhwZrNH1SsONyXx7k3kpMbsoLA1aMDMbaPoB40HbnBiYAhsoV5jSbDYmg0lPoQrDo64OEHGgyLEjoPJUK54Lw6R68MrEPRE7uoqCBa2UV275MPkWYyVGTy0lOaNgZN5acxfRY29Kgtu1g673hefOafMBqMsm7dVRx15DNYsOKE1nYxnlMxxXZRl4hhu49ST04eCjm76Aq+svD0I0PhbI5fcA6rZ5/GYwcf5MG9N7N1ZL39/pb91wJw4vDZzA0X9SiX08PoIeWCkzfbz3PnCvv3K5e+ehiALVsiPnXFftata7B3n+GuO5u89W1zepXdGU0gAQgsX3YmixYd1+0u4pmL4gflzUQEqKedl0LmTQlEiXlHGNtsbdC06RBjvZU4nXiNZJ5OVwb0WC9Ls07NKBsoJ5E68y9lcgPGZBKDiaWBY4dP4Zj6SXb7oeZ+bj9wPc+c/WJmBbOzjk3Nci7psSCe6yc09lzpdjVhFm2jTuSNiiNtTP4FCaw0oDZdF2OjVdzZQocksoOs4rThjtvisNgvXXUEa46vMXdIGRJhzx7l9FO3sWeP8oqXD3P6qXVe+MJZLFsasPaEGobBtQv7DKMov72sZdGhXcyqz2XpwgVgQPvYLtL90xZkndgu4mND2OX1NvygvBmIAHVJB19p1muvikmMyxDRSIwr1JqVFXIRE0EmN4zRpsnS7Ux8kxur5fwkF0GkZNFH6sgNrtYcRc7LHGUSQyPTpod0iPOGXx5HaESRLRRENStEIFucRiSrUMLQHkeMQRNJJl4XwclnmbZOW7cCXD2acm06kGwpzRC1BXWQKyyaDBFx+mmx2b7xtTu54HlD/MuVS5g7Wzj5aVsB+MzHl/DaS+OWRJRk1qBEjpbfS7uwRdYU24VWVBZTaRdWieuBXQAMYaz0NCxiK4gA6b4MdZhFQw1EZeE5vFiwMODd75vPNV87yA++P8raEx/njNPjzuirv7KUCy8Y7nEOPYcDbVWIA8RAVBaCUCdrWaS9g0PpZ2AUJdS0aRtlcoOYzGsk79Gk5OPpcycuz49rRFrwGqn2JrP4+Mzzw51+PHKkB9VCR6Y7IiqZD0o168iMTNZrGjrHNMbG9MfnKuS76nrSPJdds93oJJ37Fzg75zx2sntvp6h2YutD5/lc/o65/O475/Gxj+3jr/9mL3fd3eCVLxvmhRfMzskNgWYti17aRfEe2m2DZhfj2QRMqV3UJZOeQsnKhQDprgxVuNbDgYGoLADqyRQYBmMfYmQnpY6bp+nqaKFqtpQmWXM2FGcVNKcp3FXcl9954WNpgCydKwiiZPeiNp1sN8V8JhPMqWTatIgdzyFRlC1KUMv0cTVZdE6uMFKdks68QDR3v9OCIMQUno9TmNvt8d9vvnMe737nAuoSICIEhaFDabSotws4nOwCoE5g5emiXUwemTEd3CKylxZVm6ouaOc4A1NZeA5fRCTT3T2e6WKGtCxUdT6AiHwAeBz4F+L23ZuB+e0eZyAqC0HiVcGA2IdIPY7AdmQ2wMoN7gyUrqcYSOa5xPuN7XCDLDa8JRXSjTvrQK5PVDXvvZkSTy7Kx9Bbz1FNzotM11RWky1go1GEpF5jTm7Ie7Fq03kZwZ3WoXht45GbA0jKo4kCTE56yDzL7JnEckO8/5CI4zVmLYpYbsi8yMjKMN4uvF10uXUx8+KIL1bVZzufPyEiPwb+rp0ft3X3RORv29nWKwQIJSCUwGqTRX0y1TWDpFkbMvbPpSosslMJomz30sAZQ1YApFqzmuxzUVKw3ycFQuFzcVv+XJqdS012/FTyKOxedcmd3IpiwWq3o2Pu9Zhnk2jV8WRwOmYkbiiSFAyS2ED2rx/toipwahDsos2gMEsv7KIrKLEM1c5f/xCJyJtFJBSRQETejDtCdRzavXMvKdn2snZPMhFE5BIReUBE1ovIe8f/hcfj8UwfacU43t+4xxE5RkS+KyLrROReEfntkn1ERP4xKQ/vEpFnOt+1W1a+Cfh5YFvy9/pkW1u0lKFE5HLgHcBaEbnL+Wo+8IN2T9IpIhICHyOupDYDt4jINap6X9VvUo8hctyfEOnIj3CjKqYFJWv2t0Mxbt4ex1R8Dp1tmkXStJU3nVZNNvUO296f8ukbxnRwe7so+eztoit07z40gd9X1dtFZD5wm4hcXyjvXgacmPw9G/gE8OxOykpVfQS4dKKZHK/P4ovAN4G/Adwaa6+qPjnRk7bBucD6ZDFyROTLxBdZWVmkhBJg7Ivj6szOPq4m2mLq4kyz7oJVaOH/0n20NK3Flx4ySaHV4Uzm2kiufHBWRas459iDFf6fBO3cV/e55J5X2b6J3NQKbxf5fbxd9A+quhXYmqT3isg6YCX58u5S4PMadx79SEQWichRwGrGKStF5J9oHQ31W+3ks2Vloap7gD3AG5OTHgkMA/NEZJ6q/qydk0yAlcAm5/Nm4trUIiKXAZcBHLtyIPrpPR7PDKIDX2GpiNzqfL5CVa8oPabIauAs4MeFr8rKxJUV23NlJXArXaCtUlZEXgV8BDga2A4cB6wDTu1GJspOWbIt92iSm30FwNlnDk9jw9jj8Rz2KJ1M97FDVc8ebycRmQd8DfgdVX2q+HVFLtopKz9XOM9cVd0/Xn6KtNtG+0vgPOBBVV0DvIgp7LMgrh2PcT6vAra088NIDek/F/eTu96F0epbYFTs36RpZy2UeM7mMWmRACk0pyUQ+1d5uIp9csdyo4NajVVoJ/9t0s59NRrYZxORLYlaEpSDwRCpccJkx+LtYvx9Dke7mBTa5l8biEiduKL4gqpeVbJLVZnYdlkpIs8RkfuIHX1E5EwR+Xh7OWy/smio6k4gEJFAVb8LPKPdk0yAW4ATRWSNiAwBbwCumcLzeTweT0d0MRpKgM8A61T1IxW7XQP8UhIVdR6wJ+nr6KSs/D/AxcBOAFW9E3h+u9fbrti/O2ki3Qh8QUS2Y5fX6j6q2hSRdwHXEYduXKmq97b6TdFjhHjG0U58iogueYvtIrT22ooEQfn+EuQjX8o69USyqRzaypt0xVtsF6NC2QqHlfuTzShbz203ucgXbxfeLuLtpvsRUd0Tvy8A3gLcLSJ3JNv+CDgWQFU/CVwLvBxYDxwAfjX5rqOyUlU3FWY7aHucRbuVxaXACPC7xEPEFwIfaPckE0FVryW+QePvSzZSNzcVtYP7wkca2NmBIgJnpqAMU2GgnRYaZbuXHiIge+EDyV7seIa1sftLkEW0BAY1gfOVOPuUvCCBK2sE2fFFStuaVZfcya2oum9l99l9JhGBXcnOqBCVFIqROqGtanKFWT/aRUf3c4bZRad1aq/soit0qbJQ1ZsYJ3dJFNQ7K75rt6zcJCLnA5q0Qn6LRJJqh7Yqi0JnyOcqd/R4PJ7DgE5HqvcJvwF8lDiCajPwbSoqoDLGG5RXNVuhEFd2bc1WONUoaqecjrsw43RDDaPpFNXqdoZlXknkdKgZDTA5CUPsb10PSNtxnySftstWOx5abjYAEdR2XEreq0u9pjBATJjkIVsLW02Qi5e3pw2y40gYQuh4pTlvNTmOuOm8l5jmW4udmW3cCvd+uffS9R4NWcelUSFyvH0jWYdmOoRglGyG2FDEkZsCjI5tWXu7cE7r7aI7zLDFj1R1B7EyNCHGG2fR9oyEvaaRGIJBaZBKD1mExCgBDY3fnAixL79xm7YaOMY4RTq1q0dL/oUU++I5s6iGASQFgUQGTfRlCbMhSSJark0Hkq1bEDiRM4UCQm0BQZYeE4HThWsvkH/5JSuoCQrPJ96noSGj6Up2zpjshjP9eIDm5n9KpSdvF3i76PLMxDOlZSEif6Cqf1c1OK8rg/I8Ho/HU8EMqSyAPySeWfZhYNdEDzIQlYXmvMasI3NUlYaVEiTviaTepBO3bchLDylGJbfei3Pi8vyUSQ0y9vvcdtezDCWLUCnIDbbdb4z1DlVM/qSpbOFGyYRBhdwQxOez2wv5rrqeNM9l12w3Oknn/pmc/JBJDJm8k8lBoxpmnr9kz6eBZFOFA+lc2UEhcqlsDe5e2EXxHtptg2YX49kE9IVdTJqZ1WexTUSOI46gumiiBxmQyiLWoSGWF1I9uoEwmjRbRwlpaHy5kaODNjTMCgsT5o001X4dg9U2px22LxSMfcHcF895Ia2U4GrTYZgtmRmG2eplYLeLSL60yr3wmTyBIz1kaUfycAoRDVsUEM61tfX+qVh92uCkC/0CjURWiUKxBYFxIpQaWrNyQ+Bo06gylJwqQDGlMlT/2IV9UgNoF3aqqD63i64wcyqLTwDfAtaSn/pDiK9ibTsHGZzKInlykcaFAUBDhYZ9+QPrNTa0lnmQOc06sPvnvMluCLP2BRM0XXgmzHRhDR1tOgisBynOi41qth230KkqFPLatPVKwyytQWALI4KskFIRJ+0cswu3ItdBbJ+V+/KHNh17kLGZjhJRT55PSJBNJCfY6w8lv8RD5O3COZe3i27SYq7JvkJV/wn4JxH5hKpePtHjlERPezwej2fQmExFAQPSsjCqjNhmazbHz4ija45onUYyj/+IqVuvpKEhDROnm5rJDU0T0nSkB3WkB+usFed+yenMWTpWKNImvaI2/FGz0MPAkRtCA7VUgw4gyWt84iTdbGaRMe5ym1AeXlmrZZ5oLp1p07E3maaza4jzlqRF8n0WZenCvYkXX3OkmyTd1ICmyWSFZurhmxqNIPP2R0w8FrceNhnROJ1f4yCLgAkK39mw2F7bRcl9Gki7qLKNlD6xi64wc2SorjAQlYUitsMybram6dAaUUNDa1yjzvZDpm4LjoZx5IaiNp2mO82cMKaAyL9skqVrSdoEsW4CUAuQ9IWv1chmWallEkMUlbcRwzDfkVlLHncQoDVHbkjTtTg/zeYIUSTArPgcMlRdQbRJet/cAta9xw0NaCSjjRtByKG0IJCIusbXPGLqhEH8+o/kJnPICoI6xspNLn1lFxUVR7/bRTFvY4Idig5DG/TaLibMzOrg7goDUVl4uocxETd9/4O5hXVOOv+XWXjCGT3MVX+w8/6d/OgfbmP3hl3MW3sEJ3/w58F2oXoOO3xlMfMwwP7EI8yHQtYYSZrqIzpkJYZD6niNjiTR0NA2f+OmsCs3xOfSXDhLC3JSgxM1EkByurgJnzb1Q0GjrNkvteQkJrBqQ84vkijzIIOgsiMz8yCzqBethVmnZi3ryNRAeGzrzYThLM590XvZ+vjtDM9ZwrqbP8eRu5/L8qc/n9rwkvIInlbk5IYsRLJpAivpNE1Y+UwOJc92SCNGKuZ3Sj3IUZTQmSYw7RCdrF08fs+TfOut/8lxrzqFtb/1Yu7+82/w5C0bGT771OS62rcLLdgG9L9d2Ois0JGeQnIjuHMBYX1uF13BVxYzD4MwkujLxZGdIxp7fiOmbiWGUa3Z5mxDQ6uJunp00wTZBHMmwDgTsmmutCwn91VAJgcIOf1XA7XpVFPWmqDpVAK1wL53ivMOupEupnqkrhsuqbXQHtNKDGEmczyxcx3rH/gGZ55zGeGs2Rx14vMAOOMV72bbxh9z91f/kqWnXsDRL/5f/Oy6LzD/5NOZe3qLFocj16QYExA5lXBaIbs6ddMpFA6ZOnWJR+ePaJ2cCJ3gFgqBGkInTMVONteBXTz8451c8dbvcNl1r+THn17H4tNWcMOf/YDTf/1sjnvzeexrDDH3uCN46t6tDD3z9OxyB9QucPsmAnKSlD2xey2tb0Ff2MVkEWZONFS38NFQHksUHaI+NJdFS/Jh17MXLueYc38OgB33/oD9Wx9h9/23senrn+1BLqeenZsOAHDFxd9g8+1PsO7qhzjt157Fqb9yFgB77nmMJ258iN0/fghzqNHLrHp6RdJn0Y31LGYKg9Gy0MB6hwZh1I1osXJD3XZkup2XI6bOoSjdP8jJUKlXmhtd2unDL5Ohyry0UDCJVxcYRdOoFJM7lP2tREE2KCtQoGzGOLHRLeqM/tVagKbSQ5hFumza8iPmLzo29irHeJNQmz2f5sG9PPyVfwSgvmBx+/fBOruZpBdpUW6I04eiGjXJOjLTVmCIydwbg52Wui4RQxJ3dgZo6TNqxy6u/odHWfa0xezbE3usZ7/tVM5+66mM1OewrzFM08D+J0cIF81DwoCRTTtpPDVCfelQZ3ZRJUP1qV2YWpCNrSixi9y1dNqH3GO7mBQDVBG0w2BUFggHdBaQTDzmjOwcceSmAybeZ8TZfsjUrDGOmlq5DKVBaZ+FFj2HdFCVMz2CHdnqhhimEkNo54JDQrGKgYYBJRIsiCDN5KUVk03HUNUcDijIDU6kS90tIAJ2797I/n3bOO/8y+ICKsxr6CrCGW/+C5587C4e+eZnWXLG+ay45HV24s34Op1Ma6aGSEGbdtchcOWG0URKbARNDiXputQJ0va+c52RBAwnhVxE01b+YYUr145d/PenNgIbAbjovc/kpNc9HQVGmzVGG/C9X/sqex7awZp3vghNCuTaEQuyENAWdqHO4DVX1+93u4jPm/VZlNlFmh5zbTbdv3YxKXxl4Tkc2fSz74Mawtqsyn32P/EzHvnmZ1nz+ncwb9UJ7S+xNUN40TtP4PGNI1zyR2cQzZ3PQecCNTLseWgHAI9f81MAzvzaHzDai4x6+oJBkpjaYSAqC0XYb7KWhXGmB0ilh0OmbiNgDpm69VAOmZpNN02YeTSFDjcbA24k1zFXmSen40/BiaF35YZCZEnqmebOkeuxyzxCI9ajEnXcNbJ9VCT7eeBICeHYjswndz3M8mPOKcgNbocr7Nm8jllLljN31QnxNQXjdGSm90IFY7JOTTdwwL3fqQR4yNTsCmd1ibKOyQB7zbOChj1mXUJnDYPyaKh27OJZbzyBj736Bh7fcJAFpy7O2YXW6zzvM7/Aw1+/n+03PMiZ//dtMGsIHZ2AXaTPHPreLiDpWG9hFyTXkdqEvc7x7kUP7aIr+Mpi5hGpcMAMJen8XDJu9ESaPhAN2YJg1NRsM3fUhIymk5apU1kYZ8ZKN0Syylgc7TYtHOxL5aRNCJI8AWOc44VZiRJg0CSKJXCkCnEKBdooFOJ5hVINWrJ0PcCEwtLlp7Lj8btZE7423r+WSSGaRLrs27KeRSc+I18gOLLDGNz75A64MlmhYGcRNSE1k0g7pmZlg0OSb7+4C940kpV96hLZgiAUY2cdhWy+oXbsQucNsX/HCA9+bxtnnHIsBw/Cvi17ac5Wtv70Z2z/wQZ23/s4p3z4LQytOoKRRgd2USbPzAC7SPNjbaHELtx06eDNIn1gF5NGuxcNJSJXAq8EtqvqaSXfv4ds0aIacAqwTFWfFJFHgL3ELkNTVc/uTq7GMhCVhWfymKiBabYWVfY9vgHqNZaqUlj0faC45TP3cetn16FRvnBZ+8vnsva3L6ZRm9ujnHn6iu7VPZ8F/hn4fOlpVP8e+HsAEXkV8Luq+qSzy0XJKnhTykBUFobAylCmMK1DOr9PQ8Oc9HQwilsih6Iao1EqN2TN39i7yTwdtyMuH86S5aNqSgx1OzJzcerZfEA5ucEJG4cAacZfGNFsiodI7A9KOz1tPrKomrIpJEwoaF3Yse0elh3zzIqolzh94mt+k4eu+kf2PHgHC592Vul1Zjcjf0PcAIG0xRY5affej0Y1Z02C/MXNCprJbwPqSTr2IJP9K9y9du3irT98Mw/91waaQZ2Fp61k3qqFHGyG1l5GmnWi5gTsosQ2ZoJd2H1a2IXNf4vrzG4GfWUXk6FbfRaqeqOIrG5z9zcCX+rOmTtjQCqLTIZyF0hx5YamhjZE1u2niGWoMEmHmcE6fRaqZIOhlLY8ijERL6XadNakz71E6kQCSiZHixEk8XYlyCSGKqNVt1AI3BeeLF0Ttj9+NwDHnPFyTCp/hBQKMphz9GqGFi1lz/0/ZcHTzxpbQVTeDOw9i3X3OO3e48gE9jkEolabDpyLc9fCbkhEXR25IV3PoFVl0YZdNMMaa37uFEaiOiNRjVETR0911S4cW+hnu8hkKGlpF/HJCtczQ+xiUrRfWSwVEXcdiStU9YpOTycic4BLgHcVcvFtEVHgUxM5brsMRGXhmRz1WfMAGJq9YIz9b771Wg7seow5K45h/+OPolGToy75+enPpMfTT7TpHCTs6FJfwquAHxQkqAtUdYuIHAlcLyL3q+qNXTjXGAaisjAasK85bD+nHZlNDbMZK515nw6ZGiNRMsWD06k9GoXZqlwm63AzpjDOImHsGIu0AzFrjWviAZYOvgoFtZJBlnYX1Qkgm+4hUjtjgxjJr47mOk5OoIwbB29nC3Xm9zGh8NDtX2H2/COhnkXGNKKDqARsu+s7AMxecQxHnP4c5qx9GjJ3lr227DyMIRde784BZAMHDFFyjxsmtN6x60EWSaWkugQ00+mqnY7M9Pcp+SVWe2AX6fWnrQA33ed2YWw0lNPiCMm3RArBDpU24dyUfrCLySJUt9ymkDdQkKBUdUvy/3YRuRo4F/CVRRVuNFQ8vXEWGZEWBO7o7Jz0FIWMRulU1CENJ50LnU3nhmrXoyjIULmmexpNou7IcLVvV4BmBYNkURcSQDoKTk22f2u5IflQoTWbZJnMZWvPiQuIJG93/Muf2OMMLzmK5Re8LN6/Fi+lYK+xY7khwCSFWWSyJTMlAkklhkhLI1eMCrV00KSENDSOiAmpWZmhqkDoK7uokJ76zS60lqVTuzChZPksk6TcPovx6AO7mAzTWVmIyELgBcAvOtvmAoGq7k3SLwU+MFV5GIjKwjM5Du3bSW3WnNLvlj/zxaw4/+XdjlD3eGY+XaosRORLwIXEfRubgT+DeGEOVf1ksttrgG+r6n7np8uBq5PIxBrwRVX9VndyNZaBqCwMwsEonaI8cGLrhfyU41lHdhphMRLVHIkhi4Ayjtzgdr6pcbUD8m3tEq9Kg0RucJrrrjdm49RDVzEQ600hknmQxont1syzaRX1knW0V3iQiQVs+Mm/s+Tp51mvMWXp6RdkkS7kvUfbuZlecy7qJYvKQbPZUlXVpt17HAVZBEwaJ59i17zOeZDGpkPJPM5WHdw9tYuCTUD/20U2N5TT6im2hgp2McYm3P/thfaPXUyK7kVDvbGNfT5LHGLrbtsAnNmdXIzPQFQWqsL+ZjZNRVooNE1YWKoxsOkyiaEZBTSjwKazFdHIol6ccD+pCPlQ52VMCwQ7LCHIbMxoMtdbkufUsA0QJGsYGMGGRUok9sdiCs3gXP9JPi+QvPhpOswKAvfl1hAO7H2CdV/5EAAnvP63CRcsjAuKYtRLmk6PWXYrNCsIJO3PcQpYo2Lvt0i5xAAwlEws4hYKNTEESU1XC6Lcb6v6LHppF+59sque9rldGEduMk6fhbt9jF0411kZEdUHdjFplMNuuo9g/F26j4j8uYg8JiJ3JH8vd757n4isF5EHROTiXuTvcGR44XIAdj38UySsseaStzJnxXE9ztXMobn7KfZ+78e9zoZnOtE2/waEXrYs/kFV/7e7QUSeTtzjfypwNPDfInKSqracs85oYKNYDJJrnpZ6kE4HWmQyr7FhQiciI8Ak29U4bl+HHdw2Gsp2HGZeozidmmne4/uAdVtiiUHs/lkMfd6bLGNM3HsuAibefmDPVgBWnf8aTCgse8YLePy269j4rSs54+SP2N/mpnjIHbONe+HeM8V2Cpso+7kxSqPkWO4zrAWGWjIFa01MFnNvwixd8XCm2i5EAqReH982Slqd0H92EbcsnO0lrYkyu8gNTB2PPrCLyXC4LX7UbzLUpcCXVfUQsFFE1hOHgt3c6kcGONCsO59LCgVnvhl30E/DHQDkyA3uXDVqKAy+ytL58FnnfyetATa0EXUKCC007dLfRGCSF1iM5rTpXKSLc+5cOKL7cpUU7GmhoKrc+x//m9lLVzL/2JPY/+RmHvz3uIJY+cLX2agnV0PPzWdUuE63gIinn3YylIZzGokjdoivURLZI73vKekzrGs2eV6khmY6T1BgcgO0OpGhpsIuwvnzmHfes2wBOMYuSu5Tv9pFnK5YPjWktV1UvANufvrBLrqBl6Gmj3eJyF0icqWILE62rQQ2OftsTraNQUQuE5FbReTW0d0HpzqvA0k6v1Nj/1MAtqI4/tXv4IjTz+9ZvjyevqddCWqAKpQpa1mIyH8DK0q++mPgE8AHiW/lB4EPA2+lvPFaeruTYe1XACw4ebmOJFM2uJ6kqmQRMBVr/MaDrLLOS7tP5ERqmMxNUiNteRRuc9yVG1Byi5e5Swe78kSQCG8auFEvmkkMVZ2YLsUOTadTM83P8a/8DR7+xidZ97UP233nHncCpiA95CJ12unEdLPhdARL4LTSJL7P2R1w8msjeoSaZs8wDOKbEWk2i6jI+N5jX9lFiQzVb3ZBUWJKWz3j2EXLgIdiNvrALibFAFUE7TBllYWqvrid/UTk/wLfSD5uBo5xvl4FbBnvGEbFzu+TWwTebaqafAHhTofsprOCIMgXCk5oYsupqCF+C9yV0YLCb0p+b3AiY8SOsUJM9kKKkS4VCll63uqTOPmX3sfo3l3suOsH7H34bvbvfIzhFStL5wAqzgeUSQ0lmSi7ZgPqDCxLp9lWMXaENUBgI2aMLejDIBvZGwZBNljLHR1cUSj03C6ckNfcc+hTu6gMka2YG6p0UF7lqMCSa+6RXUyUVpc3qPSkz0JEjlLVrcnH1wD3JOlrgC+KyEeIO7hPBH7SgyweVsxatIxZi5YRHTrA3ofv5uF/zVoZJ/za+6gvXdbD3Hk8/Ymd6fcwoVcd3H8nIs8g9iseAX4dQFXvFZGvAPcBTeCd40VCQeyRHGpkl+LOoeN6kOo0W21Mt7PWthpxIl2cZnGxfZ0bfOXko6Rzz/UKIfa+7KCp4oWYbN/Ma8x2dGPotXCA8Toy3e1jBtalrYyTTmfZ7kuIRvYzsn0Ljb27qC1ePLZT2zl22W2xGcx1ZDrb06QBTTJuCOyqZhGBMzdSFnVigmxQWqialxlwvUjnHDkvvYd24fxESp5PP9tFccCgjdQqs4uy29GHdjFpBqw/oh16Ulmo6ltafPdXwF91djyxg6nSz8V0LorFiYZBcaJbApIIvPyIXIPdXow2qcR9UYL8Zjt4yz2UKwc4L4KK88IH+cKvk3wU8+RWbNmI4pBlz3tpvI+jUxfT7qjjdkMkbWFmnJc2EKfgVTvNtKjJNHLNwh+NEYIgTecLhXZkhn6xC5uNGWIXbno8u+g0dLYf7GKieBnK4/F4POPjK4uZhyo0mq1bFupud71GzUe35MZTlHmTOQ+uwn0qem1FbzL1plxJwvXStZC2ESD587fpvI39QZseZOn2gqTSMjPqXIQjPVhvHLWdmqrZ1Nqimdwgos6UEmJnJhXoSsti2uyi7F4Nol2MZxN25yQXfWAXE8W3LGYgitB0ZCjni0IBkfzvaNMozmRmeYmhdPCdq7O2IvmBJtEv4r4fSVbFlMsN7gAtMXmJQdx8jz2dzeLY/BS2p5dWiGDJvfyu3FBynLiAaONmaJYpUc0PynJiRNPCIlcQilOIBCa39nc2qV7hUp08adnN6LFdqBMpN2h2kVU2/W0XXcFXFh6Px+NpidN/dLgwGJWFuoN4Ct6kjcvGabZLadr1DtWNXTfuPtlxxjgqJc15SWUo1wtzOibtb4zzGy3sU3Lu3HFo4eSMvRX5TtOqdDD+9qJsMcaTK5NJFDunkaLWtRYnrZI9QxHN7osJ7DnceZLS/Uov2qGnduHcS1diGji7KPmu3+1iIrTbqB4kBqKyUAXTDEq2u28EOW06w2nyquRkCMnp1E66qq+icNj4f82NttUg80hyEnTobHdf9pzUkWW9U0MtSg2lUSuFwms8/ToefNjmyW1MqtgZ88QImqTVSF6nT6WOnIyShUjGeXMKjg4Kgp7bRXoNg2gXFRJQ5cn7yC46pnLBkMFkICoLj8fjmW58y2Imotmgqfhz2U5CaVO94DVmHW7OPkW5oeV5yMXAq8RN6XRRmdw0DUr+3GWSRC6vhU3tGGuZY+V6kIV9xpUkHO+xRNXJY93dQl5TbzzQTHoQ7EAsxJlnyZUzlKxz2DmgIm16smUbp8cuXCnJpgfQLsa1CTd//WIXE6Eo/R0GDEZlAZmhQflDLEa/6NjvcgWBFgoCZ3tbHkUhO0ImPVTp3GUDq9zCZew1dZaP0kNIYZ+ygiCnP2tpwdEyC4ojBziZNtnLLM58WrkCAnefsfms/FyGt4vSfJQeYqJ2USZjVWWhX+xignSrg1tErgReCWxX1dNKvr8Q+DqwMdl0lap+IPnuEuCjxNNQflpVP9SdXI1lcCoLj8fjmUa6GA31WeCfgc+32Of7qvrK3PlFQuBjwEuIJ2G9RUSuUdX7upYzh8GoLBQ0KnEhih6W44qVzs5ZjJUfT4YoHqekcVP0m8TpvlS3GV5x3EpZITtMe7TwuErj7Md4ljp2uytbFPbPeb7q3oei9uBkosSbdPOnTj7yJxt7uJb02C7c7A6cXZS2OvK7961ddIJCtzq4VfVGEVk9gZ+eC6xX1Q0AIvJl4gXkfGXRCimrLFxaFBDj6cBtSwxVFE+VHCz+r6SAKJ6/1bZJ5seeu+y7kr6JdqWnylMXpYeSd9wtIIDy8MduFQTeLkrzY89d9t3hYBdt0MGzXyoitzqfr0jW4+mE54jIncRLNrxbVe+lfLG4Z3d43LYZmMrC4/F4ppX2K4sdqnr2JM50O3Ccqu4TkZcD/0G8fENZ1dgt12EMg1FZqOQ7uHPfVf1mnH2K0kSnj8BtItPiqboDiJxzlE5T4R5wIlR5Xe70GGX7FDzISo+zFa6UYjcVBJkK6aZybunJeJE9tIuqy8yd1ttFtsN02kWbCB21LCaFqj7lpK8VkY+LyFImuFjcRBmMygKyAVRVVBYOecsqNYCK31YaiysdJOksXLJiPzcrWn5w6VSPrqAyiqaFfl3623EKBzd8OH+Qkp+o5Pez90/L73+FxN0x02UXuWebpQfOLtqoMGaEXYyH6rQtfiQiK4Btqqoici5xMPVOYDdwooisAR4D3gC8aaryMTiVhcfj8UwnXaorRORLwIXEfRubgT8D6gCq+kngdcDlItIEDgJv0Hj63aaIvAu4jjh09sqkL2NKGIzKQqmWoUqobBG0+fA7iqcveH1jfpp6llVN7TIPs5tUqRqT8DLH/EQrjqdj9xtzTHccRMV96RbTZhcl1+ntwj3J2P3GHHMa7aKKbslQqvrGcb7/Z+LQ2rLvrgWu7U5OWjMYlQWdyQRdP0/LHxTyUdH0bvslnGA+JlSgjPebDo5ZKT2UMRHZrxN6bReuI1HcRv67GWcXHR6vr+yiExTwa3B7PB6PZ1wOr7pigCqLXnuLLQ9U+DxeOEwbTIn00Iouna94T7tyHdP80nq7cDiM7cJPJDhDmVEPbrpf6D5mRj23qcbbhWUm2MV0RUP1CwNTWXg8Hs+0MZExNjMcX1l4PB5Ph8SD8g6v2sJXFh6PxzMR/BrcM4/RTZt3bPztdz86RYdfCuyYomNPFJ+n9vB5ap9+zNdU5em4bhzEtyxmIKq6bKqOLSK3TnISsK7j89QePk/t04/56sc8WXyfhcfj8XjGZ/rmhuoXfGXh8Xg8E8HLUJ4CnS5SMh34PLWHz1P79GO++jFPMUo3l1WdEYgeZrWjx+PxTJYF81bqs8+8vK19//uH77+tb/teOsC3LDwej2ciHGZ+tq8sPB6PZwKIObx0qKDXGegXROT1InKviBgROdvZvlpEDorIHcnfJ53vniUid4vIehH5R5GqdR67m6fku/cl531ARC6erjyV5PHPReQx5/68fLw8Tgcickly3vUi8t7pPHchH48kz+MOEbk12bZERK4XkYeS/xdPcR6uFJHtInKPs60yD9Px3Cry1Je2VIoSD8pr529A8JVFxj3Aa4EbS757WFWfkfz9hrP9E8BlxIunnwhcMh15EpGnEy+heGpyzo+LSDhNeSrjH5z7c20beZxSkvN8DHgZ8HTgjUl+esVFyb1JK/z3At9R1ROB7ySfp5LPMtYOSvMwjc+tLE/QZ7ZUhaCItvc3KPjKIkFV16nqA+3uLyJHAQtU9eZkicPPA6+epjxdCnxZVQ+p6kZgPXDudOSpA0rzOE3nPhdYr6obVHUU+HKSn37hUuBzSfpzTPEzUtUbgSfbzMO0PLeKPFXRS1uqRrW9vwHBVxbtsUZEfioi3xOR5yXbVgKbnX02J9umg5XAppJz9ypP7xKRuxJpIZUzqvI4HfTy3EUU+LaI3CYilyXblqvqVoDk/yN7kK+qPPT63vWbLVVzmFUWh1UHt4j8N7Ci5Ks/VtWvV/xsK3Csqu4UkWcB/yEip1K++kDHljHBPFWduyt5GnOyFnkklr0+mJzng8CHgbdOVV7apJfnLnKBqm4RkSOB60Xk/h7lo116ee/60ZbKSfssuoCIXAm8EtiuqqeVfP9m4A+Tj/uAy1X1zuS7R4C9QAQ0pzJE97CqLFT1xRP4zSHgUJK+TUQeBk4i9m5WObuuArZMR56Scx9Tcu6u5KlIu3kUkf8LfGOcPE4HvTx3DlXdkvy/XUSuJpZPtonIUaq6NZEOt/cga1V56Nm9U9VtabqPbKmSLkZDfRb4Z2LZuIyNwAtUdZeIvIx4sOKzne8vUtUpnwTSy1DjICLL0s40EVlL3Gm8IWm67xWR85KIo18CqloC3eYa4A0iMktE1iR5+kkv8pQUNCmvIe6Ur8zjVObF4RbgRBFZIyJDxJ2j10zTuS0iMldE5qdp4KXE9+ca4JeT3X6Z6bMbl6o89Oy59aktVdCmBNWGDDVe/42q/lBVdyUff0TeIZw2DquWRStE5DXAPwHLgP8SkTtU9WLg+cAHRKRJ3NT7DVVNH+zlxF7BbOCbyd+U50lV7xWRrwD3AU3gnaoaTUeeSvg7EXkGccP8EeDXAcbJ45Siqk0ReRdwHRACV6rqvdNx7gLLgauT6OUa8EVV/ZaI3AJ8RUTeBvwMeP1UZkJEvgRcCCwVkc3AnwEfKsvDdD23ijxd2G+2VInSSX/E0jRsOuEKVZ3oVCZvI/9Op31iCnxqEscdFz/dh8fj8XTIwtlH6XPWvrWtfa+776/Hne5DRFYD3yjrs3D2uQj4OPBcVd2ZbDva7RMDfjNpqXQdL0N5PB7PBJjOcRYicgbwaeDStKKAfJ8YkPaJTQm+svB4PJ6JME2hsyJyLHAV8BZVfdDZXtUnNiX4PguPx+PpFFWIuhMNVdF/U49Po58E/hQ4gnjkOmQhsqV9Yl3JVAm+svB4PJ6J0CWJSVXfOM73bwfeXrJ9A3BmVzLRBr6y8Hg8nolwmAUH+crCM+MRkX2qOq/X+fAcRijg1+D2eDweT2sUdIDmH28DHw3lGRgk5u9F5B6J15D4hWT7hSJyg4h8VUTuF5EvJCPcPZ6JocQd3O38DQi+ZeEZJF4LPIO4028pcIuIpAOUziJeD2EL8APgAuCmHuTRMygcZn0WvmXhGSSeC3xJVaNkUrrvAeck3/1EVTerqgHuAFb3JouegcFPUe7xzFhaSUuHnHSEt33PpBisiqAdfMvCM0jcCPyCiIQisox4Esgez07qGUgUMKa9vwHBe1eeQeJq4DnAncSv8x+o6uMi8rTeZsszkBxmLQtfWXhmPOkYi2Td8fckf+73NwA3OJ/fNY3Z8wwk3ZvuY6bgKwuPx+PpFAU9zMZZ+MrC4/F4JoIfwe3xeDyecfF9Fh6Px+NpiepARTq1g68sPB6PZyL4loXH4/F4WqNoFPU6E9OKryw8Ho+nU/wU5R6Px+NpCx866/F4PJ5WKKC+ZeHxeDyelujht/iRryw8Ho9nAhxuHdyih1n4l8fj8UwWEfkW8QJb7bBDVS+ZyvxMB76y8Hg8Hs+4+PUsPB6PxzMuvrLweDwez7j4ysLj8Xg84+IrC4/H4/GMi68sPB6PxzMu/x/5DOefHVlzMQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Display the global field and countries' outline.\n", "fig, ax = plt.subplots()\n", "ds.field.plot(ax=ax, x=\"lon\", y=\"lat\")\n", "regs.plot(ax=ax, edgecolor=\"k\", facecolor=\"none\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Compute the field average over each country\n", "\n", "`xesmf.SpatialAverager` is a class designed to average an `xarray.DataArray`\n", "over a list of polygons. It behaves similarly to `xesmf.Regridder`, but has\n", "options to deal specifically with polygon outputs. It uses the `conservative`\n", "regridding, and can store and reuse weights.\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "xESMF SpatialAverager \n", "Weight filename: spatialavg_90x180_6.nc \n", "Reuse pre-computed weights? False \n", "Input grid shape: (90, 180) \n", "Output list length: 6 " ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "savg = xe.SpatialAverager(ds, regs.geometry, geom_dim_name=\"country\")\n", "savg" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "When called, the `SpatialAverager` instance returns a `DataArray` of averages\n", "over the `geom` dimension, here countries. `lon` and `lat` coordinates are the\n", "centroids each polygon.\n" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
<xarray.DataArray 'field' (country: 6)>\n",
       "array([1.58448717, 1.49195654, 2.48197742, 2.5740491 , 2.89305772,\n",
       "       1.26271837])\n",
       "Coordinates:\n",
       "    lon      (country) float64 -65.18 103.8 25.09 -3.646 -3.541 -102.5\n",
       "    lat      (country) float64 -35.39 36.56 -29.0 40.23 17.35 23.95\n",
       "  * country  (country) object 'Argentina' 'China' ... 'Mali' 'Mexico'\n",
       "Attributes:\n",
       "    regrid_method:  conservative
" ], "text/plain": [ "\n", "array([1.58448717, 1.49195654, 2.48197742, 2.5740491 , 2.89305772,\n", " 1.26271837])\n", "Coordinates:\n", " lon (country) float64 -65.18 103.8 25.09 -3.646 -3.541 -102.5\n", " lat (country) float64 -35.39 36.56 -29.0 40.23 17.35 23.95\n", " * country (country) object 'Argentina' 'China' ... 'Mali' 'Mexico'\n", "Attributes:\n", " regrid_method: conservative" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "out = savg(ds.field)\n", "out = out.assign_coords(country=xr.DataArray(regs[\"name\"], dims=(\"country\",)))\n", "out" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As the order of the polygons is conserved in the output, we can easily include\n", "the results back into our `geopandas` dataframe.\n" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
idnamegeometryfield_avg
5032ArgentinaMULTIPOLYGON (((-67.19287 -22.82225, -67.02727...1.584487
9156ChinaMULTIPOLYGON (((77.79858 35.49614, 77.66178 35...1.491957
37710South AfricaMULTIPOLYGON (((19.98200 -24.75230, 20.10800 -...2.481977
67724SpainMULTIPOLYGON (((-5.34065 35.84736, -5.37665 35...2.574049
98466MaliPOLYGON ((-12.26352 14.77561, -12.13752 14.784...2.893058
155484MexicoMULTIPOLYGON (((-97.13797 25.96581, -97.16677 ...1.262718
\n", "
" ], "text/plain": [ " id name geometry \\\n", "5 032 Argentina MULTIPOLYGON (((-67.19287 -22.82225, -67.02727... \n", "9 156 China MULTIPOLYGON (((77.79858 35.49614, 77.66178 35... \n", "37 710 South Africa MULTIPOLYGON (((19.98200 -24.75230, 20.10800 -... \n", "67 724 Spain MULTIPOLYGON (((-5.34065 35.84736, -5.37665 35... \n", "98 466 Mali POLYGON ((-12.26352 14.77561, -12.13752 14.784... \n", "155 484 Mexico MULTIPOLYGON (((-97.13797 25.96581, -97.16677 ... \n", "\n", " field_avg \n", "5 1.584487 \n", "9 1.491957 \n", "37 2.481977 \n", "67 2.574049 \n", "98 2.893058 \n", "155 1.262718 " ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "regs[\"field_avg\"] = out.values\n", "regs" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqoAAAE7CAYAAAAVV6MRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAABmmklEQVR4nO3dd3wc5Z348c9XXVa1bLnKtlwxGAMGY3rvJUBISCDlSCA4kHJJfiRHSC6XXC5cCndwOUgCvoSDXEISEpqT0AkldAy42xg33Itc1Sxppef3x67EWNrZfWZnRrNafd+vl14ezT478+zs7teP5jvfZ8QYg1JKKaWUUtkmL+oOKKWUUkoplYwOVJVSSimlVFbSgapSSimllMpKOlBVSimllFJZSQeqSimllFIqK+lAVSmllFJKZaWCqDuglFJKKaWSO//8801DQ4Pn57311ltPGmPOD6FL/UoHqkoppZRSWaqhoYEFCxZ4fp6IDA+hO/0uJwaq5eXlZtiwYVF3QymllFI5ZMOGDQ3GmNqo+zGY5cRAddiwYdx8881RdyNrdHV1Rd0FlSXy8vQydJWaxgvVTeNFXzfccMP7UfcBYDDfRVQ/lUoppZRSKivlxBlVpZRSSqlcpWdUlVJKKaWUyjJ6RlUppZRSKovpGVWllFJKKaWyjJ5RVUoppZTKUsYYPaOqlFJKKaVUtsmZM6oiEnUXrA3Ev4wGYp8HmoH0Ge42EPs80AzE795A7PNAMxC/ewOxz9liMH+n9IyqUkoppZTKSpEOVEXkayKyTESWisjvRKRERGpE5GkReS/x79Ao+6iUUkopFaXu61S9/OSKyFL/IjIW+EfgMGNMq4g8AFwJHAY8a4z5kYh8E/gmcFOabXm69Vt/voHJ9mWT/oi6j/3x3GzbV3+lpdz2Y/P6+jN15mdfA6WfNvrrM+62H40X2bevbPh8a7zI/n2pYESd+i8ASkWkABgCbAEuBe5LPH4fcFk0XVNKKaWUit5gPqMa2UDVGLMZ+A9gA7AV2GeMeQoYaYzZmmizFRgRVR+VUkoppVR0okz9DyV+9nQisBf4o4h8ysPz5wJzAWpqajyl/r3q6ury1N75l0y6NINNui+Mv4xstul1v2G3D5vX/nhJIXl9P23SfWGksGy26XW/XtuH+V0OW5ixond7J40X/SvMWNG7vcaL1AZyvPAimz7//S3Kd/hsYJ0xZqcxpgN4CDgR2C4iowES/+5I9mRjzDxjzGxjzOyKiop+67RSSimllOofUQ5UNwDHi8gQif8JdRawApgPXJ1oczXwaET9U0oppZRSEYos9W+MeV1E/gS8DcSAd4B5QDnwgIhcS3wwe4XN9oI4/e92aj0/P99TeyevqcBk2w4qref2XK/rw2jjp33YvFbdJmsfVArO7bMQRsWt1/VONt/HsFOCYUv2/kcRK3pvX+NFdIKIFb3ba7yw346f9ip7RXpnKmPMd4Hv9lrdRvzsqlJKKaXUoJZrVfxeDY6rkJVSSiml1IAT6RnVoHid8N/JT9WtTYqtOxXo9rgz3Wezba9pPT/puzBSfAMlreen6ra7fVBpPa+fBT+VuDbr3b5rYezX63bClumMHl5iRao2Gi8y206YgogVvdc7abzI7Ll+2mcjPaOqlFJKKaVUlsmJM6pKKaWUUrlqMJ9RzYmBqoi4Vtt64acyN9PUlTMtYlP1G1QVb1CpvDAqg/1UP4ch04rUoNJWNoKq4nXbjvMY2Lwur689jKrfsKX7LAcdK0DjRTLZFC/8fI41XqRu47Z9m36qgS0nBqpKKaWUUrlqMJ9R1WtUlVJKKaVUVsqZM6oFBfYvxU/KzKYyM9n2bSpDvab1guKWmvOT1gsj7Rn2X5Ru763ztXhJUdmkv4Kqvg2DW/+9puy8pgFt1rv1MwzpvotBx4rez9V4Yb/ea5tMBR0rQONFuu24tbdZ79bPgWQwn1HNmYGqUkoppZSyIyLjgF8Do4AuYJ4x5qe92nwD+GTi1wLgUKDWGLNbRNYDjUAnEDPGzA6jnzpQVUoppZTKUiHemSoG3GiMeVtEKoC3RORpY8xyx75vBW4FEJEPAV8zxux2bOMMY0xDGJ3rlhMDVZuqf+eb7Hbq3ybNZJOiSpa285qa81Pd25/LbhOQRzVBuFdBTTqdrI3bZ8ImtRUGm3Sc1/Sdn+c6+Un3eW1jI9N4EUSsSPVcNxovUuuveOGnwl3jReplp2yLFwOVMWYrsDWx3CgiK4CxwHKXp1wF/K6futdjYF6soZRSSik1SHSfVfXy44WI1AOzgNddHh8CnA886OwW8JSIvCUiczN7ZenlxBlVpZRSSil1kOEissDx+zxjzLzejUSknPgA9KvGmP0u2/oQ8HKvtP9JxpgtIjICeFpEVhpjXgys9wk5MVAVkbRV/14njraptPWS3nJL8XV2dqbtl3OfNu3d+EnT2aTswqgGdhNV6t9Las/teTYpPq+VqX5SWF5TcM7LbMJO8bltx00YqX833d+JoGMFaLxItmzTNzdRpP69Xgag8SKz7TtFFS+yWEO6AicRKSQ+SP2tMeahFE2vpFfa3xizJfHvDhF5GJgDBD5Q1dS/UkoppVQWCyP1L/GR+q+AFcaY21K0qwJOAx51rCtLFGAhImXAucBSny8zqZw4o6qUUkoplatCqvo/Cfg0sEREFibWfQsYn9jnXYl1HwaeMsY0O547Eng4cVa6ALjfGPNEGJ3MiYGqTerfya0q1iZV57Ydt4mek1X9O5edaRG3NJ1NStCN19SZWz+9pvJsLqPw0z4MQU9G7SeF5eSWwrJJD7rxmpZ0S9+5pR9t0no2r8trBXAY0sWLoGNFqvYaLzJrH7SgY0Wq9hov+i7b9CGqeDGQGGNeAtJ+AIwx9wL39lq3FjgylI71khMDVaWUUkqpXJRJFX8u0YGqUkoppVQW04FqDvBS9e9MS/hJS7mlKJwpuXQfLpu0XhiVzH7Sd0GtD6rq1w+vKS23tG2yz4Kzj25pLq8TtbtdnuLnMxJU+s7tNUZV9euHl3jRX7Gi9/Y1XqRednstfiQ7tkHEit591HiRel/ZFi9UuHJmoKqUUkoplYsG8xnVSK8qFpFqEfmTiKwUkRUicoKI1IjI0yLyXuLfoVH2USmllFIqSmFMTzVQRH1G9afAE8aYj4pIETCE+NQIzxpjfiQi3wS+CdyUaiMiQmFhIeC9stVrla7bc91Sb8kqc93SDX5SRX6W/aTy/FT9eq3c9VO17LWy1U+Vbvfrsknr+emj1/Snn3SZs71Nis+tjU1az0+lrxu3NkHEizBjRe82bvvVeJH6uW4yjRdBx4rebTRepF4fVbyIUi4NPL2K7IyqiFQCpxKfbBZjTLsxZi9wKXBfotl9wGVR9E8ppZRSSkUrytT/JGAn8L8i8o6I/FLidzcYaYzZCpD4d0SEfVRKKaWUikwmaf9cOgMbZeq/ADga+LIx5nUR+SnxNL8VEZkLzAUYOXLkQSkCL9xSS25pGLcq3XQpPOf2YrFY2v24rffDT7rP7XU71/upAPaT4vMjjLRddxuv1cheZ6MI4xjYVBs71ztn27BJ6/mp9LWZQD1syT6nQceK3tvUeJF6X1HHiyBiRe8+arxIvZ2BEi+ClEsDT6+iPKO6CdhkjHk98fufiA9ct4vIaIDEvzuSPdkYM88YM9sYM7uqqqpfOqyUUkoppfpPZANVY8w2YKOIHJJYdRawHJgPXJ1YdzXwaATdU0oppZTKCpr6j86Xgd8mKv7XAp8lPnh+QESuBTYAV6TbiIh4mvDfbb1Nhaxbqsatkrc75eN83K0y0a2K18lrCsPtdXut0LVJ2dmk8txSgjYpPqf+TP3bpJOSpaJsUnBe7wduc3mKDT+v223ZLWUXVHVvf07g7SVeBB0rerfReJF6v1HHiyBiRe8+arzo2yab40V/yKWBp1eRDlSNMQuB2UkeOqufu6KUUkoplZUG80A10gn/lVJKKaWUchN16j8QztS/1786bNJJbik2tzSfs72zYjfZ85zcqjjDSGfapO9slt1SfG7rvR5vr6k8t/VeUz9uaSy3VFqy9W7vpxu3e3e7pbDCSGe6pdfc0nFubZyX4tik8vxUEtvwM+G/W/vu9yvMWJGqjxovUq93ey02673EiyBiBWi8SPfcbIgXKho5MVBVSimllMpFuVYc5ZUOVJVSSqkByhjDvn37aGlpYffu3ZSVlTFp0iQ9K5hjdKCaA5JV/duke2xSZjYpLedznSm87n7ZTNrttj23tJFNOjOoSbtt0nTO1+in0jfsit6gKv3dKlW7++PsV7pZKcC9GtzPe+jk9fV5reJ1S9/ZTPJtc79um346BVXx7iVeBBEreq/XeNF3Oep4EXSs6L1sEy/279/P4088zer3ViF5+RQVD6Gzq4u21kZMVyfDa0cybeokLrroIkRE40WSZadsH9jrQFUppZRSA8bjTzzJ8mVLmX7spxk2ZuZB0+I17n6fPTve5aWXnuGdd97hO9/5TsS9VSpzWvWvlFJKDTAzD58BgOTl9zljWDmsngmHnsf46efS3NzMnj17ouqmCohO+D/AichB1YHdbN4om2pWm0mEne2TpfBsqjvd0hk26Ts3XqtmvU7U7Uzf2TzXZgLvbE79u70vyd5Tt365pbacy26XAdhUD9sIe6Jut7SezbJNf2y+B0Gl/t3adB//oGMFaLxI1iab4kXQsSJVv5LFix07dvCHP/yBQ2Z/gmGjZ7j2ffz0c9iw8inuuOMOvve97yXtmw2NF9HLpYGnV3pGVSmllBogurq6uPfeexk5YQ61dbPStj/6zK/T0tLC2rVr+6F3SgUvJ86oKqWUUrnqb8+9AKaLqqoq3nhjAZ2dndRNO8PquUMqR1JWOYr58+dzzDHHUFdXx6RJk5JmIVV2yrVUvlc5MVAVEQoLCwHvlbtuVbw2qTSbSsJ0sxE4g4VbKjeMCbz93LvbLX3nterXa3Wvn7Se13S/1/RWsvfUT4rRrerfT7rfZr9e03fOz7fNpN1BpfWCuqd3EPEizFjRe58aL1Iv27xvra2tvP7GW+zbt4c5x85m1KhR9JbucxR0rOi93undd1exZMkyRow/lo07Ghg/8wqqa6cgYv/9P/S4z7B59d949vlXaWlsoKpqKJdffhmHH3649Ta6DeZ4ESUdqCqllFI5bvnylSxfsZz8wgpampp59tlnOf300xkxYkTPyY5scuDAAR5/4nFmnHAdVcMn0xlrRyTP0yAVoKRsGJOPvAIAY7p4f/lj/O///i833ngjEyZMIBaLsW7dOtasWYNBuOjCC/SMa5bRgapSSimVw2KxGG+/vYBDj/ts/IxkXj47NizglVcfo72thTFjx3HKySdazWHaX9avX09XZyfrljwI5NPavAsQ6qaeRt20sz0PWAFE8qifcTG7tizm3nvvY+q0Q3jt1ZepqKplSNV4tm94i5NPOpHa2trAX4/KLiIyDvg1MAroAuYZY37aq83pwKPAusSqh4wx3088dj7wUyAf+KUx5kdh9DN7vpE+eK3693oPa69pg2TcUnZuVbxeU9Ju1bF+Jup2S8d5XfaT1vMzsbeboCbrTleZbXNGwmb/zmNmM1m42/Hw8/nyeu9ur+v9pPiimvA/06p/Gxov7PeVKl68/MprtLe3U11VSXt7OwDl1WPJLygCYFT9cYyqP44DzbtY8tIvePPNBcyefQxOyY5n0LGi9/puU6ZMobi4mLKyMvbt28eoUaPYtWsXjz76KPsaVjPjxBsyTlMffvIX2LV1GdubuphxwucYOvIQlr38C0aNHkNFRYXGiywT0hnVGHCjMeZtEakA3hKRp40xy3u1+7sx5mLnChHJB34GnANsAt4UkflJnutbTgxUlVJKKac9e/awdes2Jhx2IY2NO2g7sIcZJ36O4tKqPm1LyoZxxClf4s2nbmHmzMMpLi6OoMd9iQj19fW89fY7LFq0iKbG/RQWlgCwd+c6ujo7egbdXhWXVjNm0kk9vzfv38aenWv5wQ9+kJWXQQx2YQxUjTFbga2J5UYRWQGMBWwGm3OA1caYtQAi8nvgUsvneqIDVaWUUjklFouxbPlKRk6Yw+iJJ1g9p/3AfsDu9qX9obGxkfb2dhYtXsz7G3cyevK5VA6byJCKERgTP4ucSeq/t+5tLXo+nvEdNmyY722qgUdE6oFZwOtJHj5BRBYBW4CvG2OWER/QbnS02QQcF0bfsuMb6ZPXqn+bVJpNSssmVdC9fbcJn23ShkFVd9uk9bym3cJO8dm8b05+qv5t0nluqah0KTy3fdp8trym+L3yeu9u5+sLO/XntW9hVP2nixdBx4re29F4Ybevffv28dJLL4Pk09y0D4BR9aPdX2AvnZ3xSwNeePEVSkuLyZM8KirKmDixvuf/Fy+pfy+xIhaL8Xpi2qmJ9eNpbGzizTcXUFxaQXl1HdPnfJbS8uGOffp/j9ta97F+2V/YuWkhpWVDAcMPfvCDtM8bzPEiKj6mpxouIgscv88zxszr3UhEyoEHga8aY/b3evhtYIIxpklELgQeAaYCyQ5SKNcn5MRAVSml1ODW1tZOc3MTI8cfy+GnXkxby25Kyoanf2JCde0Ujjn7Jpr3bSHW0YIxXWzdsZIVK59gXN1Yxo0bz4gR4RQYrVz5Lhs372T42Fm8tXgdeXmFHH7yFyivHhvK/gD2717Pzk0L+cY3vsG+ffuora3Vs6lZLMOBaoMxZnaqBiJSSHyQ+ltjzENJ9rvfsfyYiPxcRIYTP4M6ztG0jvgZ18DpQFUpNSC1t7fzta99jeeee468vDw+97nP0dbWRllZGaeddhpPPPEE06dPp66ujsmTJzN06NCsO1OigjNiRC114+rJKxtKYdEQCouGeN5Gafnwg85cjp54Ii2N21n+2r1s2Pgyl3zoQ6FM25Sfn0fz/h0cetzRjDvk7MC3n0xF9Tjy8vKoq6ujrq6uX/apsovEA+KvgBXGmNtc2owCthtjjIjMIX5H013AXmCqiEwENgNXAp8Io585M1BNl+5yS/G7pbFsqhydE1mn4/bXkDOF5QyANvd697pfm/SdW3Wvn2XncXK777efib2d6234qdh1u0+7l792bT5bNlXffi4D8JO6dEvN2SwHdX/vvXv3ctxxH1wONebIYby59zmadx/gQEOM+c89xKYV2+jt5ptvZurUqbS2tnLcccdRXV3dpw14ixdBx4re+3HSeNF3vXN5966dtGxcz/hDzun74jI0pGIks8+5ieWvzuPV115nUuJSgOLiYkpKShARurq6KCwsdP3eusWKWCzGosVLWLVqDeOnn03xkOrA+u2mM9bOkpd+gUgeXV1dGGPS/gE30ONFULNxRCmkqv+TgE8DS0RkYWLdt4DxiX3eBXwUuEFEYkArcKWJdyYmIl8CniQ+PdU9iWtXA5czA1Wl1OBRUlLSszz3yYspGlJAfnE+BUXJz3bt3dTEgvve5d6H/4ftq3fTFYsPIgqLCulo7+DHP/kxH7n8I/3SdxWOWCxGS0szU2d9LJTtT5v9D6x4/R5WvLeFzs4O2g800n6gGcnrnoDfUFZWzmGHHsLQoUOpqKjoeW57ezsbNmxkZ8MuWloOJAbdhn1791BRM44jz7iRouLyUPrdLdbeSmvTTha9eEfPuo985KMDarA2mIVU9f8Sya81dba5E7jT5bHHgMcC71gvkQ9UE3NxLQA2G2MuFpEa4A9APbAe+JgxZk90PVRKZZvS0lLKy8tpampi3nl/6Vl/zaMXMHRCRZ/21XXlnP3tD+bHNF2G9a9sY9mf1/Huk5vYsjmUS6tUP2hubmbVe2vYtGkTxaWV1NbNCmU/BYUlzDz5CwetM6YLjAHJo7PjALu3LWPVurdpemcxRUVFNDXupap6KM1NzVSPmEpV7UyGjawhP78QY2By5chQB6jGGF5+9J/Iyyugqyt+RnpIWRk3/r//R1lZWdbMcKDS0ztTResrwAqgMvH7N4FnjTE/EpFvJn6/KdUG3Cb8d/Jy725wnxzb7a9PZ1oqWZrMLX0cxn2KvabpbG56YLPslr7zc69vm9kAwpjw3yYVFUSlv9fLHNxSvzYzInidtDuo9W7pO5tUYap9LVq0iLa2NrZs2cL8Pz/Kf//0Djpa06fY25s7eOSLr7BnfSPDh9dy++23c8kll3i+mUKyCvAgYkXvbWq8OHi5q6uLWCzG4iXLiMU62bFjG7V1xzDzlMsYUjEy6WsMi0hez/mogqJSRoyfzYjxszFdnTTv30pr0066OmMMGzOTgsKS1BvLQFvrPtYseghjuph+7KfIL/hg/tf9u9fz7pu/AaCrK8bQmuFcfNGFTJ9+SM8sBt2p/+SvLbfiRbrlbOaj6j8nRDpQFZE64CLgFuD/JVZfCpyeWL4PeJ40A1Wl1OBUXFzMpEmTWLxoMYddPIER04embG+6DHec9DAAL7/8MmPGjKGtrY2GhgbWrFnDCSfYzbmp+pcxhq1bt/Lqq68CkJeXT3n1WEbWn8CMKSOpGDouzRb6l+TlU15dR3l1eEVKzfu28M5zt/f8vmbRQ0w75iq6Ojt4983/Y9e2FdTU1PCF679GZWUlRUWZ3RhAqahFfUb1v4B/Apy5upGJuyVgjNkqIiOSPVFE5gJzAcaMGRNyN5VS2epvf/sbzz//AqffeGTatq/dtZKjZh/Jb+79LWVlZfzbLd/nnl/+L0PKSmlpbmXp0qWUlZX1Q6+VF8uXL2fVqjWMm3YWFTXjaW3excjxsykoLI26a5F46ZFvICLMnHkE7723igMHDtDRfoCXHvkGAOPG1/OJL3+ZmpqaiHuqgqJnVCMgIhcDO4wxb4nI6V6fn5i0dh7AEUccYbpTBF7v1+2W+ndb71a967bN7rSEWxWxWzojjPSEn0m7/VTu+rkkwG05jAn/vVakpgscNukmm8tN3PYZVODyeimEzUTdXi+jcKvotUnldd/uctKp6f9gXfv8Vs469jwef/xx7v7lXVxy8SUAXH/99Xzm6s9SWRm/AslLvAg6VvRer/EC1q5bT1n1WMZPPwfJC356qIGg5wxqYrL/f/iHf6C0tJQlSxYDsGf7cmpra7nkkkuora1FRDRepOnzQJnwH3SgGpWTgEsSdzooASpF5DfAdhEZnTibOhrYEWEflVJZrrQ0flYtvyj1VEwH9rezbVUDT+55ghNPPJE88vnqV77G1f/wGYYOTX3JgIrW6FGjWLduLa8//q8cd+G/ZuVAImwtjTsQEYYPH8bxxx9PVVUVryQuhQC44YYbaG1tZcSIpElINcAN5oFqMPfay4Ax5mZjTJ0xpp74RLF/M8Z8CpgPXJ1odjXwaERdVEoNAA0NDQB0HOh0bbN9+W5+duojlJWX8daCt7nkkkt49plnAXSQOgAcddSRXHjhhcQ6Wtm/a13U3YnE8LFHMOHQ89m7r5k///nPvPPOOyx4800AysvLWbZsmX6Wc1h3QZWXn1wR9TWqyfwIeEBErgU2AFeke4Jb1b/XyavdKrOdnCkHZ7oqXUrYa+Wu20TdbmcS3Prr9VIIr9W6XtN6zuWOjo6k673c37t3/224HXO3fTnfWy/TuXi977RN1b/XSx7cuH2ObCbM9nqP7qAm8HbrQ1tbGwD3fvgJPvPQ+QybVEnD6n3s29KM6TTsWrOPBfe8x4knncjn536+p+K5t0zjRdCxovc2w4gXbrI5XixfsQKAtUseZebJN4RSRZ/NRPKom3YmVbVTWPTCHfz9739n7JTT2Lz6BZqamnjxxReZNGkSo0aN6nmO13jh9lnOpXgxEGcAyLWBp1dZMVA1xjxPvLofY8wu4Kwo+6OUGjg+8pGPcM4553DYYYdx7+VPUDWygvamDmbOmonp6sJ0FfODH9zCJ64M5e5+qp9MmjiRkuJiNmzYyHtv/45Dj/ts1F2KRMXQ8Rx34ffYuvZlNqx8umd9fn6+nlHNYTpQVUqpAayyspKnnnqKb9x8I7H2Lv76+l97zpxm85kSZa+oqIj6+noqKip48803efvZnzBm8mmUlA2junZK1N3rV50dbT2D1CuvvJLRo0dH3CMVNh2oDnB+Uv82lf5u/9G5pfCSTcQeVHrCD6/37vYzabfb+vb29qTrvc4A4NZPG27H2S395DUN373NoNL9Qd3owI2fybxtjqXXSt9MvyszZ87kyb/G//N+6aWX2Lp1K+Xl5ZQOKeW0U09Le9/7TONF0LGi9+vTeHHw+srKSk488UR27tzJhpV/IdbRxfEX/1viNqaDw4Knf0hxcTHXXHNNn+p+jRf2fRsoqX/QgapSSuWUk046qec/J2NM1v8npLwpLCxkzJgxjB49mkWLlvLaX/6Z2ed+m8LiwTEHbmXNBJr3bWbFync57NDpUXdHqVANnj9BlVKDxkA6U6IyJyJMmTKRzs4O8vIHz3mXI079EjNOnMvzz/3Ndb5elVu06j8HBJH6d5sBwI2Xybr9pCS8pvhsXrefKt7u+0O3tbXR1NREW1sb+fn5VFR8cIMxm3S/c71NStAt9W9zGYDNMXRL/TurTW0uM0iW+ndLW9kcbz/vZ7J+pepbUGk9P+k7m3t3h5HKyzReBB0rei/nQrzoZlPp7zVexC/9KmDJS7/gqNO/2uf15arKYRMpKCzh+eef5/zzz+9Zr/Ei9fYHYuo/1waeXuXMQFX1rwULFrBixQrKSqrp6GynvaOFcWMnMKZuFCNHjoy6e0qpQWLNmnV0dsbo6mxP3zhH7NqylBVv3AfAqlWrqK+vZ/p0vQRA5SYdqKqMTJgwgRUrVjBr/FVs37eC97Y/z8bN77Nx8/sU5BcytHoYo8bUUllZ2TNf6u7du1m9ejXFRaWAMLZuFGPGjMnqv2SVUtmtoqKMPXt2MePEz0fdlX5TUlF70O9PPfU09fX1lJQMrrllBxM9ozrAOav+vU5Y7Xbq3ybF60wnuKXwuvvlTFXZpDxsUio2bNI9bhPpp0rfDRs2jOHDalmy6VGmjjoDiFfj7t+/n+HlUxlZdig7Nq1gY8daSOyqsyvGpBEnM27Y0bS272PZ5vk0NjZTXz/+oD47bwTgdcJ/P2kst315rZztfh/d3me3mwnYpEu9pu/ceL0Uwk/lrp90n9eJvb3euzuIeBF0rIDcixfplt1mALGZJaS8vBwgPl3VlNOYMP3c5C8uh5RVjOTky26lq7ODBc/cSnvrHnbt2sXo0aM1Xnjog9d4ESUdqCqVgZNPOYlHHnmEzq74wPKEE05g1ar32LuzkVHVhzKq+lDX55YUVnL85M+xbNOfWf3eOqZMnQh8cC1ONgYKpVT2KSoq4sgjj6SxsZG1K59m65qXiHW0Mmz0DKYd8wnyC4qi7mIotq1/g9UL/wjEB1aPPfY41157TcS9UmEZzANVrfpXnjgv6i4rK6OubhzvbX+O0qIq/v7i3+nq6mRP00Z27l+ddlsF+cUcOvYi9u/fx+5du2lpaWHBggUsWrSIhQsXHnQ2RQ1OLS0tNDc3R90NleVEhMrKSmbOnMmwmmqqq4fSvHcdr/7l22x49xnaDzRmvO2OtmZi7S3sa1jD/l3rg+u0TxU1H2SijDEaL3OcVv3ngGQprqDuZ+y2TZt0XrKqf7c0hJNbWs8Pr1XlydovXryYnTsaOPmUkwCYPfsYVr37HstXLAPyeP/99wHhnff/yOmHfoWigiEp+1RUUMqsCVewbNNfaW3fx4Thc5g66gyWbJzPmtXrqJ84vuf1+0nHO9mkrpyV/jbbT5aWctu2TerUa5W1H279dGvjNdWWbjkWi3HX3Xdx7TXX9qwfMmQInZ2dHHHkEUyZMoWnn3o6svvcZxovMo0VqZYHWrwA98+78xIft7S+zSUBvS8hGDq0CqgiFovR1NTEhhVPsmHFkxx+0uepHDaRvLzkn/HeYh0HaNi8iNUL/3TQ+iNP+zJlVWOttxOWsspRnHTpj3n50ZsA+PCHL+szQMnFeJFq2esMAG59zja5NvD0KmcGqqp/TJs2jeHDh/f83tnZybbtWxk7diybN28GJHGnlC5eevfnHDnhowwrr0+5zWEVkzjlkC/RdGAH5SUjEBGOHH85f3/3ThoaGqitrU35fDWwdXZ28sN//yFr167lD7//A2eedSYf/chHefHvL9LU2MTdd90ddRfVACQiVFRUUFFRwZ49e1n6cvxzVDPqMGrrZjF0xDT27VpPQWEJDZsXsWPjAoaPPYpRE45D8vJZ9uovKSwQRo8eTWlpKSLCmjVrWPTCHQBMOeoKRtXPifIlsm7pYz3LOtuKylU6UFWelJSUMGzYsJ7fRYTm5hZ2794NwAlTr+X1Nf/L1KlTWLNmLQvW/oaj66+itnJyyu2KCBWlIw/6ffbET/LSqrsoLCyiuroqnBekIrd48WKG1w6nvaOdc84/h7Vr1/Ddf/0uH7viCv75n/+ZSZMmRd1FNcBVV1cxZEgpO3bsYPe25ezetvygx4eUVTCktIT9O5ax/f03AKisqmL4sGEHnYGbMmUKra2tNDU1sXrhHykqraRmZDTTQh1o3s2WNS8AcPnll0fSB9V/9IzqACciaVP/buk2r5X+XtN53dv3OuGw12W3PvqZwNsmhZSXl8c555zNCy+8SGPjfl5975eMrDyUVatWMH36dFauXMnb63/HmOojmTn+Q0n762ZIcQ0zxl7Equ1PU1Y2xLXy2E/q323yda+XFiRLM3mt9Ld5H8KYwNvPsp+JvbuXZ8+ezUMPPsS0adOAePo2Pz+fgoICX5N2BzXhf6bxItNYkWr9QI8XblX/zssA3C4JsJkBwO0yg87OTkSEkSNH9vTTGENXVxe7du2mZmhVzzHv6hqKSDwzFIvF+sSL4uJiSkpKaGlpYe/2lZEMVPfuXN1zhvjUU09l5MiRPa831+NF72Wvsw0MxAn/YXAPVLWYSvlWUFDAWWedyVFHzQJg+/4VAKxcuZITTjiBkuIStuxdxKZd73je9uihM8Hks3HjZvbu3UtbW1ugfVfRy8/P7xmkAhQXFx90jbBSYeie1nDEiNo+A5x0gxYRoaSkhLbWvSH3si/T1dkzSB0zZgxTp07t9z6o/jeYi6l0oKoCM2rUSM4//3wmT/4gzb9o0SIOmX4IAMs2/5W2jiZP28yTPE455AtMGn4qpq2czZs3B1YUoJRSmWhpaaGxsZHy6nH9vu/O2Adnlp23TlW5bTAPVHPmtEX3X8Ruk3P7Seu6tfdyz+5sS0N4TevZpPK6U3B79uwBYGTF4WxvXEpnZyf5+fl0dnby/Ir/4qwZ36Agv9i6rwX5xUwYPoc8KaCtazednZ2BpbF6X8aQbjtu23R+FrqPg9dKf6/p/jD4SeWFPVF3qu/Qnj17uOzD59Pe1gLAySefyXe/ewsVFRV92jsFHS+CiBW9l3M9XvReDmpmAD/fp1TxwhjDli1bGDn+WMYdclbS54dl5+bFvPvm//X8nmxQovGi77LX7022XQYQ1sBTRMYBvwZGAV3APGPMT3u1+SRwU+LXJuAGY8yixGPrgUagE4gZY2YH3kn0jKoKwahRowDY3rgUgKVLlzJjxgyOOeYYAJ5ddivtsRbP221oXE1JaW5O3q38ERGWLX2Xf7+1kR/8uJFdex5izpwjePbZZ6PumsoxPYNnERo2L+q3/bY07uDdN38DQG1tLZ/97Gf7bd8qeiGdUY0BNxpjDgWOB74oIof1arMOOM0YcwTwb8C8Xo+fYYw5KqxBKuhAVYVgzJgxzJkzh9GjR/es27RxM4WFhcycOZPqqqG8te63xLq8TVC9Y/+qg86YKNWtqqqKwsICxo3LZ+YRhfzoP0r4/r/H+OIXP8HnP/+ZnrP8SvkVP1MnbH//DVa++ZuUAwJjDHu2r6StdZ+nfXR2dhz0+4JnfsLbz94KGD784Q9z7rm5f5tYdbAwBqrGmK3GmLcTy43ACmBsrzavGGO6A+hrQF3ALy2tnEn9d5+qd0vfua13S/062aS60qU6vFY7OtlM+O3GJg1kk35OVU2bbLmgoICJEydSU1PDsmXL6OzsoqOjAxFh3Pg6NmzYzGvv/YoTpl5Hfp79x7D7EgKbtF539Xh3BW83t+PvNbXiTEs5U5Td23Q7Nm7HzyaN6uQnxed2DJxs0mh+Um1+0nfJ1o+tq2XLllYqq+KPnXp6MX99upD/+NEzzD52Jrff9nMuvfTSg15j0PEiiFjh9vpStXEaiPHCbQYA53Nt0v1+bqbhlOo9EhHq6sbS1tbGzp07Wfry3cw8+frEfmJsW/86sfYW9u58j86OAzTv3wpAcWk1U4++kura1FP1vfTIN1wfu+yyyyguLu5z+ZPGC/v9DpR0f38TkXpgFvB6imbXAo87fjfAUyJigLuNMb3PtgYisjOqIjJORJ4TkRUiskxEvpJYXyMiT4vIe4l/h0bVR+VfZWUl06ZNY9Lk+p51IsK4cWNo7djL1r1LQ9t3S3OL3lZwEKkbO5atWzoPWldensf3flDCf90hfOvb1/HJT36E7du3R9RDlUuKi4spKCxiX8MaOtqaaW1u4JX5N7N28SNsWPkU+3eto7N9LxMnTmTKlCm0te5ly5q/A7Dt/TeIdRw4aHtNezfx3jt/SrYrZs2axeWXX05JSUnor0tlpwzPqA4XkQWOn7nJti0i5cCDwFeNMftd2pxBfKB6k2P1ScaYo4ELiF82cGqgLzohyjOq3ddGvC0iFcBbIvI08BngWWPMj0Tkm8A3OfjAqAFEEneHcerq6mL16tVUVVWyvuFVRlfPID+v0Gp7e/bsobS0NG07YwzFJfEzDzt27OhZX11drcE+R9XV1bN1y4qkjx17XBF/eaKQn97+KnPmHMEPf3g7V1111aA/i6L8qawoZ/fu3bz++Pd61k2YMIGCgoKDplgrLi6mrq6OTZuW9ZwxXf3OH3tt7YPPYlVVFfUTJzJ2zBiNVwpjTKaz3TSku3ZURAqJD1J/a4x5yKXNEcAvgQuMMbsc/dqS+HeHiDwMzAFezKSjqUQ2UDXGbAW2JpYbRaT72ohLgdMTze4DnifNQLV3ate53rG/pM91pgTc0nNubbykH2zSEGFX9NpUhvqp4rVp053WKy0tpbCwkD17dtDZFbMaqE4ZeRqrt79wUGqw97667dq1izwpYFT1YYweVcOm3Qtpbd9Lc3MLhYWFSZ/rltJycr4XydL9zm16PU5+Jvn3I6hJu72m6fxsM1mb8eOnsmXzB7eU7K2kVLjpWyVccFEHX/r8P/LWW29x2223Jd2ml3gRdKxItc3BGi96r3dL99tU/Xv9Dz/Ze9q9rqioiFGjRtHS0kJ7ewcHDrT29Kd3rKioqGDy5Mns3r076TXT1dVVTJgwgSFDhlBZWdnz2jRe+N9mVLNlZDuJH4BfASuMMbe5tBkPPAR82hizyrG+DMhLjN/KgHOB74fRz6y4RrXXtREjE4NYjDFbRWRElH1T4Rg1ahSbN28GsJ6qavLIU9jVtI7mpl0MKStNG2ROO/TLFBWUAVCYX8ryzU+Qlze4A1Muq6sbx9Jl6f/gOeLIQv79J118+YZf8dZbr3Djjd/iwgsv1JsMqIwUFBT0DCyd160nU1RUxIQJE5gwYUL8NqyrVwMwefJk6ur6vUZFDSAhTTV2EvBpYImILEys+xYwPrHPu4B/AYYBP0/8n9s9DdVI4OHEugLgfmPME2F0MvLI3PvaCNu/cBLXWswFGDeu/yddVuml+2JVVVXR3NzKjv3vMrJyutVft0dPvJKX3v05sVjsoDOjTh0dHRTml1KY/8ElAss3x6//Lisr8/AK1EBSV1fH1q128ePkU4v569MFPPLgOuZ+/lPc8oP4pQBDhgwJuZdKQWNjI2vWrOn5fc2aNTpQVSmFMVA1xryE85qT5G0+B3wuyfq1wJGBdyqJSAeqLtdGbBeR0YmzqaOBHcmem6gumwdwzDHHmGQT/mfQn55lm3tzu1UJJ9uOn3RJGPykjWwrd/fv209jUxPDhw9Lut+CggKqq6tY9P6DlBRWcvzUaykuSD2QLMgron748by79Zmee25396+tra3nFquHjb0QkfgxXLLh0YP27/Z+OjlTiM73xfkak03y3/3au19fsue5Vfp7Ta+GwU+a2U+KL4hU9/jx49m6pQOwu955zJh8vvDlIWzcAF/60pf4l3/5Flu2JA03SXXvN+hYkapNLseLbm7pe7fLBry291r1n+7yoExixbbtO6mtO4r6wy7izaduoba21uo4abzIbDteU//ZeElASGdUB4Qoq/7dro2YD1ydWL4aeLT3c1V26ezs5P333++zvri4mOqqqpTPraysoL6+noJiw6vv/Q+rt73A9n0rU34p62uPj+83JjQ1trBv3z7279/fM0gFGF4xiZb2Pazb8Qpb9i4BsjP4qOCMHTuWrVub6eryFtC//++lvPhqLa2trSH1TKmDlZeXsXPTQt586hZAs4Iqte5iKq8/uSLKM6pu10b8CHhARK4FNgBXRNM9ZSs/P59x48b1mQqquKSYYooPGkAmIyJUV1fR1tbGjtZFvL/rAJWlYzhs7IWUFdckfc55R/wzxhh27n8PQxe1FVPJy8sn1tnOCyt+ysur7kIkj1hnB6WlQ8jPz9NrEHNcaWkpFRWl7GroonZE6msFnQoLhZpheQziExaqn9UOH0ZFeRmtra09d/JTKpXBfEY1yqr/VNdGeL6BcrLqeq9vrE3Frs3k38nO3AV1H2+btIXX9FBQVbw2y6km5y4sLKSwsJCSkmKam3fyyqq76TKdHDb2AmorplJSVNnntdZWTj3oODS17STW1UZBQQElJUVA/JarIpK0GtfmM+Lsp/O5ztfifE+729tM+B9UFa+fScz9pOzcUs42qbkw0npj60ayefMeTwNVgLYDps+20unuW9CxoncbjRep29hsP6gJ/9PFCy+xorCwkPLy8p7K/gMHDtDa2kptba3rzCQaL+zbBPX9yBY6UFUqi4gI5eVllJUNYe/evaxteIHlmx+nuKCcGXUXMaxiEvtaNvPGml/3POecmd8iT/JoPtAAxP9j6OzsTFuFq3LLETOP4eE/PcVRs7w97+232pk4cUw4nVIqjY6ODtra2qiurg79emM1MOlAVaksJCI9076UlBbR3NzC0s3zae9opSC/mMLCIjo64pcbvLLqbuprT2DllicpKirCGKMBP8d845/+kfmPPsYZZ5zGypWreeONN7nkkjNxJmbWrdvICSe730KyvQ2KimFXQxfPP9dG2wFYs6aL+/+vmdtu+9d+eiVKHUxEGDp0KEVFRVF3RWUpHajmgKBP2wed9hgo6QavVbw26/1MBO5cX1oar/AvLimgq6sLEaG4uIjOzk5aD+xh2aa/kJ9fQH5+QZ/teK2Qdqvedbs8JFkazutxCmNy7qAEVcWb6fLWrVv52Z2/ZPqhJTTsfoRjj8/jhJPLmXbI2336OvvYg+fl7eoyvP5aO3fc3sXrr+0DYEhZIeeecybDh4+gpLCYW289nLlzr8/o+xhGilTjReplm3jhNWXu5CVeBBErCgsLyc/P73m+xotol7NNdzHVYJUzA1U1eIhInwBTXFJMZ6wLJLsCtgrG3XffBcBfn6r29LzXXm3jC3ObaGrs5Pbbb+fZZ65l7dq1DBs2jBEjDr6XSLb9Z6+UUt0Gc3zS3KjKCSJCQWG+XpOao04++RQA1q2LpWkZ9/STB5g8bhuf/NgevvKPN7NmzVpuuOELFBcXc+ihh/YZpCqlVDYzxnj+yRU5cUa19xk2tzbdbE6h20wE7CWl4SfFEFRKwk8Vr5+0nk0KzuvE3n4qWN1eo9sE/W6fnXT995q+y4Z7d3ttH0YaO1mbs846ixEjavnwRbtYuHwEXV2G9es6mTS5bwgzxnD95/YCcP31n+M73/mOr9eYLl4EHStSLdvQeBEXZrwIOlb03o/Gi9TLXrdvcylHNqf+Qc+oKqVUVsvLy+O66+bS2NjFOafv4bhZuzjn9AZ+fW/zQe1iMcMv7oyvW7hwIXfeeVfW/weklFKp6IT/Sik1ANxwww20trayefMGjjvuRBYvXsy/fudeLry4hOHD45d8PPdsG//5kyb+9Kc/cfjhh0fcY6WUCsZgPqOaMwPVdGdNnG+yzcTXbtsOumLQT5rDhttr8poe8pP681p967WNn8mr3SpzbS4VSNdPr6k5Jz/tbfhJo3kV1PentraWH/7whz3tlyxZwr333ssf7m/l8itKuf3WJhYv6mDUqFouu+wyX33ONF6EXV2s8cK+fZjxIuhYYdtHjRf22/dz2UA20oGqUkoNMDNnzuSzn/0Ut936G267taln/f/9320R9koppYKnA9UckOlfQ25nTvycrUj23Gz/ay2ZsM8O2BQVeC2ssOH1Nn5ucyOme11+jt9AEeZ3xqb9L37xS2666Z/Zv38/48ePp729nbFjx6Z8rh/J4kV/vu5sNpjiRdCxItV6jRfBbieXvnODRc4MVJVSg09+fj5Tpkzp+V3/s1FK5Zpcm27KKx2oKqWUUkplsVyq4vcqZwaq3WdSbIogkj3P9rlBpCLCLohwE3ahhJ95/fqzOMKtEMJtvzYpvGTLXvsYdkGEm/5Mr3ndr9fte91mEPEiG9KWNvv1SuNF3+d2rw86VmTSR40X/rcf9ncoSHpGVSmllFJKZSUdqCqllFJKqayj16jmGLdUvp/nBjEHm9dUQtipB6/HxqZ9UCk+r9v3uh23bfrZb7LX6HXbNsIOVmF8TsNI69nwc+mPl+0HHSts2/hp75XGC7v2YaT4vezfb3uvBlu8iJJeo6qUUkoppbLSQD2jKiKNgGvnjTGV6bahA1WllFJKqSw1kFP/xpgKABH5PrAN+D9AgE8CFTbbGPQDVT+pP5tt5qKwJ532un2vbWwqwJ19yM/PT7t9L/sfDCmcML4D2fC9CjpeZMNrClsuxoswY0Wq7eeqXI0XQRqoA1WH84wxxzl+/4WIvA78JN0Tk9+iRymllFJKZYXus6pefrJMp4h8UkTyRSRPRD4JdKZ9FpYDVRH5sc26IInI+SLyroisFpFvhrkvpZRSSqls1dXV5fknHREZJyLPicgKEVkmIl9J0kZE5L8TY7HFInK04zEv47RPAB8Dtid+rkisS8s29X8OcFOvdRckWRcIEckHfpbY7ybgTRGZb4xZ7nE7SddHneL3WgHcn/xUxDq5fUnc1tvsK6gUoltaLdk9vVPty+tr9NKvbPtrOKjPbNgpvjD6FnWKX+NF6v3abDPoeBF1rOi9L40XmW0z275PbkI8QxoDbjTGvC0iFcBbIvJ0r7HWBcDUxM9xwC+A47yO04wx64FLM+lkyoGqiNwAfAGYJCKLHQ9VAC9nskNLc4DVxpi1iX78nvgL9DRQVUoppZQa6MIYqBpjtgJbE8uNIrICGMvBY61LgV+beAdeE5FqERkN1GMxThORO0hd9f+P6fqZ7ozq/cDjwA8B52ndRmPM7nQb92EssNHx+ybiI3mllFJKqUGlH+bErQdmAa/3eijZeGysy/pk47QFfvuWcqBqjNkH7AOuAhCREUAJUC4i5caYDX474CLZ+fiD3iURmQvMBRg/fnxI3VBKKaWUilaGl4kMFxHnQHGeMWZe70YiUg48CHzVGLO/98NJtmtSrD94hTH39dpXmTGmOW3PHWyLqT4kIu8B64AXgPXEz7SGZRMwzvF7HbDF2cAYM88YM9sYM7u2tjbEriillFJKRSOTiv/EGdiG7nFS4ifZILWQ+CD1t8aYh5Ls3m08lnac1ms/J4jIcmBF4vcjReTnNq/fdnqqHwDHA6uMMROBswj3GtU3gakiMlFEioArgfkh7k8ppZRSatCQeDXZr4AVxpjbXJrNB/4hUf1/PLAvcW2r13HafwHnAbsAjDGLgFNt+mlb9d9hjNmVmPsqzxjzXJjTUxljYiLyJeBJIB+4xxizLIPtJF0OSqbbHAzVmm5VsW7rbfblfK6fY+jcl1t/bPrm9TV62U+2VaMG9Zntz++h12MY5vcyqGOm8cJ+X1HEi/6KFb33pfHC/zaz7Rj2FtJ3/yTg08ASEVmYWPctYHxin3cBjwEXAquBFuCzicc8j9OMMRt7HWereVRtB6p7E9cwvAj8VkR2EJ/WIDTGmMeIHyCllFJKqUErpEH7SyS/1tTZxgBfdHnMyzhto4icCJjEGdh/JHEZQDq2f85dCrQCXwOeANYAH7J8rlJKKaWUylAO3JnqeuID3rHEr289CpcBcG9WZ1R7VWjd59pwAMqmSwIGCme6Kox0iXP7lnfX6Fn2k9bxekmAF2FuOxuF/b2KKk2XTZcEDBS5GC/C/j5rvAh2m9me1rcx0GOFMaYB+GQmz0034X8jySdqlfh+TWUmO1VKKaWUUull6RlSKyLyT8aYn7hN/O97wn9jTIWP/imllFJKKZ8G6kAVuAn4CfFLRvdksgHbYqoBI4wKQ5v16fbrtV/9cBeKwNu7VaG6pf68Vq36qXINY7/pXqPXbdsIO4UVxufU5ntis95PFb/bc4OOF0HEikz6pfHCfvs2Mt1vprHCto8aL+zXhxEvojSAB6rbRWQC8dkCzshkAzk3UFVKKaWUyiUDeKD6C+JF+JM4+HaqQvxSgEnpNqADVaWUUkqpLGWMyfQWqpEzxtwB3CEivzDG3JDJNnJmoJouDWeznG57tuvTbT+ovnjlls7wmuYKKmVnsx23SbvdKn29TgRus99MX4vXFJ+f98GPID7Tftu7rXer0Lb5TthUd7vtN9Pvs9f1/XlcvdJ44d7PoGOF7XacNF70Xd+f8aK/ZWOfvMh0kAo5NFBVSimllMpFA32g6ocOVJVSSimlspgOVHNA0Cn2oFIR3csD8UNmk16zSaPbpLry8/N7lp2pObdUnls/bdp4Tfc7++bldfk5fgNFmN+Z3sthp+YyfS1hv+6BYjDFi6BjRar1Gi8y2062p/K9yoXXkKmcGagqpZRSSuUaYwbuhP9B0IGqUkoppVQW04FqDuivybeDTtWFXdEbduVuGGkst/SdTYrPayrP62UAXl6Ln4nLw67o7c9UdFDfHz+pXBtBxIswY0Um29R4kbpNpvEi6FiRar3Gi8y2H3a86G8Dqa9By5mBqlJKKaVULhrMA9W89E2UUkoppZTqfzlxRtXtQmOvKYSgZgZwppy8VAa7CWMC76AmoHZWudpUv9osO5/rfO0FBR98XJ3HONnx7s1rSi6I/tscp6AmQPfDz+fRy3cg1bLN9t14ncDdTaYxIuhYkWrZhsaLuDDjRdCxovc+NV6kXrbZvpug4kV/G8xnVHNioKqUUkoplYuMGbi3UA2CDlSVUkoppbKYnlHNAWGm/v0sd/8VFEZ1ZBhs0kZe01V+0nfO9TZ99lrFa5PC85ra696+n7RetqWkwkjfhV0h71UQqf8gYoVt+2wwWONF0LGi93M1XvhfzjW5/NrSyZmBqlJKKaVULhrMA9VIqv5F5FYRWSkii0XkYRGpdjx2s4isFpF3ReS8KPqnlFJKKZUNugvGvf7kiqjOqD4N3GyMiYnIj4GbgZtE5DDgSmAGMAZ4RkSmGWM6020wzApcm4rRdJWkNtsIKs3htQLU6wT4bqkum2W3Sly3Y+nGLX0Xdurf7bUka+P2PK8pQT/LTmGkqN3eK6/fK7f+u03absO5L5tj4rbey2U9QcQK2+1ovOi73k3Q8SLoWJFqWeNF3/VRxYso5dLA06tIzqgaY54yxsQSv74G1CWWLwV+b4xpM8asA1YDc6Loo1JKKaVUNhjMZ1SzYcL/a4DHE8tjgY2OxzYl1imllFJKDUqDeaAaWupfRJ4BRiV56NvGmEcTbb4NxIDfdj8tSfukR1tE5gJzAcaNG5c2/RNUBa7XlJyXqn+v6Sw//KSNbNJVhYWFPcudnZ1Jl53tna/XmSaLxWIk4+y/12Nlc69tZx9sLg9wtu9+7TbHyU96NQx+0sxevg/g/rps7tHuldfnenntQccKr/vv/dwwaLyI635dQceK3u01XvRtk83xoj+E0ScRuQe4GNhhjDk8yePfAD6Z+LUAOBSoNcbsFpH1QCPQCcSMMbMD76Bjx6Ewxpyd6nERuZr4ATrLfPAObALGOZrVAVtctj8PmAdwzDHHZN+nSimHAwcOUFxcnJXXPimllMpeIZ4hvRe4E/i1y35vBW4FEJEPAV8zxux2NDnDGNMQRsecIimmEpHzgZuA04wxLY6H5gP3i8htxIuppgJvRNBFpQKzY+dOVixfTnFJKWUVlRTl57Fr925EhDnHzmbMmDFRd1FloT179nDHz+9g957dzDx0JpWVlZxzzjkHnZFTSg0OYQxUjTEviki9ZfOrgN8F3gkLUUW8O4Fi4OnEGabXjDHXG2OWicgDwHLilwR80VhW/Cd7E20qev2k5NxSe8naB5Xi8yPsKl6bNs5Ul031s5PbMbZJ93h9jTYVu25pO+f9wGOxGCuWL2fKqVdTUFJG4/Y1bF70JIec/Xk6Wht56aX7Offcc6mpqcmKe3eHkbKzed/cUnxu23RLafqp3M00XgQdKwDu/Nmd3H333Yw/fRI1R47gL2ueZNsbm8nLy+Pcc8/t017jRV9hxougY0WqbYZR9e+HxotoRHkLVREZApwPfMmx2gBPiYgB7k5kuUMRyUDVGDMlxWO3ALf0Y3eUCk1TUxPFQyqpmXAEABUjJlFeW0/V6GlAPDg++7cHmH3M0UyfPj3rg6UKV1tbG7NmzaKzs5MLfvsRqibV9Dz27u+X8MLLL/YMVJVSg0eGf4AOF5EFjt/nZTig/BDwcq+0/0nGmC0iMoL4SceVxpgXM+lkOppDUiokDQ0NrHt/A8MmH9uzTkR6BqkAwyYeTeGQKha9fD/vLFzEpEmTOOboWRQVFUXRZRWxBx98kM7OTi6d/wlKa8sOeqzujIk8cdef2L9/P5WVlRH1UCnV33xco9oQUJHTlfRK+xtjtiT+3SEiDxOfSlQHqqkkS8MlezzVsk36wWvavru915Sgk9eUl5NN6sfrfai9TuDtVtHr9Vo7tzSjTdrTa+rS2Te3ZWf73qm9PXv2sHTpUiad9AlqJhyZ4lVB5cjJzPzwP9PWtIstCx/nwYceYfLkSRx15BFUVFS4pq38pPVsKsa9Tl6f6XcDvH9ObZ5rk9YLIl7YxIoXXnqB1197g2/c+PWDZgvo6urq6ds/3vgVnn78Kc78xcV9BqkAZSPLGT51JKtWreLoo4/WeJFGmPEiyFiRbtnt+NlcdmHzmmxovIheVDMRiEgVcBrwKce6MiDPGNOYWD4X+H5YfciZgapS2aSlpYXismqGTzrGqr2IUFIxnEmnfJrm3ZvYvPBx3n3gAa695pqeNkuXLuPZZ5/hE5/4BHV1dSm2prLNnXfeySt/f4XaYcP51Kc+xZIlS/jQhz7E5OlTuO4zn6O0tJTn//Ycl/71kxRXlSTdxt7Vu2nZ3ZSV/4kqpQYeEfkdcDrxSwQ2Ad8FCgGMMXclmn0YeMoY0+x46kjg4UQsKgDuN8Y8EVY/daCqVAiKioooLOl7VsxGWU0dezevBA7+y37V6tUAvPLa63zsozpQHUjG1MXvW/LC31/gsssu4/9940YOv/YYSkcM4ecP3w2dhqO/caLrIHXfuj38be6faW9tp6qqqj+7rpTKAiFV/V9l0eZe4tNYOdetBVKnCgOUMwPVZOkIP9WJNsvOtFS69lFV+rsJo4rXmepym7TbmfayeY1u93F2btNPetMmxeaWznO+lt6pwPb2dlr27cAY4/kM2P6t7wFw1VVXISKICG+88SYb31/P0PFHsG/vltDOqvXnd8ZPdbLXKmebz0imr93LZUIvPPcCRxxxBOPPmMRRnz6C/OICJl8yPW3f1j+xmrGjx3Lm2Wdy2523c9nFl3LKSad4en1+aLyI6+5/kLEi2XIYVf9hGMzxor9lY5/6S84MVJXKFrFYjPc3bGRo3YyMnr9/W3yg+tDDjzJkSBl7du/seWzPhsUcfYzd5QQquxxy1UxqDq2luKqEUXO83Rl64gVTeempp/jVvF8BcPwxx4XRRaVUltKBqlIqMKveW015bT2TT/l0Rs+vm3UhdbMu5MD+Blr3baNo+xq2r4gXU1511VWMHz8+yO6qgP3mN7/h61//Ovfddx9nn302xhi2bNpMxUlVTDhnsvV2Ygdi5BfnIyJU1ldzwUNX0NXRxZNXPMgJx58Q4itQSmUTH1X/OSEnBqrGmLSpfyebykMvaX2bbfpJ/Xu9PMBrmsMtfed1YvygJup265vzPfFzT2e37dtM3O2cNipZaq+1tZXt27Zy1Ee/66lPyZRUDgcR1r70G0aPGcOVH/84xcXFaV+HTYV7UMteq3vdUrM26Tvn9t0+j34qdzONF71jxf0PxGdxufrqqznzvDN5+823ORBr49RrzrPqR3tTO69/53k2v/Y+406ZyHHfP52CkgJiLR08ePZ9VFRVUF9fT0dHR9K+a7yICzNeBBErei/bVP3bHG8/76eTxovsMpgHqsnnsVBKZWT79u3UjJ9JUan/eS7bW/ez5vlfcvTRR3Plxz8eQO9Uf/j6V2+kZtxwzvv15bQc38Xka2dw3v2XU3NobdrndrbFeOazjxJbf4D58+ez/Y3NvP/UakyX4bm5f6Vm2DB+8+vfZP1/qkqpYHWfVfXykyty4oyqUtniwIEDHGg74Hs7HQeaWPXUz6gfN4qTTtQ070By2mmnUZpfTFd7J/XnT/X0XMnPY/+mvZxx6aV87Rtfo72tnXFnTmLFrxexe20DL774ous8mUqp3JVLA0+vcmag2p0W8DqBt5+0vpdlm8sNnPxUqbrxWjFqk76zqdB1u4zCjTOl6dyv23b8pDdt7q/tfF2pJvfu7Oxk+/btjJ99adI+eLFj1SuUFMAZp59mVfXrh9t75dbGz/fBTyrPz/26g5rw3/a7PffaufzPXfdw0n+f4+nsZ15BHif9+9k8+q1HqT9vCh/5+Wdo33eANQ+s4Le//S0ikrTPGi/67jfMeOE3ViRbdpsZwOtlAF5vBODVYI4XURrMA1X901ypgDQ0NAAw8pCTfG+rdc9mhg8flpUBU6V39aevprytlNV/WOb5uePOmMiVr17H8d87g73v7eKF6x/nus9+jsmT7QuxlFK5I5O0fy4NbHWgqlRAuu+/Hmtv9b2t/dtWM336Ib63o6JRUFDAf//nf/Pe/csz3sb2t7bw2k1/45tfvYkrP3ZlgL1TSg00g3mgmhOpf2NM0hSEnwmIg5oBoPu5Nvt0S9+5PdeGTaWn18ml3apTnX1zm8zbaz9tKndtjo9NmsztddlW9VZUVFBQWETDmjcZPeOMpP2w0RXrINZ+gJqaGvLy8lzv7+1nYm+bVKjz2Lu9z27tbe6zbfNZtvkOB3VMMo0Xbsfg2s9dy/6d+zz1AaCzo5OV9y1i/UOr+Mm//4RZs2YlvcxA40Xf5f6KF35jRbr2XmcAcDveGi/6Lkd1cwS/cmng6VVODFSVyhbDhg1j24oXfA1Ud61/m9IhZQwZMiTAnqn+NqQ0/v51xbrIKzg4eWW6Emc9ugz5hfm0N7Wz6fl1tDe2sekv65g4agL3//p+Ro0addA1mEqpwUkHqkqpQAiGwuIyX9voONCM6eo66GyLGngeffRRpk+fzgOn/IoLfvdRmjbuZ+fibex4cyu73/3gbmNDqobQsq8FgEs/chnXfPWTHH/88Vl/hkcp1X90oJoDgkj9e03P26xPlrKz2b/XNk5BTeDtNc3ktZ9OsVgs6fbd7gHuJ73pltazSb25pe2Kioro7Oxk27ZtTDrpE57601tZzVg2HWjl6Wf+xumnnUJJSUnS/gY1gbefz6Cf+3XbpHht0rdRpf6717u97oKCAm666SZ+/OMf8/hVfyIvL48jjjqCS0+6mCtuv4KioiJqampYu3YtlZWVVFdX93zWYrFYz/Y1XvQVdbzwGyuSrfeT7u/PCf81XvS/XLvm1KucGagqFbXt27cDMHzSMb62UzXmEGZ99Huseu5/eP2NNzn7rDOD6J6KwHXXXcd1113Hli1bqKyspLS09KBBVmdnJ9OmTQPQFL9SytVgHqhq1b9SATDGsGPHjsC2V1hawZDKWt5+awH79+8PbLsqGiNHjqS0tDTqbiilBiit+h/ggqr695OiSFeNazO7gJ9UmI2g7tftlsqzOd5e++M2WbjX+4G77ddmYm23dJ5zgu6dO3eyceNGDj3vS576kkrd7A+zc91CNm/eTHV1dZ++B5WqCuq74VbR65a+85pydLYfCKn/TGNFqu1ovEi93/6KF35ihVu6322913S/zWUAfmi8UP0tJwaqSmWLihETA9vWng2LKa+s1onelVJqkMulM6ReRZr6F5Gvi4gRkeGOdTeLyGoReVdEzouyf0rZGjp0KCJCV2csfeMUjDG0Ne5i6/Ln2bjgUa74yIcPKqZSSik1+GjqPwIiMg44B9jgWHcYcCUwAxgDPCMi04wxKcv9nKl/tzcnqMsAnIUQbim8ZMteU3lulalur8/tvsY29zIOI33np7LW7Vi6VXQ62aQN3dbbVPK6peq2bN2KMYb25j2UVNYm3Vdv7S37WPjg96kaPZWi0nI6WvbSvGcb+Xl51NfXc9VVV1JTU2OVnvKT5nby8zl1q8p167/zu+Q2FZdzv27vv5/jE0S8CDpW9G4fRrxwo/EiLt1x8BMrvE7y7+cyAK/fB7f/RzReRCvXBp5eRZn6vx34J+BRx7pLgd8bY9qAdSKyGpgDvBpB/5SyYozhrQULGHf0xdaDVID8wviZ0gO7N3D0iSdSXV3P0KFDKSsrO2gKG6WUUoObDlT7mYhcAmw2xizq9ZfLWOA1x++bEuuSbWMuMBdg7NikTZTqF42NjQCMmHZi2ramqwtJ/KV/oHEnSB5z5sxh+vTpBxVcKKWUUt10oBoCEXkGGJXkoW8D3wLOTfa0JOuSvjvGmHnAPIAjjjjCONMCLu3TLttUJNqk8JK18VoNHFUVr1vVrNdqXZs+u6XdbY6rn+Pj9Tikq/BtaWkhv6CY/MLilPtt2bOFpX/5T5A8xs26kAP7doAxDB06lPz8fKsZCKKq4rW5L7efyl3n99ftXuX9Wbmbaeo/iFjRe73Gi759iDpeZBorei8HdRmAxou+7XOt0l8HqiEwxpydbL2IzAQmAt1nU+uAt0VkDvEzqOMczeuALWH1UamglNUk+5vsYF2dMfLy8zluznE07FpDdVkRZ378Y9TW2l8uoJRSanAxxni+jjuX9Hvq3xizBBjR/buIrAdmG2MaRGQ+cL+I3Ea8mGoq8EZ/91EpLzZu2sT+nZvStisur6Grs5NXX30FgLlz5+q1qEoppdLSM6pZwhizTEQeAJYDMeCL6Sr+E89Lex9gr6l/m2W3qt50k4J7nWnAz8TYNve1d0vT+LmPtxu3/jhfr9d7gwdx7+7e622Wu9N2be0xhk9Mf9vUguKyg34XEc9VvG5pLrf2TkHNfOG23i3F51a5a9M3r5W7UaT+g44VQa7XeBEXdLzINFakWu+1ul/jRd/lqOJFfxjMA9XIb6FqjKk3xjQ4fr/FGDPZGHOIMebxKPumlI2y0mIa1i5I2273+4sAmDVrFtdcc43Vf9RKKaVUsnlS0/2kIyL3iMgOEVnq8vjpIrJPRBYmfv7F8dj5ifnuV4vINwN8qX1EPlBVaqCbMGECAHs2LEnZbs3f/4+ammGccMIJOom/Ukopa2EMVIF7gfPTtPm7MeaoxM/3AUQkH/gZcAFwGHBVYh78UGRV6t+PZBW4Tl7TZH7u2Z1s2U/Vr5+0lU3KwyZt5Kda18mm0tNrOsnJJr3ptt7mkoBk6bmdO+MJgS1LnmLo+Jl99rHjvdfYuuQZACorK6wqem365SeFZXOM3aq4bSbqtuFWuWvznfSaxvQz4X+6eBF0rEi1rPGi7/oo4kWmsaL3+qCWNV6k7ltQ8SIqHgaeXrf7oojUZ/DUOcBqY8xaABH5PfF58JcH2L0eekZVKZ+mTp0CwNDxR/V5bNuKF1n/2h/JN23MmDGD889P98erUkopdbCQzqjaOEFEFonI4yIyI7FuLLDR0cZ1zvsg5MwZVaWi0traCsCoQ0/t89i2Zc9y7rnnMnny5P7ulhqktm3bxsiRI7PurJBSKnMZDjyHi4izgGKeic9Bb+ttYIIxpklELgQeIT4bk/Wc90HIiYGqMYaOjo6k693aJ1v2WvXvZbJu5368Tkod1PxpXtN6bhW0Xit33fZlUwHq9ZING15TmukqbTdv3hz/vaDvnaVGHXY6L7z4NzZu2sL4cXVMnTrFqgo6jPSdG6+fQZvJvG36ZlO5a1PZ7Pb98Hp8Mo0XQceK3vvxEi8WL17Mvb++l+98+zuUl5cnfT22NF6498fPzCFeLwnQeJF6X1HFi/6W4UC1wRgz28c+9zuWHxORn4vIcPp5zntN/SvlQ0dHB8uWLWPEtBOSPj7qsDOYfMZ1rFi+lCeffKKfe6cGm+nTp3P93Ot9D1KVUtnDmPiE/15//BKRUZIYwSduypQH7ALeBKaKyEQRKQKuBOb73qGLnDijqlRUdu3aBcC4oz+U9HFjDC179OZqqn8UFRUxceLEqLuhlBoAROR3wOnELxHYBHwXKAQwxtwFfBS4QURiQCtwpYmf2o2JyJeAJ4F84B5jzLKw+pkTA1VjPpjw32v6ziZVFERqz6aK1+uk3TYpDJtlm0my3aos3bhV6LqlYPwcByc/Vf826U3nclNTEy+/8hpD62aQV9D3DlPGdLHwT/9KaXEhF110EVOmxIuu/KT4vL63Tn4+XzapPD+pdps++Elv+qn6T9e3oGNF7/UaL/ruK+p44TVWBHUZgMYL+z6EES+iFFLV/1VpHr8TuNPlsceAxwLvVBI5MVBVKgqPzp/PgdZWZl90U5/Atu7VB9i3dRUdB5q45urP67ypSimlMhZUrcpApANVpTJUVlZOQfkImna+T+WoKT3rTVcnO1e/DsBRs46hsLBvkZVSSillI6x5VAeKnBioeq36d2vj9R7cXu7x7Xavb69pPT/cUlRu6RK3VJ7Xfbml9dz64LVa108Vr9t6m3Teaaeewl/++hgrn/4FR33kXygaUgXAxnfi2ZDPf/6DM6nO4+f1/t5eJ6v2ymv6zrlft3txO1+j2768pu+cbWxee1BV/25tklX9BxErei9rvEjdhyjihddY4XbphE36XuNF6uWo4kV/G8xnVLXqX6kMjRgxgk9/6pMA7HzvNfZtXcX+bavZtvx5gKR/PCmllFJeRTjhf+R0oKqUD4WFhZx99tl07FjKvuWPsXvpn3seu+eeeyLsmVJKqVwxmAeqOZP6T3Y/4Uwn8O7Na7VpsrSd1/RdGJNV26TvvN6v2y1lY1Nl6XUC72xI/SdbP336dA477LCD1j/22GOsWrWKWCxGcXHxQcfJTzVwUGk9r5W7Nmk9tzS287XbfKZsPqdhVO4GES+CiBW9lzVepO5Dtqb+/cwGoPEi9X6zIV70J2OMpv6VUsHqLqD6458epLm5OeLeKKWUGsj0jKpSKlBVVVXUTj6W9gON/PJX9zBz5kzG1Y1l2LBhbN26lWnTplFWVhZ1N5VSSg0Ag/mMas4MVLtTB17/inBLY9isdy67pee6++UnZec1teV1kmebib1tjqszZeOW7vFzj26noCfwTrXeayotLy+PyspKdr7yBHn5BXR1xlj4ztssfOdtioqKaG9vZ/WaNVz58Y8HNlG3TWrLT9Ws2325bdJ3zvVunymblJ3XKl4nrxP+u0nW56BjRe/1Gi/isjVe+I0V3Wwm/9d4kboPYceLqOTaGVKvcmagqlQ2mTFjBocffjidnZ00NzfT1NTEkCFDGDp0KAUFBQcFV6WUUioVHagqpQKXl5dHXl4e1dXV1NTUHLReKaWUsqWp/wHOGOM6kbCzTbr1NsvOFJVNOq972ebe3W7bDuovKZtKVbd92Qyu3Lbvdpy8TiLuFMYxcVtvk8JLlobzWulrs2yTsvPK5rPu9dIJ5/fRZhJutxSfW9ow7MpdL/Ei6FjRu43Gi77bzKZ4EUSs6N1e40Xf7WdzvOgPg/mMamSndkTkyyLyrogsE5GfONbfLCKrE4+dF1X/lFJKKaWi1j09ldefXBHJGVUROQO4FDjCGNMmIiMS6w8DrgRmAGOAZ0RkmjGm7ySpSimllFKDwGA+oxpV6v8G4EfGmDYAY8yOxPpLgd8n1q8TkdXAHODVVBuzSf2nem43rxW7bs9NlqqzmSS5Pye3d0uvuKXs3CpP3dI9NtWXXquWnfoz9W+TSkvWxmZybpvUn1sbm9fhJojJ7Xuvd/sc2aTp3NYHNWF5UNJV/QcRK1K113jRdznqeBFErOi9XuNF3zYDMV4EaTAPVKNK/U8DThGR10XkBRE5NrF+LLDR0W5TYp1SSimllBpkQjujKiLPAKOSPPTtxH6HAscDxwIPiMgkINmfPUn/jBCRucBcgJEjRwbRZaWUUkqprNJ9jepgFdpA1RhztttjInID8JCJn8t+Q0S6gOHEz6COczStA7a4bH8eMA9g2rRppjulZHN63CZl5jXdlq7C1+tk1TbpLK+Cmszb7Qtjk8qzmcDb5v7hbrxO4O3WJohLAoKq0PU6UbdXXifqdq53q9Z1275bGtOtittr2tiNzaTmbtJ9BoOOFam2o/Gib5so4kXQsaJ3G40XcdkWL6Kkqf/+9whwJoCITAOKgAZgPnCliBSLyERgKvBGRH1USimllIqcVv33v3uAe0RkKdAOXJ04u7pMRB4AlgMx4Ita8a+UUkqpwWwwn1GNZKBqjGkHPuXy2C3ALR63R0dHR9o26dZ7rSr1ktrzc+9uryk+r2kgm8m5ve7LJmXnVkkcRjWzG5u0XabH02ulr9eKXq8pPpuUndcqXq/pOz/Vul7TrX54iRdBx4re29F40bdNNsWLIGJFkOs1XqRedsrGFH8yeo2qUkoppZTKWnpGVSmllFJKZSUdqOYALxP+u51C95pO8pKSs2mbDffudmNzb22bqk+nqCbtduOnejdZe68pqWy+d7dNFbeTTRW32+fO7bPj9nptPr9+pIsXYcwiovGir2yKF0HHilTtNV70XXbKtngRBk39K6WUUkqprKVnVJVSSimlVFYK44yqiNwDXAzsMMYcnuTxTwI3JX5tAm4wxixKPLYeaAQ6gZgxZnbgHUzIiYGqMSZt6t/mrxGv6SQvlbZeq36DSme5pVT8pEVsjoFNZajbfv1M1O2V18mfvcwS4KcC2Gv6zmtaz+vnyyat5/b+u723Nulem5Sz28TefniJF1HPKGLbBxsaL1JL9/nKNFZksqzxoq+o4kV/COmM6r3AncCvXR5fB5xmjNkjIhcQv8nScY7HzzDGNITRMaecGKgqpZRSSuUiY0woA1VjzIsiUp/i8Vccv75G/G6h/U4HqkoppZRSWSwLiqmuBR53/G6Ap0TEAHeb+G3tQ5ETA9WgUv826Y1MZwzwmuKz2adXXtMcbukkr6lLp6DSlVGl/t3We5kIPKj7dQeVtnKr0Hbjljqz+bzYpJDd2nitHvYj089gELGi93M1Xti3cRNF6t9tvdc2Gi/6yrZ40R8y/AwPF5EFjt/nZTKgFJEziA9UT3asPskYs0VERgBPi8hKY8yLmXQynZwYqCqllFJK5aoMB6oNfoucROQI4JfABcaYXY7+bEn8u0NEHgbmAKEMVAfehGJKKaWUUoNE9zWqXn/8EpHxwEPAp40xqxzry0SkonsZOBdY6nuHLnLijKox5qD0Qjo2KQGvaSab9FaydUGl/pzCTme4pZO8pj1t1ntt40dQlwF081OJG0aKz0slO3ifGNtris+mbzbvSVQT/ncLM1b0Xq/xIvV6r20yFXSsAI0XybYzEONFWML4PIvI74DTiV8isAn4LlCY2N9dwL8Aw4CfJ45t9zRUI4GHE+sKgPuNMU8E3sGEnBioKqWUUkope8aYq9I8/jngc0nWrwWODKtfvelAVSmllFIqi2VB1X9kcmagmq7q34bXVIfNc7vX+5lRIGxe01g2aZeg0knZwCZVlK7q3229nzZhs0nr2VTo2kywnekx7t3P/pTus+wnna3xIvVznbIpXvj5HGu8iMvVeOFX2Je+ZbOcGagqpZRSSuWasCb8HyhyYqDqtZiq93ODaJ9uO26Pey2m8COMv7S9ngmxeS1eb48YBq/HyssZVZs2Xs+i+BFG0YTXYgev73NUZ4289DOMQiCNF31FHS+CiBW229F4Yd8Hm+0MJDpQVUoppZRSWUkHqkoppZRSKivpQHWAM8YEcoF0UJcBOIU5B6MNm1vT2VzgbrNNr/t147V92IJK7XnZnp+CCxtePwtei3i8FlPYfK69tg9bEJcBOIU9v7MNjRf+BB0r/LbReGHfPpsFNcYZqHJioKqUUkoplasG+mDbj0hu0SAiR4nIayKyUEQWiMgcx2M3i8hqEXlXRM6Lon9KKaWUUtki2S1S0/3kiqjOqP4E+FdjzOMicmHi99NF5DDgSmAGMAZ4RkSmGWPSlvSHeVrc67aDTgmGUcUbVLrMa1rPKapbOgbFSz+jSgna8JNec2MzF6IbmzRgtgozVti213iRfcL4/nttr/FiYMu11+NFVDe9NUBlYrkK2JJYvhT4vTGmzRizDlgNzEnyfKWUUkqpQUHPqPa/rwJPish/EB8sn5hYPxZ4zdFuU2KdUkoppdSgk2sDT69CG6iKyDPAqCQPfRs4C/iaMeZBEfkY8CvgbCDZOf+k746IzAXmAtTU1HhKufXnG57pvsLuo58UnJOflKCfytMw9FcKMajUaRhs3s+wK2693ryjP49Jtu5H40XqfQWtPy830HiRWrbFi7Bo1X8IjDFnuz0mIr8GvpL49Y/ALxPLm4BxjqZ1fHBZQO/tzwPmAdTX1w/ePzWUUkopldMG8xnVqK5R3QKcllg+E3gvsTwfuFJEikVkIjAVeCOC/imllFJKZQW9RrX/XQf8VEQKgAMkUvjGmGUi8gCwHIgBX7Sp+B9ok+Fm2wfIa8rO63ODMlBTNray7fV5Tdl5fW5Qsu37FLRse30aL7JDtr0+jRfhybWBp1eRDFSNMS8Bx7g8dgtwS//2SCmllFJKZRu9M5VSSimlVBbTM6o5YDC/if0h29JMKjvo904lo/FCJaPxInOD+djlzEBVKaWUUioXDaQ6nKDpQFUppZRSKktpMZVSSimllMpaOlBVSimllFJZSQeqSimllFIqK+lAVSmllFJKZSUdqCqllFJKqawz0O6+GTQdqCqllFJKZbHBfEZVcuHFi8hO4P2o++HRcKAh6k4MAHqc7OmxsqPHyZ4eKzt6nOwMxOM0wRhTG2UHRMRkchMNY8xbxpjZIXSpX+XEGdWoP0SZEJEFufABCpseJ3t6rOzocbKnx8qOHic7epwylwsnFTOVEwNVpZRSSqkc9STxs9FeDbSz10npQFUppZRSKksZY86Pug9Ryou6A4PYvKg7MEDocbKnx8qOHid7eqzs6HGyo8dJeZYTxVRKKaWUUir36BlVpZRSSimVlXSgGjIRuUJElolIl4jMdqyvF5FWEVmY+LnL8dgxIrJERFaLyH9LJvNSDEBuxyrx2M2J4/GuiJznWD8oj1U3EfmeiGx2fI4udDyW9JgNZiJyfuJ4rBaRb0bdn2wiIusT36WFIrIgsa5GRJ4WkfcS/w6Nup9REJF7RGSHiCx1rHM9NoP1u+dynDRGKV90oBq+pcDlwItJHltjjDkq8XO9Y/0vgLnA1MTPYLmQOumxEpHDgCuBGcSPxc9FJD/x8GA9Vk63Oz5Hj0HaYzYoJV7/z4ALgMOAqxLHSX3gjMTnqPsPxW8CzxpjpgLPJn4fjO6lb2xJemwG+XfvXpLHYI1RKmM6UA2ZMWaFMeZd2/YiMhqoNMa8auIXEP8auCys/mWTFMfqUuD3xpg2Y8w6YDUwZzAfKwtJj1nEfYraHGC1MWatMaYd+D3x46TcXQrcl1i+j0H6/TLGvAjs7rXa7dgM2u+ey3FyM2iPk/JGB6rRmigi74jICyJySmLdWGCTo82mxLrBbCyw0fF79zHRYxX3JRFZnEi7dacf3Y7ZYKbHJDUDPCUib4nI3MS6kcaYrQCJf0dE1rvs43Zs9HPWl8YolTGdRzUAIvIMMCrJQ982xjzq8rStwHhjzC4ROQZ4RERmAMmuscyZqRkyPFZuxySnj1W3VMeM+KUP/0b8df8b8J/ANQySY+ORHpPUTjLGbBGREcDTIrIy6g4NUPo5O5jGKOWLDlQDYIw5O4PntAFtieW3RGQNMI34X5V1jqZ1wJYg+pkNMjlWxI/JOMfv3cckp49VN9tjJiL/A/wl8avbMRvM9JikYIzZkvh3h4g8TDwNu11ERhtjtiYutdkRaSezi9ux0c+ZgzFme/eyxiiVCU39R0REarsvHBeRScQLgdYmUkiNInJ8ooL9HwC3M42DxXzgShEpFpGJxI/VG3qseq5p7vZh4gVp4HLM+rt/WeZNYKqITBSRIuKFHPMj7lNWEJEyEanoXgbOJf5Zmg9cnWh2NYPs+5WG27HR756Dxijll55RDZmIfBi4A6gF/ioiC40x5wGnAt8XkRjQCVxvjOm+CP0G4tWTpcDjiZ+c53asjDHLROQBYDkQA75ojOlMPG1QHiuHn4jIUcRTZuuBzwOkOWaDkjEmJiJfIn7f7HzgHmPMsoi7lS1GAg8nZncrAO43xjwhIm8CD4jItcAG4IoI+xgZEfkdcDowXEQ2Ad8FfkSSYzOYv3sux+l0jVHKD70zlVJKKaWUykqa+ldKKaWUUllJB6pKKaWUUior6UBVKaWUUkplJR2oKqWUUkqprKQDVaWUUkoplZV0oKqUykki0hR1H5RSSvmjA1WllFJKKZWVdKCqlMppEneriCwVkSUi8vHE+tNF5HkR+ZOIrBSR3ybucKaUUipL6J2plFK57nLgKOBIYDjwpoi8mHhsFjCD+D3GXwZOAl6KoI9KKaWS0DOqSqlcdzLwO2NMpzFmO/ACcGzisTeMMZuMMV3AQqA+mi4qpZRKRgeqSqlclyqd3+ZY7kSzTEoplVV0oKqUynUvAh8XkXwRqQVOBd6IuE9KKaUs6NkDpVSuexg4AVgEGOCfjDHbRGR6tN1SSimVjhhjou6DUkoppZRSfWjqXymllFJKZSUdqCqllFJKqaykA1WllFJKKZWVdKCqlFJKKaWykg5UlVJKKaVUVtKBqlJKKaWUyko6UFVKKaWUUllJB6pKKaWUUior/X9xzfvCTbrFXAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(figsize=(12, 5))\n", "ds.field.plot(ax=ax, x=\"lon\", y=\"lat\", cmap=\"Greys_r\")\n", "handles = regs.plot(\n", " column=\"field_avg\", ax=ax, edgecolor=\"k\", vmin=1, vmax=3, cmap=\"viridis\"\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Extract the weight mask from the averager\n", "\n", "The weights are stored in a sparse matrix structure in\n", "`SpatialAverager.weights`. The sparse matrix can be converted to a full\n", "DataArray, but note that this will increase memory usage proportional to the\n", "number of polygons.\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA2MAAAKACAYAAAARlqCRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAABUG0lEQVR4nO39eZxlZXnv/X++PdF0AzJjCyioHRWJAxJxSDwqmqAxgZPHAeOAhsij0TgkMSGJxwy/kxyPyTHqo5F0HIIJDogYiIc4EYkxRmQQZWgIBBlaGhpkBqXpquv3x1qlRVHVvat27bV3VX3evvZr77XWvfZ97Y19d1/7Wve9UlVIkiRJkrq1bNgBSJIkSdJSZDImSZIkSUNgMiZJkiRJQ2AyJkmSJElDYDImSZIkSUNgMiZJkiRJQ2AypiUlyVuTrBlCv+9L8v0knfyZS3JMkkMmbf9pkud10bckDcowxvAkT03ytSRXJLk8yYeTrEnyx0l+Z4ZzvtFljJIWLpMxLTVvBab9izzJ8kF02CZg/x24HnjWDG1WzHO3xwA/Tsaq6p1V9ZV57kOSuvZWOhzDk+wHfAb4vap6DPA44AvArts7r6qeMd+xSFqcTMY0UpK8Osl3k3wnyd+3+x6R5Ox2/9lJHt7u/7skL5507t3t87OTnJPktPZXzFPSeDPwMOCrSb46cU5bNToXeEeSz016v+cnOX0ePtZzgEuADwEvn/T+f5xkQ5IvAR9Psk+SLye5MMnfJLk2yd5t21cm+VaSi9pjyyfF/2ft9/XNJPsleQbwy8BftO0fNfm7SnJNkj9p+7k4yWPb/U9N8o0k326fHzMPn13SErIIx/A3AidX1X8AVOO0qrqpPX5IG+vVbXw9fZb22DuTnJfkkvbvgvQZq6QFyGRMIyPJ44E/BJ5bVU8E3tIe+gDw8ap6AnAK8P4e3u7JNL+gHgI8EnhmVb0fuAF4TlU9p223Frikqo4A/hR4XJJ92mOvBT42TZx/1SY5Ux8nzhDLy4FPAp8DXpRk5aRjTwGOrqpfBf4I+JeqOqxtO/EPlscBL2s/w5OAMeAVk+L/Zvt9fQ14XVV9AzgTeHtVPamq/muamG5p+/kQMHGZzeXAs6rqycA7gT+f4fNI0oMs0jH8UOCC7cT5WOAXgKcCfzRlfJ/xs7T7P1BVP1NVhwI7Ay/aTj+SFqn5vjRK6sdzgdOq6haAqrq13f904Ffa138PvLuH9/pWVW0CSHIRcBDw9WnajQGfbfur9pfcVyb5WNvvq6eeUFVv6/HzkGQV8ELgbVV1V/vr7c8D/7dtcmZV/bB9/bM0lzNSVV9Iclu7/0iapO289ofTnYEt7bGtwOfb1xcAz+8xtIlfiy/gJ9/tQ4CTk6wHCpjuHxWSNJNFN4b34P9W1X3AfUm2APsBm6a0memzPCfJ79JcdrkncCnwT/MYm6QFwGRMoyQ0ScCOTLTZRlvdbS/vWDWpzX2TXo8x8//Xf1RVY5O2P0bzl+GPgM9U1bYHBZn8Fc2lh1N9qqreNWXfUTRJzsVtIrUGuJefJGP3TH7rGWIMzWUyvz/NsfurauL72N7nnGri+5l8zv8P+GpV/fckBwHn9PhekgSLcwy/lObHsDNm6L+XOB/UJslq4K+Bw6vq+iR/DKyeoQ9Ji5iXKWqUnA28NMleAEn2bPd/Azi2ff0KfvLr6DU0f0kCHE1vlZy72M7E66q6geYymHcAfzdDm7e1l/9NfUz9SxyaSxR/vaoOqqqDgIOBn8/0q4F9HXgpQJKfB/Zo958NvDjJvu2xPZM8op/POYOHAN9vX79mludK0mIcwz8AHJfkiIkdaebwPrSHWLdnIvG6JckuwIu311jS4mUyppFRVZcCfwb8a5LvAO9pD70ZeG2S7wKv4ifzEP4W+G9JvgUcwQOrTDPZAPzzxOTvGZwCXF9Vl83hY/xYm3D9Aj+pglFV99D8Q+SXpjnlT2gStQuBFwCbgbvaON4BfKn9Dr4MrNtB958C3t4uxvGoHkN+N/C/kvw7MJCVJSUtXottDAdoF+o4FvjLNEvbbwR+Drizz/e9nebzXwz8I3Bef5FKWqjykyucJAEk+QDw7ar6SMf97gSMVdW2JE8HPtQu2CFJ6tGwxnBJmguTMWmSJBfQ/Dr7/HZSdpd9rwdOpalYbwV+o6r8tVSSejTMMVyS5sJkTJIkSZKGwDljkiRJkjQEJmOSJEmSNAQL4j5jq7JTrWbtsMOQpEXrLm67par26bW947IkDdZsx+WFKMkhv/QLay/9py/es2LKPQOXjAWRjK1mLUfkyGGHIUmL1lfqtGtn095xWZIGa7bj8kL0smN2uXTzjWN84qSHbqO5cfyS42WKkiRJkjqV5JBbbh3jk3/zUP7PX99GkiV5j9MFURmTJEmStHi87JhdLn39q3fnofuu4OeetjO//Rt7LMnqmJUxSZIkSZ2ZqIo96+k7A/D2N+6xZKtjVsYkSZIkdWaiKjZhKVfHrIxJkiRJ6sTUqtiEpVodszImSZIkqRNTq2ITlmp1zMqYJEmSpIGbqSo2YSlWx6yMSZIkSRq4mapiE5ZidczKmCRJkqSB2lFVbMJSq45ZGZMkSZI0UDuqik1YatUxK2OSJEmSBqbXqtiEpVQdszImSZIkaWB6rYpNWErVMStjkiRJkgZitlWxCUulOmZlTJIkSdJAzLYqNmGpVMcGWhlL8rYklya5JMknk6xOsmeSLye5sn3eY5AxSJIkSereXKtiE5ZCdWxgyViS/YE3A4dX1aHAcuBY4ETg7KpaD5zdbkuSJElaRF52zC6XvuOte835/Inq2CdOeui2eQxrpAx6ztgKYOckK4A1wA3A0cDJ7fGTgWMGHIMkSZKkDvVbFZuw2KtjA0vGqur7wF8C1wGbgTuq6kvAflW1uW2zGdh3uvOTnJDk/CTn3899gwpTktQjx2VJUq/6rYpNWOzVsUFeprgHTRXsYOBhwNokr+z1/KraUFWHV9XhK9lpUGFKknrkuCxJ6sV8VcUmLObq2CAvU3we8L2qurmq7gdOB54B3JRkHUD7vGWAMUiSJEnq0HxVxSYs5urYIJOx64CnJVmTJMCRwEbgTOC4ts1xwBkDjEGSJElSR+a7KjZhsVbHBnafsao6N8lpwIXANuDbwAZgF+DUJMfTJGwvGVQMkiRJkroz1/uK7chive/YQFdTrKo/qqrHVtWhVfWqqrqvqn5QVUdW1fr2+dZBxiBJkiRp8AZVFZuwGKtjg17aXpIkSdISMN9zxaaa7dyxJEcluSLJVUkedG/jNN7fHv9uksMmHbsmycVJLkpy/nx+jslMxiRJkiT1ZdBVsQm9Vsfa4x8EXgAcArw8ySFTmr0AWN8+TgA+NOX4c6rqSVV1+PxE/2AmY5IkSZL69T8GWRWbMFEdY8frTjwVuKqqrq6qrcCnaG67NdnRwMer8U1g94lV37tiMiZJkiSpX48ddFVswot/aReA1yc5f9LjhCnN9geun7S9qd3Xa5sCvpTkgmnee94MbDVFSZIkSUvGeFcdLV8egI1V9YbtNJtuxcWaRZtnVtUNSfYFvpzk8qr62uyj3T4rY5IkSZIWm03AgZO2DwBu6LVNVU08bwE+R3PZ47wzGZMkSZLUt/EO/9eD84D1SQ5Osgo4FjhzSpszgVe3qyo+DbijqjYnWZtkV4Aka4GfBy6Zv2/qJ7xMUZIkSdKiUlXbkrwJ+CKwHPhoVV2a5PXt8ZOAs4AXAlcB9wKvbU/fD/hcEmjypU9U1RcGEafJmCRJkqS+jVU308bGa+rUr+lV1Vk0CdfkfSdNel3AG6c572rgif1F2RsvU5QkSZKkIbAyJkmSJKlv4w9arHBh99MFkzFJkiQtGMsf8hAAsnonAOr++wEYu+2OpkFHl8pJ88FkTJIkSVLfelzlsG/V3S3NBs5kTJIkSSMvK1Y2z7vu0uxY1tyvNzutAmD5eHPp2tjtt3cemzRXJmOSJEmS+jbW4yqH/ep1NcWFwGRMkiRJPVm+664AjN19T7Ojg/lZK/bbt3mxsqmMTVTEpsoua5v2q5p2Y7feDkBtu79t0C4i7pwyjRCTMUmSJM1o2Zo1QLNgxo8Tnj2aRTTGb7m1eb733oH1P377nU3/O69unnfbZYZA2yStXdhj+X77PODw2I03AVBjAwhSmiOTMUmSJEl9625p+8XDZEySJEkPsnz33QHI2jXtjmUParOsrZSlvYRw/J6mQvbjSwNnKcuXN++7Zs2Pl66nrYjN2rZtTUx33tXENGZJTKPHZEySJElS38a86fOsmYxJkiTpQSaWjJ+uIvZjbfVqooq1vF12vm5vbsD842pUu+z8+H0/emAfE8vVT1TE2koba3eee+Dbmj7H77q7ef7hD+f+XtKAmYxJkiRJ6lt3c8asjEmSJGkRmKhOLZuock3cVHnlHP6ZOHEj5j13b56nHm4rZEzcJ2p7VbeZ/LCtrk3MJbtvK/CTVRenVt+kUWYyJkmSJKlv3vR59kzGJEmSlqL2JsjL99272V6xfPB9/viGzdPfuHlaY81C5uO33gZAbW0qYRPzzMbbbWkhMhmTJEmS1Leu7v/lfcYkSZK0sFVbcWpXPly2957DjGZm7byyZbs/BIBq54yN3XHn0EKS5ovJmCRJkqS+eZ+x2TMZkyRJWsLG770XgLrhfqDjOWRzMHb77cMOQZo3c1hPtHdJdk9yWpLLk2xM8vQkeyb5cpIr2+c9BhmDJEmSJI2igSZjwPuAL1TVY4EnAhuBE4Gzq2o9cHa7LUmSpCGqbfc3j3t/SN37w2GH82ArV8DKFSzfdVeW77rrsKPRNMaqu8diMbBkLMluwLOAjwBU1daquh04Gji5bXYycMygYpAkSZKkUTXIOWOPBG4GPpbkicAFwFuA/apqM0BVbU6y73QnJzkBOAFgNWsGGKYkqReOy9LSMLFaYXZe3exYOVpLDGTNzs2Lu+4abiB6kK6WnF9EhbGBXqa4AjgM+FBVPRm4h1lcklhVG6rq8Ko6fCU7DSpGSVKPHJclSZpfg/ypYxOwqarObbdPo0nGbkqyrq2KrQO2DDAGSZIkzcL4fU1lbHzzjQAsW7Wqed79IbC6/SFmfEptYmwMgLr7HgAy0W6iujaflg16yQPN1RhZVP10YWD/b66qG4Hrkzym3XUkcBlwJnBcu+844IxBxSBJkiRJo2rQFwH/JnBKklXA1cBraRLAU5McD1wHvGTAMUiSJGmOxrdubZ633AyZ8jt+zTBLqJ3PtWLffZrt1fNwaXNbjRu7+Qf9v5cGYmrBdFAW05yxgSZjVXURcPg0h44cZL+SJEmSNOpGa3kcSZIkja6ZKmFTZPny5sXKlfPX9R13Ns/b7p+399T8cs7Y7DkDUpIkSZKGwMqYJEmS5lVNrK54ZzN3LA/ZrTlw333N/q1NdStr23sWrlj+wDcYG4cfNW3H2xUaJ1Z51OjqqmI1bmVMkiRJktQPK2OSJEkaiLF2VcVlbSXsQdWttnI2cS+zrGz+aTp2z709z0+TFjKTMUmSJA3UjJcYtgnXxE2iq7000URsYRqvji5T7KifLniZoiRJkiQNgZUxSZIkDdXYHXcMOwTNg+4W8Fg8rIxJkiRJ0hBYGZMkSZLUt7GO6jzji6ietHg+iSRJkiQtIFbGJEmSJPXN1RRnz8qYJEmSJA2BlTFJkiRJfetuNUUrY5IkSZKkPlgZkyRJktS3seqmzjPmnDFJkiRJUj+sjEmSJEnqW1f3/6pFVE9aPJ9EkiRJkhYQkzFJkiRJGgIvU5QkSZLUt+6Wtl88rIxJkiRJ0hBYGZMkSZLUt86Wtl9E9aQdfpIkb0qyRxfBSJIkSdJS0Utl7KHAeUkuBD4KfLGqarBhSZIkSVpIxjuaM1Yd9dOFHVbGquodwHrgI8BrgCuT/HmSRw04NkmSJElatHqaM1ZVleRG4EZgG7AHcFqSL1fV7w4yQEmSJEmjr6u5XF2t2tiFHSZjSd4MHAfcAnwYeHtV3Z9kGXAlYDImSZIkSbPUS2Vsb+BXqurayTurajzJiwYTliRJkqSFpKvVFMc76qcLO0zGquqd2zm2cX7DkSRJkqSlYeD3GUuyHDgf+H5VvSjJnsCngYOAa4CXVtVtg45DkiRJ0uCMdzRnrKt+utDFJ3kLMLmCdiJwdlWtB85utyVJkiRpSRloMpbkAOAXaRb+mHA0cHL7+mTgmEHGIEmSJEmjaNCXKb6XZrXFXSft26+qNgNU1eYk+053YpITgBMAVrNmwGFKknbEcVmStD1j1c2S8+PVSTedGFhlrF1pcUtVXTCX86tqQ1UdXlWHr2SneY5OkjRbjsuSJM2vQVbGngn8cpIXAquB3ZL8A3BTknVtVWwdsGWAMUiSJEnqQFc3fXYBjx5U1e9X1QFVdRBwLPAvVfVK4Eyam0jTPp8xqBgkSZIkaVQNfGn7abwLODXJ8cB1wEuGEIMkSZKkedTVzZjHO5qb1oVOkrGqOgc4p339A+DILvqVJEmSpFE1jMqYJEmSpEXGOWOzt3g+iSRJkiQtIFbGJEmSJPWtq/uMddVPF6yMSZIkSdIQWBmTJEmS1Leu5nIVVsYkSZIkSX0wGZMkSZKkIfAyRc2rZYc9HoC7D9ql2ZGmjLzrV68AYOzW24YSlyRJkgZrrKObPne1hH4XFs8nkSRJkqQFxMqY5sWKffYG4LbH7ApALXvgxMrbXvBYANbesBWAlf9+8QOOj9+/rXlR44MMU5IkSQMy3tHCGi7gIUmSJEnqi5UxzY+H7AY8uCI21X17rgRgp4fuC8D4PrsDsOySq5rt+340oAAlSZI0SJ3NGfOmz5IkSZKkflgZ07y48wn79NRu2+rml4wfPm4dAKu+ciEANTY2mMAkSZLUia5WOezq5tJdWDyfRJIkSZIWECtj6svy3XcHYGx1b3l9qnne6ZxmNcVxK2KSJEmLwnhHc7m66qcLVsYkSZIkaQisjKkvde+97Yu25JXt/1Kx8p7mPmKumihJkrS4dDdnzMqYJEmSJKkPVsbUl+y8c/uit18oxlY27VbstRcA237wg4HEJUmSJI06kzFJkiRJfRvv6KbPXfXTBZMx9aUedUDz3OOlu2M7NQ3HH75fs8PKmCRJkpYokzH1Zdsuq7Z7fPnWZmGPXf/zTgCW3dUs+DF27fWDDUySJEmdGutoYQ0X8JAkSZIk9cXKmPryw313mnb/2hu3ArBqyz0AjF18OQDj3YQlSZKkjnU1l6u86bMkSZIkqR9WxjQ3afL4Wj794eX3bmuO/+c1HQUkSZKkYepqzlhX/XTBypgkSZIkDYGVMc1NNbO/Vt3dPP9o96ZEtmbL/QAsv+waAMbu+1H3sUmSJKlznc0ZW0T1pIF9kiQHJvlqko1JLk3ylnb/nkm+nOTK9nmPQcUgSZIkSaNqkJWxbcBvV9WFSXYFLkjyZeA1wNlV9a4kJwInAr83wDg0AMtWNfcX27rLA/P5ZfeNATB+z72dxyRJkqThGeuoMjbmaoo7VlWbq+rC9vVdwEZgf+Bo4OS22cnAMYOKQZIkSZJGVSdzxpIcBDwZOBfYr6o2Q5OwJdm3ixg0zx73aADGVzzwl4lVW+4GYNu2+zsPSZIkSVpIBp6MJdkF+Czw1qq6M+mtrJjkBOAEgNWsGVyAkqSeOC5LkrZnvKMl58ul7XuTZCVNInZKVZ3e7r4pybr2+Dpgy3TnVtWGqjq8qg5fyU6DDFNzsOyOe1h2xz2kIAUPueJuHnLF3Wy7/Eq2XX7lsMOTNACOy5KkhSTJUUmuSHJVu1bF1ONJ8v72+HeTHDbl+PIk307y+UHFOMjVFAN8BNhYVe+ZdOhM4Lj29XHAGYOKQZIkSVI3xmpZZ48dSbIc+CDwAuAQ4OVJDpnS7AXA+vZxAvChKcffQrPuxcAM8jLFZwKvAi5OclG77w+AdwGnJjkeuA54yQBj0IBsu+ZaAHZtn5fvtRcAY8ub+43V2NhwApMkSZLgqcBVVXU1QJJP0SwkeNmkNkcDH6+qAr6ZZPck69p1LQ4AfhH4M+C3BhXkwJKxqvo6zHhB55GD6leSJElS98Y7WnK+CoDHJTl/0u4NVbVh0vb+wPWTtjcBR0x5q+na7A9sBt4L/C6w67wEPYNOVlPU4ldbtwKw/OBHADB+3abmud0vSZIkzaONVfWG7RyfLjOsXtokeRGwpaouSPLsOcbXE5MxSZIkSX0bG+zagLPtZxNw4KTtA4AbemzzYuCXk7wQWA3sluQfquqVcw56Bt18Y1r0xu66i7G77mLbVVez7aqrGd+61aqYJEmShuU8YH2Sg5OsAo6lWUhwsjOBV7erKj4NuKOqNlfV71fVAVV1UHvevwwiEQMrY5IkSZLmQXdzxnbcT1VtS/Im4IvAcuCjVXVpkte3x08CzgJeCFwF3Au8dmBBz8BkTJIkSdKiU1Vn0SRck/edNOl1AW/cwXucA5wzgPAAkzFJkiRJ82C8oxlQ4zMu2L7wOGdMkiRJkobAZEySJEmShsDLFCVJkiT1bayjBTy8TFGSJEmS1BcrY5IkSZL61tXS9l310wUrY5IkSZI0BFbGJEmSJPVtvLqp8/Ry0+eFwsqYJEmSJA2BlTFJkiRJfRvraJXDrvrpgpUxSZIkSRoCK2OSJEmS+tbVKofOGZMkSZIk9cXKmCRJkqS+dbWa4vgiqictnk8iSZIkSQuIyZgkSZIkDYGXKUqSJEnq23hHS85XJ710w8qYJEmSJA2BlTFJkiRJfRvraMn5rvrpgpUxSZIkSRoCK2OSJEmS+tbV0vbVUT9dWDyfRJIkSZIWECtjkiRJkvo23tFcrq766YKVMUmSJEkaAitjkiRJkvrW1X3GxjvppRtDqYwlOSrJFUmuSnLiMGKQJEmSpGHqvDKWZDnwQeD5wCbgvCRnVtVlXcciSZIkaX50NZernDPWl6cCV1XV1VW1FfgUcPQQ4pAkSZKkoRlGMrY/cP2k7U3tvgdIckKS85Ocfz/3dRacJGl6jsuSJM2vYSzgMV1dsR60o2oDsAEgyc1fqdPuAW4ZcGxzsTejGReMbmyjGheMbmyjGheMbmyjGheMZmyP2FEDx+V5MaqxjWpcMLqxjWpcMLqxjWpcMJqx7XBcHjVd3fR5fBEtCD+MZGwTcOCk7QOAG7Z3QlXtk+T8qjp8oJHNwajGBaMb26jGBaMb26jGBaMb26jGBaMdW68cl+dmVGMb1bhgdGMb1bhgdGMb1bhgtGPT4jaMZOw8YH2Sg4HvA8cCvzqEOCRJkiTNk+4W8Oikm050noxV1bYkbwK+CCwHPlpVl3YdhyRJkiQN01Bu+lxVZwFnzfK0DYOIZR6MalwwurGNalwwurGNalwwurGNalww2rHNxqh+jlGNC0Y3tlGNC0Y3tlGNC0Y3tlGNC0Y7tgWju5s+L56l7VOLqc4nSZIkqXNJLviVf3/DYV30detlN3HO604/qare0EV/gzSUypgkSZKkxcWbPs/e4lkXUpIkSZIWECtjkiRJkvrWVWWsq366YGVMkiRJkobAypgkSZKkvlkZmz0rY5IkSZI0BFbGJEmSJPWts9UUO+mlG1bGJEmSJGkITMYkSZIkaQi8TFGSJElS38bpaAGPjvrpgpUxSZIkSRoCkzEtSknemmRNh/2tSXJKkouTXJLk60l2meN7vSbJwyZtX5Nk7x7PfV+S7ydZNmnfzyW5NMlFSXae5pxvzCVOSeraAh/b/zTJ8+Y7RmmUjFc6eZRL20sj763AtH9hJ1k+gP7eAtxUVT9dVYcCxwP3z/G9XgM8bEeNpmoTsP8OXA88a9KhVwB/WVVPqqofTmq/HKCqnjHHOCWpa29lgY7tVfXOqvrKvEYnacEzGdNQJHl1ku8m+U6Sv2/3PSLJ2e3+s5M8vN3/d0lePOncu9vnZyc5J8lpSS5vf71MkjfTJDNfTfLViXPaXyXPBd6R5HOT3u/5SU7v8yOtA74/sVFVV1TVfe37/1b7i+olSd7a7jsoySWTYvidJH/cfs7DgVOmVLJ+M8mF7a+zj50hhucAlwAfAl7evu+vAy8F3tl+P89O8tUknwAunvhuJsXxu20f30nyrnbf65Kc1+77bJe/SktaWJbK2N6O4ZcnObn9XKdNjI1J3tmOmZck2ZAkUz9vmise/qSHcV1aULqqjHnTZ6kPSR4P/CHw3Kp6Is0vjwAfAD5eVU8ATgHe38PbPZnml9JDgEcCz6yq9wM3AM+pque07dYCl1TVEcCfAo9Lsk977LXAx6aJ86/ahGjq48Rp4vgo8HtJ/iPJ/0yyvn2Pp7TvfwTwNOB1SZ4804epqtOA84FXTKlk3VJVh9EkWr8zw+kvBz4JfA54UZKVVfVh4Ezg7VX1irbdU4E/rKpDpnzeFwDHAEe0/13e3R46vap+pt23keaXYUl6gKU0trceA2xoP9edwG9MfN52zDwU2Bl40QyfsZdxXdIiZzKmYXgucFpV3QJQVbe2+58OfKJ9/ffAz/bwXt+qqk1VNQ5cBBw0Q7sx4LNtf9W+/yuT7N72+89TT6iqt7UJ0dTHu6ZpexHNPxj+AtgTOC/J49rP8Lmquqeq7gZOB36uh8811cSvuxdM9xmTrAJeCPxjVd0JnAv8/Azv9a2q+t40+58HfKyq7m0/08R/l0OT/FuSi2kueXz8HOKXtPgtpbEd4Pqq+vf29T9M+lzPSXJuO2Y+l5nHzO2O69JC5Jyx2XNpew1D6O3m6RNtttH+cNBe7rFqUpv7Jr0eY+b/T/+oqsYmbX8M+CfgR8Bnqmrbg4JM/orm0r+pPjXDX9oTydbpScZpkqOxqe2mfqbW6hnaTZj4nDN9xqOAhwAXt1fErAHuBf7vNG3vmaGPmf67/B1wTFV9J8lrgGfvIFZJS9NSGts/y4M/ayVZDfw1cHhVXZ/kj5l5fN/RuC5pCbAypmE4G3hpkr0AkuzZ7v8GcGz7+hXA19vX1wBPaV8fDazsoY+7gF1nOlhVN9Bc7vIOmmRjujY9/3qa5JlJ9mhfr6K5tOZa4GvAMWlW5FpLs8DGvwE3Afsm2SvJTjzwMpbtxj6DlwO/XlUHVdVBwMHAz89yfteXgF+bNO9h4r/LrsDmJCtp/rtI0nSW0tgO8PAkT29fv7z9XBOJ1y1pVl18MdIS0tmcsWF/0HnkLzHqXFVdmuTPgH9NMgZ8m2YFwTcDH03yduBmmuv9Af4WOCPJt2j+sp+psjPZBuCfk2yeNLdgqlOAfarqsrl/mh97FPCh9tfdZTQVqc9WVSX5O+BbbbsPV9W3oVnmmOZywu8Bl096r78DTkryQ5rLbLarTZ5+Afh/J/ZV1T1Jvg78Uq8foKq+kORJwPlJtgJnAX8A/I82zmtpFv2YbaIoaQlYSmM78AiaObTHJfkb4ErgQ1V1b5K/pRkrrwHOm4cYJC1iaS6xlpaeJB8Avl1VHxl2LJKk+dHF2J7kIODz7SIdkoAkF/y3r/z2YV30deflm/n2mz5xUlW9oYv+BsnKmJakJBfQ/Ar728OORZI0PxzbJS00JmNakqrqKTtuJUlaSLoa26vqGsCqmKS+mYxJkiRJ6ts43Sw5Xx310wVXU5QkSZKkIVgQlbFV2alWs3bYYUjSonUXt91SVfv02t5xWZIGa7bj8igY7+hmzF3104UFkYytZi1H5MhhhyFJi9ZX6rRrd9zqJxyXJWmwZjsua2FaEMmYJEmSpNFWHVWsuuqnC84ZkyRJkqQhsDImSZIkqW9dzeWq6qSbTlgZkyRJkqQhsDImSZIkqW+dzRnzPmOSJEmSpH5YGZMkSZLUN+8zNntWxiRJkiRpCEzGJEmSJGkIvExRkiRJUt+6WnLepe0lSZIkSX2xMiZJkiSpb+MdLTnv0vaSJEmSpL4MNBlL8rYklya5JMknk6xOsmeSLye5sn3eY5AxSJIkSRq8qnT0GPYnnT8DS8aS7A+8GTi8qg4FlgPHAicCZ1fVeuDsdluSJEmSlpRBX6a4Atg5yQpgDXADcDRwcnv8ZOCYAccgSZIkacDGK508yps+71hVfR/4S+A6YDNwR1V9Cdivqja3bTYD+w4qBkmSJEkaVYO8THEPmirYwcDDgLVJXjmL809Icn6S8+/nvkGFKUnqkeOyJGl7qjp6DPuDzqNBXqb4POB7VXVzVd0PnA48A7gpyTqA9nnLdCdX1YaqOryqDl/JTgMMU5LUC8dlSZLm1yDvM3Yd8LQka4AfAkcC5wP3AMcB72qfzxhgDJIkSZI60NVcrsU0Z2xgyVhVnZvkNOBCYBvwbWADsAtwapLjaRK2lwwqBkmSJEkaVYOsjFFVfwT80ZTd99FUySRJkiRpyRpoMiZJkiRpafAyxdkb9H3GJEmSJEnTsDImSZIkqW/jVsZmzcqYJEmSJA2BlTFJkiRJfauO7sbcVT9dsDImSZIkSUNgZUySJElS3zpbTbGTXrphZUySJEmShsDKmCRJkqS+eZ+x2bMyJkmSJElDYGVMkiRJUt+6msvlnDFJkiRJUl9MxiRJkiRpCLxMUZIkSVLfOltYwwU8JEmSJEn9sDImSZIkqX+u4DFrVsYkSZIkaQisjEmSJEnqW2c3fcY5Y5IkSZI0spIcleSKJFclOXGa40ny/vb4d5Mc1u5fneRbSb6T5NIkfzKLPpcl2a3X9iZjkiRJkvpW1d1jR5IsBz4IvAA4BHh5kkOmNHsBsL59nAB8qN1/H/Dcqnoi8CTgqCRP205fn0iyW5K1wGXAFUne3st3ZjImSZIkabF5KnBVVV1dVVuBTwFHT2lzNPDxanwT2D3Junb77rbNyvaxvRTwkKq6EzgGOAt4OPCqXoJ0zpgkSVqw/uuvng7AVS87adrjx1//TAC+9m8/DcDBZ/wIgPzbtzuITlpaOpsz1ls/+wPXT9reBBzRQ5v9gc1tZe0C4NHAB6vq3O30tTLJSppk7ANVdX/S23dhZUySJEnSQvO4JOdPepww5fh02dDU6taMbapqrKqeBBwAPDXJoduJ5W+Aa4C1wNeSPAK4o4fPYGVMkiQtXOt/7wIAHrnnrwFw9fM/+oDjHznw3wH4x6MvAuCDn38pAMkyqPGOopSWiI4qY21KtbGq3rCdVpuAAydtHwDcMNs2VXV7knOAo4BLZujrn6rq/RMbSa4Dfm07sf2YlTFJkiRJi815wPokBydZBRwLnDmlzZnAq9tVFZ8G3FFVm5Psk2R3gCQ7A88DLt9OX5+dvFFVRTNHbYesjEmSpAVrfOtWAN7zjE9vt92Gw58CwIpHN+1r+XJqm5UxabGqqm1J3gR8EVgOfLSqLk3y+vb4STSLbbwQuAq4F3hte/o64OR23tgy4NSq+vzUPpI8Fng88JAkvzLp0G7A6l7iNBmTJEmS1Ldelpyfn456bFZ1Fk3CNXnfSZNeF/DGac77LvDkHrp4DPAiYHfglybtvwt4XS8xmoxJkqQF77e+fiwAx/zCRwC4eOsPAXjD298KwNo7vtk0vKCnOfWStENVdQZwRpKnV9V/zOU9TMYkSZIk9a+rytjouSrJHwAHMSm/qqodLuJhMiZJkha89a89H4Bf4IkP2L+Wbw4jHElLyxnAvwFfAcZmc+JAk7F2FZIPA4fS5Mq/BlwBfJomc7wGeGlV3TbIOCRJkiQN1ojd9LlLa6rq9+Zy4qCXtn8f8IWqeizwRGAjcCJwdlWtB85utyVJkiRpIfp8khfO5cSBVcaS7AY8C3gNQFVtBbYmORp4dtvsZOAcYE6ZpCRJkqQRMWKrKQ5akrtoognwB0nuA+5vt6uqdtvRewzyMsVHAjcDH0vyROAC4C3AflW1mSbCzUn2ne7kJCcAJwCsZs0Aw5Qk9cJxWZKkn6iqXft9j0EmYyuAw4DfrKpzk7yPWVySWFUbgA0Au2XPEcl/JWnpclyWJG3PUp0zluSwaXbfAVxbVdu2d+4gk7FNwKaqOrfdPo0mGbspybq2KrYO2DLAGCRJ0gK24tGPBKBubdb6Sprp7tt+8IOhxSRJU/w1TRHq4nb7p4HvAHsleX1VfWmmEwe2gEdV3Qhcn+Qx7a4jgcuAM4Hj2n3H0SwFKUmSJGkhqw4fo+Ua4MlV9ZSqegrwJOAS4HnAu7d34qDvM/abwClJVgFXA6+lSQBPTXI8cB3wkgHHIEmSFpqJCth/XQPA8rXNPMVtd989rIgkaSaPrapLJzaq6rIkT66qq5PtX1I50GSsqi4CDp/m0JGD7FeSJC1wNf6AzTGTMEmj64okHwI+1W6/DPjPJDvRrK44o0FXxiRJkiQtCV0trDFaC3jQ3MrrN4C30gT3deB3aBKx52zvRJMxSZIkSZqjqvoh8H/ax1TbLeubjEmSJEnq3+gtrDFQSU6tqpcmuZhpPn1VPWFH72EyJkmSJEmz95b2+UVzfYOBLW0vSZIkaQlZYkvbV9Xm9vnadtf69vUW4NZe3sNkTJIkSZLmKMnrgNOAv2l3HQD8Yy/nepmiJEmSpP5VR6scdtVP794IPBU4F6Cqrkyyby8nWhmTJEmSpLm7r6q2TmwkWUGPF1NaGZMkSZLUt+poLldX/czCvyb5A2DnJM+nuefYP/VyopUxSZIkSZq7E4GbgYuB/xc4C3hHLydaGZMkSZLUv64qVqNXGXs2cEpV/e1sTzQZkyRJkqS5ew1wUpIfAP/WPr5eVbft6ESTMUmSJEn9G71VDjtRVa8GSPIw4MXAB4GH0UOuZTImSZIkSXOU5JXAzwE/DdwCfICmOrZDJmOSJEmSNHfvBf4LOAn4alVd0+uJJmOSJEmS+paOFtboqp9eVdXeSR4PPAv4syTrgSuq6lU7Otel7SVJkiRpjpLsBjwceARwEPAQYLyXc62MSZIkSerf0l3a/uuTHh+oqk29nmgyJkmSJElzVFVPmOu5JmOSJEmS+tfV0vaLaAl954xJkiRJ0hBYGZMkSZLUv9GbyzXyTMYkSZIkaY6SvH+a3XcA51fVGds7d4eXKSZ5U5I95hqcJEmSpCWgOnyMltXAk4Ar28cTgD2B45O8d3sn9lIZeyhwXpILgY8CX6yq0fsKJEmSJKl7jwaeW1XbAJJ8CPgS8Hzg4u2duMPKWFW9A1gPfAR4DXBlkj9P8qg+g5YkSZK0WCzdytj+wNpJ22uBh1XVGHDf9k7sac5YVVWSG4EbgW3AHsBpSb5cVb87t5glSZIkacF7N3BRknOAAM8C/jzJWuAr2ztxh8lYkjcDxwG3AB8G3l5V9ydZRnNNpMmYJEmSpCWpqj6S5CzgqTTJ2B9U1Q3t4bdv79xeKmN7A79SVddO6XQ8yYvmErAkSZKkRWZp3/R5GXAzTX716CSPrqqv7eikHSZjVfXO7RzbOKsQJUmSJGkRSfK/gZcBlwLj7e4C+k/G+pVkOXA+8P2qelGSPYFPAwcB1wAvrarbBh2HJEmSpMFJRwtrdNXPLBwDPKaqtrtYx3R2uJriPHgLMLmCdiJwdlWtB85utyVJkiRpIboaWDmXEwdaGUtyAPCLwJ8Bv9XuPhp4dvv6ZOAc4PcGGYckSZKkARu9ilVX7qVZTfFsJi1lX1Vv3tGJg75M8b00qy3uOmnfflW1GaCqNifZd7oTk5wAnACwmjUDDlOStCOOy5IkTevM9jFrA0vG2pUWt1TVBUmePdvzq2oDsAFgt+y5dPNsSRoRjsuSJD1YVZ0813MHWRl7JvDLSV4IrAZ2S/IPwE1J1rVVsXXAlgHGIEmSJEnzLsmpVfXSJBczzUWaVfWEHb3HwJKxqvp94PcB2srY71TVK5P8Bc1NpN/VPp8xqBgkSZIkdWMJrqb4lvZ5zvdeHvjS9tN4F3BqkuOB64CXDCEGSZIkSZqzSetgXDvX9+gkGauqc2hWTaSqfgAc2UW/kiRJkjpSWVz99CjJXTz4MsU7aO61/NtVdfVM5w6jMiZJkiRJi8V7gBuATwABjgUeClwBfJSf3NbrQbq46bMkSZIkLVZHVdXfVNVdVXVnu/rwC6vq08Ae2zvRZEySJElS/6rDx2gZT/LSJMvax0snHdtutCZjkiRJkjR3rwBeRXPLrpva169MsjPwpu2d6JwxSZIkSf3rqmI1YpWxdoGOX5rh8Ne3d67JmCRJkiTNUZLVwPHA44HVE/ur6td2dK6XKUqSJEnqW6qbxwj6e5rVE38B+FfgAOCuXk60MiYtUMsf/1MAXP72XQFYtWknAB7xjm8MLSZJkqQl6NFV9ZIkR1fVyUk+AXyxlxNNxiRJkiT1b4nOGQPub59vT3IocCNwUC8nmoxJoyrNVcQrfuqRANSq5o9rLWv2X/vO5vlxe20G4Lb91nQdoSRJkmBDkj2AdwBnArsA/6OXE03GJEmSJPVvCVbGkiwD7qyq24CvAY+czfkmY9KI2/ibewGw7yNvBWCvNXcCD65933LRvgA85LCdGL/w0q7CkyRJWrKqajzJm4BT53K+yZgkSZKkvnW10uEIrqj45SS/A3wauGdiZ1XduqMTTcakEXPrrz8DgLsPaLYfd+g107a7/PqHArDqv5rbWWw7+D4All1/E+ODDVGSJEk/MXE/sTdO2lf0cMmiyZgkSZIkzVFVHTzXc03GpBGzz3/cAsCu771vu+12e8i9APzwkG0APO6dPwRg2823DDA6SZKkGVSGHcFQJFkD/Bbw8Ko6Icl64DFV9fkdnWsyJo2YscuuAmDzHY8H4JF7/wCA8XaAu29b88f2jusfAsDj/upmALZddXWncUqSJAmAjwEXAM9otzcBnwF2mIwtG2BQkiRJkpaK6vAxWh5VVe+mvflzVf0Q6KlMaGVMGjXVLL9x8J80N3MfW7M3APn25QDs1N4Mev191wKwrev4JEmSNNnWJDvTpolJHgVsf75Jy2RMkiRJUt+W8NL2fwx8ATgwySnAM4HX9HKiyZg0osYuvvwB26M37kiSJKmqvpTkAuBpNJcnvqWqelpRzWRMkiRJUv+6+uV4xH6hTnIm8EngzKq6Z0ftJ3MBD0mSJEmau/8D/BxwWZLPJHlxktW9nGhlTJIkSVLfluqcsar6V+BfkywHngu8DvgosNuOzjUZkyRJkqQ+tKsp/hLwMuAw4ORezjMZkyRJktS/EatYdSXJp4EjaFZU/CBwTlV7r6IdMBmTJEmSpLn7GPCrVTUGkOSZSX61qt64oxNNxiRJkiRpjqrqC0melOTlNJcpfg84vZdzTcYkSZIk9W+JLW2f5KeAY4GXAz8APg2kqp7T63sMbGn7JAcm+WqSjUkuTfKWdv+eSb6c5Mr2eY9BxSBJkiRJA3I5cCTwS1X1s1X1/wFjs3mDQd5nbBvw21X1OJq7Ub8xySHAicDZVbUeOLvdljQIWdY8JEmSBizVzWNUKmPA/wPcCHw1yd8mORLIbN5gYP9Kq6rNVXVh+/ouYCOwP3A0P1nq8WTgmEHFIEmSJEmDUFWfq6qXAY8FzgHeBuyX5ENJfr6X9+jkJ/MkBwFPBs4F9quqzdAkbMC+M5xzQpLzk5x/P/d1Eaa0+NR485DmgeOyJEkPVlX3VNUpVfUi4ADgInq8+m/gyViSXYDPAm+tqjt7Pa+qNlTV4VV1+Ep2GlyAkqSeOC5LkrR9VXVrVf1NVT23l/YDTcaSrKRJxE6pqonlHW9Ksq49vg7YMsgYJEmSJHWgOnwsEoNcTTHAR4CNVfWeSYfOBI5rXx8HnDGoGCRJkiRpVA3yPmPPBF4FXJzkonbfHwDvAk5NcjxwHfCSAcYgSZIkqQPpqGLVVT9dGFgyVlVfZ+alHY8cVL+SJEmStBAMsjImSZIkaalYRBWrrng3WEmSJEkaApMxSZIkSRoCL1OUJEmS1L+uLlNcRJdDWhmTJEmSpCGwMiZJkiSpby5tP3tWxiRJkiRpCKyMSZIkSeqfc8ZmzcqYJEmSJA2BlTFJkiRJfXPO2OxZGZMkSZKkIbAyJkmSJKl/i6hi1RUrY5IkSZI0BCZjkiRJkvpXHT56kOSoJFckuSrJidMcT5L3t8e/m+Swdv+BSb6aZGOSS5O8ZS5fRy9MxiRJkiQtKkmWAx8EXgAcArw8ySFTmr0AWN8+TgA+1O7fBvx2VT0OeBrwxmnOnRcmY5IkSZIWm6cCV1XV1VW1FfgUcPSUNkcDH6/GN4Hdk6yrqs1VdSFAVd0FbAT2H0SQLuAhSZIkqW8dL23/uCTnT9q9oao2TNreH7h+0vYm4IgpbzVdm/2BzT/uKzkIeDJwbr9xT8dkTJIkSdJCs7Gq3rCd45lm39R0cbttkuwCfBZ4a1XdOfsQd8xkTJIkSVL/ulravrd+NgEHTto+ALih1zZJVtIkYqdU1elzDXVHnDMmSZIkabE5D1if5OAkq4BjgTOntDkTeHW7quLTgDuqanOSAB+hqb69Z5BBWhmTJEmS1L8RuulzVW1L8ibgi8By4KNVdWmS17fHTwLOAl4IXAXcC7y2Pf2ZwKuAi5Nc1O77g6o6a77jNBmTJEmStOi0ydNZU/adNOl1AW+c5ryvM/18snlnMiZJkiSpb12tpjhKFbh+OWdMkiRJkobAypgkSZKk/lkZmzUrY5IkSZI0BFbGJEmSJPWtqzljnc1N64CVMUmSJEkaAitjkiRJkvrnnLFZszImSZIkSUMwlGQsyVFJrkhyVZIThxGDJEmSJA1T58lYkuXAB4EXAIcAL09ySNdxSJIkSZpH1dFjERlGZeypwFVVdXVVbQU+BRw9hDgkSZIkaWiGsYDH/sD1k7Y3AUdMbZTkBOAEgNWs6SYySdKMHJclSduTRdZPF4aRjE33/T2o4FhVG4ANAElu/kqddg9wy4Bjm4u9Gc24YHRjG9W4YHRjG9W4YHRjG9W4YDRje8SOGjguz4tRjW1U44LRjW1U44LRjW1U44LRjG2H47IWvmEkY5uAAydtHwDcsL0TqmqfJOdX1eEDjWwORjUuGN3YRjUuGN3YRjUuGN3YRjUuGO3YeuW4PDejGtuoxgWjG9uoxgWjG9uoxgWjHduC4tL2szaMOWPnAeuTHJxkFXAscOYQ4pAkSZKkoem8MlZV25K8CfgisBz4aFVd2nUckiRJkuZPOqpYddVPF4ZxmSJVdRZw1ixP2zCIWObBqMYFoxvbqMYFoxvbqMYFoxvbqMYFox3bbIzq5xjVuGB0YxvVuGB0YxvVuGB0YxvVuGC0Y9MilqpFlFpKkiRJ6lySC57w5vcc1kVf9954LVed+r6TquoNXfQ3SMOYMyZJkiRJS95QLlOUJEmStMi4muKsWRmTJEmSpCEwGZMkSZKkIfAyRUmSJEl962xp+2666YSVMUmSJEkaAitjkiRJkvrnAh6zZmVMkiRJkobAypgkSZKkvnU1Z8zKmCRJkiSpL1bGJEmSJPXPytisWRmTepTkrUnWdNjfs5NUkuMn7Xtyu+93dnDuH0+0SfKnSZ436Hglab4tpHF3O+/54SSHzF+UkhYTkzGpd28Fpv1HQZLlA+rzYuBlk7aPBb4zmzeoqndW1VfmNSpJ6sZbWYDj7mRV9etVdVnfUUkLQKq7x2JhMqZFI8mrk3w3yXeS/H277xFJzm73n53k4e3+v0vy4knn3t0+PzvJOUlOS3J5klPSeDPwMOCrSb46cU5bdToXeEeSz016v+cnOX0ePtZ1wOok+yUJcBTwz5P6eV2S89rP/NnpfkGe+lklab4s0XH3UUm+kOSCJP+W5LFJVrRj8bPbNv8ryZ+1r89Jcnj7+qgkF7bf19ntvj2T/GP7fX0zyRPm4TNIWiCcM6ZFIcnjgT8EnllVtyTZsz30AeDjVXVykl8D3g8cs4O3ezLweOAG4N/b93x/kt8CnlNVt7Tt1gKXVNU727+wNybZp6puBl4LfGyaOP8KeM40fX6qqt41QzynAS8Bvg1cCNw36djpVfW37Xv/T+B44P/bweeTpL4t4XF3A/D6qroyyRHAX1fVc5O8BjitTSKPAo6YEsc+wN8Cz6qq7036vv4E+HZVHZPkucDHgSfNEJc02hZRxaorJmNaLJ4LnDbxF3ZV3drufzrwK+3rvwfe3cN7fauqNgEkuQg4CPj6NO3GgM+2/VX7q/Ark3ys7ffVU0+oqrf1+HkmOxX4NPBY4JPAMyYdO7RNwnYHdgG+OIf3l6S5WHLjbpJd2tefaXJBAHZq+7m0jeefgKdX1dYp7/k04GtV9b22/cT39bPA/9Pu+5ckeyV5SFXdMYe4JS0wJmNaLEJvv8dMtNlGe5lu++vqqkltJv8COsbMf05+VFVjk7Y/RvOX8I+Az1TVtgcFOYdfaKvqxiT3A88H3sIDk7G/A46pqu+0v8o+e4ZYJWm+LcVxdxlwe1U9aYb4fhq4HdhvmmMzfV+ZZp/1BWmJcM6YFouzgZcm2Quaa/Db/d+gmXwN8Ap+8kvrNcBT2tdHAyt76OMuYNeZDlbVDTSX2LyDJkmars3bqupJ0zxmulRmwjuB35vyjxDaeDYnWUnz+SSpK0tu3K2qO4HvJXkJNEllkie2r38F2At4FvD+JLtPeb//AP5bkoPb9hPf19dox+92ztktbT/SwlMdPhYJK2NaFNrLQ/4M+NckYzTX+b8GeDPw0SRvBybmFEBz3f4ZSb5F8w+Ke3roZgPwz0k2V9V0v7ICnALsM98rZ1XVN2Y49D+Ac4FraVYAm/EfLZI0n5bwuPsK4ENJ3kGTUH4qyfeBdwFHVtX1ST4AvA84btL73ZzkBOD0JMuALTSVtz8GPpbku8C9k8+RtPilahGlltKQtX8Bf7uqPjLsWCRpKXDclUZDkgsOe917Duuir3u2XMsVZ7zvpKp6Qxf9DZKVMWmeJLmA5pfe3x52LJK0FDjuSlroTMakeVJVT9lxK0nSfHHclUZMVxfcLaIL+1zAQ5IkSZKGwMqYJEmSpL6lo7UoprsfxEK1IJKxVdmpVrN22GFI0qJ1F7fdUlX79NrecVmSBmu247IWpgWRjK1mLUfkyGGHIUmL1lfqtGtn095xWZIGa7bj8khwztisOWdMkiRJkoZgQVTGJEmSJI22WBmbNStjkiRJkjQEJmOSJEmSNARepihJkiSpf16mOGtWxiRJkiRpCKyMSZIkSepbVwt4dLZQSAesjEmSJEnSEFgZkyRJktS/RVSx6oqVMUmSJEkaAitjkiRJkvrmnLHZszImSZIkSUNgZUySJElS/7zP2KwNtDKW5G1JLk1ySZJPJlmdZM8kX05yZfu8xyBjkCRJkqRRNLBkLMn+wJuBw6vqUGA5cCxwInB2Va0Hzm63JUmSJC1gqe4ei8Wg54ytAHZOsgJYA9wAHA2c3B4/GThmwDFIkiRJ0sgZWDJWVd8H/hK4DtgM3FFVXwL2q6rNbZvNwL7TnZ/khCTnJzn/fu4bVJiSpB45LkuSNL8GeZniHjRVsIOBhwFrk7yy1/OrakNVHV5Vh69kp0GFKUnqkeOyJGm7qrp7LBKDvEzxecD3qurmqrofOB14BnBTknUA7fOWAcYgSZIkSSNpkMnYdcDTkqxJEuBIYCNwJnBc2+Y44IwBxiBJkiSpA50t4DHsDzqPBnafsao6N8lpwIXANuDbwAZgF+DUJMfTJGwvGVQMkiRJkjSqBnrT56r6I+CPpuy+j6ZKJkmSJGmx8KbPszbope0lSZIkSdMYaGVMkiRJ0tKQ8Y466qqfDlgZkyRJkqQhsDImSZIkqX+LaC5XV6yMSZIkSdIQWBmTJEmS1Ld0VBnrqp8uWBmTJEmSpCEwGZMkSZKkIfAyRUmSJEn9q46uH+yqnw5YGZMkSZKkIbAyJkmSJKlvLuAxe1bGJEmSJGkIrIxJkiRJ6l9XFSsrY5IkSZKkflgZkyRJktS3zuaMddNNJ6yMSZIkSdIQWBmTJEmS1D/vMzZrVsYkSZIkaQisjEmSJEnqm/cZmz0rY5IkSZI0BCZjkiRJkjQEXqYoSZIkqX/e9HnWrIxJkiRJ0hBYGZMkSZLUt84W1rAyJkmSJEnqh5UxSZIkSf0b7+qmz9100wUrY5IkSZI0BFbGJEmSJPXPOWOzZmVMkiRJkobAypgkSZKkvnW1mmJnqzZ2YKCVsSS7JzktyeVJNiZ5epI9k3w5yZXt8x6DjEGSJEmSRtGgL1N8H/CFqnos8ERgI3AicHZVrQfObrclSRpJK9Y9lBXrHsryQx/D8kMfQ1asJCtWDjssSRo9Vd09FomBJWNJdgOeBXwEoKq2VtXtwNHAyW2zk4FjBhWDJEmSJI2qQc4ZeyRwM/CxJE8ELgDeAuxXVZsBqmpzkn2nOznJCcAJAKtZM8AwJUm9WOrj8i0/sycAKw9trq7fukv7e2b7A+2P9mqeD/znWwEYu/jyTuOTpGHrbM5YN910YpDJ2ArgMOA3q+rcJO9jFpckVtUGYAPAbtlz8dQiJWmBWqrj8rbNNwKw58X7NNsPWQXA7Y/eqdne+YHtrzmmSdrWPO0ZAKy6e5xdP/nNLkKVJC0wg5wztgnYVFXnttun0SRnNyVZB9A+bxlgDJIkSZI0kgaWjFXVjcD1SR7T7joSuAw4Eziu3XcccMagYpAkab4s/94NLP/eDdyzbhX3rFvFQZ/ZwkGf2cKqO2DVHT9pN7a6edx1UPP4waHL4OlPbB6StJhVh48eJDkqyRVJrkryoCv00nh/e/y7SQ6bdOyjSbYkuWS2X8NsDPo+Y78JnJJkFXA18FqaBPDUJMcD1wEvGXAMkiRJkpaQJMuBDwLPp7li77wkZ1bVZZOavQBY3z6OAD7UPgP8HfAB4OODjHOgyVhVXQQcPs2hIwfZrzRoKx51MAA3PW8dAPue8l0Axu6+e2gxSRqsbT/4AQC7f+bbANzz/CcA8PDTmzll1xz70Kbdzg8+95oXrQXgUbf/FABjG/9zoLFK0jCkoyXne+znqcBVVXU1QJJP0azqPjkZOxr4eFUV8M32HsnrqmpzVX0tyUHzHPqDDPo+Y5IkSZI03x6X5PxJjxOmHN8fuH7S9qZ232zbDNSgL1OUFqVauxqAx7/mUgAu++WHAbDbe3cDYMWXzx9OYJIGbvy+HwGw5utNdauWNYssH/SxZv/VrzuoaTfNfaGz5QeDD1CShmW8o36awtjGqnrDdlpNtwL+1JJaL20GysqYJEmSpMVmE3DgpO0DgBvm0GagrIxJc3Hd5gdsHrLnTc2LP22eL3hKc3+hA9/TVMjGt27tLjZJnRi7/fYpe24DYO33DwKalRSn2vaY9u/8bzywQrZin70BuOmYZk7ZXpfc0xz4j+/MQ6SS1I0RmzN2HrA+ycHA94FjgV+d0uZM4E3tfLIjgDuqajMdsjImSZIkaVGpqm3Am4AvAhuBU6vq0iSvT/L6ttlZNCu+XwX8LfAbE+cn+STwH8BjkmxqV4Kfd1bGpDmY+EX8lvuauWJ773TPA44/5ehmLtk3nvB4ANa/6RrGbr2tuwAlDc2+X78FgB/t1VS77t/1J8euP3INAAfd1lTAJn7dve3JewHww32bdiuu3QLAtoFHK0nzqKvZVj32U1Vn0SRck/edNOl1AW+c4dyXzz3A3lkZkyRJkqQhsDIm9eG/tjS/fO994D3THn/GwVcDcMtO09x4SNKiNHEPsf3b52VPOZQtP9OstHp3O2Xsmhc3Y8dYszAre1/cLEH28L+8EIBt7YqNkrSgdDRnrLN+OmBlTJIkSZKGwMqYtB3LVq0C4D//4jAAxteOAZBlzS8yj9qrWT3x9q1N5es/z3kkAGvbRVF3+16ziuLKzd53TFqqxi+4hL0vaF7vt2szgWzs7qaavuywQ5o2F17WPFdXN+mRJI0CkzFJkiRJfcuILeCxEJiMSdsxcX+w2v1+AB59QLPC2cPW3AnAxbesA2Dl77bzPTZ+o+sQJS0gY3fd9YDt8QsuGVIkkqRRYDImSZIkqX+LaGGNrpiMST1Yf1wz4SNP/WkAzv2dPQA48K+bP0JjGy8cTmCSJElasEzGJEmSJPUtHa1B1FU/XTAZk2ahvnUxAI946ZADkSRJ0oJnMiZJkiSpf970eda86bMkSZIkDYGVMUmSJEn98z5js2ZlTJIkSZKGwMqYJEmSpL6lo7lcWUSlMStjkiRJkjQEJmOSJEmSNAQ7TMaSvCnJHl0EI0mSJGmBqurusUj0Uhl7KHBeklOTHJUkgw5KkiRJkha7HSZjVfUOYD3wEeA1wJVJ/jzJowYcmyRJkqSFYrzDxyLR05yxqirgxvaxDdgDOC3JuwcYmyRJkiQtWjtc2j7Jm4HjgFuADwNvr6r7kywDrgR+d7AhSpIkSRp1nS1tv4jmjPVyn7G9gV+pqmsn76yq8SQvGkxYkiRJkrS47TAZq6p3bufYxh2dn2Q5cD7w/ap6UZI9gU8DBwHXAC+tqtt6DViSJEnSCFpEFauudHGfsbcAk5O2E4Gzq2o9cHa7LUmSJElLykCTsSQHAL9IM9dswtHAye3rk4FjBhmDJEmSpA54n7FZG3Rl7L00C3xMXoByv6raDNA+7zvgGCRJkiRp5AwsGWsX99hSVRfM8fwTkpyf5Pz7uW+eo5MkzZbjsiRpu7zP2Kz1spriXD0T+OUkLwRWA7sl+QfgpiTrqmpzknXAlulOrqoNwAaA3bLn4qlFStIC5bgsSdL8GlhlrKp+v6oOqKqDgGOBf6mqVwJn0ty3jPb5jEHFIEmSJEmjapCVsZm8Czg1yfHAdcBLhhCDJEmSpHnU3c2YF8/FGZ0kY1V1DnBO+/oHwJFd9CtJkiRJo2oYlTFJkiRJi01XlTGXtpckSZIk9cPKmCRJkqT+WRmbNStjkiRJkjQEVsYkSZIk9a+zylg33XTBypgkSZIkDYGVMUmSJEn9G++oHytjkiRJkqR+WBmTJEmS1Ld0NGesq366YGVMkiRJkobAZEySJEmShsDLFCVJkiT1z5s+z5qVMUmSJEkaAitjkiRJkvo33lHFqqt+OmBlTJIkSZKGwMqYJEmSpP4torlcXbEyJkmSJElDYGVMkiRJUv9cTXHWrIxJkiRJ0hBYGZMkSZLUPytjs2ZlTJIkSZKGwMqYJEmSpP51df8vK2OSJEmSpH6YjEmSJEnSEHiZoiRJkqT+1XhH/XiZoiRJkiSpD1bGJEmSJPXPpe1nzcqYJEmSJA2BlTFJkiRJ/etqafuu+umAlTFJkiRJGoKBJWNJDkzy1SQbk1ya5C3t/j2TfDnJle3zHoOKQZIkSVJHqrp5YGWsF9uA366qxwFPA96Y5BDgRODsqloPnN1uS5IkSdKSMrA5Y1W1Gdjcvr4ryUZgf+Bo4Nlts5OBc4DfG1QckiRJkjrgaoqz1skCHkkOAp4MnAvs1yZqVNXmJPvOcM4JwAkAq1nTRZiSpO1wXJYkaX4NPBlLsgvwWeCtVXVnkp7Oq6oNwAaA3bLn4kl/JWmBclyWJG1XZ5WxbrrpwkBXU0yykiYRO6WqTm9335RkXXt8HbBlkDFIkiRJ0iga5GqKAT4CbKyq90w6dCZwXPv6OOCMQcUgSZIkSaNqkJcpPhN4FXBxkovafX8AvAs4NcnxwHXASwYYgyRJkqQujI8vrn46MMjVFL8OzDRB7MhB9StJkiRJC0EnqylKkiRJWuQ6W3J+8azgMdAFPCRJkiRJ07MyJkmSJKl/3vR51qyMSZIkSdIQWBmTJEmS1L/xjipWi2cxRStjkiRJkjQMVsYkSZIk9a2qm5JVV/10wcqYJEmSJA2BlTFJkiRJ/etqzpirKUqSJEmS+mEyJkmSJElD4GWKkiRJkvrX2U2fu+mmC1bGJEmSJGkIrIxJkiRJ6t94R0vOu7S9JEmSJKkfVsYkSZIk9a+zOWOLZ9KYlTFJkiRJGgIrY5IkSZL6Vl3NGevq5tIdsDImSZIkSUNgZUySJElS/5wzNmtWxiRJkiRpCKyMSZIkSepfV3O5rIxJkiRJkvphZUySJElS/6qj1RS76qcDVsYkSZIkaQhMxiRJkiRpCLxMUZIkSVLfqqMFPBbRVYpWxiRJkiRpGKyMSZIkSepfZyWrxVMaszImSZIkSUMwlGQsyVFJrkhyVZIThxGDJEmSpPlT49XZY7HoPBlLshz4IPAC4BDg5UkO6ToOSZIkSYvXjgpAaby/Pf7dJIf1eu58GUZl7KnAVVV1dVVtBT4FHD2EOCRJkiTNlxrv6LHjyliPBaAXAOvbxwnAh2Zx7rwYxgIe+wPXT9reBBwxtVGSE2i+FFazppvIJEkzclyWJG3Hmo11QScdbeU+gPt30OzHBSCAJBMFoMsmtTka+HhVFfDNJLsnWQcc1MO582IYyVim2feg9LaqNgAbAJLc/JU67R7glgHHNhd7M5pxwejGNqpxwejGNqpxwejGNqpxwWjG9ogdNXBcnhejGtuoxgWjG9uoxgWjG9uoxgWjGdsOx+URc9T3+d7eHfb3s0nOn7S9of17akIvBaDp2uzf47nzYhjJ2CbgwEnbBwA3bO+EqtonyflVdfhAI5uDUY0LRje2UY0LRje2UY0LRje2UY0LRju2Xjkuz82oxjaqccHoxjaqccHoxjaqccFox7ZQVNW1wLUddnkB8L7tHO+lADRTm56KR/NhGHPGzgPWJzk4ySrgWODMIcQhSZIkaXHqpQA0U5tZF4/mqvNkrKq2AW8CvghsBE6tqku7jkOSJEnSotVLAehM4NXtqopPA+6oqs09njsvhnGZIlV1FnDWLE/bsOMmQzGqccHoxjaqccHoxjaqccHoxjaqccFoxzYbo/o5RjUuGN3YRjUuGN3YRjUuGN3YRjUuGO3YNAdVtS3JRAFoOfDRqro0yevb4yfR5CMvBK4C7gVeu71zBxFnqoelISVJkiRJ82sYc8YkSZIkackzGZMkSZKkIRj5ZCzJUUmuSHJVkhOHHMuBSb6aZGOSS5O8pd2/Z5IvJ7myfd5jSPEtT/LtJJ8fsbh2T3Jaksvb7+7poxBbkre1/x0vSfLJJKuHFVeSjybZkuSSSftmjCXJ77d/Jq5I8gtDiO0v2v+e303yuSS7dx3bdHFNOvY7SSrJ3pP2DfU7a/f/Ztv/pUnePYzY5oPj8qzic1yeXVyOy3OLa+hj8kyxTTo2tHF5sY/JWuCqamQfNBPm/gt4JLAK+A5wyBDjWQcc1r7eFfhP4BDg3cCJ7f4Tgf89pPh+C/gE8Pl2e1TiOhn49fb1KmD3YcdGczO/7wE7t9unAq8ZVlzAs4DDgEsm7Zs2lvb/c98BdgIObv+MLO84tp8HVrSv//cwYpsurnb/gTQTbq8F9h6h7+w5wFeAndrtfYcR2zx8Nsfl2cXnuNx7TI7Lc49r6GPyTLG1+4c6Li/mMdnHwn+MemXsqcBVVXV1VW0FPgUcPaxgqmpzVV3Yvr6LZmn+/duYTm6bnQwc03VsSQ4AfhH48KTdoxDXbjSD4EcAqmprVd0+CrHRrCa6c5IVwBqa+0cMJa6q+hpw65TdM8VyNPCpqrqvqr5HswLQU7uMraq+VM1tKgC+SXP/jU5jm+E7A/gr4Hd54M0Zh/6dAW8A3lVV97VttgwjtnnguNwjx+U5cVyeQ1yjMCbPFFtrqOPyIh+TtcCNejK2P3D9pO1N7b6hS3IQ8GTgXGC/au5JQPu87xBCei/NQDc+ad8oxPVI4GbgY+2lOh9OsnbYsVXV94G/BK4DNtPcV+JLw45ripliGbU/F78G/HP7eqixJfll4PtV9Z0ph0bhO/sp4OeSnJvkX5P8zAjFNhsjG6/jcs8cl+duIYzLIzMmw0iPy4tlTNYCN+rJWKbZN/S1+JPsAnwWeGtV3TkC8bwI2FJVFww7lmmsoLk04ENV9WTgHppLO4aqvc7/aJpLEB4GrE3yyuFG1bOR+XOR5A+BbcApE7umadZJbEnWAH8IvHO6w9Ps6/o7WwHsATwNeDtwapIwGrHNxkjG67g8K47L828k/lyM0pjcxjPK4/JiGZO1wI16MraJ5jrjCQfQXLIwNElW0vyFf0pVnd7uvinJuvb4OmDLTOcPyDOBX05yDc0lQ89N8g8jEBc0/w03VdW57fZpNP8IGHZszwO+V1U3V9X9wOnAM0YgrslmimUk/lwkOQ54EfCKqpr4i2qYsT2K5h9x32n/LBwAXJjkoUOOa8Im4PRqfIumWrL3iMQ2GyMXr+PyrDkuz93IjssjOCbDaI/Li2VM1gI36snYecD6JAcnWQUcC5w5rGDaX0w+AmysqvdMOnQmcFz7+jjgjC7jqqrfr6oDquogmu/oX6rqlcOOq43tRuD6JI9pdx0JXDYCsV0HPC3Jmva/65E0c02GHddkM8VyJnBskp2SHAysB77VZWBJjgJ+D/jlqrp30qGhxVZVF1fVvlV1UPtnYRPNwg43DjOuSf4ReC5Akp+iWTThlhGJbTYcl3vguDwnjstzNIpjMoz8uPyPLI4xWQtdDWnlkF4fwAtpVsf6L+APhxzLz9KUqr8LXNQ+XgjsBZwNXNk+7znEGJ/NT1btGom4gCcB57ff2z/SXBYw9NiAPwEuBy4B/p5m5aShxAV8kmaOxP00f1kdv71YaC77+C/gCuAFQ4jtKppr6if+HJzUdWzTxTXl+DW0q3aNyHe2CviH9v9vFwLPHUZs8/T5HJdnF6Pjcu9xOS7PLa6hj8kzxTbl+FDG5cU+JvtY2I9UeRmsJEmSJHVt1C9TlCRJkqRFyWRMkiRJkobAZEySJEmShsBkTJIkSZKGwGRMkiRJkobAZEySJEmShsBkTJIkSZKGwGRMmkdJfibJd5OsTrI2yaVJDh12XJK0VDkuSxpl3vRZmmdJ/iewGtgZ2FRV/2vIIUnSkua4LGlUmYxJ8yzJKuA84EfAM6pqbMghSdKS5rgsaVR5maI0//YEdgF2pfklVpI0XI7LkkaSlTFpniU5E/gUcDCwrqreNOSQJGlJc1yWNKpWDDsAaTFJ8mpgW1V9Isly4BtJnltV/zLs2CRpKXJcljTKrIxJkiRJ0hA4Z0ySJEmShsBkTJIkSZKGwGRMkiRJkobAZEySJEmShsBkTJIkSZKGwGRMkiRJkobAZEySJEmShuD/D04s6CAaMePGAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Convert sparse matrix to numpy array, it has size : (n_in, n_out)\n", "# So reshape to the same shape as ds + polygons\n", "w = xr.DataArray(\n", " savg.weights.toarray().reshape(regs.geometry.size, *ds.lon.shape),\n", " dims=(\"country\", *ds.lon.dims),\n", " coords=dict(country=out.country, **ds.lon.coords),\n", ")\n", "plt.subplots_adjust(top=0.9)\n", "facets = w.plot(col=\"country\", col_wrap=2, aspect=2, vmin=0, vmax=0.05)\n", "facets.cbar.set_label(\"Averaging weights\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This also allows to quickly check that the weights are indeed normalized, that\n", "the sum of each mask is 1.\n" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "pycharm": { "name": "#%%\n" } }, "outputs": [ { "data": { "text/plain": [ "array([1., 1., 1., 1., 1., 1.])" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "w.sum(dim=[\"y\", \"x\"]).values" ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }