diff --git a/.gitignore b/.gitignore index 4a4fa10..883508f 100644 --- a/.gitignore +++ b/.gitignore @@ -99,6 +99,7 @@ tests/out/*.ipynb !tests/localize-data.ipynb /isce2_topsapp.egg-info/ !tests/prepare-for-delivery.ipynb +!tests/solid_earth_tides.ipynb # Jsons in Test Directory tests/**/*.json diff --git a/CHANGELOG.md b/CHANGELOG.md index 14f9668..bd4111e 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -6,6 +6,11 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/), and this project adheres to [PEP 440](https://www.python.org/dev/peps/pep-0440/) and uses [Semantic Versioning](https://semver.org/spec/v2.0.0.html). +## [0.2.3] + +### Added +* Added support to compute and embed solid earth tide correction layers into GUNW products (see PR #91) + ## [0.2.2] ### Added diff --git a/README.md b/README.md index 80b3048..bd24ddf 100644 --- a/README.md +++ b/README.md @@ -64,6 +64,17 @@ isce2_topsapp --reference-scenes S1B_IW_SLC__1SDV_20171117T145926_20171117T14595 ``` Not including `--esd-coherence-threshold` means no ESD correction will be applied. The ESD threshold refers to a coherence value and therefore must be in $[0, 1]$. +### Apply Solid Earth Tide Correction + +``` +isce2_topsapp --reference-scenes S1B_IW_SLC__1SDV_20210723T014947_20210723T015014_027915_0354B4_B3A9 \ + --secondary-scenes S1B_IW_SLC__1SDV_20210711T014922_20210711T014949_027740_034F80_859D \ + S1B_IW_SLC__1SDV_20210711T014947_20210711T015013_027740_034F80_D404 \ + S1B_IW_SLC__1SDV_20210711T015011_20210711T015038_027740_034F80_376C \ + --compute-solid-earth-tide True \ + > topsapp_img_set.out 2> topsapp_img_set.err +``` + ### Using "fixed frames" (experimental) Sentinel-1 Frames are not constant over passes. We generate fixed frames [here](https://github.com/ACCESS-Cloud-Based-InSAR/s1-frame-generation) and enumerate interferograms using this [repo](https://github.com/ACCESS-Cloud-Based-InSAR/s1-frame-enumerator). This is highly experimental. We then ensure ISCE processes only over the frame. The key is overlap. We provide some examples of the additional options (you will need to run this in *two* separate directories because ISCE2 outputs are organized with respect to the working directory of the processing). For one frame over CA: diff --git a/environment.yml b/environment.yml index d39604c..7d463f6 100644 --- a/environment.yml +++ b/environment.yml @@ -29,6 +29,7 @@ dependencies: - notebook - numpy<1.24 - pandas + - pysolid - papermill - pytest - pytest-cov @@ -39,6 +40,7 @@ dependencies: - scipy<1.10 - setuptools - setuptools_scm + - scipy<1.10 - shapely - tqdm - dem_stitcher>=2.3.1 diff --git a/isce2_topsapp/__main__.py b/isce2_topsapp/__main__.py index fb56d1f..f0cf264 100644 --- a/isce2_topsapp/__main__.py +++ b/isce2_topsapp/__main__.py @@ -8,6 +8,7 @@ from pathlib import Path from typing import Optional +import h5py from isce2_topsapp import (BurstParams, aws, download_aux_cal, download_bursts, download_dem_for_isce2, download_orbits, @@ -15,6 +16,7 @@ package_gunw_product, prepare_for_delivery, topsapp_processing) from isce2_topsapp.json_encoder import MetadataEncoder +from isce2_topsapp.solid_earth_tides import update_gunw_with_solid_earth_tide def localize_data(reference_scenes: list, @@ -114,6 +116,7 @@ def gunw_slc(): parser.add_argument('--secondary-scenes', type=str.split, nargs='+', required=True) parser.add_argument('--estimate-ionosphere-delay', type=true_false_string_argument, default=False) parser.add_argument('--frame-id', type=int, default=-1) + parser.add_argument('--compute-solid-earth-tide', type=true_false_string_argument, default=False) parser.add_argument('--esd-coherence-threshold', type=float, default=-1.) args = parser.parse_args() @@ -162,9 +165,15 @@ def gunw_slc(): reference_properties=ref_properties, secondary_properties=sec_properties, extent=extent, - additional_2d_layers=additional_2d_layers + additional_2d_layers=additional_2d_layers, ) + if args.compute_solid_earth_tide: + nc_path = update_gunw_with_solid_earth_tide(nc_path) + # Update to 1c + with h5py.File(nc_path, mode='a') as file: + file.attrs.modify('version', '1c') + # Move final product to current working directory final_directory = prepare_for_delivery(nc_path, loc_data) diff --git a/isce2_topsapp/packaging_utils/additional_layers.py b/isce2_topsapp/packaging_utils/additional_layers.py index 8dc0cbc..43ede48 100644 --- a/isce2_topsapp/packaging_utils/additional_layers.py +++ b/isce2_topsapp/packaging_utils/additional_layers.py @@ -24,9 +24,12 @@ def add_2d_layer(layer_name: str, gunw_netcdf_path: Path) -> Path: # The layers generally already exist within the file with h5py.File(gunw_netcdf_path, mode='a') as file: if dst_group in file: - for key in file[dst_group].keys(): - del file[dst_group][key] - del file[dst_group] + # Delete the variable to be written to + if dst_variable in file[dst_group]: + del file[dst_group][dst_variable] + # Delete the group if there are no variables left + if len(file[dst_group].keys()) == 0: + del file[dst_group] ds = xr.open_dataset(layer_data['input_relative_path'], engine='rasterio') diff --git a/isce2_topsapp/solid_earth_tides.py b/isce2_topsapp/solid_earth_tides.py new file mode 100644 index 0000000..5c4dd73 --- /dev/null +++ b/isce2_topsapp/solid_earth_tides.py @@ -0,0 +1,167 @@ +# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +# +# Authors: Simran Sangha, Charles Marshak, Brett Buzzanga, David Bekaert, +# Marin Govorcin, Zhang Yunjun +# Copyright 2023, by the California Institute of Technology. ALL RIGHTS +# RESERVED. United States Government Sponsorship acknowledged. +# +# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +import datetime as dt +from pathlib import Path +from typing import Tuple + +import h5py +import numpy as np +import pysolid +import xarray as xr + + +def compute_solid_earth_tide_from_gunw(gunw_path: str) -> xr.Dataset: + """ Read GUNW and compute/export differential SET """ + + gunw_basename = Path(gunw_path).name + dates = gunw_basename.split('-')[6].split('_') + ref_time = gunw_basename.split('-')[7] + + # get GUNW attributes + (wavelength, + z_meta, + inc_angle, + az_angle, + set_attrs) = get_gunw_attrs(gunw_path) + + # compute differential SET ENU + # Use reference time twice (assuming exact same pass time) + tide_e_ref, tide_n_ref, tide_u_ref = compute_enu_se_tide(dates[0], + ref_time, + set_attrs) + tide_e_sec, tide_n_sec, tide_u_sec = compute_enu_se_tide(dates[1], + ref_time, + set_attrs) + tide_e = tide_e_sec - tide_e_ref + tide_n = tide_n_sec - tide_n_ref + tide_u = tide_u_sec - tide_u_ref + + # compute SET LOS estimate for each height level + tide_los = np.zeros(inc_angle.shape) + for i in range(len(z_meta)): + tide_los[i] = compute_los_solid_earth_tide(tide_e, + tide_n, + tide_u, + inc_angle[i], + az_angle[i], + wavelength) + + solidtide_corr_ds = export_se_tides_to_dataset(gunw_path, tide_los) + + return solidtide_corr_ds + + +def get_gunw_attrs(gunw_path: str) -> Tuple[float, np.array, np.array, + np.array, xr.Dataset, dict]: + """ Access necessary GUNW attributes to compute SET """ + + group = 'science/radarMetaData' + with xr.open_dataset(gunw_path, group=group) as ds: + wavelength = ds['wavelength'].item() + + group = 'science/grids/imagingGeometry' + with xr.open_dataset(gunw_path, group=group, engine='rasterio') as ds: + z_meta = ds.heightsMeta.data + lat_meta = ds.y.data + lon_meta = ds.x.data + # convert angles to rad + inc_angle = np.deg2rad(ds.incidenceAngle.data) + az_angle = np.deg2rad(ds.azimuthAngle.data-90) + solidtide_atr = { + 'LENGTH': len(lat_meta), + 'WIDTH': len(lon_meta), + 'X_FIRST': lon_meta[0], + 'Y_FIRST': lat_meta[0], + 'X_STEP': lon_meta[1] - lon_meta[0], + 'Y_STEP': lat_meta[1] - lat_meta[0], + } + + return wavelength, z_meta, inc_angle, az_angle, solidtide_atr + + +def compute_enu_se_tide(acq_date: str, + acq_time: str, + solidtide_atr: dict) -> np.array: + """ Compute SET in ENU """ + # source: + # https://github.com/insarlab/ and + # PySolid/blob/main/docs/plot_grid_SET.ipynb + dt_obj = dt.datetime.strptime(acq_date + '-' + acq_time, + "%Y%m%d-%H%M%S") + tide_e, tide_n, tide_u = pysolid.calc_solid_earth_tides_grid(dt_obj, + solidtide_atr, + display=False, + verbose=False) + + return tide_e, tide_n, tide_u + + +def compute_los_solid_earth_tide(tide_e: np.array, + tide_n: np.array, + tide_u: np.array, + inc_angle: np.array, + az_angle: np.array, + wavelength: float) -> np.array: + """ Compute SET in LOS """ + # source: + # https://github.com/insarlab/ + # PySolid/blob/main/docs/plot_grid_SET.ipynb + + # project ENU to radar line-of-sight (LOS) + # with positive for motion towards satellite + tide_e_slant = tide_e * np.sin(inc_angle) * np.sin(az_angle) * -1 + tide_n_slant = tide_n * np.sin(inc_angle) * np.cos(az_angle) + tide_u_slant = tide_u * np.cos(inc_angle) + tide_los = tide_e_slant + tide_n_slant + tide_u_slant + + # Convert m to rad + tide_los = tide_los / (wavelength / (4*np.pi)) + + return tide_los + + +def export_se_tides_to_dataset(gunw_path: str, + tide_los: np.array, + lyr_name='solidEarthTide') -> xr.Dataset: + """ Create SET array and populate with metadata leveraging + the same geodata from the imaging geometry in the gunw""" + + # obtain affine transformation and coordinate metadata + group = 'science/grids/imagingGeometry' + solidtide_corr_ds = xr.open_dataset(gunw_path, + group=group) + + solidtide_corr_ds = solidtide_corr_ds.rename({'longitudeMeta': 'longitude', + 'latitudeMeta': 'latitude', + 'heightsMeta': 'height'}) + attrs = {'description': 'Solid Earth tide', + 'units': 'radians', + 'long_name': lyr_name, + 'standard_name': lyr_name} + + dim_order = ['height', 'latitude', 'longitude'] + solidtide_corr_ds[lyr_name] = (dim_order, tide_los) + solidtide_corr_ds[lyr_name].attrs.update(attrs) + solidtide_corr_ds = solidtide_corr_ds[[lyr_name]] + solidtide_corr_ds = solidtide_corr_ds.astype(np.float32) + solidtide_corr_ds.rio.write_crs('epsg:4326', inplace=True) + solidtide_corr_ds['solidEarthTide'].rio.write_nodata(0, inplace=True) + return solidtide_corr_ds + + +def update_gunw_with_solid_earth_tide(gunw_path: Path) -> Path: + se_tide_group = '/science/grids/corrections/external/tides' + with h5py.File(gunw_path, 'a') as file: + if se_tide_group in file: + del file[se_tide_group] + solid_earth_tide_ds = compute_solid_earth_tide_from_gunw(gunw_path) + solid_earth_tide_ds.to_netcdf(gunw_path, + mode='a', + group=se_tide_group) + return gunw_path diff --git a/tests/solid_earth_tides.ipynb b/tests/solid_earth_tides.ipynb new file mode 100644 index 0000000..4b61f5d --- /dev/null +++ b/tests/solid_earth_tides.ipynb @@ -0,0 +1,732 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "source": [ + "from pathlib import Path\n", + "import requests\n", + "from isce2_topsapp.solid_earth_tides import compute_solid_earth_tide_from_gunw, update_gunw_with_solid_earth_tide" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:03.484032Z", + "start_time": "2023-02-07T02:10:00.556938Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 2, + "source": [ + "url = 'https://grfn.asf.alaska.edu/door/download/S1-GUNW-A-R-064-tops-20210723_20210711-015001-35393N_33512N-PP-6267-v2_0_4.nc'\n", + "gunw_path = Path(url.split('/')[-1])\n", + "\n", + "resp = requests.get(url)\n", + "with open(gunw_path, 'wb') as file:\n", + " file.write(resp.content)" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:13.866460Z", + "start_time": "2023-02-07T02:10:03.486567Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 3, + "source": [ + "set_ds = compute_solid_earth_tide_from_gunw(gunw_path)\n", + "set_ds" + ], + "outputs": [ + { + "output_type": "stream", + "name": "stderr", + "text": [ + "/Users/ssangha/Downloads/snap_setup/stable_oct5_2021/envs/topsapp_env/lib/python3.9/site-packages/pysolid/grid.py:93: UserWarning: Input line 1 contained no data and will not be counted towards `max_rows=1062`. This differs from the behaviour in NumPy <=1.22 which counted lines rather than rows. If desired, the previous behaviour can be achieved by using `itertools.islice`.\n", + "Please see the 1.23 release notes for an example on how to do this. If you wish to ignore this warning, use `warnings.filterwarnings`. This warning is expected to be removed in the future and is given only once per `loadtxt` call.\n", + " fc = np.loadtxt(txt_file,\n" + ] + }, + { + "output_type": "execute_result", + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
<xarray.Dataset>\n",
+       "Dimensions:         (height: 4, latitude: 26, longitude: 37)\n",
+       "Coordinates:\n",
+       "  * longitude       (longitude) float64 -119.2 -119.1 -119.0 ... -115.7 -115.6\n",
+       "  * latitude        (latitude) float64 35.7 35.6 35.5 35.4 ... 33.4 33.3 33.2\n",
+       "  * height          (height) float64 -1.5e+03 0.0 3e+03 9e+03\n",
+       "    spatial_ref     int64 0\n",
+       "Data variables:\n",
+       "    solidEarthTide  (height, latitude, longitude) float32 17.8 17.7 ... 14.47
" + ], + "text/plain": [ + "\n", + "Dimensions: (height: 4, latitude: 26, longitude: 37)\n", + "Coordinates:\n", + " * longitude (longitude) float64 -119.2 -119.1 -119.0 ... -115.7 -115.6\n", + " * latitude (latitude) float64 35.7 35.6 35.5 35.4 ... 33.4 33.3 33.2\n", + " * height (height) float64 -1.5e+03 0.0 3e+03 9e+03\n", + " spatial_ref int64 0\n", + "Data variables:\n", + " solidEarthTide (height, latitude, longitude) float32 17.8 17.7 ... 14.47" + ] + }, + "metadata": {}, + "execution_count": 3 + } + ], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:14.199042Z", + "start_time": "2023-02-07T02:10:13.871876Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 4, + "source": [ + "update_gunw_with_solid_earth_tide(gunw_path)" + ], + "outputs": [ + { + "output_type": "stream", + "name": "stderr", + "text": [ + "/Users/ssangha/Downloads/snap_setup/stable_oct5_2021/envs/topsapp_env/lib/python3.9/site-packages/pysolid/grid.py:93: UserWarning: Input line 1 contained no data and will not be counted towards `max_rows=1062`. This differs from the behaviour in NumPy <=1.22 which counted lines rather than rows. If desired, the previous behaviour can be achieved by using `itertools.islice`.\n", + "Please see the 1.23 release notes for an example on how to do this. If you wish to ignore this warning, use `warnings.filterwarnings`. This warning is expected to be removed in the future and is given only once per `loadtxt` call.\n", + " fc = np.loadtxt(txt_file,\n" + ] + }, + { + "output_type": "execute_result", + "data": { + "text/plain": [ + "PosixPath('S1-GUNW-A-R-064-tops-20210723_20210711-015001-35393N_33512N-PP-6267-v2_0_4.nc')" + ] + }, + "metadata": {}, + "execution_count": 4 + } + ], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:14.391162Z", + "start_time": "2023-02-07T02:10:14.200716Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 5, + "source": [ + "import geopandas as gpd\n", + "df_world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))\n", + "xmin, ymin, xmax, ymax = set_ds.rio.bounds()\n", + "\n", + "s = set_ds['solidEarthTide'].plot.imshow(col='height', \n", + " col_wrap=4, \n", + " extent=[xmin, ymin, xmax, ymax])\n", + "[df_world.boundary.plot(ax=ax, color='black') for ax in s.axes[0]]\n", + "[ax.set_ylim(ymin, ymax) for ax in s.axes[0]]\n", + "[ax.set_xlim(xmin, xmax) for ax in s.axes[0]]" + ], + "outputs": [ + { + "output_type": "stream", + "name": "stderr", + "text": [ + "/var/folders/13/168512yn66d8nmvt7g6jsybh0000gq/T/ipykernel_59692/183130835.py:8: DeprecationWarning: self.axes is deprecated since 2022.11 in order to align with matplotlibs plt.subplots, use self.axs instead.\n", + " [df_world.boundary.plot(ax=ax, color='black') for ax in s.axes[0]]\n", + "/var/folders/13/168512yn66d8nmvt7g6jsybh0000gq/T/ipykernel_59692/183130835.py:9: DeprecationWarning: self.axes is deprecated since 2022.11 in order to align with matplotlibs plt.subplots, use self.axs instead.\n", + " [ax.set_ylim(ymin, ymax) for ax in s.axes[0]]\n", + "/var/folders/13/168512yn66d8nmvt7g6jsybh0000gq/T/ipykernel_59692/183130835.py:10: DeprecationWarning: self.axes is deprecated since 2022.11 in order to align with matplotlibs plt.subplots, use self.axs instead.\n", + " [ax.set_xlim(xmin, xmax) for ax in s.axes[0]]\n" + ] + }, + { + "output_type": "execute_result", + "data": { + "text/plain": [ + "[(-119.25, -115.55000000000001),\n", + " (-119.25, -115.55000000000001),\n", + " (-119.25, -115.55000000000001),\n", + " (-119.25, -115.55000000000001)]" + ] + }, + "metadata": {}, + "execution_count": 5 + }, + { + "output_type": "display_data", + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAABIAAAAD8CAYAAAAPB9/zAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB//klEQVR4nO3dd3gU1f4G8Hd2N733QiChCFIEQlGwASogoiJYULkCYvenCBZQLBQVBAvoRcGK4PXaBRtSRBBBQVoUBUEgJAHSey+78/sjZC8heyY7k8m2vB+f80imnHNmk++cycnMfCVZlmUQEREREREREZHHMji7A0RERERERERE1Lo4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AURERERERERE5OE4AUREREREREREpMLWrVtxzTXXID4+HpIkYc2aNY3Wl5WV4YEHHkBCQgL8/PzQvXt3LFu2zDmdPY0TQEREREREREREKpSXl6NPnz5YunSpzfXTp0/HunXr8J///AcHDx7E9OnT8eCDD+Krr75ycE//R5JlWXZa60REREREREREraSqqgo1NTV2bevt7Q1fX1/VbUiShNWrV+O6666zLuvVqxfGjx+Pp59+2rqsf//+uOqqq/Dss8+qbkMPJqe0SkRERERERETUiqqqqtAxMRBZOWa7to+NjcXvv//eaBLIx8cHPj4+qtu++OKL8fXXX2PKlCmIj4/Hli1bcPjwYbz66quq69ILJ4CIiIiIiIiIyOPU1NQgK8eM1D2JCA5SfgNOSakFHfunISYmptHy2bNnY86cOarbfu2113DXXXchISEBJpMJBoMB77zzDi6++GLVdemFE0BERERERERE5LECAuuLEvPpl+NkZGQgODjYulzL3T9A/QTQjh078PXXXyMxMRFbt27F/fffj7i4OFxxxRWa6mwpTgARERERERERkceqgxl1UH79cR0sAIDg4OBGE0BaVFZWYtasWVi9ejVGjx4NAOjduzdSUlLw0ksvOW0CiFnA3MzQoUMxbdq0FtUxZ84c9O3b1+HtEpFj8XxB5NoYo0Seh3FN5JrMsmxX0UttbS1qa2thMDSecjEajbBYLLq1oxYngNqgRx99FJs2bdK9XkmSsGbNGt3rtcfzzz+PCy+8EP7+/ggNDbVrn8mTJ0OSpEZl0KBBrdtRAG+88QY6duwIX19f9O/fHz///HOj9XPmzMG5556LgIAAhIWF4YorrsDOnTtbvV9Etrj7+aK6uhoPPvggIiMjERAQgGuvvRYnTpxodr/m4pTIVbh7jN5zzz3o3Lkz/Pz8EBUVhTFjxuDvv/9u1TZlWcacOXMQHx8PPz8/DB06FH/99ZfT+0XUwN3j+ujRoxg7diyioqIQHByMm266CdnZ2a3aJuOammOBbFdRo6ysDCkpKUhJSQEApKamIiUlBenp6QgODsaQIUPw2GOPYcuWLUhNTcX777+PVatWYezYsa1whPbhBFAbFBgYiIiICGd3Q1c1NTW48cYbcd9996na78orr0RmZqa1rF27tkX9mDNnDiZPnixc/8knn2DatGl48sknsW/fPlxyySUYNWoU0tPTrdt07doVS5cuxf79+7Ft2zYkJSVhxIgRyM3NbVHfiLRw9/PFtGnTsHr1anz88cfYtm0bysrKcPXVV8NsFmeCsCdOiVyFu8do//79sWLFChw8eBDr16+HLMsYMWKEYow2Z/LkyYov61y0aBFeeeUVLF26FLt27UJsbCyGDx+O0tLSVu0Xkb3cOa7Ly8sxYsQISJKEH3/8Edu3b0dNTQ2uueaaFt31wLimlrJAhrmZonYCaPfu3UhOTkZycjIA4OGHH0ZycjKeeeYZAMDHH3+MgQMHYsKECejRowdeeOEFPP/887j33nt1Pz67yeRWhgwZIj/44IPyY489JoeFhckxMTHy7NmzG21TVFQk33XXXXJUVJQcFBQkDxs2TE5JSbGunz17ttynTx/r17W1tfKDDz4oh4SEyOHh4fKMGTPkiRMnymPGjLG73cTERBmAtSQmJrbOB9CMFStWyCEhIXZtO2nSpEbHaEtzn+XZZs+eLU+aNEm4/vzzz5fvvffeRsvOPfdc+fHHHxfuU1xcLAOQf/jhB8W+Ep2trZ8vioqKZC8vL/njjz+2Ljt58qRsMBjkdevWCffTEqdEWrT1GLXl999/lwHIR44csS47ceKEfNNNN8mhoaFyeHi4fO2118qpqanCOiZNmtTkc2xgsVjk2NhY+YUXXrAuq6qqkkNCQuTly5er6heRLW09rtevXy8bDAa5uLjYuqygoEAGIG/cuNG6jHFNjtLwu9ThgzFy5ok4xXL4YIwMoNHPr6fhHUBuaOXKlQgICMDOnTuxaNEizJs3Dxs3bgRQf/vj6NGjkZWVhbVr12LPnj3o168fLr/8chQUFNisb+HChfjwww+xYsUKbN++HSUlJTZvD1Vqd9euXQCAFStWIDMz0/q1LT179kRgYKCw9OzZs4WfkP22bNmC6OhodO3aFXfddRdycnKs67R8lkpqamqwZ88ejBgxotHyESNG4JdffhHu89ZbbyEkJAR9+vRR3SZRWz5f7NmzB7W1tY1iLj4+Hr169VKMObVxStQSbTlGz1ZeXo4VK1agY8eOaN++PQCgoqICw4YNQ2BgILZu3Ypt27YhMDAQV155JWpqauyuu0FqaiqysrIaxbiPjw+GDBkijHFb/SJS0pbjurq6GpIkNcqa5OvrC4PBgG3btgFgXJNzWOwsno5ZwNxQ7969MXv2bADAOeecg6VLl2LTpk0YPnw4Nm/ejP379yMnJ8d64n3ppZewZs0afP7557j77rub1Pfvf/8bTzzxhPVZxKVLl9p8FEqp3aioKABAaGgoYmNjFfu/du1a1NbWCtd7eXnZ8Sm03KhRo3DjjTciMTERqampePrpp3HZZZdhz5498PHx0fRZKsnLy4PZbEZMTEyj5TExMcjKymq07Ntvv8XNN9+MiooKxMXFYePGjYiMjGzZAVOb1JbPF1lZWfD29kZYWFij5bZiroGaOCXSQ1uO0QZvvPEGZsyYgfLycpx77rnYuHEjvL29AdTfPm8wGPDOO+9AkiQA9b/AhoaGYsuWLU0ma5vTEMe2YjwtLc3ufhEpactxPWjQIAQEBGDmzJmYP38+ZFnGzJkzYbFYkJmZCYBxTc7R8JhXc9t4Ok4AuaHevXs3+jouLs5658qePXtQVlbW5LnhyspKHD16tEldxcXFyM7Oxvnnn29dZjQa0b9//ybP6Sq1q0ZiYqLqfRrce++9+M9//mP9uqysTHNd48ePt/67V69eGDBgABITE/Hdd99h3Lhxdn2WP//8M0aNGmVdV1NTA1mW8fnnn1uXzZo1C7NmzbJ+3TDQNZBlucmyYcOGISUlBXl5eXj77bdx0003YefOnYiOjtZ8vNQ2teXzhYitmDubPXFKpAfGKDBhwgQMHz4cmZmZeOmll3DTTTdh+/bt8PX1xZ49e3DkyBEEBQU12qeqqsr6GXz44Ye45557rOsa7kB46aWXrMvefPNNTJgwwfq1PTGu1C8iJW05rqOiovDZZ5/hvvvuw2uvvQaDwYBbbrkF/fr1g9FoBADGNTmFWa4vzW3j6TgB5IbOnnWXJMk6AFgsFsTFxWHLli1N9lPKjmXrhKmmXTV69uzZZDb+TImJiU3e2t9g3rx5ePTRR1W3aY+4uDgkJibin3/+AWDfZzlgwADrW98B4LXXXsPJkyexcOFC67Lw8HAAQGRkJIxGY5O7CHJycpr8xSIgIABdunRBly5dMGjQIJxzzjl499138cQTT+hwpNSWtOXzRWxsLGpqalBYWNjoLqCcnBxceOGFNvdRE6dEemjLMdogJCQEISEhOOecczBo0CCEhYVh9erVuOWWW2CxWNC/f398+OGHTfZruKPh2muvxQUXXGBdPnPmTLRr1w5Tp061LmuI34Y7H7KyshAXF2ddbyvGlfpFpKStx/WIESNw9OhR5OXlwWQyWe866tixIwAwrskp6iChFsp/zKtrZr0n4ASQh+nXrx+ysrJgMpmQlJTU7PYhISGIiYnBb7/9hksuuQQAYDabsW/fPvTt21dV215eXna9Rb8lt5VGR0e32l0w+fn5yMjIsA4c9nyWfn5+6NKli/Xr8PBwlJSUNFrWwNvbG/3798fGjRsbpf7buHEjxowZo9g3WZZRXV2t4aiIxDz9fNG/f394eXlh48aNuOmmmwAAmZmZ+PPPP7Fo0SKb+7QkTon05ukxKnLmmNevXz988skniI6ORnBwsM3tg4KCGt1JEBQUhPDwcJtjcceOHREbG4uNGzdas7bU1NTgp59+avTHm+b6RaRVW4rrhtcX/Pjjj8jJycG1114LgHFNzmGR60tz23g6TgB5mCuuuAKDBw/Gddddh4ULF6Jbt244deoU1q5di+uuuw4DBgxoss+DDz6IBQsWoEuXLjj33HPx73//G4WFhaofd0hKSsKmTZtw0UUXwcfHp8l7Nxq0xiMd6enpKCgoQHp6Osxms/WunC5duiAwMBAAcO6552LBggUYO3YsysrKMGfOHFx//fWIi4vD8ePHMWvWLERGRlp/6dPyWTbn4Ycfxm233YYBAwZg8ODBeOutt5Cenm5NBVheXo7nn38e1157LeLi4pCfn4833ngDJ06cwI033qjPh0V0mqefL0JCQnDHHXfgkUceQUREBMLDw/Hoo4/ivPPOwxVXXGHd7vLLL8fYsWPxwAMPAGg+TokcxdNj9NixY/jkk08wYsQIREVFWe+g9fPzw1VXXQWg/nGNF198EWPGjMG8efOQkJCA9PR0fPnll3jssceQkJCgqk1JkjBt2jTMnz8f55xzDs455xzMnz8f/v7+uPXWW+3uF5FWnh7XQP37fLp3746oqCj8+uuveOihhzB9+nR069YNAOOanMMMCeZm7vBpbr0n4ASQh5EkCWvXrsWTTz6JKVOmIDc3F7Gxsbj00kuFjy/MnDkTWVlZmDhxIoxGI+6++26MHDnS+pyuvV5++WU8/PDDePvtt9GuXTscP35chyOyzzPPPIOVK1dav26Y/d+8eTOGDh0KADh06BCKi4sB1D87vX//fqxatQpFRUWIi4vDsGHD8Mknn1j/2qDls2zO+PHjkZ+fj3nz5iEzMxO9evXC2rVrrQOt0WjE33//jZUrVyIvLw8REREYOHAgfv75Z4dmR6O2oS2cLxYvXgyTyYSbbroJlZWVuPzyy/H+++836m/DbeoNmotTIkfx9Bj19fXFzz//jCVLlqCwsBAxMTG49NJL8csvv1jv9vX398fWrVsxc+ZMjBs3DqWlpWjXrh0uv/xy4Z0DzZkxYwYqKytx//33o7CwEBdccAE2bNhgHf/t6ReRVp4e10D9NfcTTzyBgoICJCUl4cknn8T06dOt6xnX5AycAKonybYeIKU2zWKxoHv37rjpppvw7LPPOrs7ROTCeL4gcm2MUSLPw7gmsl9JSQlCQkLw05/tEBhkUNy2rNSCIb1Oori4WPNkpKvjHUCEtLQ0bNiwAUOGDEF1dTWWLl2K1NRU6y2TREQNeL4gcm2MUSLPw7gmajkzDDBDeQKo+TdkuT/lT4DaBIPBgPfffx8DBw7ERRddhP379+OHH35A9+7dnd01InIxPF8QuTbGKJHnYVwTtZwsS7A0U2SZj4AREREREREREbmdhkfANuxPREAzj4CVl1ow4rw0PgJGREREREREROSOzLIBZrmZR8DawK0xnAAiIiIiIiIiIo9VCwNqoZw5r9ZBfXEmTgARERERERERkcey7w4gz78FyOMngCwWC06dOoWgoCBIkue/1InIHciyjNLSUsTHx8NgUP8uesY1kethXBN5HsY1kedpaVy7KwskWKB8HmpuvSfw+AmgU6dOoX379s7uBhHZkJGRgYSEBNX7Ma6JXBfjmsjzMK6JPI/WuHZXFjvSwFvg/DuAwsPDVW0vSRL27t2LxMREu7b3+AmgoKAgAMDmnVEIDGz8DS+2eNvcp8jsJ6yv0BJoc3mxYJ/CugCby0sU2iiu9bW5vKzO9vISwfYVtV7CNipqbR97ZY3tfWqqbT8vaa4VP0dpqRGsq7EdeIZawfJq8UysoUawvNb2PqLtjdXCJmCssX0iENVlqhZsr/BQqbHaos/yqjphG4Zqs83lUrXtjklVgg+lSnDgAOTKSpvLLZWN66qTa/Fz3RprfKqlJa4LzbZjEQCKLP42lzsirkvqbO9TVutjc3l5je3jA8RxXVVr+1Qvius6wfYAIAviVxjXouVKcS2IFUON7X2MovOA+EdVGNeic4H4PCC+UDBVic4FtuPXVGU7RkWxCwAGQcxLVYK4rrZ9gLJCXKOqyuZiS0XjeK+Ta/Gz+esWx/WmHdEIOCuuSyy246HYIo65AovtOC0y2473ojrby0UxCgDFtYL4Ndvub2mN7eVlguUAUFVne1wWjde1wrgWj9eyYLyWBPErCWNRHNeSME7VxbXSeG2oVRfXJg1xbRTEtVEQ18ZKUVwrjdcq41owXivFtVxZYXO5rfF6m/xNi+N6446YJnFdLIjrAoXxWhTzhWbb1+ei+BXFOwCUCsbfEsF1uCh+ywRjMgBUCa7Rq2psj7+1guVm0bU2IByXRXEtikXRckA8lovGcWFcK43XoutqYV2i84BCXKu+DtcQ18LxWnCCqrZ9gHKl7TG5fp3gOvysNupQi21Yqzmu3VWtbEKt3Mw7gFwgDXxRURGWLFmCkJCQZreVZRn3338/zGbxteLZPH4CqOF208BAAwLPSvtWZ7F9Aqw1i38wqgXrasy2P8pKwYVbtVk8OeMtGDC8VC43KQw8RsEvkUaT7X4ZjLaPT1a4oIRRsM4oGHhEyxVuGRbdtWg0CAYkQVVKpwKjYCZYNH9sFDw7qtiGRTDAmFUuNykMPHWCgUfwuUuiA1SYOJcl2ycfi2S7v1pvB9cS16LYBYAad4prL/EvisK4FkzoiOLaILjQBABZ8PMijGtBkCrGtejUIdjHKKhKFO+AOK6NgutD8fYKE0CCVBIGWTABZBJcUCoM6gaj4BdClXGt+Ei8IH4tkuCc0sK4DrAR12bReG0Rx3WlIH6rBctFEy3VguWAQlzXqYtfxbgW/KIoGq/NwrhWmAASjNfinyPBGKsQdKLxWnQuEMW10lhqkFTGtWi8VohroyCujYK4NgrjWmm8bv24lgXxa3O8llsnrkXjdZXieC0Yl0XxLohfby1xrXJcVr4OVxnXJsF1uOhaGxCOy6KfI1H8KsW1aFwWXp9ruQ4XxKPox9skimuLQlyrvg4XxLUgdgHAIDrPCq+3bX9YsmBMrl9n53X46Y+irT2WaZYlmJuZ4GluvaPcfPPNiI6OtmvbBx98UFXdHj8BRERERERERERtl9mOR8DMLvAImEUwISlSWlqqantOABERERERERGRx7LIBliayQJmaQNZwNrOa7+JiIiIiIiIqM2phQE1slGx1LrY9MjKlSvx3XffWb+eMWMGQkNDceGFFyItLU1TnW3mDiAj5CbvcBC+00Hh2UqDYJ3wuXPB9kaI2zAK6hK1IVquRLSPJFwuqEjpMUmV+wgfudTShtrtFeqRRQcv+KyE2yvdUijsl6htQQsKz/KK1gn30PJcsOg58bOX6/R8rZq4NijEnCiuhe2q3L5+H1G/1MW1KEaV1qmNa6U2hGvUxpYjzh0KRPEgqz2fKsWcQVCXsG11yxXbF8Wi6DygdBzCts++SNLnoskoyU3iRTj2Ko2lgnWi5WpjUalf4u31G8eFcSpcrroJYTyIfl40/d1Uz3jXqS5NbYhoCQu146+9Y2+jdfa930mSJSiEmd0MkJvEmJbrcNFYqrYupetwEdXXzgp1qR6vhfWI23Cn8Vop5tT+fqDr7xOqf89Qf6KVBXEqrElDXFM9CwywNJsFzLU+w/nz52PZsmUAgF9//RVLly7FkiVL8O2332L69On48ssvVdfZZiaAiIiIiIiIiKjtMcsGmJt5BKy59Y6WkZGBLl26AADWrFmDG264AXfffTcuuugiDB06VFOdrnWEREREREREREQ6skCyq7iSwMBA5OfnAwA2bNiAK664AgDg6+uLyspKTXXyDiAiIiIiIiIi8lg1sglGWXn6o8bF3gE9fPhw3HnnnUhOTsbhw4cxevRoAMBff/2FpKQkTXXyDiAiIiIiIiIi8lgWWbKruJLXX38dgwcPRm5uLr744gtEREQAAPbs2YNbbrlFU528A4iIiIiIiIiIPJYFBpjd7CXQoaGhWLp0aZPlc+fO1Vxnm5kAktD0didRxgAlwiwDoqwiKrOGAUqZSPTLXKI+y4DK7GCKlSnsY4OWbBxqJ291zfihKcOBXhmJ1G1ev4+6bEGK3/QmWYFal624FlHMKqIyhkSxqJRVRG2castIpO7nRZxFSGkndcs1ZeMQUF2XrsehIUb1ylKoQLfsfkpx7cLZ/ZQ4NWunu2X3E54L1GWiVM6oqW4fp547lOoSDDri7EYaMnTpOS4LmxCd05yf3U8p3lVfb+s4luq1XInq+NUyBukYc6pP+3qO1yq315YxVGVdWrLxOiDem+4jaUzb6N4ssgGWZl7y3Nx6ZygqKsJvv/2GnJwcWCz/O9dJkoTbbrtNdX1tZgKIiIiIiIiIiNoeMySYm5lRbG69o33zzTeYMGECysvLERQU1GjCUOsEkOtNcRERERERERER6aRWNqBWNjZTXGt65JFHHsGUKVNQWlqKoqIiFBYWWktBQYGmOnkHEBERERERERF5LHd8BOzkyZOYOnUq/P39davTtY6QiIiIiIiIiEhHZtlgV3ElI0eOxO7du3Wtk3cAEREREREREZHHkiHB0sw7fmQXewfQ6NGj8dhjj+HAgQM477zz4OXl1Wj9tddeq7pOTgARERERERERkceqtRhhsBib2Uac+dMZ7rrrLgDAvHnzmqyTJAlms1l1nW1mAsgo1ZczCVPBKqWfFKWsVJnKUikltSiVpdrttaSfFKeoVV2VcmpKVdsrpVRUt4umVLBOTF2rNiW1cipnURsqO2ZQf2vk2Sku9ZpbtxXXwpjTEteCdLNq08bXr1OXulZcj4Z00cLt1S1XrEy4ve0+CX+GNbShOk2rEl1TT6tLHS+Oa4VG1Ma1pvSxthtprbiW0PSwtIylalPHq41FQH2KaT1TUouI00Ur7aSqCX3HUrXbK6akFv3c6xhzws9X3QEqngNFlFLH2yBMLw1AVllXSxkh2x2TBsHYCyiljtdvLNUrpbxyG7aXq0737oj06Yr7aBjjbW6vfp36FO1KbQhWir9R4spE1F4+65LuvWH52Y0b2mQaeDMMMDfzjWhuvaNZWmFCqs1MABERERERERFR22ORJVia+etgc+s9ASeAiIiIiIiIiMhjWWCApZk7fJpb7wzl5eX46aefkJ6ejpqamkbrpk6dqro+TgARERERERERkceqtRhgsChP8NQ2s97R9u3bh6uuugoVFRUoLy9HeHg48vLy4O/vj+joaE0TQK51hEREREREREREOpJlAyzNFNnF0sBPnz4d11xzDQoKCuDn54cdO3YgLS0N/fv3x0svvaSpTtc6QiIiIiIiIiIiHZkh2VVcSUpKCh555BEYjUYYjUZUV1ejffv2WLRoEWbNmqWpzjbzCJgB9s92KWUCUcpMoLYu8T6ibCf6ZBtR2keYRUhtVgIlajNoKLWhJYOHSnplGFLcXm3GAk3ZjVRmZBBkPlDKKiLMlnB25rBWnF0XZ+NQH4uiDEPCLCSa2lDXXy1xrTojoGLMqa1LfT2aMvOo3d6ZWY9UH4eOmUBElLL7OTpbkI3sfsJtmd3PxvbqlitXJlquX3Y/h4ylKrfX8g5Q1ecIDRmJhOOvnln/zo53nV6Iqld2PxHxtbP6rGGicVlt28Jr52bW2d5efT3ijJOic4r6sVfPn3shnerS9dwhGDKV2lAdv8IMZErjNe/tUGKRm3/Js8XFsqN5eXlZf3ZiYmKQnp6O7t27IyQkBOnp6ZrqbDMTQERERERERETU9tTJRhhkY7PbuJLk5GTs3r0bXbt2xbBhw/DMM88gLy8PH3zwAc477zxNdXKakIiIiIiIiIg8llmW7CquZP78+YiLiwMAPPvss4iIiMB9992HnJwcvPXWW5rq5B1AREREREREROSxGl703Nw2rmTAgAHWf0dFRWHt2rUtrtO1jpCIiIiIiIiISEcWSLDIzRSVL4TaunUrrrnmGsTHx0OSJKxZs6bRekmSbJYXX3xRxyNTh3cAEREREREREZHHktH8BI+scgKovLwcffr0we23347rr7++yfrMzMxGX3///fe44447bG7boF+/fti0aRPCwsKQnJysmIBn7969qvoLtPEJIGH2AcXMAKLMPCqziihlLhHUpbSPLYrZB1TVpD4LSf0+tpc7ImOA2jY0Pe7piP6q3F45w4HtlaqzjSjRsk8L2MrupyWuRfGrOuuflgxdGjIPiYizAqnM+qfUhvBnT9RfDRm9HJDxQ3UMOfHcoWdci7KKaIr3VsruZyuutWT3E47LzO6nbp1ObajN7ufczHviVeJMZxoyo4moDSVN47XtRs4+F7TmqC4er5WukdWN8eJMfQoxp3Jc1jKO65eNV9iEgzLMqluuJa5FMSS87tByHa4hq5dNWrJmqs7aqb4N6ax9JFmCyktcj1BnMUKyNPMS6GbWn23UqFEYNWqUcH1sbGyjr7/66isMGzYMnTp1Eu4zZswY+Pj4AACuu+46Vf2xR5ueACIiIiIiIiIiz9bwmFdz2wBASUlJo+U+Pj7WSRmtsrOz8d1332HlypWK282ePdvmv/XCdwARERERERERkceynH4ErLkCAO3bt0dISIi1LFiwoMXtr1y5EkFBQRg3blyL62oJ3gFERERERERERB5LzR1AGRkZCA4Oti5v6d0/APDee+9hwoQJ8PX1VdwuLCzM7kfzCwoKVPeDE0BERERERERE5LHqLAZIFuUHoOpOrw8ODm40AdRSP//8Mw4dOoRPPvmk2W2XLFli/Xd+fj6ee+45jBw5EoMHDwYA/Prrr1i/fj2efvppTX3hBBAREREREREReSw1dwDp7d1330X//v3Rp0+fZredNGmS9d/XX3895s2bhwceeMC6bOrUqVi6dCl++OEHTJ8+XXVf2swEkPF0OZPorf1K2QeE9avMoKElq4jaTCD6Zh+wXY9SpjHVWUU0ZcBR14SeWUX0yiKk3IbKA3TZbCNn76PPyVVNXCvWI8z+o+4coRTXSlnIbFEbo1r2UYxfEbVZN3TMKqIptlS2oW+2IHVVacso5oBzhIOz+9miJbufSFvN7qfYhtrsfpp+jlQu10DtuKwto6Zeyx2RLUhhgNeSragFjFJ9adwF/cZr4fYargnEmcP0GccBbVk4bdejtFLdcj3HWEdk/FV97azYhrpznTgzmUIjDhivRY8NyU2y/rXN1wDLgB1p4NUpKyvDkSNHrF+npqYiJSUF4eHh6NChA4D6F0p/9tlnePnll1XWDqxfvx4LFy5ssnzkyJF4/PHHVdcHtNXvPhERERERERG1CQ13ADVX1Ni9ezeSk5ORnJwMAHj44YeRnJyMZ555xrrNxx9/DFmWccstt6juc0REBFavXt1k+Zo1axAREaG6PqAN3QFERERERERERG1PncUA2PkOIHsNHToUsqx839Ddd9+Nu+++W1W9DebOnYs77rgDW7Zssb4DaMeOHVi3bh3eeecdTXVyAoiIiIiIiIiIPJYz3wGk1eTJk9G9e3e89tpr+PLLLyHLMnr06IHt27fjggsu0FQnJ4CIiIiIiIiIyGPJsgS5mQme5tY7wwUXXIAPP/xQt/qc+g6gZcuWoXfv3tY0a4MHD8b3339vXT958mRIktSoDBo0yIk9JiIiIiIiIiJ3YoFkV3FVlZWVKCkpaVS0cOodQAkJCXjhhRfQpUsXAMDKlSsxZswY7Nu3Dz179gQAXHnllVixYoV1H29vb6f0lYiIiIiIiIjcjzs+AlZRUYEZM2bg008/RX5+fpP1ZrNZdZ1OnQC65pprGn39/PPPY9myZdixY4d1AsjHxwexsbEtbssoSTCelRpPbep2QH0qWmEqSS1pZYVpKbWkuBS1obIixTTwKperrQcKqfqE/VKfGlL1eUBLWlmVdQm/UQr39InaF3ZLlH5S6YekSZrJhn3OWi7rc/OhrbgWbqsQc2rTxKpNBVu/j7q00GqXK7dte7kwU6qm9PDqlmuLB9v90pby1TZNaaFF1KaoVfuNUqI2xBRT1woqa7KPPhdNxtPlTEpjprAeQZyKxn7h+K4Q72rT0IviV1vqdnV1afruOCCVs9p00bqOpY5Ib622HkB16mnRuUOUElpRK43XBjQ9LYliTmmMNQhiThSnorqU4loUp6JzhyOuCbR8KxWv0VVsr5zaXF0TmsZYR1yHq91H07lDEKcqx36luJZV/yLXtpgtBkjNvOTZrPIl0K3tsccew+bNm/HGG29g4sSJeP3113Hy5Em8+eabeOGFFzTV6TJHaDab8fHHH6O8vNz6hmsA2LJlC6Kjo9G1a1fcddddyMnJcWIviYiIiIiIiMidNLwDqLniSr755hu88cYbuOGGG2AymXDJJZfgqaeewvz58zW/F8jpL4Hev38/Bg8ejKqqKgQGBmL16tXo0aMHAGDUqFG48cYbkZiYiNTUVDz99NO47LLLsGfPHvj4+Nisr7q6GtXV1davtT4bR0Sug3FN5HkY10Seh3FNRK5KtuMRMFebACooKEDHjh0BAMHBwSgoKAAAXHzxxbjvvvs01en0O4C6deuGlJQU7NixA/fddx8mTZqEAwcOAADGjx+P0aNHo1evXrjmmmvw/fff4/Dhw/juu++E9S1YsAAhISHW0r59e0cdChG1EsY1kedhXBN5HsY1EbkqGYAsN1Oc3cmzdOrUCcePHwcA9OjRA59++imA+juDQkNDNdXp9Akgb29vdOnSBQMGDMCCBQvQp08fvPrqqza3jYuLQ2JiIv755x9hfU888QSKi4utJSMjo7W6TkQOwrgm8jyMayLPw7gmIldllg12FVdy++234/fffwdQf35944034OPjg+nTp+Oxxx7TVKfTHwE7myzLjW4dPVN+fj4yMjIQFxcn3N/Hx0f4eBgRuSfGNZHnYVwTeR7GNRG5KossQXKzLGDTp0+3/nvYsGH4+++/sXv3bnTu3Bl9+vTRVKdTJ4BmzZqFUaNGoX379igtLcXHH3+MLVu2YN26dSgrK8OcOXNw/fXXIy4uDsePH8esWbMQGRmJsWPHqm7LVvYB4bYasg+IiDKEaMk+IO6T+mxBqrOH6JihRFa5XFJ61b4D3trvkGxmemYPEVGfCkq/NhxIzywdwvgVZf3TMbuf2nqaW2eLMK61ZPcTEMa1UhYLdU1oilHV5whHZCQS0ZLdT21ca4ldZvdzSHY/tVk7NWX3E54LRNvbXqxnJh9N47gj2hDWJcrcpb4u1eOvlvFaU1op7Wxl9xNvq5B5T5hdVzReC+JHw4Meaq+3ldoQ/uipvQ7XlLVTVJfo5CFuQ7csnFriWtfrcNHK1s8+6pB4P3vAcLFJDkdpeMyruW1cRW1tLUaMGIE333wTXbt2BQB06NABHTp0aFG9Tp0Ays7Oxm233YbMzEyEhISgd+/eWLduHYYPH47Kykrs378fq1atQlFREeLi4jBs2DB88sknCAoKcma3iYiIiIiIiMhN2JPly5VeAu3l5YU///xT8Y+mWjh1Aujdd98VrvPz88P69esd2BsiIiIiIiIi8jRmiwGwKN+tbG5mvaNNnDgR7777Ll544QXd6nS5dwAREREREREREenF3R4BA4Camhq888472LhxIwYMGICAgIBG61955RXVdXICiIiIiIiIiIg8Vv0EUHOPgDmoM3b6888/0a9fPwDA4cOHG63T+mgYJ4CIiIiIiIiIyGO52zuAAGDz5s2619mmJ4BEmQEUs4qIsg8Il6ufRhRnDtMvq4iIaB9xVgLVTah/m78D2tCSuUSc8UNLRgZBXQZ1mRrEWQzUt63b9mg6Q62Y2U0FW9n9tMS1KHuI2qxAStn9xPuoyx6iZ3Y/cT1K69Rl8RMGkKaYU9hH7fYOyBYk2keULEvX7H4ibSi7nzgrELP7NV2nqgl9s/vpeO4QxpwTzx2q07UptC+pfTWF4slcUFmTfVov/kWxpTaDHyA+R4i3V8g0pjLjr9prZyWqs/Qq1mV7uTgWVS7Xuo9KauNa1+ycqjONKewgil89s/EK4rrpdXjb5I5p4FuDa73liIiIiIiIiIhIT7KdxcnGjRuHkpISu7efMGECcnJy7N5e8wTQzz//jH/9618YPHgwTp48CQD44IMPsG3bNq1VEhERERERERHp6/QjYEpF2y1k+vrqq6+Qm5uLkpKSZktxcTG++eYblJWV2V2/pkfAvvjiC9x2222YMGEC9u3bh+rqagBAaWkp5s+fj7Vr12qploiIiIiIiIhIV+6SBUyWZXTt2rXV6tc0AfTcc89h+fLlmDhxIj7++GPr8gsvvBDz5s3TrXNERERERERERC3hLi+B1vLi53bt2tm9raYJoEOHDuHSSy9tsjw4OBhFRUVaqiQiIiIiIiIi0p1skSBbmpkAama9IwwZMqRV69c0ARQXF4cjR44gKSmp0fJt27ahU6dOevRLd0ZIMJ71znNdsw+o3EdthhBAOROJ2u1VZxURLdcxq4ieGbpEtGTyUd22yqxhp9eqa0NEMTuKyso0ZP4RZYCRz65Lp9l1W3Gt5e1tomwgorjWkvlHbYYhPbP4qc/upyGunRhberYhjhNR9jOlLDuiFczup0Sv7H7C+pndz+66VGf3U+yA7cXOHJe1tKFbf7W8hVNttiAtGfwMZ3VMlLZQJaMkwXhWP9Vm7gKUzgUqs/Rqyu6nXzZe9dn9RMuVxmudri2V6BTXjsjc5dTsYIDqcbbJtXNDNWfH6JlcIGunS7PnJc8u8AhYa9N0Vr/nnnvw0EMPYefOnZAkCadOncKHH36IRx99FPfff7/efSQiIiIiIiIi0qS5F0Db84iYJ9B0B9CMGTNQXFyMYcOGoaqqCpdeeil8fHzw6KOP4oEHHtC7j0RERERERERE2rWBO3yao2kCCACef/55PPnkkzhw4AAsFgt69OiBwMBAPftGRERERERERNQi7vIOoNameQIIAPz9/TFgwAC9+kJEREREREREpDMJzb/AyfUmgOrq6rBlyxYcPXoUt956K4KCgnDq1CkEBwdrugHH7gmgcePG2V3pl19+qbojRERERERERES6c8OXQKelpeHKK69Eeno6qqurMXz4cAQFBWHRokWoqqrC8uXLVddp90ugQ0JCrCU4OBibNm3C7t27rev37NmDTZs2ISQkRHUniIiIiIiIiIhahWxncSEPPfQQBgwYgMLCQvj5+VmXjx07Fps2bdJUp913AK1YscL675kzZ+Kmm27C8uXLYTQaAQBmsxn3338/goODNXWktUmSBIOWVLpnEaafFKaRVpdeGlCfYlpTSnmHpJUVrVDVhDiNJRTSI+uVHr6Zdara1rENPVNZij5DUUp3bd90xxHNaivFnNpU0qK4VqI25auW1LUiwjTSqmvSENcOSMfq1DTSCnTrl+K5Q13qeElDiljRuaBJilqdMmcYIcF41kFr+bkXj8tqUzmL412YelpDuncRtSmm1Y7j9Tu18nId29byYyaOE/3SZIuypWsKC7VxqmXsdXC6aAPs/6uzUswZVI6/eo7Xaq+3taSBV/1tUUwDr89ybeO17X5pum7X6xrZAdfhep6fhFVpivezok90wvJw7vgOoG3btmH79u3w9vZutDwxMREnT57UVKem7/57772HRx991Dr5AwBGoxEPP/ww3nvvPU0dISIiIiIiIiLSnRveAWSxWGA2m5ssP3HiBIKCgjTVqWkCqK6uDgcPHmyy/ODBg7BY1M+0ExERERERERG1Clmyr7iQ4cOHY8mSJdavJUlCWVkZZs+ejauuukpTnZqygN1+++2YMmUKjhw5gkGDBgEAduzYgRdeeAG33367po4QEREREREREelNkpWfmmzYxpUsXrwYw4YNQ48ePVBVVYVbb70V//zzDyIjI/HRRx9pqlPTBNBLL72E2NhYLF68GJmZmQCAuLg4zJgxA4888oimjhARERERERER6c4i1ZfmtnEh8fHxSElJwUcffYS9e/fCYrHgjjvuwIQJExq9FFoNTRNABoMBM2bMwIwZM1BSUgIALvvyZyIiIiIiIiJqw9wwDTwA+Pn5YcqUKZgyZYou9WmaADqTO0/8GAUTfMpv7VeZfUBDhi61+xg1ZCtQm6FES+YS4T10et5bp2eWAWEbKpdroTJDl+rjVqL2TWBa2pDObkSf7AO2svtpiWsRYRYhldmFlPZRn1VEIa5VZxoTZRFSVY02WrL7CesS1SPexRHZ/cRtq8xIpET1cbhndj8RZvezsVzl9sqNqFuuGFcOyOSjV4Y9xXOQcJXomkddpj7FdarPjS3P+ifperHTmDCzroZ4EGfKVZe9V4nwelvP+NUzu5+wcdE1pKhthZ9VXToETRm6dM2gq+VcYLMeDWOp6JJY0/W2647XLsFNJoC+/vpru7e99tprVdevaQKoY8eOiieDY8eOaamWiIiIiIiIiEhfbjIBdN111zX6WpIkyLLcZBkAmxnCmqNpAmjatGmNvq6trcW+ffuwbt06PPbYY1qqJCIiIiIiIiLSnz1ZvlwgC9iZWdV/+OEHzJw5E/Pnz8fgwYMhSRJ++eUXPPXUU5g/f76m+jVNAD300EM2l7/++uvYvXu3po4QEREREREREelNstSX5rZxJdOmTcPy5ctx8cUXW5eNHDkS/v7+uPvuu3Hw4EHVderzIo7TRo0ahS+++ELPKomIiIiIiIiI2pSjR48iJCSkyfKQkBAcP35cU526TgB9/vnnCA8P17NKIiIiIiIiIiLNJNS/A12xOLuTZxk4cCCmTZuGzMxM67KsrCw88sgjOP/88zXVqekRsOTk5EYvgZZlGVlZWcjNzcUbb7yhqSOtzQgJxrO+pXrOfhkE2UOEy/XM0KVyeyVqs4doyUogen+4MAOOprfgq1yugdqsBJoeKdUrUxHgkOwhMOg6p9wsW3GthTgzj34ZusT7iDKUqI9rvc4FinGtV3Y/xZ9VHbPcqW3fAVlF1G7P7H7M7qeuLvXZ/XQLLS3Z/XSMOU2Zw1TSrQ0dMwhKBvVZw0TJXOSz69LpfRi2r8N1jF+VdSnFnNqMgHpm0BXHr5a4Foylwh1ULldYp2ssqq3LAVl6IYo5BaLPRFiTpqydgnG4yT6uNs3hIG7yDqAzvffeexg7diwSExPRoUMHAEB6ejq6du2KNWvWaKpT0wTQmDFjGg0cBoMBUVFRGDp0KM4991xNHSEiIiIiIiIi0p3ldGluGxfSpUsX/PHHH9i4cSP+/vtvyLKMHj164IorrlDMyq5E0wTQnDlzNDVGRERERERERORIDY95NbeNq5EkCSNGjMCIESN0qU/TBJDRaERmZiaio6MbLc/Pz0d0dLSmfPRERERERERERLqTofAM5BnbqLB161a8+OKL2LNnDzIzM7F69Wpcd911jbY5ePAgZs6ciZ9++gkWiwU9e/bEp59+an2k62yvvfYa7r77bvj6+uK1115TbH/q1KnqOgyNE0CybPuTqa6uhre3t5YqiYiIiIiIiIj01woTQOXl5ejTpw9uv/12XH/99U3WHz16FBdffDHuuOMOzJ07FyEhITh48CB8fX2FdS5evBgTJkyAr68vFi9eLNxOkqTWnwBqmIGSJAnvvPMOAgMDrevMZjO2bt3KdwARERERERERkcuQLBIki/J7c5pbf7ZRo0Zh1KhRwvVPPvkkrrrqKixatMi6rFOnTop1pqamYuvWrbjwwguRmpqqqj/2UDUB1DADJcsyli9fDqPRaF3n7e2NpKQkLF++XN8eOoFixg/Bm6FEWQm0tKG2bRGl7ANqs3dpyQKm+rVUOmYfENGSlUD1y+Cd2YYC8bGrq6xJhpAzqxKtOHsflSdXNUR5iLTEtbgNDdm+1GYo0XCOEMapyu0V2xDWZXu5MLuf0g+xEzN+6JnFzyHZgpjdT8idsvsp1uVO2f0UO6Ayu58jMvnoeE4RHYf4HKhA9TlCSxY/x2a6sZXdT0R5vFb3eaod3wEtWXdF5xqFrJ1q41e0XDGuVS7XQqfY0jPxkiOu9YXbaxku1caihgxkTcZx2bHjustQcQdQSUlJo8U+Pj7w8fFR1ZzFYsF3332HGTNmYOTIkdi3bx86duyIJ554osljYmcbNmyYzVfu6EHVBFDDDNSwYcPw5ZdfIiwsTPcOERERERERERHpRc1LoNu3b99o+ezZs1UnwsrJyUFZWRleeOEFPPfcc1i4cCHWrVuHcePGYfPmzRgyZIhwX9Erd/Sg6R1Amzdv1rsfRERERERERET6U3EHUEZGBoKDg62L1d79A9TfAQQAY8aMwfTp0wEAffv2xS+//ILly5crTgAB0JzmvTl2TwA9/PDDePbZZxEQEICHH35YcdtXXnmlxR0jIiIiIiIiImoxC9DsU92n1wcHBzeaANIiMjISJpMJPXr0aLS8e/fu2LZtW7P7P/300/D391fcRsu8i90TQPv27UNtbS0AYO/eva02I0VEREREREREpJtWyAKmxNvbGwMHDsShQ4caLT98+DASExOb3X///v2KGda1zsfYPQF05mNfW7Zs0dQYEREREREREZEjqXkHkL3Kyspw5MgR69epqalISUlBeHg4OnTogMceewzjx4/HpZdeimHDhmHdunX45ptv7JpPWb16tfNfAt1gypQpePXVVxEUFNRoeXl5OR588EG89957unROT4bT/53JKNxWv6k/LRk/1GYUU5tJQGkf9VlFlNbplFVEy1v7nZnhQMc2dMtKoNi2KHOJ7eWKs83CbC6SfduppCautRDFoiiuFTOXiOrSMbuf2n20ZPcTEsa7aHv1TYi4anY/1W0zu5+Qp2f30xLXemb3E7ahNrufnneHaxrP1G2uKeufTucOYfYzJWoT92jKDnZ2I62XLcgo6J6WeBCP16KxV/14rTYjoNL2asdlLYTffjfL0KU+hlRmHFRqw6Du+6HYhtrsXZrOEXxCx9F2796NYcOGWb9ueFXOpEmT8P7772Ps2LFYvnw5FixYgKlTp6Jbt2744osvcPHFFyvW25pPW2k6q69cuRKVlZVNlldWVmLVqlUt7hQRERERERERkS5kO4sKQ4cOhSzLTcr7779v3WbKlCn4559/UFlZiZSUFIwZM6b5rrZiFjBVE0AlJSUoLi6GLMsoLS1FSUmJtRQWFmLt2rWtcpuSHnbtrXZ2F4hIZwf/rnF2F4hIZ6lpdc7uAhHpLCvL7OwuEFEbJ8n1L4FWLK0376LKihUrEBIS0ip1q3oELDQ0FJIkQZIkdO3atcl6SZIwd+5c3Tqnp+tvy8W8WaF48J4gvsCayEOMHp+Dl58Px8SbAxjXRB5i9PU5eP2VcIwZrZz5gojcx1Vjc/H2G+EYdqmvs7tCRG2Vg18C3RKTJk2y/vvw4cPYsmULcnJyrKnlGzzzzDOq61Y1AbR582bIsozLLrsMX3zxBcLDw63rvL29kZiYiPj4eNWdcASzGXjy2SL88ls1lr0SgbDQ1nummYgco7oaeODRAmzfUY3FC8IQ4M+4JnJ3pWXAxLsLcO8dNXj2qRB4e3Nyl8jdFRTKuOHWfDw6LQgzpgfBKHoBEBFRK2mNl0C3trfffhv33XcfIiMjERsb2+gP3pIktf4E0JAhQwDUv926ffv2MBjc55etBc+EYvYLRfhufSUuuTITq96MROfzOPgQubPHpwdj0asl+Ojzcuz7owYfvBWJuM7O7hURtcS9dwZg+TvlWP5uGXbtrcaKZREIaefsXhFRS9xykx8++rQSLy4uxW+7avDm0jD4RDi7V0TUprjRHUANnnvuOTz//POYOXOmbnVqmsFJTEyEwWBARUUF/v77b/zxxx+NiiuaeGsgNn0di46JJqRlmDH8umyseL+8VV+wRESt68F7gvHdp9GIjTHi78O1GDIqC1982fQF9UTkPp58NASfvB+B0FAJe/bV4tIrs/HDxipnd4uIWmDBvFC8+e8wBPhL+GlbNYaMyMGOX/l+TiJynGbf/3O6uJLCwkLceOONutapaQIoNzcXV199NYKCgtCzZ08kJyc3Kq6q73ne2Pp9LK4Z5YeaGuDpp0tx3/3FKC21/ztthKyqiBgki7hAtl0kURHVJd5H3C/RPrBZNJFUFg2EL3VvuPfvrCJLEBZhGwr7qC1qPytZkmwWRZJku4gYBEXPNnRw8WBfbF8fi6EX+6CiUsZDDxVjxsxiVFY2/jkXxpXKmNYa1+J9lGJbnyJu23YRfRvri2yzqKYQ88I4UXmOUIw5nc5DmtpQSdu5w/Y30CHnjhYywICrhgdg+/o4DEj2RlGRjDumFGHB86Ww1MowAtaihVGyCIpsuyidD0R1wXYRHrNC/KqNOdH2igUqf1S1xJBOMafl3KE6fhToFu+Kx68y5rTEqMFgu7QSIyTcPC4Am9dGo3s3E7JzLJhwcyHe+HcZYJHtvswQXzvbLqrr0TBei84diseh+lpfdH7Q8O0XXAur3l7h+lnXa3q1bejYtm4xCoVrd73iXWGfhnf4nlnapFbIAtbabrzxRmzYsEHXOlU9AtZg2rRpKCwsxI4dOzBs2DCsXr0a2dnZeO655/Dyyy/r2kG9hYYY8OHbkXjjnVI89VwRvv22Cn/+WYu33gxFjx5ezu4eEWkQHWXEmv9GY+GSYrywuAT//W8lUlJqsXx5KDp11HSaIyIn65BgwvovY/DM80V4/Z1SLF9ejj17avD666GIj9M6/UNEztTtHC9s+jYaj84qwn8/q8CLi8rw2281WPJqKMLD3efVEkTkftzlHUCvvfaa9d9dunTB008/jR07duC8886Dl1fj+YqpU6eqrl/Tb0Y//vgjvvrqKwwcOBAGgwGJiYkYPnw4goODsWDBAowePVpLtQ4jSRL+765gdE824P77i3D8uBnXXJuPuXODceMNfgDngYjcjtEoYdYjoTivvxEPTi3GgQN1uOqqfCxaGIxRo3y13y5ARE7j7S3hhblh6D3QiEcfLcauXbW48so8LF4ciiGXemu8j5mInCnA34BlS8KRfIEJTz1Zgp+21OCqK/Pw6muhGDDQS/vdGkREStzkHUCLFy9u9HVgYCB++ukn/PTTT42WS5LkuAmg8vJyREdHAwDCw8ORm5uLrl274rzzzsPevXu1VOkU/ft5Y933kZg2vRibNlXj8cdL8PjjJQgOlhARZUREpOF0Of3vKAN8w+sQHmlEeJQR4ZFG+Pjy6pPIlVx6qQ/Wr4vA/z1QhJ07a3H//xUDKEZIqNQolq3/jjTAL8KM8EgjwiLr49rbh1efRK5k9FW+6NnDhHvvK8Kff9Zh0qRCSBIQGmZA+FnjdMO/fc6I67AII7yYTYzIpdw03h/n9fbC/fcW4dgxM266sQAGAxAaXh/HDbEdHmlARFT9/30jzAiLNCEs0ojQcCNMXoxrIrKPPe/4cYV3AKWmprZq/ZomgLp164ZDhw4hKSkJffv2xZtvvomkpCQsX74ccXFxevexVYWFGbDivVAsf7Mcr75ajvJyGSUlMkpK6pB6tPn9/QMlhEcaERzhjdBII0IjTQiLNCEkwoTQSBMMYTJCIr0QHOEF3wBD233mksiBYmON+OTjcLz8chneersc1dVAcZGM4qI6HDtia4/iRl8FBhkQFmlEcKTXGTFtRNjpuJbCZYREeCEk0gs+/ry1iMgRkpJMWLM6As8+V4r//rcCtbVAYYEFhQUWHD1c1+z+QSH1cR0UUR/XoafjuaGYwmUEn45rb/5xh8ghunf3wjffReDpp0qwZnUVLBagIM+Cgjz7fgsLDm2Ia2+ERpoQEmlC6OkYDzl9HR4c6Y3gSC94eTOuido0N7kD6Ezz5s3Do48+Cn9//0bLKysr8eKLL7Z+GvgG06ZNQ2ZmJgBg9uzZGDlyJP7zn//A29sbK1eutLueZcuWYdmyZTh+/DgAoGfPnnjmmWcwatSoJtvec889eOutt7B48WJMmzZNS7eFDAYJ998XiPvuDUBRkYxjOTLy8yynixn5uf/7d06ujII8MwrzLKipkVFRJqOirA443vzFp7evAcGRXgiJMMEvwg/BEV4IivBGUKR3/f9P/1sO8YJvkImTRUQtYDJJmDkzCI89Fvi/uM49HdN5FuTlWZCfW//vnDwZhXlmFOSZUVcLlJVaUFZqAVJrm23Hx9+AkAgvBEd6wT/c96x49rL+2xzqDZ8AxjVRS/j6Snj+uWDMmxuEwkILjuRIKLCO0+b/jd25ZmtcF+abYa4DSostKC22AEebj2vfAIP1jzd+EX7/i+mzxmxLmBd8/PmeMaKWCAw0YPGSULz4koyCAguO5kgoOH3dXXBGbBecvg4vzDOjqMAMixkoKbKgpMgCHGk+rv2CjAiO8EJwpDf8zhivA62x7YPACG+YQ3zg5ce4JvI07vIOoDPNnTsX9957b5MJoIqKCsydO9dxE0ATJkyw/js5ORnHjx/H33//jQ4dOiAyMtLuehISEvDCCy+gS5cuAICVK1dizJgx2LdvH3r27Gndbs2aNdi5cyfi4+O1dBcAYJQkGM/6xevsryFJiAwHTKEAutqup8jiDQCQZRnlpfWTQQV5ZqTneKEor65+UMqrqy/5dcjPtaAkvxbVlRbUVFmQd6IaeSeqAZQr9tc3yIQJ/x6IpP7h1mVK2XxsMShmK1JXlyjjiGL2H7URJNheVsx2ou6XaWHWj+YydehBoR7V2UhE22tIzSZqW1iTlgkMyaD8tUb2xLXRKCEqAvAKA9DNdj0FZl8A9XFdWmKpnwzKNeNErgmFDfGcV4eigvr4LsgzozivDjVVFlRXWJBTUY2cjGoAZYr9DQj3xuS3L0Bs1+Bmj00pE4l4H9uxojbeFYnq0hJbwjZsLxYehfD49IsHrRmDbBK2oSW2VO4jCr3mMo7Ys9yBk5tGo4TISCNM4QaILmWKLT4AAItFRmmxxTpen8gxWcfowrw6FOfXoSjPjPw8M0ryalFbI6Oq3IKq8mpkpzUf18Exvrhj1YUIjf/fxZmesahXXcrfYlFc63fuUDv+OiLm9Dxvabq+UF2XymsehWsC4Zqz97E4Lq5NJgnR0Ub4RIrvrC06I65LCi2n/yhrRnqO6XQs18d28elxuyDPjJL8WphrZVSWmlFZakb28SoAJYp9iUgMwO0fXAL/UO9m+y263tZz7HXIdbiwHqV1omt3dReXTo13De3rG++2d5I0XNOLfi9ysTkN53HDO4BkWbb5ff39998RHh5uY4/m2T0B9PDDD9td6SuvvGLXdtdcc02jr59//nksW7YMO3bssE4AnTx5Eg888ADWr1/vUi+XliQJgcESAoMN6NDJCx3MgTa3KzLXXxBWVdRfWBbn16E4rxaZuRJK82pQml+LkvwalObXoDSv/t9VZWZUldZh2/tHG00AEVHrkiQJwSFGBIcYkdgZ6Gjxt7ldkdkfslz/S2Jxfi2K82pRkleLk7nG07F8Oqbz62O8NK8G1RVmlBfU4Nf/pGLsvD4OPjKitstgkBASZkRImBEdzwE6C+I63xwIWZZRWWa2xnRxfi2ycqT/xfJZ8V1TaUFJdhV2fZKG4dO7O/jIiNoug0FCaIQRoRFGoBuQpHAdLssyKkrM1pguOX0dXpb/v1guOyOu66otyE8rx+9fZ2DwxM4OPjIiajVuNAEUFhYGSZIgSRK6du3aaBLIbDajrKwM9957r6a67Z4A2rdvn13baX28wWw247PPPkN5eTkGDx4MALBYLLjtttvw2GOPNbojSEl1dTWqq6utX5eUKM/wO4qvvxG+HYyI7lD/dcPEkC0H/5axdOxW/PNzLkrzqhAU6eugXhK5JleMa0mS4BdohF+gEbGJ9THaSXABCgApv1bivSk78Oe6TFz1eE8+NkJtnqvGtX+QCf5BJsR19AOgPF7//G0xPntsL1K+PoHLHuwGo4nvGKG2zVXjOiDEhIAQE+I618d1gWC8lmUZmz7Mwffz9yNldToG3daJj24TeQh3egRsyZIlkGUZU6ZMwdy5cxESEmJd5+3tjaSkJOuciVp2/wayefNmTQ00Z//+/Rg8eDCqqqoQGBiI1atXo0ePHgCAhQsXwmQyqUpvtmDBAsydO7dV+uooMV2C0L5PKDJ+L8Lv35zExbfzrw/UtnlCXCcNCEd4B38UpFfgrw2Z6Hdde2d3icipPCGuz708Fv5h3ijLq8bRX3LR9dIYZ3eJyKncPa4lScJ5VyVg48t/IfdoKU79VYR2vcKc3S0i0oE7TQBNmjQJdXX17xi+4oorkJCQoFvdTv9TVbdu3ZCSkoIdO3bgvvvuw6RJk3DgwAHs2bMHr776Kt5//31VM+9PPPEEiouLrSUjI6MVe996Gn453LvmBGTZRX4SiZzEE+JakiRrXO9bc8LJvSFyPk+Ia5OXAb1HtwMA7Fvjfv0n0psnxLVvsBfOvbw+q/HvjGsizyHbWVyEyWTC/fffD7PZrGu9Tp8A8vb2RpcuXTBgwAAsWLAAffr0wauvvoqff/4ZOTk56NChA0wmE0wmE9LS0vDII48gKSlJWJ+Pjw+Cg4MbFXfU68o4ePkakHusDCf+KHJ2d4icylPiuu817SAZgON7CpCfrvwieCJP5ylxnXx6YvfQ5myUF9Y4uTdEzuUpcd33uvp3Nuz//gRqq/T95YuInMhNJn8aXHDBBXa/isdeLvcSClmWUV1djdtuuw1XXHFFo3UjR47Ebbfdhttvv111vYbT/zVeph8DbGfs0ZLJxyhZEBBkRK/hsdj3zSns+yoDSX1DFLbXL6uIYjYBlduLE8cIMgaoarmhMpXLNVD9pn89sx7plZUAUP8Dr+WZdw0ZC1yBKE5FcS1ibGb7kFg/dBkchX+252LfmhO4Ymo3XbP7iajNHqL2PKC0jyiLn2J2P5VB5JDsfg7IKsLsfvXsytppXa7f1Vpz8WuLATLiuwUhvnswTh0swZ/fncDgf3UUby861+iaHUx1VWKOGP+EbavMLqRj+4qxq3pcFmTl0XStIKpLkEVIQ12txdZ1uCjXl5ZxziiILdE5wthMGx3Pj0RInB+KMytx6MdM9LoqQfU5QimuheOvqhaauQ4X7mN7uXhc1vFnRcdrZPWxqL4N8fnGAecnezNw2uPsgUHXtGvuQ7LUl+a2cSX3338/HnnkEZw4cQL9+/dHQEBAo/W9e/dWXadTJ4BmzZqFUaNGoX379igtLcXHH3+MLVu2YN26dYiIiEBERESj7b28vBAbG4tu3QS5nD1Mv7EJ2PfNKfyx9hRGz+gOeDm7R0TUUsnXJdRPAH19Apf9X1dnd4eIdJB8XXucOvgX9qzOwKAJSc7uDhG1kGSQ0Ofa9tj65mHsW5OOXlfp9/4NInIOd3oHUIPx48cDQKN3IkuSZE0Pr+XxMKdOAGVnZ+O2225DZmYmQkJC0Lt3b6xbtw7Dhw93ZrdcRscB4QhL8EPhiUr89UMWEke55220RPQ/3S+LgV+IF0qyq3B0Rx6iB3VwdpeIqIV6j47HupcOIvtwKTIPlsCva0DzOxGRS+szpn4CKHVnHoozKyBFMSsvkVtzozTwDVJTU3Wv06kTQO+++66q7Y8fP946HXFRBoOEfmPaYdPrR7BnzUkkjuLdAkTuzuRtRO+r4rHzozTsXZOBKzkBROT2/EO80f2yGPy5PhN712Tgohlxzu4SEbVQWEIAkgZG4viuPPz+dQb63hHu7C4RUQu44x1AiYmJutfp9JdAk7J+YxIgScCxnfkoOsmXxhJ5goZsYH//mI3KYr40lsgT9BtbH9d/fHcKddV8aSyRJ+h7erxOWZMB2eJivxkSkToWO4sLOnDgANatW4evv/66UdHC5V4CTY2Fxfuh8wUROLIjH399fRwX3dfT2V0iohaK6x6M2G7ByDpUgr+/T0fyzV2c3SUiaqHOgyIRHONb/3jnT6fQbUR7Z3eJiFqo+xVxWDt/P4pOVuDE3ly0HxDt7C4RkUbueAfQsWPHMHbsWOzfv9/67h+g/j1AADS9A4h3ANlghKy+SBabRcQgycJytn5j618899fXaTBYzDDC0qjo0YY9+6ipR5GkY3FmGwKyJAkKhEV8HJLNIhtgs2jpFwyCorJPilkJ1NSvgsHmf7BZlOJXdbuSRVDEsXLmecFkkNH/unYAgL++Pt4kpo2wwCjJNotyv/SJU0mSFYq6b724ER2LBsJY1LENUV1qzwOK5w7RD7ywTxq+gWrPES5CGKewXURsxqeN8ddglJB8bcN4napu7IUsLirjWjl+bReIirAR29vLCkUYD9AxE6/KmFM9JsNR5w61Y6ygaBmvdTnBtx611+ci4nG88XW8r78B510ZDwA48JXtuNYy9qq93hbGrhbCmEfrFw3Uxq+WuBZywPHpNo5rife2prkU8C6YCv6hhx5Cx44dkZ2dDX9/f/z111/YunUrBgwYgC1btmiqkxNAbqDn5THwDTKhJLMCabtynd0dItJB36vjYTRJyDpQhJzDxc7uDhHpIPm6+gmgtF+zUZpd4eTeEJEekk8/BvbPphOoLqt1cm+ISCtJlu0qruTXX3/FvHnzEBUVBYPBAIPBgIsvvhgLFixolBlMDU4AuQEvXyN6j6p/oeT+NWlO7g0R6SEgzBvnDqu/lXz/V4xrIk8Q0SEAif3CIVuAv75hXBN5goQ+oYjsGIi6KjMObchwdneISCPJYl9xJWazGYGBgQCAyMhInDp1CkD9y6EPHTqkqU5OALmJ/qcfAzv840lUlfClsUSeoP/puwX+WpsOc62LjThEpEnDY9t/fpVmfVafiNyXJEnWu/v++uq4cztDRNq54SNgvXr1wh9//AEAuOCCC7Bo0SJs374d8+bNQ6dOnTTVyQkgN5HQKwQRnYJQV23B3xtOOLs7RKSDcy6KRGCULyoLa3D050xnd4eIdNBzRBy8/IwoyijDyZR8Z3eHiHTQ55oESEYJmX/koyC1xNndISINlF55Z8/r75zhqaeegsVS/0fi5557Dmlpabjkkkuwdu1avPbaa5rq5ASQm5AkCb2vSwIA/MHHwIg8gtFkQM+rOwBgXBN5Ch9/E7qNrH9nyJ9rjju3M0Ski6AoX3S8KBYA8OfXx53bGSLSxg3vABo5ciTGjRsHAOjUqRMOHDiAvLw85OTk4LLLLtNUJyeA3EjP0e3r//rwZyHeu2kTflj0O/7ZcgrVpXwkjMhdnTcmEQBwdGsWVt76IzYv3o9j27NQW8EXTRK5q15jkgAAf31zHP+d+CN+/vefSNuZjbqqOud2jIg063ltEgBg738O45M7NuPXN//CiT25MNeoT8NMRI7njncA2RIeHg5JkpCTk6Npf5PO/aGzaE0xbUtQpD/639wJuz88itzDxcg9XIw9/z0KybATUeeGI35ADNoNjEFc32h4+Sl/aw06TW8qJRZUnZ5SVJlCI5rTOqqpR0O/1G4vC1M0qv0MdUz1qKUuN001KYpTo45vghPFdURSEHpc1R4H1mYg60ARsg4U4beV/8BglBDdMwLtBsai3cAYxJwXBZOPUVNad9E+WlPE214h2kF1EwqN67RcQxuiGJW1XCmI6jLoeNWhNk27ltht0oY+32zD6f8aL7NNyxgrbFdDvIvOEe36RiBpcAyO/5qNU38U4NQfBdj57t8weBkQ2zsK8QNiED8wFtE9I2D0MjbTL9vHqHaMVdpe9O0X7aPpU1cZp8Jx2ZnnFD3bUBujEH8mwpqaS/fuQEZJgvGsNs/++n/L9YtrI9THtegaudOl8YjtFY6sPwtwcm8eTu7Nw684AKOPEbF9oxE/IBZxA2MR2S0cBpNB0xgr7JMwRsX7aE4Tr4LaOHXVuNbr94lmfjFSWZeWcVkwWkpnL2+j94DIdrzk2UUmgPz9/ZGWloaoqCgAwJVXXokVK1YgLq4+MVR2djbi4+NhNqufgOYEkJsZPqM3LrqzG9J25yHtt1yk/ZaLgrQy5BzIR86BfKSsOgCDyYDoXhGISE5AzIA4RPaMhtGH32oiV3XN/IEYOq0X0nflIm1XLtJ35aL4ZAWy/shD1h952PPunzB6GxDTOwoR/RIQ0z8eET2iYDC10QGcyMVJkoQbll2C4pPlSN+Vi4zdOUjfmYOy3Cqc2pONU3uygTf/gMnXiNjkaIQnt0d0v3iEdYuAwci4JnJFRi8Dbll1GYoyypGxKwcZv+UgY3cOKgqqcXJnJk7urH+Xn1eAV31c92uP6P7tENI5HJKGCT8i0pks15fmtnEBVVVVjRJJbN++HZWVlY220ZpogrMCbsg/3AfdR7RD9xHtAAAnTgAnd2fh1O5snNiVhbKsCmSl5CIrJRd/rdgHo7cRkb1jENM/DjH94xHePcrJR0BEZwuK9kPP0R3Qc3T9O4FOpJlxcnc2Tu7Owsld2ajIq8Sp3dk4tTsb+7EHJj8TovrEIqZ/PGIGxCP0nHAnHwERnS2kXQDOaxeA865LgizLOHGsFqd2ZeHknmyc2pWFqqJqnPg1Eyd+Pf2LY6A3ovrGIfr0eB3SiXFN5EokSUJYh0CEdQhE7+s7QZZlpB+uxqndWTi1OwuZe7JRU1qDjG0nkbHtJADAO8QHUcnxiO4fj+j+7RCUGOLkoyBqm+x5xMsdHgFrIGm8k5MTQB4gMMYf3UZ3QrfR9QNR6ckynNidjeO/5SF7zylU5Vcie/cpZO8+BWAPTP5eCO8dj8h+7RCRnICQLhGQ+BdHIpcS3C4Qwe0C0X1MZ8iyjKK0UpzcnYXjO/OQvfcUaoqrkbnjBDJ31GcF9AryRnifdojol4DIfgkI7BiueWAgIv1JkoTQxGCEJgajxw1dIVtkFBwrwqld2Uj7LRc5+zJRW1aDU9vScGpb/UvhfUJ9EZ6cgIh+CYhIboeA9qGMayIXIkkSwjqHIqxzKHqOPxcWswUF/xTi1K4spO/KRW5KJmqKq3FySypObkkFAPhG+CEsuQPCkxMQ3i8B/nGcECJyCHte8uxGE0BacQLIw0iShOCEIPRICEL86PMgyzJK0oqQvTsT2XtOIWdvJmpKqpGzIw05O+ovML2CfBDRtx0i+7WDX++O8E+M4AUmkQuRJAlhScEISwpGwpg+kC0yio4VIHt3JnL2nKr/xbG0BtnbUpG9rf4C0zvUDxHJ9RNC/r07wq8df3EkciWSQUJElzBEdAlD++uTYTFbUHQ4v36s3nMKub9nobqoCpmbjyBz8xEAgE9kACL7JSCiXzv49uoE31j+4kjkSgxGAyLPjUDkuRFIvNkPljoLCg/mImfPSeTsOYW8/dmoyq9E5g+HkPnDIQCAb0wQIvq1R3hyAnx6dYJPZJCTj4LIM0lmG69DsrGNK5AkqdF1+9lftwQngDycJEkISQpDSFIYut7Qo/4XxyMFSP8tB3l7T6Dg91OoLa1G1s/HkPXzMQA/wyvMH2F92yO0b3uEJreHX3yosw+DiM4gGSSEdYlAWJcInHtzr/oLzEN5SN+VWx/Xf2SipqjyjF8ct8A7MhBhye0R2rcDwpLbwzcm2NmHQURnMBgNCO8ehfDuUej+rz4w15pRcDAXJ3bVj9dFf2WhOq8cJzccwskNhwD8CN/YEIQmt0dI3w4I7dMePpGBzj4MIjqDwWRAxHkxiDgvBt0n94O5ug75f+YgY1cOCvZmoPhgNqqyS3Hy+wM4+f0BAIBvQhhC+3ZASN9EhPRuD69QfycfBZFncKdHwGRZRteuXa2TPmVlZUhOTobh9Iu+tb7/B2hDE0AGSDCc9Wp2o+BV7UoZsoTZglTeL6aUlUBtXaLMJTaXG4GIbmHw6hiPzuP7wlJnQfHhHOTtPWn9xbG2sAI5mw8hZ3P9XyZ8ooMQMrAzEqYMhSnAx+5+ibOK2F2FHY0IMpQIG9evaUdkMtC1DZXZjYQ0vMjw7BlrvWawWzuu1W6vKduIypHGVlwbvICoXpHw6doe50zoB3ONGUUHs5G39wTy951EwV/ZqMkrQ/bGg8jeeBAA4BsXgtDBXZEw6RIYfLwa1afrfULult1PLSdn91Mdv8zupzq7n1KMisZf0blAdB6y1YbB24CYPjHw65GIcyYNhLm6DoV/ZiF/74n6CaGD2ajKKkbW98XI+v5PAIB/h3CEXnQu4m69CAaTcnYxuzgi855ebUO/MVPLuUPteUtxe7VP5euaAdSg/LWTCa95VY6/jhivbbbrY0J0/3j49+4I3HEB6ipqULg/E/n7TqBg7wkUH85B1YlCZJ0oRNa3vwMA/DtGInToeYgZN0jdy6T1ytrpzOt2pfYdcL7R8zpcXJe6ymSFnwHhmrP3sbjnuN5ibvQS6BUrVrRa3W1mAohsM5gMCOsRi7AesTjnX/1RWmFA8cEsFO7LQOG+DJQczER1TilyvktB6V8n0e25G+HNW1OJXJrR24iIPvGI6BMP3A4UlRpQ/NcpFO3LQFFKBkoOZaEqsxhZX+5C2d+ncM6c6+EV7OfsbhORAqOPCZH9ExDZPwHdABQXA8X7T6IoJR2F+zJQdiQHFekFqEj/BeX/ZKHzrOtg9PN2dreJSIHJ3xtRFyQi6oJEAEBhvoyS/SdQlJKOopR0VKTmnS6bUXk0G4kPXw2DF399I9LCne4AmjRpUqvVzTMINWLwNiGsTwLC+iQAkwfDXFmLwpQMHHzpB1Qez8WB6R+g67M3wj+JmcSI3IXR1wvh/RMR3r/+ArOuvBoFu47j0OIfUHbgJA4+/B90e+4m+PB9IkRuw+TvjYgLOiLigo4AgNqSSuT9chT/vPYjincfw98z/4uuc2+EV1iAk3tKRPbyCvJFxIVdEHFhFwBATWE5cn86hGPLt6Bw6wHUFpWj81PXwxjg6+SeErkfyVJfmtvG07nWfZ3kcox+Xogc3Ak9Ft8G34Rw1OSW4uCjH6Lkj3Rnd42INDIF+CB6aDd0f2kCvCODUHWiAAemf4DyI9nO7hoRaeQV7Ie4K3uh2wu3whTsh4p/snDw4Q9QdaLA2V0jIo28wwLQ7rp+6DJ3PAx+3ij7Iw2HZnyAmrwSZ3eNyP00PALWXHGysLAwhIeH21W04AQQ2cUnNgTdX/kXAnu0g7msGoee/BT5Ww86u1tE1AL+SVHosfhf8EuKQm1hOQ4+9l8U7011dreIqAUCz41H91dug09sKKqzinDwkQ9QdvCks7tFRC0QnNwRXRf+C6awAFQdz8WhR1aiMi3X2d0icisNj4A1V5xtyZIlWLx4MRYvXoynnnoKADBy5EjMmTMHc+bMwciRIwEATz/9tKb6OQFEdvMK9sO5C8Yj7MJzINeacXTB18j7eqezu0VELeAdFYzuL09AUJ8OsFTW4PDTn6Nw8x/O7hYRtYBvu3B0f+U2+J8Ti7qSShx64iOU/nbY2d0iohbw7xyLbi9Pgk9COGrzSnH4sQ9Q/udxZ3eLyH3IdhYnmzRpkrVs374d8+bNw0cffYSpU6di6tSp+OijjzBv3jz89NNPmurnBFArM0oWm8URdRklWVgMksVmETFIMgySDJOvCV2fGoOYa5IBGch6dwOyVmyEJFus2zQU1RwwDStLtgsUinAftZppx2bR6fiazWhi87glm0X5GCXbxU0ZIdsuGuJaVJeI2hhV2ufs2Dw7Rk0BPuj27I0IH9IdstmCk0u+Qt7n2yCh6b6ib7Gu33rBuUB1/Grgdm3oFb8GSVxETUtSk6KHhux+ZxajoBgAYVEbc8J4h0VYhMfQTMyp+zxkm6U5XmEBOHfhrQgZ0AmW6jqkL/gUhev22P6RkWSbRTUtY6mOsaW2X5rGOWEbtk+CsgE2iyYqY1Tx5Ky2rhYy2PxPXeyqzdgJaBtL9breFsWuUvyKYrGh+MaG4NyXbkNA93Ywl1chY+5/ULL9T53iVzT2ioswttD6v1crXfOqLerbtn3ucEQbijzsOlxvklm2q7iS9evX48orr2yyfOTIkfjhhx801ckJIFJNMhqQ9H9XoP2USwEAeWt2IP3l1bDU1jm5Z0SklcHbhM4zr0HsDecDALI+2IxTb66DbG4Db8Mj8lBGP290mX09IoefB1hknFq2FtkfboHsAu84ICJtTMH+6Dr/FoQO7gq5zoyTL3+B/K93OLtbRK7PTe4AOlNERARWr17dZPmaNWsQERGhqU5OAJEmkiSh3fhBaD/9WsBoQPHPB5A65yOYy6qc3TUi0kgySOhw5zDE3zkCkID8tXuQtvALWKprnd01ItLIYDIiafpViBp/CQAg99OfcfLf30CuMzu5Z0SklcHHC51mjUXYVQMBGch+bwOyVmyAbHGx316JXIgE4Y1u/ysq69y6dSuuueYaxMfHQ5IkrFmzptH6yZMnN7lbetCgQXbXP3fuXDz++OMYPXo0nnvuOTz33HO4+uqr8cQTT2Du3Lkqe1uPE0DUImHDeqPjMzfD4OeN8v1pOPrESmYmIHJzkdeejw4zrofkZUTJjkM49vSHqCupcHa3iEgjSZIQc+tQxP/faMAgoWjT70h7/hOYK2uc3TUi0kgyGhB715WInng5AKDgqx04ufhL3pFPJNIKWcDKy8vRp08fLF26VLjNlVdeiczMTGtZu3at3fVPnjwZv/zyC0JDQ/Hll1/iiy++QEhICLZv347Jkyer6msDk6a9iM4QlNwJnRdMROrcj1GVloujM1ag45xbgdh4Z3eNiDQKvag7TCH+SHv+M1T8fQJHZ65Exzm3AGFRzu4aEWkUPqIfTGGByFj0Bcr2HkXqk6uQ9MzNgH+os7tGRBpIkoTIcRfBFB6EU//+GiU//4W6wjK0f3w84BXo7O4RuRR7Xi+r9vVZo0aNwqhRoxS38fHxQWxsrLqKz3DBBRfgww8/1Lz/2XgHEOnCr1MsuiyaDJ+ECNTmleLozJWo+Ou4s7tFRC0Q2CsRnRdOgldkMKpP5uPIjPdRdSzT2d0iohYIHtgVHZ+bCGOwP6qOZuLojBWoOZXn7G4RUQuEDu2NDk/fAoOfNyr+TMPxWe+jLr/Y2d0icimSRbarAEBJSUmjUl1drbndLVu2IDo6Gl27dsVdd92FnJwcxe1LSkoa/VupaNFmJoCMkqFJMUiSzWKUICyizB5qswyIMgwYhVl89M02IjyOFmQr8IkNQZdFk+DfPQHm8iqcfHYVSn/Z3zQrgdpsX809qGmztH5GsdbO6KXUhpbMALplO9GS8qmVMhLoFdd6UYprtfuItxdn91PLnnOHb4codHlxMnwTo1FXWIaMZ95DecqRppU1+1A1s/s5NLufsC5mFdEru5CWusSZ+vQb++3J6OXfrR06LZwMr9gw1GYXIeOpd1H1T4b92f10HOfUxq9iPOjVL8X2W38sVd2GlhhVkx3MheiVnVN8Pa+UKVCfbIBa2hY5M8aDkjshaf4kmMICUZ2eg5NPv42aEzk2vs06Zf1T7JjtomW8dkxcOyBLob0nWC3X1E66Dnc7FjsLgPbt2yMkJMRaFixYoKnJUaNG4cMPP8SPP/6Il19+Gbt27cJll12mOKEUFhZmnSQKDQ1FWFhYk9KwXAs+Aka6MgX5ofOztyL95a9Q/OshZC7+HHWFpQi/erCzu0ZEGnlFBKPzCxNxfP5nKN+fhpMLPkTsfWMQPLSvs7tGRBr5xEeg88LJSHv2Y1QeycSJOSsR9/CNCBzQzdldIyKN/DrFouPC25E277+oOZGPk8+8g7jHboFfj47O7hqR00myDKmZd/w0rM/IyEBwcLB1uY+Pj6Y2x48fb/13r169MGDAACQmJuK7777DuHHjbO7z448/Ijw8HACwefNmTe0q4QQQ6c7g44XEmeOQvvwHFK3bhdwV61GXX4Ko24ZDMrSZm86IPIoxwBcd59yCtMXfoXTbfmQtXY3a/BKEj7sEUlv9SxKRmzOFBiLpuYlIW/glKvb9g1OLPkL0naMROmKgs7tGRBp5x4Si4wu3I+3ZT1F1KB0nn1uFmAeuR9CFvZzdNSLnsifN++n1wcHBjSaA9BIXF4fExET8888/wm2GDBli89964QQQtQrJaED0nVfBFBmCvP/8gMKvf0VdfiliH7wO/LEjck8GLxNip46DKSIYhV9tR/5Hm1BXUILoKVcBMDq7e0SkgdHPG+1m3oLst75ByY/7kPPWt6grKEXE+GHO7hoRaWQK8kP805OQ/drnKP/tILJf/QzmwhKEjr7Q2V0jcpoz3/GjtE1rys/PR0ZGBuLi4oTb/PHHH3bX17t3b9V94G/i1GokSULE2IthCg9C1utfoXT7n6grKkPU9AkwBvg5u3tEpIFkMCDqthEwhQcj9/11KF6/C3WFpYj8v/EweHs5u3tEpIFkMiLmvjEwRQSj4LOfUPD5T6jLL0H47WMhmTi5S+SODN5eiH14PPJWrEXx+t+Qt3Id6vJLEHLTKN6RT22TPWneVaaBLysrw5Ej/3s3ZmpqKlJSUhAeHo7w8HDMmTMH119/PeLi4nD8+HHMmjULkZGRGDt2rLDOvn37QpIkyM09riZJMJvNqvoLcAKIHCBkSB+YQgNx6sVPUPnXcWTOfQexMyfCFBHi7K4RkUZhowfBFBGMrFe/QPlvf6Ou8D3EPHYbjEH+zu4aEWkgSRIix18GU3gwct7+FiWb96G2oAzRD90Cg6+2dx8QkXNJBgMip4yGKSIE+f/diKJvf0FtXhki77kekhd/DaS2RbLUl+a2UWP37t0YNux/d8w+/PDDAIBJkyZh2bJl2L9/P1atWoWioiLExcVh2LBh+OSTTxAUFCSsMzU1VV0nVGLkk0ME9OmM9s/ejhPPfYjajGyceuZNxD4+Cd7tY5zdNSLSKGhQDxiDA3Bq4X9R/U8GMue8hZiZE+EVHe7srhGRRqHDB8AUFoTMVz5D5e//IHPeO4iZMRGmUPHFKhG5LkmSEHbdJTCFByN72WqU7/gD5uJSRE3jHfnUxrTCHUBDhw5VvFNn/fr1quoDgMTERNX7qMH7/1qZARabRVtd6lK0i7ZX2kdEjyyIvh3jkLjgTnjFR8FcUILMOW+j6kDq/+rTkMVUNQekmVRKUau2OOQ4hKlgYbsocdN00eKUrILUzBriWq8Uz5qOQ+V5QJQi1lYJ6NkBHZ6fAmNECGpP5SHzmbdQnXpK+QMHoJwmHoLS+inlWz0NrUIbWtJLC/fR6zzghLSyRsnQpBgkyWYxShAW9e3aTiOtVNTWJd5eFha9zh1qUk8HDuiGhLmTYQjyR83xU8ic/SZqM/P+t4GWWHRA/KoeY3WMa13b0OscYSulu5NSuxsgNSlGQRFdfmj5pUWUNl6xr2pjTsO1tvAavRXTxgdd2gfxT9wGydcHVQdTkfXs26grKNZ2TdkaY6OAbqnj9SRqQ0PMCY+vDV2HO0rDO4CaK67m6NGjePDBB3HFFVdg+PDhmDp1Ko4ePaq5Pk4AkUN5RYcifu5d8OnWAZaKKmQueB9lv+53dreIqAV82kcjft498O4QC3NxGTLnvYOKP8TZDYjI9fmdk4D4uffCFBOOutxCnJr9Jqr+SXd2t4ioBfx7d0bs03fCGBqE2hPZyJzzJmoysp3dLSLHaLgDqLniQtavX48ePXrgt99+Q+/evdGrVy/s3LkTPXv2xMaNGzXVyQkgcjhjoD/iZt0O//N7AHVm5Lz2CYrWbnd2t4ioBUzhwYibfSd8e3aCXFWD7EUfoPSnvc7uFhG1gFdsBOLn3APvTu1gKatA1vPvoXzPQWd3i4hawCcpHnFz7jl9R34xMue9haqDrfvOESKXIAOwNFNca/4Hjz/+OKZPn46dO3filVdeweLFi7Fz505MmzYNM2fO1FQnJ4DIKQzeXoh56GYEjxgEACj44HsUfPQdZIu2x+OIyPkM/r6IfXwiAi7qDZgtyFv+JYq/2dxsFgMicl3GkEDEPXUH/Pp2hVxTi5xXPkTp5p3O7hYRtYApKgyxs++GT9dEyBVVyHphBcp/sz/1NJE7kmTZruJKDh48iDvuuKPJ8ilTpuDAgQOa6uQEEDmNZDAgYvJohN8yEgBQum4b8pZ/Arm2zsk9IyKtJJMJUfffgJBrLgEAFH+xAYWrvuLkLpEbM/j6IOaRfyFwaH9AllG4cg2KvtjAyV0iN2YM9EfME7fDf2BPoM6M/OX/RcmGn53dLaLWY5EBi6WZ4lrjWlRUFFJSUposT0lJQXR0tKY6mQWMnEqSJIReewmM4UHIXb4aFTv/QHZxGaKn/gsGZiYgckuSwYDwW0fCGB6MglVrUbZ5J8xFpYi4dzwMPt7O7h4RaSAZjYi8ayxM4SEo+vJHlHyzGebCEoRPHgvJZHR294hIA4O3F6Km3oyCD75D6YYdKPr4O5gLihF601WQDLxPgDyMBc2/ENzF/l5511134e6778axY8dw4YUXQpIkbNu2DQsXLsQjjzyiqc42PQFkFPwEKJ3umsscYC+ljEFq2zBq+EkVZQ6QVGYUUNpeTV1BF/eFFBCC3Nf+g+q/jyFr/puIfuR2mMJDmslco3K5FoK6FLN0uVEb6uvRUFGTiwjPzEigJcOf2qx8mrL4ieJddU0KbdhYFnrlYBiCQpH35ieo3HcAOYveQfT0STAGBWg7k6qNdx0PUBSLSjGqOn4dcT7Tcj3vphlEROOcUjYf9W2oq0tt1h6lfbScC0SE4/VZyyUJCLvxMhhCQlGwcg3Kt+2BuagEkQ9MgMHXR2PjttsWZrdyxLio5bpDJT3Hd9FnpaUJ6ay6zv7a2UTXyM68dhZvr9+5Rss4bu91uGQwIHzi1TCGhKHos+9RumEbzEUliLjjJkheJsVGhD/HThyvRZrNrqt2H7Xb6zWfpiUmJYPy122EPY94udojYE8//TSCgoLw8ssv44knngAAxMfHY86cOZg6daqmOtvmd59ckl/PLoiddbc1M0HWs8tQcyLL2d0iohYIGNgLMTPugCHADzVHM5D13HLU5hY4u1tE1AKBQwYiauptkLy9UPXnP8he8BbMRaXO7hYRaSRJEoJHDUHEXeMBoxEVv/2BnMXvwVJR6eyuEenHDbOASZKE6dOn48SJEyguLkZxcTFOnDiBhx56SPMEPSeAyKV4J8Yj9un7YIqrz0yQ9fybqDp0zNndIqIW8O2ahNin7oUxIhR1WXn1k7tpJ53dLSJqAb++5yL68btgCApAbdopZD2/DLXZuc7uFhG1QMDgZERPmwzJ1wfVfx9D9oLlqCsscna3iPThhhNAlZWVqKioAAAEBQWhoKAAS5YswYYNGzTXyQkgcjmmyDDEPnUvfM6pz0yQs+QdlO/63dndIqIW8IqPRuzT98GrfSwsxWXIfmk5Kv865OxuEVEL+HRqj5in7oUpOhzm3EJkL3wd1cfSnN0tImoB357nIObxe2AMCULtyWxkL3wdNSd5Rz55ALNsX3EhY8aMwapVqwAARUVFOP/88/Hyyy9jzJgxWLZsmaY6OQFELskY6I/oGXfAf8DpzARv/RclG5mZgMidmcKCETvrHvj26Ay5uga5/34fZb/scXa3iKgFvGIiEfPUffDu2A6W8grkvPIWKn7/y9ndIqIW8O4Qj5gn76+/I7+wGNkvvoGqw0ed3S2iFnHHNPB79+7FJZfUZ9b9/PPPERsbi7S0NKxatQqvvfaapjo9/iXQDSlKS8qavoitTJCWuFQhXXGZxfazduWCfSrMZpvLKy3iVOdVZtvrqutqBcttZ9+oqa0RtlFrtn0cdTWi5bbrMdfZPj4AMNfY/kws1bb3sdQ2PY7wKeMg+QegfOtvKPr0W9Tl5CHkupHWzARSjWAOU3AcwuXV4mcoZdsfO2TRx1stWF6jcEIR1VUr2EewXK4V/+zKdbbXyYLvoUHwcyhaDgCSRXDwlsYHWHf6a60phJXiukRTXNte58y4rhEtrxXP29fW2V4nqEoc1wo/R+YaQfzW2B5OmsS1BETefzPyV3yDyj1/oOD9T1GXk4+gEZdYn2UWxnWtjnEtijnBcvH2+sW1rGtci+JXFO+CHxIAkkVwIIxrDXFtO05qbIx/DWrrVI7Xgp8js8LPqrnW9nEI47qmcX8lLxMiH5qI/Dc/R/XBf5D3xkqE3ng1Ai86/3/biOJUz7gW/RiLxmXR9hriWhS/EMSvQbQ9AIjGZeF4LVouvg509HhdaiOuRfGrLa5tfwYOieta28tF19oAUFtjO+ZF47VoXFaM65qWxbUhwA9R06Yg/82PUHMsHTlL3kbYhOsR0O886zZSte3xWhjvwrFUw3it8npbWA8ASbSPyutwg2BMBqAQ1+qutyXF8dr2hyLrHNduy55HvFzsM6moqEBQUBAAYMOGDRg3bhwMBgMGDRqEtDRtd9x6/ARQaWn9SwkT+x13bkdIF2VbdqBsyw5nd4N0UlpaipCQEE37AYxrT1GydhNK1m5ydjdIJ4xrAoCiz75F0WffOrsbpJOWxnXH/ul6d4kczWxB4arPULjqM2f3hHSiNa7dlkUWZp5stI0L6dKlC9asWYOxY8di/fr1mD59OgAgJycHwcHBmur0+Amg+Ph4ZGRkICgoSNObsktKStC+fXtkZGRo/pDdEY+bx92aZFlGaWkp4uPjNe3fkrhuq99joO0eO4+bce3J2uqx87gZ156qrR430HaP3d3i2m3JFkDhDkPrNi7kmWeewa233orp06fj8ssvx+DBgwHU3w2UnJysqU6PnwAyGAxISEhocT3BwcFt6kTUgMfdtjjyuFvyFwc94rqtfo+BtnvsPO7Wx7h2nrZ67Dzu1se4do62etxA2z12d4lrt+WGj4DdcMMNuPjii5GZmYk+ffpYl19++eUYO3aspjo9fgKIiIiIiIiIiNowiwzAvR4BA4DY2FjExsY2Wnb++ecLtm4eJ4CIiIiIiIiIyHPJluYf8XKxR8BaAyeAmuHj44PZs2fDx8fH2V1xKB43j9tTtaVjPVtbPXYet+cfd1s61rO11WPncXv+cbelYz1TWz1uoO0ee1s9bocz2zEB1Nw7gjyAJLe5/G9ERERERERE5OlKSkoQEhKCK+LvgcmgPMlWZ6nGD6feRHFxsce+h4p3ABERERERERGR55Jhx0ugHdITp+IEEBERERERERF5LjfMAtYaOAFERERERERERJ7LYgHAdwBxAoiIiIiIiIiIPBcngAAABmd3wFU8//zzuPDCC+Hv74/Q0FCb2zz00EPo378/fHx80LdvX5vbfPrpp+jbty/8/f2RmJiIF198sfU6rQO9jnv9+vUYNGgQgoKCEBUVheuvvx6pqamt1/EW0uO458yZA0mSmpSAgIDW7XwL6PX9lmUZL730Erp27QofHx+0b98e8+fPb72Oa8S4ZlyfjXHNuG7AuGZcM65dB+OacX02xrX7x7XLscj2FQ/HCaDTampqcOONN+K+++4TbiPLMqZMmYLx48fbXP/9999jwoQJuPfee/Hnn3/ijTfewCuvvIKlS5e2VrdbTI/jPnbsGMaMGYPLLrsMKSkpWL9+PfLy8jBu3LjW6naL6XHcjz76KDIzMxuVHj164MYbb2ytbreYHscN1A9O77zzDl566SX8/fff+Oabb3D++ee3RpdbhHHNuD4b45pxDTCuGdeMa1fDuGZcn41x7f5x7Wpk2WJX8XgyNbJixQo5JCREcZvZs2fLffr0abL8lltukW+44YZGyxYvXiwnJCTIFotFx17qryXH/dlnn8kmk0k2m83WZV9//bUsSZJcU1Ojc0/11ZLjPltKSooMQN66das+nWtFLTnuAwcOyCaTSf77779bp3OtgHEtxrhWxrh2XYxrMca1Msa162JcizGulTGuSaS4uFgGIF8eOlEeGXanYrk8dKIMQC4uLnZ2t1sN7wDSUXV1NXx9fRst8/Pzw4kTJ5CWluakXrW+AQMGwGg0YsWKFTCbzSguLsYHH3yAESNGwMvLy9ndc5h33nkHXbt2xSWXXOLsrrSqb775Bp06dcK3336Ljh07IikpCXfeeScKCgqc3bVWwbhmXDOuPQ/jmnHNuPY8jGvGNeOaFJnN9hUPxwkgHY0cORJffvklNm3aBIvFgsOHD2PJkiUAgMzMTOd2rhUlJSVhw4YNmDVrFnx8fBAaGooTJ07g448/dnbXHKa6uhoffvgh7rjjDmd3pdUdO3YMaWlp+Oyzz7Bq1Sq8//772LNnD2644QZnd61VMK4Z14xrz8O4Zlwzrj0P45pxzbgmRQ1p4JsrHs6jJ4BELwc7s+zevVu39u666y488MADuPrqq+Ht7Y1Bgwbh5ptvBgAYjUbd2mmOo487KysLd955JyZNmoRdu3bhp59+gre3N2644QbIDgwiRx/3mb788kuUlpZi4sSJrVK/Ekcft8ViQXV1NVatWoVLLrkEQ4cOxbvvvovNmzfj0KFDurUjwrhmXDOuGdctxbhmXDOuGdethXHNuPbkuHZnssViV/F0Hp0G/oEHHrCe+EWSkpJ0a0+SJCxcuBDz589HVlYWoqKisGnTJt3baY6jj/v1119HcHAwFi1aZF32n//8B+3bt8fOnTsxaNAg3dpS4ujjPtM777yDq6++GrGxsa1SvxJHH3dcXBxMJhO6du1qXda9e3cAQHp6Orp166ZbW7YwrsUY1/piXDOuWxvjWoxx3XKMa8Y147r1tbW4dmuyDKCZSdE2cAeQR08ARUZGIjIy0uHtGo1GtGvXDgDw0UcfYfDgwYiOjnZY+44+7oqKiiZ/WWn42uLAWVRnfb9TU1OxefNmfP311w5vG3D8cV900UWoq6vD0aNH0blzZwDA4cOHAQCJiYmt3j7j2jEY14xrxnXrY1w7FuOace0IjGvHYlw7Nq7dmtkCSM2846cNZAHz6AkgNdLT01FQUID09HSYzWakpKQAALp06YLAwEAAwJEjR1BWVoasrCxUVlZat+nRowe8vb2Rl5eHzz//HEOHDkVVVRVWrFiBzz77DD/99JOTjqp5ehz36NGjsXjxYsybNw+33HILSktLMWvWLCQmJiI5OdlJR6ZMj+Nu8N577yEuLg6jRo1y9GGopsdxX3HFFejXrx+mTJmCJUuWwGKx4P/+7/8wfPjwRn+NcAWMa8Y1wLgGGNeMa8Y1wLhmXLsexjXjGvDcuHY1skWGLCnf4ePIxyadxqk5yFzIpEmTGu4Ja1Q2b95s3WbIkCE2t0lNTZVlWZZzc3PlQYMGyQEBAbK/v798+eWXyzt27HDOAdlJj+OWZVn+6KOP5OTkZDkgIECOioqSr732WvngwYOOPyA76XXcZrNZTkhIkGfNmuX4g9BAr+M+efKkPG7cODkwMFCOiYmRJ0+eLOfn5zv+gJrBuGZcM67rMa4Z1w0Y14xrxrXrYFwzrj05rl1FQxr4YcZx8nDTeMUyzDhOBjw7Dbwky21hmouIiIiIiIiI2pKSkhKEhITgYulqmOCluG0darFN/hbFxcUIDg52UA8di4+AEREREREREZHH8fb2RmxsLLZlfWvX9rGxsY0eM/Q0vAOIiIiIiIiIiDxSVVUVampq7NrW29sbvr6+rdwj5+EEEBERERERERGRhzM4uwNERERERERERNS6OAFEREREREREROThOAFEREREREREROThOAFEREREREREROThOAFEqg0dOhTTpk3zmDYnT56M6667rlXqJnIHjGkiz8O4JvI8jGsiaimTsztAZI8vv/wSXl5e1q+TkpIwbdo0hw+CRKQPxjSR52FcE3kexjWRZ+EEELmF8PBwZ3eBiHTEmCbyPIxrIs/DuCbyLHwEjFqksLAQEydORFhYGPz9/TFq1Cj8888/1vXvv/8+QkNDsX79enTv3h2BgYG48sorkZmZad2mrq4OU6dORWhoKCIiIjBz5kxMmjSp0S2hZ95+OnToUKSlpWH69OmQJAmSJAEA5syZg759+zbq35IlS5CUlGT92mw24+GHH7a2NWPGDMiy3GgfWZaxaNEidOrUCX5+fujTpw8+//xzfT4wIhfHmCbyPIxrIs/DuCYiLTgBRC0yefJk7N69G19//TV+/fVXyLKMq666CrW1tdZtKioq8NJLL+GDDz7A1q1bkZ6ejkcffdS6fuHChfjwww+xYsUKbN++HSUlJVizZo2wzS+//BIJCQmYN28eMjMzGw1kzXn55Zfx3nvv4d1338W2bdtQUFCA1atXN9rmqaeewooVK7Bs2TL89ddfmD59Ov71r3/hp59+sv+DIXJTjGkiz8O4JvI8jGsi0kQmUmnIkCHyQw89JB8+fFgGIG/fvt26Li8vT/bz85M//fRTWZZlecWKFTIA+ciRI9ZtXn/9dTkmJsb6dUxMjPziiy9av66rq5M7dOggjxkzpkmbDRITE+XFixc36tfs2bPlPn36NFq2ePFiOTEx0fp1XFyc/MILL1i/rq2tlRMSEqxtlZWVyb6+vvIvv/zSqJ477rhDvuWWWxQ/FyJ3xZgm8jyMayLPw7gmopbiO4BIs4MHD8JkMuGCCy6wLouIiEC3bt1w8OBB6zJ/f3907tzZ+nVcXBxycnIAAMXFxcjOzsb5559vXW80GtG/f39YLBZd+1tcXIzMzEwMHjzYusxkMmHAgAHWW1APHDiAqqoqDB8+vNG+NTU1SE5O1rU/RK6GMU3keRjXRJ6HcU1EWnECiDSTz3pu98zlDc8EA2iUOQAAJElqsu+Z2yvVrcRgMDTZ78zbYO3RMOB99913aNeuXaN1Pj4+qvtE5E4Y00Seh3FN5HkY10SkFd8BRJr16NEDdXV12Llzp3VZfn4+Dh8+jO7du9tVR0hICGJiYvDbb79Zl5nNZuzbt09xP29vb5jN5kbLoqKikJWV1WgASklJadRWXFwcduzYYV1WV1eHPXv2NDomHx8fpKeno0uXLo1K+/bt7TomInfFmCbyPIxrIs/DuCYirXgHEGl2zjnnYMyYMbjrrrvw5ptvIigoCI8//jjatWuHMWPG2F3Pgw8+iAULFqBLly4499xz8e9//xuFhYVN/iJxpqSkJGzduhU333wzfHx8EBkZiaFDhyI3NxeLFi3CDTfcgHXr1uH7779HcHCwdb+HHnoIL7zwAs455xx0794dr7zyCoqKiqzrg4KC8Oijj2L69OmwWCy4+OKLUVJSgl9++QWBgYGYNGmSps+KyB0wpok8D+OayPMwrolIK94BRC2yYsUK9O/fH1dffTUGDx4MWZaxdu3aJrecKpk5cyZuueUWTJw4EYMHD0ZgYCBGjhwJX19f4T7z5s3D8ePH0blzZ0RFRQEAunfvjjfeeAOvv/46+vTpg99++61RpgMAeOSRRzBx4kRMnjwZgwcPRlBQEMaOHdtom2effRbPPPMMFixYgO7du2PkyJH45ptv0LFjRxWfDJF7YkwTeR7GNZHnYVwTkRaSrOVBT6JWZLFY0L17d9x000149tlnnd0dImohxjSR52FcE3kexjWR5+MjYOR0aWlp2LBhA4YMGYLq6mosXboUqampuPXWW53dNSLSgDFN5HkY10Seh3FN1PbwETByOoPBgPfffx8DBw7ERRddhP379+OHH36w+yV2RORaGNNEnodxTeR5GNdEbQ8fASMiIiIiIiIi8nC8A4iIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMNxAoiIiIiIiIiIyMP9P66320B8Ji+DAAAAAElFTkSuQmCC", + "text/plain": [ + "
" + ] + }, + "metadata": {} + } + ], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.516207Z", + "start_time": "2023-02-07T02:10:14.392886Z" + }, + "scrolled": true + } + }, + { + "cell_type": "code", + "execution_count": 6, + "source": [ + "import rasterio\n", + "from rasterio.crs import CRS\n", + "from affine import Affine\n", + "\n", + "with rasterio.open(f'netcdf:{gunw_path}:/science/grids/corrections/external/tides/solidEarthTide') as ds:\n", + " p = ds.profile" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.533079Z", + "start_time": "2023-02-07T02:10:15.518231Z" + } + } + }, + { + "cell_type": "markdown", + "source": [ + "Make sure lat/lon" + ], + "metadata": {} + }, + { + "cell_type": "code", + "execution_count": 7, + "source": [ + "assert p['crs'] == CRS.from_epsg(4326)" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.537628Z", + "start_time": "2023-02-07T02:10:15.534889Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 9, + "source": [ + "import numpy as np\n", + "\n", + "assert p['nodata'] == 0" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.541418Z", + "start_time": "2023-02-07T02:10:15.539268Z" + } + } + }, + { + "cell_type": "markdown", + "source": [ + "Make sure transform is not identity and not None" + ], + "metadata": {} + }, + { + "cell_type": "code", + "execution_count": 10, + "source": [ + "assert p['transform']" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.547207Z", + "start_time": "2023-02-07T02:10:15.544776Z" + } + } + }, + { + "cell_type": "code", + "execution_count": 11, + "source": [ + "assert p['transform'] != Affine(1, 0, 0, 0, 1, 0)" + ], + "outputs": [], + "metadata": { + "ExecuteTime": { + "end_time": "2023-02-07T02:10:15.551306Z", + "start_time": "2023-02-07T02:10:15.548818Z" + } + } + } + ], + "metadata": { + "kernelspec": { + "name": "python3", + "display_name": "Python 3.9.15 64-bit ('topsapp_env': conda)" + }, + "language_info": { + "name": "python", + "version": "3.9.16", + "mimetype": "text/x-python", + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "pygments_lexer": "ipython3", + "nbconvert_exporter": "python", + "file_extension": ".py" + }, + "interpreter": { + "hash": "83c8da8660bbe18150fdc09aeaa55e9731d119a6641346a233b76fde249768bb" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} \ No newline at end of file diff --git a/tests/test_notebooks.py b/tests/test_notebooks.py index 2dfd101..7c3c434 100644 --- a/tests/test_notebooks.py +++ b/tests/test_notebooks.py @@ -4,7 +4,8 @@ import pytest notebooks = ['localize-data.ipynb', - 'prepare-for-delivery.ipynb'] + 'prepare-for-delivery.ipynb', + 'solid_earth_tides.ipynb'] @pytest.mark.parametrize('notebook_name', notebooks)