{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The pytesmo validation framework\n",
    "\n",
    "The pytesmo validation framework takes care of iterating over datasets, spatial and temporal matching as well as\n",
    "scaling. It uses metric calculators to then calculate metrics that are returned to the user. There are several\n",
    "metrics calculators included in pytesmo but new ones can be added simply.\n",
    "\n",
    "## Overview\n",
    "\n",
    "How does the validation framework work? It makes these assumptions about the used datasets:\n",
    "\n",
    "- The dataset readers that are used have a method for reading data (usually called `read`) that can be called either by a grid point index (gpi) which can be any indicator that identifies a certain grid point or by using longitude and latitude. This means that both call signatures `read(gpi)` and `read(lon, lat)` must be valid. Please check the [pygeobase](https://github.com/TUW-GEO/pygeobase) documentation for more details on how a fully compatible dataset class should look. But a simple `read` method should do for the validation framework. This assumption can be relaxed by using the `read_ts_names` keyword in the pytesmo.validation_framework.data_manager.DataManager class.\n",
    "- **The reading method returns a pandas.DataFrame time series.**\n",
    "- Ideally the datasets classes also have a `grid` attribute that is a\n",
    "  [pygeogrids](http://pygeogrids.readthedocs.org/en/latest/) grid. This makes the calculation of lookup tables easily possible and the nearest neighbor search faster.\n",
    "\n",
    "Fortunately these assumptions are true about the dataset readers included in pytesmo. \n",
    "\n",
    "It also makes a few assumptions about how to perform a validation. For a comparison study it is often necessary to\n",
    "choose a spatial reference grid, a temporal reference and a scaling or data space reference.\n",
    "\n",
    "### Spatial reference\n",
    "The spatial reference is the one to which all the other datasets are matched spatially. Often through nearest\n",
    "neighbor search. The validation framework uses grid points of the dataset specified as the spatial reference to\n",
    "spatially match all the other datasets with nearest neighbor search. Other, more sophisticated spatial matching\n",
    "algorithms are not implemented at the moment. If you need a more complex spatial matching then a preprocessing of\n",
    "the data is the only option at the moment.\n",
    "\n",
    "### Temporal reference\n",
    "The temporal reference is the dataset to which the other dataset are temporally matched. That means that the\n",
    "nearest observation to the reference timestamps in a certain time window is chosen for each comparison dataset.\n",
    "This is by default done by the temporal matching module included in pytesmo. How many datasets should be matched to\n",
    "the reference dataset at once can be configured, we will cover how to do this later.\n",
    "\n",
    "### Data space reference\n",
    "It is often necessary to bring all the datasets into a common data space by using. Scaling is often used for that\n",
    "and pytesmo offers a choice of several scaling algorithms (e.g. CDF matching, min-max scaling, mean-std scaling,\n",
    "triple collocation based scaling). The data space reference can also be chosen independently from the other two\n",
    "references.\n",
    "\n",
    "## Data Flow\n",
    "\n",
    "After it is initialized, the validation framework works through the following steps:\n",
    "\n",
    "1. Read all the datasets for a certain job (gpi, lon, lat)\n",
    "2. Read all the masking dataset if any\n",
    "3. Mask the temporal reference dataset using the masking data\n",
    "4. Temporally match all the chosen combinations of temporal reference and other datasets\n",
    "5. Turn the temporally matched time series over to the metric calculators\n",
    "6. Get the calculated metrics from the metric calculators\n",
    "7. Put all the metrics into a dictionary by dataset combination and return them.\n",
    "\n",
    "## Masking datasets\n",
    "Masking datasets can be used if the datasets that are compared do not contain the necessary information to mask\n",
    " them. For example we might want to use modelled soil temperature data to mask our soil moisture observations\n",
    "before comparing them. To be able to do that we just need a Dataset that returns a pandas.DataFrame with one column\n",
    " of boolean data type. Everywhere where the masking dataset is `True` the data will be masked.\n",
    "\n",
    "Let's look at a first example.\n",
    "\n",
    "## Example soil moisture validation: ASCAT - ISMN\n",
    "\n",
    "This example shows how to setup the pytesmo validation framework to perform a comparison between ASCAT and ISMN data. \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "is_executing": true
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "from pathlib import Path\n",
    "from pprint import pprint\n",
    "import pytesmo.validation_framework.metric_calculators as metrics_calculators\n",
    "from datetime import datetime\n",
    "import warnings\n",
    "from ascat.read_native.cdr import AscatGriddedNcTs\n",
    "\n",
    "from ismn.interface import ISMN_Interface  # install ismn: 'pip install ismn'\n",
    "from pytesmo.validation_framework.validation import Validation\n",
    "from pytesmo.validation_framework.results_manager import netcdf_results_manager\n",
    "from pytesmo.utils import rootdir\n",
    "\n",
    "import tempfile"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You need the test data from https://github.com/TUW-GEO/pytesmo-test-data for this example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Data is stored in: /tmp/tmppeisd__e\n"
     ]
    }
   ],
   "source": [
    "from tempfile import mkdtemp\n",
    "output_folder = Path(mkdtemp())\n",
    "print('Data is stored in:', output_folder)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we initialize the data readers that we want to use. In this case the ASCAT soil moisture time series and in\n",
    "situ data from the ISMN.\n",
    "\n",
    "Initialize ASCAT reader"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "testdata_path = rootdir() / \"tests\" / \"test-data\"\n",
    "ascat_data_folder = testdata_path / \"sat\" / \"ascat\" / \"netcdf\" / \"55R22\"\n",
    "ascat_grid_fname = testdata_path / \"sat\" / \"ascat\" / \"netcdf\" / \"grid\" / \"TUW_WARP5_grid_info_2_1.nc\"\n",
    "static_layer_path = testdata_path / \"sat\" / \"h_saf\" / \"static_layer\"\n",
    "\n",
    "\n",
    "#init the AscatSsmCdr reader with the paths\n",
    "with warnings.catch_warnings():\n",
    "    warnings.filterwarnings('ignore') # some warnings are expected and ignored\n",
    "    \n",
    "    ascat_reader = AscatGriddedNcTs(\n",
    "        ascat_data_folder,\n",
    "        \"TUW_METOP_ASCAT_WARP55R22_{:04d}\",\n",
    "        grid_filename=ascat_grid_fname,\n",
    "        static_layer_path=static_layer_path\n",
    "    )\n",
    "\n",
    "ascat_reader.read_bulk = True"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Initialize ISMN reader"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Processing metadata for all ismn stations into folder /home/wpreimes/shares/home/code/pytesmo/tests/test-data/ismn/multinetwork/header_values.\n",
      "This may take a few minutes, but is only done once...\n",
      "Hint: Use `parallel=True` to speed up metadata generation for large datasets\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Files Processed: 100%|██████████| 8/8 [00:00<00:00, 65.44it/s]"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Metadata generation finished after 0 Seconds.\n",
      "Metadata and Log stored in /tmp/tmp_72gspsz\n",
      "Found existing ismn metadata in /tmp/tmp_72gspsz/header_values.csv.\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "ismn_data_folder = testdata_path / \"ismn/multinetwork/header_values\"\n",
    "meta_path = tempfile.mkdtemp()\n",
    "ismn_reader = ISMN_Interface(ismn_data_folder, meta_path=meta_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The validation is run based on jobs. A job consists of at least three lists or numpy arrays specifing the grid\n",
    "point index, its latitude and longitude. In the case of the ISMN we can use the `dataset_ids` that identify every\n",
    "time series in the downloaded ISMN data as our grid point index. We can then get longitude and latitude from the\n",
    "metadata of the dataset.\n",
    "\n",
    "**DO NOT CHANGE** the name ***jobs*** because it will be searched during the parallel processing!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Jobs (gpi, lon, lat):\n",
      "[(0, 102.1333, 33.8833)]\n"
     ]
    }
   ],
   "source": [
    "jobs = []\n",
    "\n",
    "ids = ismn_reader.get_dataset_ids(variable='soil_moisture', min_depth=0, max_depth=0.1)\n",
    "\n",
    "for idx in ids:\n",
    "    metadata = ismn_reader.read_metadata(idx)\n",
    "    jobs.append((idx, metadata['longitude'].val, metadata['latitude'].val))\n",
    "    break\n",
    "\n",
    "print(\"Jobs (gpi, lon, lat):\")\n",
    "print(jobs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For this example we only process one job / pixel / sensor.\n",
    "\n",
    "It is important here that the ISMN reader has a reading function that works by just using the `dataset_id`. In this\n",
    " way the validation framework can go through the jobs and read the correct time series."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ISMN data example:\n",
      "                     soil_moisture soil_moisture_flag soil_moisture_orig_flag\n",
      "date_time                                                                    \n",
      "2008-07-01 00:00:00            0.5                C03                       M\n",
      "2008-07-01 01:00:00            0.5                C03                       M\n",
      "2008-07-01 02:00:00            0.5                C03                       M\n",
      "2008-07-01 03:00:00            0.5                C03                       M\n",
      "2008-07-01 04:00:00            0.5                C03                       M\n"
     ]
    }
   ],
   "source": [
    "data = ismn_reader.read(ids[0])\n",
    "print('ISMN data example:')\n",
    "print(data.head())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Initialize the Validation class\n",
    "\n",
    "The Validation class is the heart of the validation framework. It contains the information about which datasets to\n",
    "read using which arguments or keywords and if they are spatially compatible. It also contains the settings about\n",
    "which metric calculators to use and how to perform the scaling into the reference data space. It is initialized in\n",
    "the following way:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "datasets = {\n",
    "    'ISMN': {\n",
    "        'class': ismn_reader,\n",
    "        'columns': ['soil_moisture']\n",
    "    },\n",
    "    'ASCAT': {\n",
    "        'class': ascat_reader,\n",
    "        'columns': ['sm'],\n",
    "        'kwargs': {'mask_frozen_prob': 80,\n",
    "                   'mask_snow_prob': 80,\n",
    "                   'mask_ssf': True}\n",
    "    }}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The datasets dictionary contains all the information about the datasets to read. The `class` is the dataset class\n",
    "to use which we have already initialized. The `columns` key describes which columns of the dataset interest us for\n",
    "validation. This a mandatory field telling the framework which other columns to ignore. In this case the columns\n",
    "`soil moisture_flag` and `soil moisture_orig_flag` will be ignored by the ISMN reader. We can also specify\n",
    "additional keywords that should be given to the `read` method of the dataset reader. In this case we want the\n",
    "ASCAT reader to mask the ASCAT soil moisture using the included frozen and snow probabilities as well as the SSF.\n",
    "There are also other keys that can be used here. Please see the documentation for explanations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/wpreimes/shares/home/code/pytesmo/src/pytesmo/validation_framework/validation.py:144: UserWarning: You are using the default temporal matcher. If you are using one of the newer metric calculators (PairwiseIntercomparisonMetrics, TripleCollocationMetrics) you should probably use `make_combined_temporal_matcher` instead. Have a look at the documentation of the metric calculators for more info.\n",
      "  warnings.warn(\n"
     ]
    }
   ],
   "source": [
    "period = [datetime(2007, 1, 1), datetime(2014, 12, 31)]\n",
    "basic_metrics = metrics_calculators.PairwiseIntercomparisonMetrics()\n",
    "\n",
    "process = Validation(\n",
    "    datasets, 'ISMN',\n",
    "    temporal_ref='ASCAT',\n",
    "    scaling='mean_std',\n",
    "    scaling_ref='ASCAT',\n",
    "    metrics_calculators={(2, 2): basic_metrics.calc_metrics},\n",
    "    period=period)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "During the initialization of the Validation class we can also tell it other things that it needs to know. In this\n",
    "case it uses the datasets we have specified earlier. The spatial reference is the `'ISMN'` dataset which is the\n",
    "second argument. The 'metrics_calculators' argument looks a little bit strange so let's look at it in more detail.\n",
    "\n",
    "It is a dictionary with a tuple as the key and a function as the value. The key tuple `(n, k)` has the following\n",
    "meaning: `n` datasets are temporally matched together and then given in sets of `k` columns to the metric\n",
    "calculator. The metric calculator then gets a DataFrame with the columns ['ref', 'k1', 'k2' ...] and so on\n",
    "depending on the value of k. The value of `(2, 2)` makes sense here since we only have two datasets and all our\n",
    "metrics also take two inputs.\n",
    "\n",
    "This can be used in more complex scenarios to e.g. have three input datasets that are all temporally matched\n",
    "together and then combinations of two input datasets are given to one metric calculator while all three datasets\n",
    "are given to another metric calculator. This could look like this:\n",
    "\n",
    "```\n",
    "{ (3 ,2): metric_calc,\n",
    "  (3, 3): triple_collocation}\n",
    "```\n",
    "\n",
    "Create the variable ***save_path*** which is a string representing the path where the results will be saved.\n",
    "**DO NOT CHANGE** the name ***save_path*** because it will be searched during the parallel processing!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{(('ASCAT', 'sm'), ('ISMN', 'soil_moisture')): {'BIAS': array([1.4210855e-14], dtype=float32),\n",
      "                                                'BIAS_ci_lower': array([-1.124606], dtype=float32),\n",
      "                                                'BIAS_ci_upper': array([1.124606], dtype=float32),\n",
      "                                                'R': array([0.5342869], dtype=float32),\n",
      "                                                'RMSD': array([10.789426], dtype=float32),\n",
      "                                                'RSS': array([41558.98], dtype=float32),\n",
      "                                                'R_ci_lower': array([0.45576647], dtype=float32),\n",
      "                                                'R_ci_upper': array([0.60455596], dtype=float32),\n",
      "                                                'gpi': array([0], dtype=int32),\n",
      "                                                'lat': array([33.8833]),\n",
      "                                                'lon': array([102.1333]),\n",
      "                                                'mse': array([116.41171], dtype=float32),\n",
      "                                                'mse_bias': array([2.019484e-28], dtype=float32),\n",
      "                                                'mse_corr': array([116.41171], dtype=float32),\n",
      "                                                'mse_var': array([3.1554436e-30], dtype=float32),\n",
      "                                                'n_obs': array([357], dtype=int32),\n",
      "                                                'p_R': array([9.664754e-28], dtype=float32),\n",
      "                                                'p_rho': array([2.471651e-28], dtype=float32),\n",
      "                                                'p_tau': array([2.1020434e-26], dtype=float32),\n",
      "                                                'rho': array([0.53934574], dtype=float32),\n",
      "                                                'rho_ci_lower': array([0.45559874], dtype=float32),\n",
      "                                                'rho_ci_upper': array([0.61362934], dtype=float32),\n",
      "                                                'status': array([0], dtype=int32),\n",
      "                                                'tau': array([0.3907124], dtype=float32),\n",
      "                                                'tau_ci_lower': array([0.35201582], dtype=float32),\n",
      "                                                'tau_ci_upper': array([0.42807564], dtype=float32),\n",
      "                                                'urmsd': array([10.789426], dtype=float32),\n",
      "                                                'urmsd_ci_lower': array([10.065898], dtype=float32),\n",
      "                                                'urmsd_ci_upper': array([11.661137], dtype=float32)}}\n"
     ]
    }
   ],
   "source": [
    "save_path = output_folder / 'ascat_ismn'\n",
    "\n",
    "for job in jobs:\n",
    "\n",
    "    results = process.calc(*job)\n",
    "    pprint(results)\n",
    "    netcdf_results_manager(results, save_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The validation is then performed by looping over all the defined jobs and storing the results.\n",
    "You can see that the results are a dictionary where the key is a tuple defining the exact combination of datasets\n",
    "and columns that were used for the calculation of the metrics. The metrics itself are a dictionary of `metric-name:\n",
    " numpy.ndarray` which also include information about the gpi, lon and lat. Since all the information contained in\n",
    "the job is given to the metric calculator they can be stored in the results.\n",
    "\n",
    "Storing of the results to disk is at the moment supported by the `netcdf_results_manager` which creates a netCDF\n",
    "file for each dataset combination and stores each metric as a variable. We can inspect the stored netCDF file which\n",
    " is named after the dictionary key:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "lon [102.1333]\n",
      "lat [33.8833]\n",
      "idx [0]\n",
      "_row_size [0]\n",
      "time []\n",
      "gpi [0]\n",
      "status [0]\n",
      "R [0.5342869]\n",
      "p_R [9.664754e-28]\n",
      "BIAS [1.4210855e-14]\n",
      "RMSD [10.789426]\n",
      "mse [116.41171]\n",
      "RSS [41558.98]\n",
      "mse_corr [116.41171]\n",
      "mse_bias [2.019484e-28]\n",
      "urmsd [10.789426]\n",
      "mse_var [3.1554436e-30]\n",
      "rho [0.53934574]\n",
      "p_rho [2.471651e-28]\n",
      "tau [0.3907124]\n",
      "p_tau [2.1020434e-26]\n",
      "BIAS_ci_lower [-1.124606]\n",
      "BIAS_ci_upper [1.124606]\n",
      "urmsd_ci_lower [10.065898]\n",
      "urmsd_ci_upper [11.661137]\n",
      "R_ci_lower [0.45576647]\n",
      "R_ci_upper [0.60455596]\n",
      "rho_ci_lower [0.45559874]\n",
      "rho_ci_upper [0.61362934]\n",
      "tau_ci_lower [0.35201582]\n",
      "tau_ci_upper [0.42807564]\n",
      "n_obs [357]\n"
     ]
    }
   ],
   "source": [
    "import netCDF4\n",
    "results_fname = Path(save_path) / \"ASCAT.sm_with_ISMN.soil_moisture.nc\"\n",
    "\n",
    "with netCDF4.Dataset(str(results_fname)) as ds:\n",
    "    for var in ds.variables:\n",
    "        print(var, ds.variables[var][:])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Metric Calculator Adapters\n",
    "\n",
    "Metric calculators compute a set of comparison metrics based on the passed input. Usually the input are two or more collocated time series. However, to validate soil moisture it is often desired to assess the quality for certain temporal subsets independently. E.g. if a varying level of agreement between two datasets is expected for two different seasons (e.g. wet/dry season, summer/winter) or if the varying level of agreement should be assessed (e.g. on a seasonal/monthly basis). In this case pytesmo provides adapters to split up the input time series before computing validation metrics. These adapters work with any (properly implemented) metric calculator. You can use one of the predefined adapters in `pytesmo.validation_framework.metric_calculators_adapters`. Currently, there are 2 options:\n",
    "\n",
    "- `SubsetsMetricsAdapter`: This adapter lets you define arbitrary temporal subsets of your time series data before metrics computation. You can compute metrics for a certain set of datetimes or datetime ranges.\n",
    "    - `MonthsMetricsAdapter`: This is a version of a SubsetsMetricsAdapter that allows splitting up the data based on their month. You can e.g. compute validation metrics for each month of a year individually, or for each season (i.e. 3 * 4 months), or any other combination of months.\n",
    "\n",
    "Let's look at an example. Here we use the base `SubsetsMetricsAdapter`. We compute the same metrics as before, but apply the metrics calculator to a set of predefined temporal subgroups. Note that groups in this example were chosen to demonstrate the feature range of the `TsDistributor` rather than to produce meaningful results.\n",
    "\n",
    "We create a SubsetsMetricsAdapter that takes the original metrics calculator and a list of named, temporal subsets (dict).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "from pytesmo.validation_framework.metric_calculators_adapters import SubsetsMetricsAdapter, TsDistributor\n",
    "from pytesmo.time_series.grouping import YearlessDatetime\n",
    "from datetime import datetime\n",
    "import shutil\n",
    "\n",
    "subsets = {\n",
    "    # The first subset does only include values in August 2010 and 2011 + July 31st of both years\n",
    "    'August2010/11': TsDistributor(dates=[datetime(year=2009, month=7, day=31), datetime(2010, 7, 31)],\n",
    "                                   date_ranges=[(datetime(2010, 7, 1, 0, 0), datetime(2010, 7, 31)),\n",
    "                                                (datetime(2011, 7, 1), datetime(2011, 7, 31))]),\n",
    "    \n",
    "    # The second subset includes all values from June to September of ANY YEAR, but not August\n",
    "    'JJS': TsDistributor(yearless_date_ranges=[\n",
    "        (YearlessDatetime(month=6, day=1, hour=0, minute=0), YearlessDatetime(7, 31)),\n",
    "        (YearlessDatetime(9, 1), YearlessDatetime(9, 30))\n",
    "    ]),\n",
    "    \n",
    "    # The third group includes all values from Feb 28th to the end of April of ANY YEAR\n",
    "    'MarchApril': TsDistributor(yearless_date_ranges=[(YearlessDatetime(2, 28), YearlessDatetime(4, 30))]),\n",
    "}\n",
    "\n",
    "adapted_intercomparison_metrics = SubsetsMetricsAdapter(\n",
    "    metrics_calculators.PairwiseIntercomparisonMetrics(calc_kendall=False, calc_spearman=False),\n",
    "    subsets,\n",
    "    group_results='join')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "The adapted metrics calculator can now be used as before. The results will contain the chosen group names accordingly, as defined in the adapter. The results manager will write them to a netcdf file. Note that this only works as `group_results='join'` was selected, otherwise the output format would different from what the results manager expects."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{(('ASCAT', 'sm'), ('ISMN', 'soil_moisture')): {'August2010/11|BIAS': array([0.02893662], dtype=float32),\n",
      "                                                'August2010/11|BIAS_ci_lower': array([0.01474692], dtype=float32),\n",
      "                                                'August2010/11|BIAS_ci_upper': array([0.04312633], dtype=float32),\n",
      "                                                'August2010/11|R': array([0.7176565], dtype=float32),\n",
      "                                                'August2010/11|RMSD': array([0.04613708], dtype=float32),\n",
      "                                                'August2010/11|RSS': array([0.05960163], dtype=float32),\n",
      "                                                'August2010/11|R_ci_lower': array([0.47057065], dtype=float32),\n",
      "                                                'August2010/11|R_ci_upper': array([0.8603755], dtype=float32),\n",
      "                                                'August2010/11|mse': array([0.00212863], dtype=float32),\n",
      "                                                'August2010/11|mse_bias': array([0.00083733], dtype=float32),\n",
      "                                                'August2010/11|mse_corr': array([0.00128012], dtype=float32),\n",
      "                                                'August2010/11|mse_var': array([1.1178111e-05], dtype=float32),\n",
      "                                                'August2010/11|n_obs': array([28], dtype=int32),\n",
      "                                                'August2010/11|p_R': array([1.718043e-05], dtype=float32),\n",
      "                                                'August2010/11|status': array([0], dtype=int32),\n",
      "                                                'August2010/11|urmsd': array([0.03593468], dtype=float32),\n",
      "                                                'August2010/11|urmsd_ci_lower': array([0.02893201], dtype=float32),\n",
      "                                                'August2010/11|urmsd_ci_upper': array([0.04980955], dtype=float32),\n",
      "                                                'JJS|BIAS': array([0.02439872], dtype=float32),\n",
      "                                                'JJS|BIAS_ci_lower': array([0.01300029], dtype=float32),\n",
      "                                                'JJS|BIAS_ci_upper': array([0.03579714], dtype=float32),\n",
      "                                                'JJS|R': array([0.43224052], dtype=float32),\n",
      "                                                'JJS|RMSD': array([0.08166534], dtype=float32),\n",
      "                                                'JJS|RSS': array([1.2204688], dtype=float32),\n",
      "                                                'JJS|R_ci_lower': array([0.3063946], dtype=float32),\n",
      "                                                'JJS|R_ci_upper': array([0.54323655], dtype=float32),\n",
      "                                                'JJS|mse': array([0.00666923], dtype=float32),\n",
      "                                                'JJS|mse_bias': array([0.0005953], dtype=float32),\n",
      "                                                'JJS|mse_corr': array([0.00590448], dtype=float32),\n",
      "                                                'JJS|mse_var': array([0.00016945], dtype=float32),\n",
      "                                                'JJS|n_obs': array([183], dtype=int32),\n",
      "                                                'JJS|p_R': array([9.96128e-10], dtype=float32),\n",
      "                                                'JJS|status': array([0], dtype=int32),\n",
      "                                                'JJS|urmsd': array([0.07793543], dtype=float32),\n",
      "                                                'JJS|urmsd_ci_lower': array([0.07087912], dtype=float32),\n",
      "                                                'JJS|urmsd_ci_upper': array([0.0870943], dtype=float32),\n",
      "                                                'MarchApril|BIAS': array([-0.0636023], dtype=float32),\n",
      "                                                'MarchApril|BIAS_ci_lower': array([-0.07745089], dtype=float32),\n",
      "                                                'MarchApril|BIAS_ci_upper': array([-0.0497537], dtype=float32),\n",
      "                                                'MarchApril|R': array([0.8920552], dtype=float32),\n",
      "                                                'MarchApril|RMSD': array([0.07466082], dtype=float32),\n",
      "                                                'MarchApril|RSS': array([0.1895241], dtype=float32),\n",
      "                                                'MarchApril|R_ci_lower': array([0.7931545], dtype=float32),\n",
      "                                                'MarchApril|R_ci_upper': array([0.94511515], dtype=float32),\n",
      "                                                'MarchApril|mse': array([0.00557424], dtype=float32),\n",
      "                                                'MarchApril|mse_bias': array([0.00404525], dtype=float32),\n",
      "                                                'MarchApril|mse_corr': array([0.00138048], dtype=float32),\n",
      "                                                'MarchApril|mse_var': array([0.0001485], dtype=float32),\n",
      "                                                'MarchApril|n_obs': array([34], dtype=int32),\n",
      "                                                'MarchApril|p_R': array([1.4273317e-12], dtype=float32),\n",
      "                                                'MarchApril|status': array([0], dtype=int32),\n",
      "                                                'MarchApril|urmsd': array([0.03910225], dtype=float32),\n",
      "                                                'MarchApril|urmsd_ci_lower': array([0.03201326], dtype=float32),\n",
      "                                                'MarchApril|urmsd_ci_upper': array([0.05224345], dtype=float32),\n",
      "                                                'gpi': array([0], dtype=int32),\n",
      "                                                'lat': array([33.8833]),\n",
      "                                                'lon': array([102.1333])}}\n",
      "Results were stored in /tmp/tmppeisd__e/ascat_ismn_adapted\n"
     ]
    }
   ],
   "source": [
    "import pandas as pd\n",
    "from pytesmo.validation_framework.temporal_matchers import make_combined_temporal_matcher\n",
    "process = Validation(\n",
    "    datasets, 'ISMN',\n",
    "    temporal_ref='ASCAT',\n",
    "    scaling='mean_std',\n",
    "    scaling_ref='ISMN',\n",
    "    metrics_calculators={(2, 2): adapted_intercomparison_metrics.calc_metrics},\n",
    "    temporal_matcher=make_combined_temporal_matcher(pd.Timedelta(6, \"H\")),\n",
    "    period=period)\n",
    "\n",
    "save_path = output_folder / 'ascat_ismn_adapted'\n",
    "\n",
    "if os.path.exists(save_path):\n",
    "    shutil.rmtree(save_path)\n",
    "\n",
    "for job in jobs:\n",
    "    results = process.calc(*job)\n",
    "    pprint(results)\n",
    "    netcdf_results_manager(results, save_path)\n",
    "\n",
    "print(f\"Results were stored in {save_path}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parallel processing\n",
    "\n",
    "The same code can be executed in parallel by defining the following `start_processing` function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "def start_processing(job):\n",
    "    try:\n",
    "        return process.calc(*job)\n",
    "    except RuntimeError:\n",
    "        return process.calc(*job)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`pytesmo.validation_framework.start_validation` can then be used to run your validation in parallel.\n",
    "Your setup code can look like this Ipython notebook without the loop over the jobs. Otherwise the validation would\n",
    "be done twice. Save it into a `.py` file e.g. `my_validation.py`.\n",
    "\n",
    "After [starting the ipyparallel cluster](http://ipyparallel.readthedocs.org/en/latest/process.html) you can then\n",
    "execute the following code:\n",
    "\n",
    "```python\n",
    "from pytesmo.validation_framework import start_validation\n",
    "\n",
    "# Note that before starting the validation you must start a controller\n",
    "# and engines, for example by using: ipcluster start -n 4\n",
    "# This command will launch a controller and 4 engines on the local machine.\n",
    "# Also, do not forget to change the setup_code path to your current setup.\n",
    "\n",
    "setup_code = \"my_validation.py\"\n",
    "start_validation(setup_code)\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Time Series Adapters\n",
    "\n",
    "Data readers extract time series from the candidate and reference data sets used in a validation run. However, often it is desired to change the data in some way after reading and before using them in the validation framework. Potential preprocessing steps include masking/filtering, conversion to anomalies and the combining multiple available variables in a dataset.\n",
    "\n",
    "## Masking / filtering\n",
    "\n",
    "Filters are used to select certain observations in a time series after reading to include/exclude from the validation run. Satellite and in situ observations often come with quality flags for their measurements. Filters are for example used to choose - based on these flags - which observations should be used and which ones should be discarded before computing validation metrics.\n",
    "\n",
    "There are 2 ways of masking input datasets.\n",
    "\n",
    "1) Directly upon loading a time series from a dataset, by removing any unwanted observations\n",
    "    For this the `SelfMaskingAdapter` and `AdvancedMaskingAdapter` are used.\n",
    "2) Indirectly after reading all datasets by using one or multiple independent masking datasets.\n",
    "    For this the `MaskingAdapter` is used. The masking datasets have the same format as the dataset dictionary and can be specified in the Validation class with the `masking_datasets` keyword.\n",
    "\n",
    "### Masking adapter\n",
    "\n",
    "To easily transform an existing dataset into a masking dataset `pytesmo` offers an adapter class that calls the\n",
    "reading method of an existing dataset and creates a masking dataset based on an operator, a given threshold, and (optionally) a column name."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                     soil_moisture\n",
      "date_time                         \n",
      "2008-07-01 00:00:00          False\n",
      "2008-07-01 01:00:00          False\n",
      "2008-07-01 02:00:00          False\n",
      "2008-07-01 03:00:00          False\n",
      "2008-07-01 04:00:00          False\n",
      "...                            ...\n",
      "2010-07-31 19:00:00          False\n",
      "2010-07-31 20:00:00          False\n",
      "2010-07-31 21:00:00          False\n",
      "2010-07-31 22:00:00          False\n",
      "2010-07-31 23:00:00          False\n",
      "\n",
      "[15927 rows x 1 columns]\n"
     ]
    }
   ],
   "source": [
    "from pytesmo.validation_framework.adapters import MaskingAdapter\n",
    "\n",
    "ds_mask = MaskingAdapter(ismn_reader,  '<', 0.2, 'soil_moisture')\n",
    "pprint(ds_mask.read(ids[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Self-masking adapter\n",
    "`pytesmo` also has a class that masks a dataset \"on-the-fly\", based on one of the columns it contains and an operator and a threshold. In contrast to the masking adapter mentioned above, the output of the self-masking adapter is the masked data, not the mask. The self-masking adapter wraps a data reader, which must have a `read_ts` or `read` method. Alternatively, to use a method with a name other than `read`/`read_ts`, use the `read_name` keyword which is available for each Adapter - it is still required that the method returns a pandas DataFrame! Calling the reading method will return the masked data - more precisely a DataFrame with only rows where the masking condition is true."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                     soil_moisture soil_moisture_flag soil_moisture_orig_flag\n",
      "date_time                                                                    \n",
      "2008-11-29 09:00:00           0.19            D01,D03                       M\n",
      "2008-11-29 10:00:00           0.19            D01,D03                       M\n",
      "2008-11-29 11:00:00           0.19            D01,D03                       M\n",
      "2008-11-30 03:00:00           0.19            D01,D03                       M\n",
      "2008-11-30 04:00:00           0.19            D01,D03                       M\n",
      "...                            ...                ...                     ...\n",
      "2010-03-16 08:00:00           0.18                D01                       M\n",
      "2010-03-16 09:00:00           0.18                D01                       M\n",
      "2010-03-16 10:00:00           0.18                D01                       M\n",
      "2010-03-16 11:00:00           0.18                D01                       M\n",
      "2010-03-16 12:00:00           0.19                D01                       M\n",
      "\n",
      "[2956 rows x 3 columns]\n"
     ]
    }
   ],
   "source": [
    "from pytesmo.validation_framework.adapters import SelfMaskingAdapter\n",
    "\n",
    "ds_mask = SelfMaskingAdapter(ismn_reader, '<', 0.2, 'soil_moisture', read_name='read')\n",
    "pprint(ds_mask.read(ids[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "While the `SelfMaskingAdapter` works only one one filtering rule, there is also the `AdvancedMaskingAdapter`. This one can take a list of conditions based on which a time series is masked."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                              sm  snow_prob  proc_flag orbit_dir\n",
      "2007-02-12 14:47:48.999982  16.0          0          0      b'A'\n",
      "2007-02-20 13:42:25.999977  35.0          0          0      b'A'\n",
      "2007-03-01 13:56:14.999987  37.0          0          0      b'A'\n",
      "2007-03-06 13:52:51.000007  35.0          0          0      b'A'\n",
      "2007-03-08 14:51:29.999990  36.0          0          0      b'A'\n",
      "...                          ...        ...        ...       ...\n",
      "2014-06-20 13:54:58.000002  80.0          0          0      b'A'\n",
      "2014-06-22 14:53:39.000015  77.0          0          0      b'A'\n",
      "2014-06-25 13:51:37.000006  83.0          0          0      b'A'\n",
      "2014-06-27 14:50:18.000019  81.0          0          0      b'A'\n",
      "2014-06-30 13:48:17.999999  83.0          0          0      b'A'\n",
      "\n",
      "[656 rows x 4 columns]\n"
     ]
    }
   ],
   "source": [
    "from pytesmo.validation_framework.adapters import AdvancedMaskingAdapter\n",
    "\n",
    "\n",
    "ds_mask = AdvancedMaskingAdapter(ascat_reader,\n",
    "                                 filter_list=[('snow_prob', np.less_equal, 10),\n",
    "                                              ('sm', '>', 0),\n",
    "                                              ('sm', '<', 100),\n",
    "                                              ('proc_flag', '==', 0),\n",
    "                                              ('orbit_dir', '==', b'A')],\n",
    "                                 read_name='read')\n",
    "pprint(ds_mask.read(1814367)[['sm', 'snow_prob', 'proc_flag', 'orbit_dir']])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## AnomalyAdapter and AnomalyClimAdapter\n",
    "\n",
    "These 2 adapters are used to transform absolute values into anomalies on-the-fly after reading. You can select one or multiple columns for which the anomalies are computed. Additional kwargs (like `timespan`) are passed to the underlying anomaly and climatology function. For more details there is a separate tutorial available."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                     soil_moisture soil_moisture_flag soil_moisture_orig_flag\n",
      "date_time                                                                    \n",
      "2008-07-01 00:00:00       0.140861                C03                       M\n",
      "2008-07-01 01:00:00       0.140861                C03                       M\n",
      "2008-07-01 02:00:00       0.140861                C03                       M\n",
      "2008-07-01 03:00:00       0.140861                C03                       M\n",
      "2008-07-01 04:00:00       0.140861                C03                       M\n",
      "...                            ...                ...                     ...\n",
      "2010-07-31 19:00:00      -0.151816                  U                       M\n",
      "2010-07-31 20:00:00      -0.151816                  U                       M\n",
      "2010-07-31 21:00:00      -0.151816                  U                       M\n",
      "2010-07-31 22:00:00      -0.151816                  U                       M\n",
      "2010-07-31 23:00:00      -0.151816                  U                       M\n",
      "\n",
      "[15927 rows x 3 columns]\n"
     ]
    },
    {
     "data": {
      "text/plain": "<Axes: title={'center': 'Anomaly (wrt. 1991-2020 avg.)'}, xlabel='date_time', ylabel='SM m3m-3'>"
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": "<Figure size 700x200 with 1 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnEAAAD2CAYAAAC5p9FxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAACOk0lEQVR4nO2dd3wT5R/HP0n3Lt0tdEFZhbI3yJA9lCFLEWTIliGiwE9lCAIOEFSWiCBLQEFQlgzZIKPsXUZpoS0t3TNNk+f3R3qXu+Sy2qRJmuf9evXV5PLc3fMk99x9n+8UEUIIKBQKhUKhUChWhdjcHaBQKBQKhUKhGA4V4igUCoVCoVCsECrEUSgUCoVCoVghVIijUCgUCoVCsUKoEEehUCgUCoVihVAhjkKhUCgUCsUKoUIchUKhUCgUihVChTgKhUKhUCgUK4QKcRQKhUKhUChWCBXiKJQK4Pvvv4dIJEL9+vXN3RWTEBERgZEjRxr1mI8fP4aTkxMuXLhg1ONySUpKwvz583H9+vVyHWfz5s0YOnQoateuDbFYjIiICI1tL126hO7du8PDwwPu7u7o1KkTzp07p9aOEILvv/8ederUgZOTE4KDgzFx4kRkZmaqtV2xYgUGDBiAyMhIiEQidOzY0aD+//vvvxg9ejTq1KkDNzc3VK1aFX379kVsbKxg+6tXr6JLly5wd3eHt7c3BgwYgCdPnvDaPHz4EDNnzkTTpk3h7e0NHx8ftG3bFn/88YfgMVNTUzFy5Ej4+fnB1dUVrVu3xvHjxw0ah6XQvn17TJ8+3dzdoNgAVIijUCqAX375BQBw584dXLx40cy9sQ5mzpyJrl27onXr1iY7R1JSEhYsWFBuIW7Lli24c+cOWrRogRo1amhsd/nyZbRv3x6FhYXYsmULtmzZgqKiInTu3FlNWJ05cyY+/PBD9O3bF/v378fs2bOxfft2dO3aFVKplNd27dq1ePbsGV5//XX4+/sb3P81a9YgPj4e06ZNw8GDB7Fy5UqkpqaiVatW+Pfff3lt79+/j44dO6K4uBi7du3CL7/8gocPH+K1115DWloa2+7IkSM4cOAA3nrrLfz+++/Ytm0batasiUGDBuGLL77gHVMikaBz5844fvw4Vq5ciX379iEwMBA9evTAqVOnDB6PuVm4cCFWr16NBw8emLsrlMoOoVAoJuXy5csEAOnduzcBQMaOHWvuLhmd8PBw8t577xnteHfv3iUAyOHDh412TC4lJSWkqKiI/W02btxYruPJZDL2de/evUl4eLhgu+7du5PAwECSn5/PbsvJySF+fn6kTZs27Lbnz58TOzs7MmXKFN7+27dvJwDITz/9pPH89erVIx06dDCo/y9fvlTblpubSwIDA0nnzp152wcNGkT8/PxIdnY2uy0+Pp44ODiQTz75hN2WlpZG5HK52nF79+5NXF1dSVFREbtt1apVBAA5f/48u00qlZLo6GjSokULg8ZiKdSvX79SznWKZUE1cRSKidmwYQMAYOnSpWjTpg127NiBgoICXpv4+HiIRCJ8++23WL58OSIjI+Hu7o7WrVvjv//+UzvmX3/9hdatW8PV1RUeHh7o2rWrmiZn/vz5EIlEuHnzJgYNGgQvLy/4+PhgxowZKCkpwYMHD9CjRw94eHggIiICX3/9NW//oqIifPTRR2jUqBG7b+vWrbFv3z6t483Ly4O3tzfGjx+v9ll8fDzs7OzwzTffaD3GmjVrEBQUhK5du7LbVq1aBbFYjNTUVHbbsmXLIBKJMHnyZHabXC5HlSpV8NFHH7HnFIlE+Prrr7Fo0SJERkbCyckJJ06cQPPmzQEAo0aNgkgkgkgkwvz587X2TQixWL9b6blz59CxY0e4urqy2zw8PNC+fXucP38eycnJAID//vsPMpkMvXr14u3fp08fAMDu3bvLdH5NBAQEqG1zd3dHdHQ0EhMT2W0lJSXYv38/3nrrLXh6erLbw8PD0alTJ/z555/sNj8/P4hEIrXjtmjRAgUFBcjIyGC3/fnnn6hduzZP62pvb493330Xly5dwosXL7T2/+jRo+jbty+qVasGZ2dnREVFYfz48Xj16hXbZu/evRCJRIIm2jVr1rBzhWH9+vWoVasWnJycEB0dje3bt2PkyJFaTeVchg8fju3btyM3N1ev9hRKWaBCHIViQgoLC/Hbb7+hefPmqF+/PkaPHo3c3Fz8/vvvgu1XrVqFo0ePYsWKFdi2bRvy8/PRq1cvZGdns222b9+Ovn37wtPTE7/99hs2bNiAzMxMdOzYEWfPnlU75uDBg9GwYUPs3r0bY8eOxXfffYcPP/wQ/fr1Q+/evfHnn3/i9ddfx6xZs7Bnzx52P4lEgoyMDMycORN79+7Fb7/9hnbt2mHAgAHYvHmzxjG7u7tj9OjR2LZtG6/fALB69Wo4Ojpi9OjRWr+3AwcOoH379jzhpEuXLiCE8B7Cx44dg4uLC44ePcpuu3LlCrKystClSxfeMb///nv8+++/+Pbbb3Ho0CGEhIRg48aNAIDPPvsMFy5cwIULF/D+++9r7Vt5KC4uhpOTk9p2ZtutW7fYdtztDA4ODmrChqnIzs7G1atXUa9ePXbb48ePUVhYiAYNGqi1b9CgAR49eoSioiKtxz1x4gT8/f15guPt27c1HhNQuCFo4/Hjx2jdujXWrFmDI0eOYO7cubh48SLatWvHmp779OmDgIAA9jfnsmnTJjRp0oQ9308//YRx48ahQYMG2LNnDz777DMsWLAAJ0+e1NoPLh07dkR+fr5B+1AoBmNuVSCFUpnZvHkzAUDWrl1LCFGYqNzd3clrr73Ga/f06VMCgMTExJCSkhJ2+6VLlwgA8ttvvxFCFGazkJAQEhMTwzOh5ebmkoCAAJ5Jbt68eQQAWbZsGe9cjRo1IgDInj172G1SqZT4+/uTAQMGaBxLSUkJkUqlZMyYMaRx48a8z1TNqY8fPyZisZh899137LbCwkLi6+tLRo0apfEchChMewDI0qVL1T6rVq0aGT16NCGEEIlEQtzc3MisWbMIAPLs2TNCCCFffvklcXBwIHl5eYQQ5Xdbo0YNUlxczDuescypXLSZUxs1akRq1arF++2kUimpXr06AUC2b99OCCHk+vXrBABZuHAhb//jx48TAMTR0VHj+ctiThVi2LBhxN7enly5coXddu7cOd71yGXx4sUEAElKStJ4zPXr1xMAZOXKlbztDg4OZPz48Wrtz58/z/te9EEulxOpVEqePXtGAJB9+/axn82YMYO4uLiQrKwsdhtjuv/hhx8IIYo5FhQURFq2bMk77rNnz4iDg4PG31aV4uJiIhKJyKxZs/TuO4ViKFQTR6GYkA0bNsDFxQVDhw4FoNBSDRo0CGfOnEFcXJxa+969e8POzo59z2gGnj17BgB48OABkpKSMHz4cJ6Wyt3dHW+99Rb+++8/NVMtY4JjqFu3LkQiEXr27Mlus7e3R1RUFHseht9//x1t27aFu7s77O3t4eDggA0bNuDevXtax129enX06dMHq1evBiEEgEKDmJ6ejg8++EDrvklJSQCETXydO3fGsWPHAADnz59HQUEBZsyYAT8/P1Ybd+zYMbRu3Rpubm68fd988004ODhoPbepmTJlCh4+fIgPPvgAL168QGJiIiZMmMB+78xv2rBhQ7Rv3x7ffPMNfv/9d2RlZeH8+fOYMGEC7OzsymQ+JYSgpKSE96eJzz//HNu2bcN3332Hpk2bqn0uZCbV9dmhQ4cwefJkDBw4EFOmTDHKMRlSU1MxYcIEhIaGstdpeHg4APCu1dGjR6OwsBA7d+5kt23cuBFOTk545513ACjmWEpKCgYPHsw7R1hYGNq2bau1H1wcHBzg7e2t0xRMoZQHKsRRKCbi0aNHOH36NHr37g1CCLKyspCVlYWBAwcCUEascvH19eW9Z8xphYWFAID09HQAQHBwsNq+ISEhkMvlaikofHx8eO8dHR3h6uoKZ2dnte1cU9iePXswePBgVK1aFVu3bsWFCxdw+fJljB49WqfJDACmTZuGuLg4VrhatWoVWrdujSZNmmjdjxmrav8AhUk1ISEBcXFxOHbsGBo3boyAgAC8/vrrOHbsGAoLC3H+/Hk1Uyog/J1VNKNHj8bSpUuxZcsWVKtWDWFhYbh79y5mzpwJAKhatSrblhGgBw8ejCpVqqBTp04YMGAAGjVqxGunL6dOnYKDgwPvLz4+Xq3dggULsGjRInz55ZdqAjdzfTLXIZeMjAyIRCJ4e3urffbPP/9gwIAB6Nq1K7Zt26YmlPn6+mo8JqB+DXORy+Xo1q0b9uzZg08++QTHjx/HpUuXWF9S5noCgHr16qF58+asSVUmk2Hr1q3o27cvew6mH4GBgWrnEtqmDWdnZ975KRRjY2/uDlAolZVffvkFhBD88ccfgrmxfv31VyxatIinedMF8xBlHOC5JCUlQSwWo0qVKmXvNIetW7ciMjISO3fu5D10JRKJXvu//vrrqF+/Pn788Ue4u7vj6tWr2Lp1q879/Pz8AIDn+M7QuXNnAApt29GjR9nAh86dO+Ozzz7D6dOnIZFIBIU4XdqcimLWrFmYPn064uLi4OHhgfDwcIwfPx5ubm48rVdAQAAOHjyI1NRUpKSkIDw8HC4uLli9ejW7EDCEpk2b4vLly7xtISEhvPcLFizA/PnzMX/+fPzvf/9TO0aNGjXg4uLC+u5xuXXrFqKiotSE73/++Qf9+vVDhw4dsHv3bjg6OqrtGxMTo/GYALTmV7x9+zZu3LiBTZs24b333mO3P3r0SLD9qFGjMGnSJNy7dw9PnjxBcnIyRo0axX7OzLGXL1+q7ZuSkqKxH0JkZmay1zOFYgqoJo5CMQEymQy//voratSogRMnTqj9ffTRR0hOTsahQ4cMOm7t2rVRtWpVbN++nTVTAkB+fj52797NRqwaA5FIBEdHR57wk5KSojM6lcvUqVNx4MABzJkzB4GBgRg0aJDOfRhh5fHjx2qfBQcHIzo6Grt370ZsbCwrxHXt2hVpaWlYvnw5PD092ahTXahqOisKJycn1K9fH+Hh4UhISMDOnTsxduxYuLi4qLUNCAhAgwYN4OXlhbVr1yI/P1+nSVoIDw8PNGvWjPfHFagWLlyI+fPn47PPPsO8efMEj2Fvb4833ngDe/bs4UVdJiQk4MSJExgwYACv/ZEjR9CvXz+0a9cOe/fuFQzqAID+/fvj/v37vByKJSUl2Lp1K1q2bKkmbHJhrk/VY69bt06w/dtvvw1nZ2ds2rQJmzZtQtWqVdGtWzf289q1ayMoKAi7du3i7ZeQkIDz589r7IcqSUlJKCoqQnR0tN77UCiGQjVxFIoJOHToEJKSkvDVV18JZs9nNFQbNmxQ81nThlgsxtdff41hw4ahT58+GD9+PCQSCb755htkZWVh6dKlRhtDnz59sGfPHkyaNAkDBw5EYmIiFi5ciODgYEF/PiHeffddzJkzB6dPn8Znn30mqIVRxdHRUWNqFUChdfvhhx/g4uLC+ihFRkYiMjISR44cwZtvvgl7e/1ubYxmadu2bahbty7c3d0REhKCkJAQbN68GaNHj8Yvv/yCESNGaD3O3bt3cffuXQAKQbegoIDVvkZHR7MP8tu3b2P37t1o1qwZnJyccOPGDSxduhQ1a9bEwoULecdcv34928esrCwcOnQIGzZswOLFi9VM0leuXGFNozk5OawGGACaN2/O+odpYtmyZZg7dy569OiB3r17q333rVq1Yl8vWLAAzZs3R58+fTB79mwUFRVh7ty58PPzY9O6AMDZs2fRr18/BAUF4X//+59aQuXo6Gg2Tcno0aOxatUqDBo0CEuXLkVAQACbLJfxgdREnTp1UKNGDcyePRuEEPj4+ODvv//mRSxz8fb2Rv/+/bFp0yZkZWVh5syZPB9DsViMBQsWYPz48Rg4cCBGjx6NrKwsLFiwAMHBwWr+iFFRUQDUNX/Md9ipUyet/adQyoUZgyoolEpLv379iKOjI0lNTdXYZujQocTe3p6kpKSwEZTffPONWjsAZN68ebxte/fuJS1btiTOzs7Ezc2NdO7cmZw7d47XholOTUtL421/7733iJubm9p5OnToQOrVq8fbtnTpUhIREUGcnJxI3bp1yfr169njctGW7HfkyJHE3t6ePH/+XNNXocaGDRuInZ2dYKTjvn37CADStWtX3vaxY8cSAOT777/nbdf23RJCyG+//Ubq1KlDHBwceN/1xo0b9Y5cZb4ToT/ub/fgwQPSvn174uPjQxwdHUlUVBT57LPP2EhaLuvWrSN169Ylrq6ubETz3r17Bc//3nvvaTy/Pv3v0KGDxv2FHhNXrlwhnTt3Jq6ursTT05P069ePPHr0SO/vBAA5ceIEr31KSgoZMWIE8fHxIc7OzqRVq1bk6NGjOvtOiCLCtGvXrsTDw4NUqVKFDBo0iCQkJAjOHUIIOXLkCNuPhw8fCh7zp59+IlFRUcTR0ZHUqlWL/PLLL6Rv376CkdlCEavDhw8nMTExevWfQikrIkI4NhkKhUIxIsXFxYiIiEC7du3UzFPaKCoqQlhYGD766CPMmjXLhD2kUPQjKysLtWrVQr9+/fDTTz9pbZuTk4OQkBB89913GDt2bAX1kGKLUJ84CoVidNLS0nD27FlMnDgRL1++xOzZsw3a39nZGQsWLMDy5cuRn59vol5SKMKkpKRgypQp2LNnD06dOoXNmzejU6dOyM3NxbRp03Tu/9133yEsLIwXMEGhmALqE0ehUIzOgQMHMGrUKAQHB2P16tU604oIMW7cOGRlZeHJkyeIiYkxQS8pFGGcnJwQHx+PSZMmISMjA66urmjVqhXWrl3Lq2ChCU9PT2zatElv30wKpaxQcyqFQqFQKBSKFULNqRQKhUKhUChWiNUJcatXr0ZkZCScnZ3RtGlTnDlzRmPbPXv2oGvXrvD394enpydat26Nf/75pwJ7S6FQKBQKhWIarMqcunPnTgwfPhyrV69G27ZtsW7dOvz888+4e/cuwsLC1NpPnz4dISEh6NSpE7y9vbFx40Z8++23uHjxIho3bqzXOeVyOZKSkuDh4WExGd8pFAqFQqFUXgghyM3NRUhIiNZayVYlxLVs2RJNmjTBmjVr2G1169ZFv379sGTJEr2OUa9ePQwZMgRz587Vq/3z588RGhpapv5SKBQKhUKhlJXExERUq1ZN4+dWEzpTXFyM2NhYtVQF3bp107sUilwuR25urtZiyhKJhFcbkpFxExMT2eziFAqFQqFQKKYiJycHoaGh8PDw0NrOaoS4V69eQSaTITAwkLc9MDBQ76LEy5YtQ35+PgYPHqyxzZIlS7BgwQK17Z6enlSIo1AoFAqFUmHocuOyusAG1QERQvTyVfvtt98wf/587Ny5EwEBARrbzZkzB9nZ2exfYmJiuftMoVAoFAqFYmysRhPn5+cHOzs7Na1bamqqmnZOlZ07d2LMmDH4/fff0aVLF61tnZyc4OTkVO7+UigUCoVCoZgSq9HEOTo6omnTpjh69Chv+9GjR9GmTRuN+/32228YOXIktm/fjt69e5u6mxQKpYJIzCjAXzeSIJdbTWwWhUKhGBWr0cQBwIwZMzB8+HA0a9YMrVu3xk8//YSEhARMmDABgMIU+uLFC2zevBmAQoAbMWIEVq5ciVatWrFaPBcXF3h5eZltHBQKpfy89vUJAIBEKsOgZjSCnGIZyGQySKVSc3eDooOHKbl4kVWITnU0u1eZEgcHB9jZ2ZX7OFYlxA0ZMgTp6en44osvkJycjPr16+PgwYMIDw8HACQnJyMhIYFtv27dOpSUlGDy5MmYPHkyu/29997Dpk2bKrr7ZUIuJ8iVlMDLxcHcXaFQLJKLTzOoEEcxO4QQpKSkICsry9xdoejBy8xC2AN4EJcLR3vzGCW9vb0RFBRUrhy0VpUnzhzk5OTAy8sL2dnZZolOjZh9AACw9t0m6FE/uMLPT7Ft8iQluJ6QhVbVfWBvZ1neF8zcGNi0Gr4d1NDMvaHYOsnJycjKykJAQABcXV1pcngL50FKDgAg0NMZ3q6OFXpuQggKCgqQmpoKb29vBAerP9v1lT2sShNny0zYehXxS6lPH6ViWfj3Xey8koiF/epjeKtwc3eHQrFIZDIZK8D5+vqauzsUPRDZFwEAHJ2c4OysCGYsKC4BCODqZHrRyMXFBYAiODMgIKDMplXLWlpTKBSLYucVRYqdtScfm7knmqG2BIq5YXzgXF1dzdwTSlmRE4JHqXl4lJYHWQUFSzHXS3l8KKkQR6FQdCK3YEmJwHL7RrEtqAnVeuHe4ypKiDPG9UKFOAqFohNLFuKoDEehUMoKc/uwVvGbCnEUCkUnlizDvcwtQlJWobm7QaFQrBELvrfpAxXiKBSKTiz5PnfuUTraLP0XhcUyc3eFQrEZ4uPjIRKJcP36dQDAyZMnIRKJzJpiZdOmTfD29i7j3lxdnCXf8fhQIY5CoejEGjIRpeVKzN0FCsVmCA0NZfO1WgpDhgzBw4cP9Wq7b9d2tKtn/RH3NMWIldAsvIq5u0CxYayhspVF++1RKJUMOzs7BAUFmbsbPFxcXNjUHeXB0DuJVCqFg4N5EvJTTZyV4EkrNlDMSEZ+sbm7oBMqwlEsCUIICopLKvzPUK35H3/8gZiYGLi4uMDX1xddunRBfn4+5HI5vvjiC1SrVg1OTk5o1KgRDh8+zO6nak41BMbsuX//ftSuXRuurq4YOHAg8vPz8euvvyIiIgJVqlTBlClTIJMp3SQyMzMxYsQIVKlSBa6urujZsyfi4uLUjstw48YNdOrUCR4eHvD09ETTpk1x5coVnDx5EnM/mozcnBz4eThDJBJhwYL5AICGoVXw1759vP56e3uzVZ6Yce/atQsdO3aEs7Mztm7dCgDYuHEj6tatC2dnZ9SpUwerV682+LsxFKqJsxKkMrm5u0CxcZ6+ykekn5u5u6ERqomjWBKFUhmi5/5T4ee9+0V3uDrq92hPTk7G22+/ja+//hr9+/dHbm4uzpw5A0IIVq5ciWXLlmHdunVo3LgxfvnlF7z55pu4c+cOatasWe5+FhQU4Pvvv8eOHTuQm5uLAQMGYMCAAfD29sbBgwfx5MkTvPXWW2jXrh2GDBkCABg5ciTi4uLw119/wdPTE7NmzUKvXr1w9+5dQU3YsGHD0LhxY/zw4yrkSwkuXo6FSGyHNm3a4JP5S7B62WLsO3kZ0SFecHF1w7OcUoFRj1vJrFmzsGzZMmzcuBFOTk5Yv3495s2bhx9//BGNGzfGtWvXMHbsWLi5ueG9994r9/elCSrEWTDFJUrB7UZilvk6QqEAyJeUmLsLWqEyHIViGMnJySgpKcGAAQPYGuQxMTEAgG+//RazZs3C0KFDAQBfffUVTpw4gRUrVmDVqlXlPrdUKsWaNWtQo0YNAMDAgQOxZcsWvHz5Eu7u7oiOjkanTp1w4sQJDBkyhBXezp07hzZt2gAAtm3bhtDQUOzduxeDBg1SO0dCQgJGT5wKuVcIXAB07FkVfu5OcHR0hLuHJ0QiEfwCAhEU5K3IDZeTrXf/p0+fjgEDBrDvFy5ciGXLlrHbIiMjcffuXaxbt44KcbZKiVwpxMVU8zJjTygUy8cagi8otoOLgx3uftHdLOfVl4YNG6Jz586IiYlB9+7d0a1bNwwcOBB2dnZISkpC27Ztee3btm2LGzduGKWfrq6urAAHAIGBgYiIiIC7uztvW2pqKgDg3r17sLe3R8uWLdnPfX19Ubt2bdy7d0/wHDNmzMCs6ZPRdMd2tGzXAd1694NPdG2j9L9Zs2bs67S0NCQmJmLMmDEYO3Ysu72kpAReXqZ9dlMhzswkZxfiblIO6oV4IcjLmfcZ15lcIqXmVErFkl3ALwVTYuHRDZbdO4qtIRKJ9DZrmgs7OzscPXoU58+fx5EjR/DDDz/g008/xdGjRwGoVxQghBitKoWq+VMkEgluk5cqMzQt0rT1af78+Wjyeh+cOX4EZ08cw5rlS7Fq/a8YO2Koxn6JRCI11wyhslhubkrXEqaP69ev5wmZAMpcE1VfaGCDmTn/KB1jfr2CPj+chVzlIcm9aK88y0SRlObBolQcu0rrpjI8SMkx6fmepOXh87238aKMiXupTxyFYjgikQht27bFggULcO3aNTg6OuL48eMICQnB2bNneW3Pnz+PunXrmqWf0dHRKCkpwcWLF9lt6enpePjwodY+RVSPwvCxk7Bu+x507tEHu37bAkAhRMo4vubM87aKrx9SUpLZ7XFxcSgoKNDat8DAQFStWhVPnjxBVFQU7y8yMrJM49UXy14m2ACujgop/VWeBJISOVwclVK7quLjw53XsebdpmU6jzFXUBTbQPVymbX7FhqFVkHtIA+TnG/oT/8hNVeCK88ycWjaa1rbxr3MVdtGZTgKxTAuXryI48ePo1u3bggICMDFixeRlpaGunXr4uOPP8a8efNQo0YNNGrUCBs3bsT169exbds2s/S1Zs2a6Nu3L8aOHYt169bBw8MDs2fPRtWqVdG3b1+19oWFhfj444/RuEMPVA0Nw8vkJNy5cQ293+wHmZwgJDQMBfl5uHj2FEK6tEVSngywd0KLNq9h7erVeK1tG8jlcsyaNUuv9CHz58/H1KlT4enpiZ49e0IikeDKlSvIzMzEjBkzTPGVAKBCnNnpEh3Ivi6WyeECpRCnqj4+dDulTOeImH0AADCxYw3M6lGnTMegUACg36pzuLewh0mOnVqarPdesm6N33WBQJ+KKlpNoVQWPD09cfr0aaxYsQI5OTkIDw/HsmXL0LNnT3Tv3h05OTn46KOPkJqaiujoaPz1119GiUwtKxs3bsS0adPQp08fFBcXo3379jh48KCgkGVnZ4f09HR8Nn0C0l+lwbuKLzr37IMZsz9DYkYBGjVriUHvjsInk0ZjXGYGJnw4CxNnzMZHny/CV/+bhvbt2yMkJAQrV65EbGyszr69//77cHV1xTfffINPPvkEbm5uiImJwfTp003wTSgREeoNrJWcnBx4eXkhOzsbnp6eRj8+IQSRcw4CAK581gV+7k7sZ5n5xWi88CivffzS3gafgxHiyro/xTZZf/oJvjyo7jBsqmvIkOt0+8UE/O/PW7xtm0e3QPta/ibpG4WijaKiIjx9+hSRkZFwdnbWvQOlwrj5PIv3voqrIzIL+HkvG1Tz5rWrFegBZwMCRMqKtutGX9mD+sSZGZFIBEc7xc+gmguO+vhQzIm3q+UmmH6Umqe2jWriKBSKLvS5S1jTnYQKcRaAg53C+UhaohLYINCWKk4pFYWTkVaiqblFOHInxahC1svcIrVtNPCHQrEMevbsCXd3d8G/xYsXm7t7urGixyz1ibMApKUPt1f5EoT5urLbhTRxzzMLEerjqradQjE2zvbCazy5nEAs1j9Ipvf3Z5GWK8HCvvUwvHWEUfrGaK+5JGRojyCjUCgVw88//4zCQuEocx8fnwrujQoCAlphseoC0HqkOCrEWQBMZYYBq8/zfYEErqNiWn6LYkLS8yTYdeU53mpSVaMm7lpiFpqGV9H7mGmlAQuH76RoFeLEIvWIbCGKS+RIyVbXxFmy+ZdCsSWqVq1q7i4A0GS5Ut8Wl5qro4XlQs2pFozQA43WUKWYkqk7ruGrw/fx3sbLrJlflfFbrpTp2NwyckL4uDmyr1UTDXMZtekSLjxJV9teRBNiU8wMdXexLJ6lq2vn9fKJq6Cf0RjXCxXiLBgicLmp+s1RKMbk3COFcHQvOUfj3e5VXrHwBzq4HJ+J7ELNwlnL6r7s6/R8ic4+qpJn4bVdKZUXJsWFrqSwlIolp0jz/UYbFfWUZa4XffLQaYKaUy0YRhPnaC9GkKczEjIKDDanWsLKkBCCQ7dTEFPVi/rzWRHcK2dSxxpYffKx1vaJGQUI9nKGvYC/GsP4LVewY1xrwc/sOX52Upnh1+2uK4mY3CnK4P0olPJiZ2cHb29vts6nq6srTa5uAZAS9QVnSbEckMm0PhslRYWwJ6ZzzyCEoKCgAKmpqfD29i5XaS4qxFkAYT6ugk7ZzEUmFnEiWA0W4srfv/Ly981kTP3tGgCap87S0eSX5uGs/YZ27O5LvL/5CjrW9semUS00truRmA25nOBSfAZqB3qgCseEyk3gWxa3gUBPmp+LYj6CgoIAgBXkKMYjX1KCzAIpPF3s4anjXsQlNVM9uCLXwQ5FJTKtz0Z5jmOF5Inz9vZmr5uyYpAQl56ejps3b6Jhw4bw8fHBq1evsGHDBkgkEgwaNKhCaqqtXr0a33zzDZKTk1GvXj2sWLECr70mXKInOTkZH330EWJjYxEXF4epU6dixYoVJu+joSwdEIN3fr6IqAB33nbmIhNBBIdS7caz9Hy04piddKEp19yJ+6k49TAN/+tVF44aohCNxX8C/ksUy8ReLGa1vdxLp3UNX9Twd8PjtHwAgKREBid75U3ul3NPAQAnH6RpPX6hVIYjd1MwYetVtKnhi+1jW7GfVfdzY31YyhLAU0L9RSlmRCQSITg4GAEBAYIF0yllp/Oyk+zr4x911Hu/9/ecVNvWNsoPV59lolBLSqKFfeujTaSfAT00HAcHh3Jp4Bj0FuIuXbqEbt26IScnB97e3jh69CgGDRoEe3t7EEKwdOlSnD17Fk2aNCl3pzSxc+dOTJ8+HatXr0bbtm2xbt069OzZE3fv3kVYWJhae4lEAn9/f3z66af47rvvTNav8sLUS1VNYCrnaOJyixT+PrN238KQ5upj1YSqVkUmJ7ATizBq02UAQISvK0a2NW2BXksw6VL0QywGUHpv4/pkNgr1xvLBjdB31TkAwNb/EjCmnfK6ERtgOvr7pqK49PnHfOGeG5gg1REEIQRN9kuxBOzs7IzycKYoeZGrFLgMqYjRIioIf157wduWWiBHaoEcOUWahbiUfLnVVN7QWwXz6aefYtCgQcjOzsb//vc/9OvXD507d8bDhw8RFxeHd955BwsXLjRlX7F8+XKMGTMG77//PurWrYsVK1YgNDQUa9asEWwfERGBlStXYsSIEfDy8tLrHBKJBDk5Obw/U+PupJSlzz16xb5mNXEiERqGKvt/63m23sdW1cSpJkR9kSWcy8eYUBnOeuAKY8zvVjdYUfIlxNuF/Wzh/ruQc4QmbTKcp7Py+h7dNpKX440r4HMjTvXVxC3sVx8Dm1YDUDY/OgqFUnkRirBPSC/Q6rcLADIremjpLcTFxsZixowZ8PDwwLRp05CUlISxY8eyn0+ePBmXL182SScBoLi4GLGxsejWrRtve7du3XD+/HmjnWfJkiXw8vJi/0JDQ412bE3U8FeaUePT89nXjAAmEgHd6ynt5vdT9BcsVa9FVRVyRVig7iSZXhCmGAc7rhBX+p/ZwmiMGfKKS9i0IVzhL0ElrD+A46vm7mwP7m115+VEwX4s2q9es1WIhtW80L+xIidVidx6zanXEjLReslx7L+ZZO6uUCiVBgcBYS3Q0xl2upKVV0Yhrri4GC4uipW4g4MDXF1d4eentBn7+voiPd10vk+vXr2CTCZDYGAgb3tgYCBSUlKMdp45c+YgOzub/UtMFH7IGBOxWIS+jUIA8DNHcx+i3EinEgPMRqqaONUAClPXZy0oLsGtF/prDinmhatRI5xFBKBewWHDmaeo9dkhzNh1HaceKn3hlhziC2Bcjd3PZ57wIpQP3haeu/qmC7EXi9moVm3pSyydiVuvIjm7CB9sVwQAZRUUo9fKM/jptPaIYIrtsPJYHNp99S+bPJuiGyEhrlgm50XCC2FNWn29hbjQ0FA8efKEfb9jxw4EBwez75OTk3lCnalQDdsmhBg1lNvJyQmenp68v4rApTQSZtcVpdDIRqeKRXB3UmpBDIncUxXSVB90phbicgpp7i5rQqicFjO9VE0QK4/HAQD2XOX7nOSrlLDhLjoKimW8VTDXvM81fWiqg6rqX2lvJ2L7/DJHgqyCsuWwMzeqc/qn009wNzkHiw/eN1OPKJbGd8ce4nlmIVadeGTurlg1sc8ydfrwWlNSfb2FuKFDh/JCp3v37s1q5gDgr7/+QosWmlMLlBc/Pz/Y2dmpad1SU1PVtHPWyKs8xerq4cs89kGkjE4FOtQKYNvqynzPRVVpp+owbixn8NUnH+HNH88iVyW5Ik2VZF2IBcypXCZ3qqHzGM84LgGA+jVWwBHy7icrTe3cZun5xYIBMarHsheLePOh2aJjVhmlqroQlZQhsKMsSEpkPE0pxfKhATz6oymoTpcvOLPwzCooxiOVklyWht5C3Lx58zB06FCNn3/66afYvn27UTolhKOjI5o2bYqjR4/yth89ehRt2rQx2Xkriqocp3FGXc7MVbFIBDuxqEwO3Ol5fNW76r5PX/EfuGXl68MPcPN5NjZfeMbbTmU464K3QuWkuGGY1FF3Mt1n6QU8v01VbS9X6HqjYQj7WvXhFDnnIBYf5JtmVR2O7cVinmavRE6wRkdSYktEVQGqy2XHGOQWSREz/wgGrbtg+pNRjIZQJZ/KzNm4V7obaaCs3xSjiWvx5XF0WX4aD19ariBXrgRh586dg0SiEBJcXV3h5ORklE5pYsaMGfj555/xyy+/4N69e/jwww+RkJCACRMmAFD4s40YMYK3z/Xr13H9+nXk5eUhLS0N169fx927d03az7LATSzIXHjMQ415rjL2fUNUvf1X84M+VPe99DTDwJ5qR5MZzBKRlMjwICWXpkDhIObJcPzrDwDcnPTLStRjxRn2tapwxjXpM0mENWmDfjr9hPde9VgikXp6ky3/8RcS1oDqGAzxey0rZ+JeobhEjthnmSY/F8V42NLtKqdIinc3XCzz/sx31TaKn1vVSUdu1PvJCqGNiZIvjyBpaspVsaFnz564fv06qlevbqz+aGXIkCFIT0/HF198geTkZNSvXx8HDx5EeHg4AIVfXkJCAm+fxo0bs69jY2Oxfft2hIeHIz4+vkL6rDece/iDlFzUCvTAimMPAShrVTI+Qzc4me11oeoDp5q6Qcjxszyo3WAsWBU3fMMlXHqage/fbow3ORohW4ar1eKa88uDqiYuT6K8JteeeozpXWpixq7reh1LyJQU6uPCe29NCwkGVc2bu57CcvnOacGTk6IRG5LhUFSsPpeZXKf6wCxE64V48Wou63JXcHXiR+JbcuR7uZ7g5tBgTJo0CfHx8ZBIJIiNjUX79u3ZzzZt2oSTJ0/y2hNC1P4sToAD0KGmP/uayRV35O5LXhvGb648FRBSsot47yUlxn3gqT6wRRYsxTFayO0XrU9zYyqE8sSV17FRVfDKl/CvuW//eYCDt/SLMD+qMicIAYK9+EJcTpH1BdOo+sTZi01bReX4vZe8dEYU68GWNHFC+dyepefj4ctcbL4Qr9P/lfmuXB3t8HYL/dOFqVoGKkIzXlZMe6eg6E2bKD+8VlMR3asa3cdwPyVX6+eqCGktlh99yHtv7GtT9Xianv/mNmFyJ7++qzpbQEh2KOu3k1Ma5KJ6HZ59xDdN/Hz2Ke894/vJcIwjuP3vz1vs69dq+rFaOG71CEM4/+gVWi4+xjuHKciXlGh1kFb93o2sIOcR+ywTY369gqWHlJGv2QXWm57F9rBcgcLYCD0nCqUydPvuNObuu4PtlxIE9uLsX/pfBBGWDGiAte/qV1EqS8WCVWLBKUfKdatYt25dpYgMtRR6xShSthRqENKyODfaawmZalGAquQXV7xGQqaidhYSAK4lZKLF4uPYd/2FwKf6kVVQzCuYbihcX8BXudaZlsIU6IpONYSfzyiEM0MXCp90r817P33ndfY1tzTXljEtWQ2WSxmLVQ//5RJe5kjw/uYrZdpfX7ouP4Uuy0/jSrywD6qqadPOhJq4BynqwqS2OpIUy4K5xVpTGoyyInTv6P39Wfb19YQsrfsrqx4p/uvrPqRaA9qSI97Ldad455134ObmZqy+2DzMg+jYPWGtADfZav/V59GHczELYQ5ll6raWdVMlC8pwcStV5GWK8G0HdfLfJ6O355Ev1XncCZOe8F1TWRyBOJGod5l7kdlg29OVQ9s0MWxGR3Y14WliwhGE9c0vIpex/BydUD80t7K4+ghYDg7CN/Kbr/Ixqd/3kJytnBKgYpK15BU6sZw4JaibuzN51m8EnvqQpzytbG11r+ce6q7EcViISDYd/0Fan56CAdK6xBXVnRF4urOc1p6Dyt9pyrEvdWkmkZLTEa+cnEvrUzm1KKiInzzzTfo1asXmjVrhiZNmvD+KGXHsVRI46Yb4eLt6sh7nysp0boa03Tz//KA6aJzdT0U5/91xygh8oxWsqxmMK5Jz9ZC9rXBj05VYIg5lStM3UhUVOpg0oIMblZNcB9V7FQEGn0ELUaLrcqMXdex7WICz3RoTpixvPnjOQz7+SLro6r6MOJOa00phbILpOi87CS+U3GR0MWj1Dy1bbag1aksyAnYBfDk7VfN2xkTo0tGO1GqMfsj9jlOPEhV+1yXJi7c1xX3F/YQPDa3upHQnLEUDA6BGj16NI4ePYqBAweiRYsWRq2WYOtE+im0mpoiZ34a0RTtvjrB21YolWlUEWuaAOvPmG4lruuBe/BWspowWh7KKn5xv7JcK3SENxVCgQ2qc3xU2whsPBfPvq/i6sBqNt0clbeUS/EZyCooZp2EXRz5t5vq/m54kqbuEqDtnlIvxBN3knLwTsswlWO5C7Z/+FJx8+WWBTMnJXLCW1y9yCpAkJczvFwceO2evlI+NPIkJfCxV58zm87H43FaPlYej8OHXWuVu18U68CWAht0jdXBTownaXmY+fsNAOBp8Ln7M/cUR3v1VD7c56eHkz38PZ3wJC2f59bkrTI/LQmDhbgDBw7g4MGDaNu2rSn6Y9MwuWskGsxH1aq4orqfG55wEvQWSGTwdBa+wLir+0961MbXhx+UuW+5RVK4O9nrFNpVhThVbWB+sQz+HurHIITgwctcRPm7C0YkaSI9Tz9/tpTsIrzMKULDUtNpYbFSUD6koX6nLeLM8y3jmyIYhjQP5Qlxm0e3xBs/Kkz73q78azEjv5gVEJqEefM+axJWRVCIY7SBvWKCcPBWClvgHlCuqLvWNcwXN19SghP3U9E4zJu3iAj2ckaySsS2KSmRyXkCE/MywIOfYzOIE3FbUFwCHzd1IU418XF5+0WxDnZffW7uLlQYusylxSUyvMzRXEtW1cqiqvBQve7bRvnhdKmLzm1OzW9LXuQYbE6tWrUqPDw8TNEXm4d5gOZKSpAvKWEFDqE2DK9UKjIQQhAx+wAiZh/g3eQntNddLkkTt19kI2b+EXxUutrRhj4Xu5CQtuHsU/RYcQYf/3HToL456kjayNBqyXH0XXUO8/+6AwC49SLLoPPYCkL5yVTl9kAPZ977mGpeWNSvPjaNaq4m5HOFehcHO14AQoi3C8J9XQXOpzhGy0hFgk5u3jfmnmtoRLFURjBq02U181NF5GPjUiSV88qOMVpK7rw+fDuZV9VCkzlV1excHiz5IUWxXXRdlTlFJSiUKi0p03Zc4++vYk6tVoV/v2HcCKqULj671Qtk5+eXnGoxljw/DBbili1bhlmzZuHZM5pby9j4uitX21eeZaJaFcVq/PM+0ex2JxUHbtU8b88zlQ7c3AhMocLmDLpW4WtOKcoYqRY6F0I1v47QpS9k/v2xtKjzn9cMi1g11Ol70/l4fLb3llq5MXOnPLEU7DlF6DXdt1wclQLHwn71AQDvtgpHx9oBam27fneafW0nFvGCFOxEIjxLV/qdrB7WBCdmdlSep1SwOXQ7BQVskIScPZYqi/vHsK+vJQhXIeAm/ASAOI6viyFJtMtKoVSG45zApZQchRaQK8RN2HoVa08pS4dpyuem26lbfyw5hQLFdtGnru/HvysX/vuuJ/EWQNwUIwDUNNrMAumPiW0w/41ojb61qlkXLAmDhbhmzZqhqKgI1atXh4eHB3x8fHh/lLLjZG+HGv4Kv7jCYhl7BdpzHlj5Er7/FjflAsDXmjC+XrqUFlGfHtKa9NfZXv/0DaorFqHnjIOdeofKmj3+oh5lw/JUvrOt/yWoaTcseKFVoXAFbCaxrmrCZm7JmuYR+kWcAoCrik+cqiwfU9WL9QsFgMYc8yvjWMxo9oSEOK6f3LXS1AOqZkpVuJq4xMwCLS2NQ+yzTN6cZeaHm6PmOfYqV9hc5OFsPC2iJWekp9gu+qxT0vP5LjXca1lVE6cKo4mr4e+OkW0j1SxdynaW+4Aw+C7w9ttv48WLF1i8eDECAwNpYIORCfZyweO0fOQUSQVrV/asH4yHL+PY96pltbgPYSbyTZ/f6Pcrz/Fuq3DBzzSlbxDiXnKOzjZCmriy5ttV9Wc6E5eG5KwiDG6uzM697IhuX0BDSrlUZrgLBtb3UuVrEYlEmNG1FlJyilA7UN214s6C7qg37x+17aqmb7FYhNqBHnhQWlzaXkW4rxnogWpVXPA8s5AVfLQJcQAwqGk1/B77HEWlixIXLcKRYizK1xWhjfJ3d+K5OTBaA21+oJpcBuy1XK8Z+cX45I8bGNQsFN3rBenslyWbiyi2S1kyBxSXyMG4vQrt7+Jgx1oE9J3zllzKz2Ah7vz587hw4QIaNmxoiv7YPMzD6fLTDMHalaoPpZXH4ngqYO7DbeVxhbCnT4oGVWFQqE/64K6iHRCaREIPn7IuBvxVNC3DN1wCoPDTqhvsCUARfq6LV3kShGhI7WJLcCN1mTxvQr/M1M41NR7DTU8/M0KAGgFurBDnKCDIKIJ2CvHwZS5aRPqwApCma9K1dH4klqYHEDomF665piIEGalMDk/OHEkr9WlVFWC5xD7LRN9GVdW222kZ2/KjD3DsXiqO3UtVi9jT1C8KxdJQnZI1A9x5LhBC7L+ZrFRICDxD570Rjdl7FJVfLgsk324Z6aNm4TkT90qtnaVgsDm1Tp06KCwUTpxJKT+MnwtPGOIIOKrOzPo69utCmxCnmv5AG/oIjEJ9LqsSTJMvG5PclRCiVwqRikr6asnkFElxgVOXNymL0eQafqzqfsJJwIO9lEER95JzWH9PR3sxfN3VTZ/MfNh1JREAICtdOWty6meurd8uKdq76hAouYJbRURoSuVyhPoonauZ6GoHLRUahNKjJGYUYOHfmvM9ZuQbVoWEXv8US0T1/v7P9PaC2n8uJ+4r88WxPnGc28XQFkq3i/EdqqvtP7y1QgDk3qssGYMlgKVLl+Kjjz7CyZMnkZ6ejpycHN4fpXw0C1f4FRZJ5YL+AI1V0jQYUi6na7TmtAw/nX6i8TNDNHHFqjnuBMagqnUjhGgNE9fGKw0pRpjv7kWWfgsOak4C4l7yyzFFl2oyVX3i9GH5kEa89yuHKt4fnPoau00qkyPYywXxS3vj4aKegsfpUV9hCmSuGV2aOG6+uOwCqU6HZG5wQEVcA4kZhcjkCFhZhcX4+cwTrT6pqr6EADDm18solnF9fxT555gAEHstQqGQdpIGNlg2qr7QqhSXyHE3KafSBWhxp+TRD9tDLBZh+ZCG6NMgGN+/3Vhwn/ucsnJs1RmVe9iDRT2wf0o7DG4WClUYIdGSTahcDDan9uihyG7cuXNn3nZCCEQiEWQy6xi4pcL4n73IKoRL6Wvu5Regkt5BNZO0tjnsJ6Dp0Adtvjeq3E/JZa8FQDg61VclQujm82yBVuWDKb6uSwvn7CBGkVReqTQRcjnBtcRM1Az00JhDUAhV4fqShjqf+tAo1FvQjFeF89vro/hqHlG6qCkN+9flEze4WSjmlJpKcoqkgsJJkVTGOjDLKlgTByj8Txn2XH2BPdAekZ0ksBBhkhgzyAkwdvMV/Hs/FWc+6cQzz35/PE6r+RvQromnmB9uugshxm25gpMP0vBl//oY1lLYt9mamPn7DfwR+xxzSzX13q4OqFkqXNUL8cKP7zTR6ALAXbgLaeIARRBh/apegvsz9wZuaUZLxmAh7sSJE7obUcoMcwGdfpiGLqUJTbkXoLMjfxVtiKlLKCeXPnCdrrkCmiZuPs8WzHHHoJpGJKOAr03jPmT1IbtACi9XYWGl58ozWvd1d7JHkbS4XEJcep4E47bEonPdAEzqGFXm4xiLw3dSMGmbIh/azfnd9BbkNEUImyp2SZ/FAXMdMH5zzO+kaV87sQg+bo7IyC9GoVQmmIYjt6gEzg52IITwVvoVpY3NLDDM1JldKMXLnCIEemo270hlcvxbakb6PfY5GlT1YlMCLT/6kCfECfmpTt95Hf0aq/vdUSyD01oqjsjkhC3YvvFcfKUQ4hg/5i/2K1wGhO5NmioVcc2gZVFManqWWGrwm8FCXIcOHUzRD0opr9X0Y18zSQy5qmBV4YYQhQaBEbS0RfOMahuBS08z2Ju9LiQlMuRLZLwHplRG1EqXAIp0B4zWiyuU6ZpEQsKToUJc9xWn8d//OvNMCUKLNDdHO+QX8zXFzM2hPCkWhv18EfdTchH7LBMDm1ZT05ZWNL+X+o8BQIP5R/B0SS+9Akc0+ZkZW4hbPrghfjr9BF+91UBnW67WtkQm16mJA5T55ZKyCgUFM4XJ0Unt2qsoIa4sQQTPMwt4QpxIxJ9b3L6LoPlBBGiek5n5xTxNKcVy0LbgOXgrmX1d2cypDIbcggI480SpidP/CJoWvRn5xWqBdJaAwT5xqubSixcv4vTp05BKrUP1aOlw82QxmaN5mjiBnG0FHNv9U4EyRhGlGjgnezusHtaE3c74PGmi36rzaPHlMV6Zr8JiYXM5V+hS84vTglQmR0I6Pz+XIX5+gDJhKq+ckZzwfI8A4L02EWr7MjfH8mjiuD4YeRVch5UQgoT0At7N+3I8P9HtTE4yTG04CAjnQNl84rQxoEk1HJ7eXqugwRDEWVUXlegnxDFO/ScfpAkmC91w9ikA9bJVFWVOvf3CcN/h4hLt5YNinyl/86ISmZqgpk/SVOoXarlou965VXsSMypn0KEmIezx4l5q26QlfF9RwDAhEADealJNbZsuv0RzobcQl5ycjHbt2sHJyQkdOnRAZmYm+vTpg9atW6Njx46oX78+kpOTdR+IohWRSMT6xTF+KtwLUChRboFEKfQM+ek/rcfnJmoN8OSvKlRXcfeSc1AiJ4jlCAVnHwmHWnN35fqhaco2zyCVydWiVZmoSG37qJ+fYA+npuCRuylovPAor42QSt7OjtHEGecBli8xjU+oTE5w4XG6WuLidaefoP03JxA55yAepeahoLhErY2+tRYrShNnCNzr9VWuRGdgAwBEhygWJ7lFJYK/a07pvFJVvj54qT11QXlgkngD4AUk6IuqpthBZfxLOD5T6XnFakLc2tPKChCarnRLNBVRFGgLVOEuQF2d9LdgWBOa7kFC1yyjGOjzwxnsv5msdX9NCFmbDFUuVBR6C3GzZs0CIQR//vkngoOD0adPH+Tk5CAxMRHPnj1DYGAgvvzyS1P21WZgJqxQcXChFcmVZ9od0Ln7cF+rBhRwSyBpUstnFWry51G2n7vvNvt64tZYrX2TyoiaUKaq7XuRVciLFBKaTOcfp/O+r2P3+CbjHvWCUDtIPTSd+a6NFdjw3sZLZd53+8UEvPfLJaw//UTtO9h47ineXv8f6s/7BzsuJbCfcxMZd1l+CvP23Snz+S1RD8O9Xvdce6GXJu6NBoq8iUUlMkENFJPHTpdgZEmozhHV5MDc/JFCvkI/n3mq8xwJGaavWEFRIpcTnHv0Co/TdC8etAkh3GoCWQVSvbSu1oYhUzMjvxjrzzzhabwNndlCc+h+imVm39BbiDt27BiWLVuGN954A6tXr8aFCxcwb948VK1aFaGhoViwYAEOHTpkyr7aLLpWEbrSc2jaXTXfzqbz8exrrlDDNTtpMqdyZT5uOSBuhM+xGe3V9muz9DjmqggeXCHtTFwa2i79F/1Xn2e3FQn04XlmgdbSKMsGNxSsi8cIA8ZKscDNz0WI4iY9fMNFrOPUwhTiakIm/vfnLZx6mIYvD95D3bmH0fiLI/hw53XkFkl5EY2z99xC3bmHETH7gNqYf9cjsbEmjFmL05hEBSjShkhKZPr5xJUKNA9TclEs8LteLS3JpaqJUzWvWhKq5lTVaFJuourbL7LVBHLuokzTAq3fqnPl6yQF5x+9wi09o+2P3XuJYT9fROdlpwSTznLRVpZQ1Q86VkPdYGvGEJeO9PxiHL6Twt/fQFWckBB3+HaKQEvzo7cQl5mZiapVFdFLPj4+cHV1RXi4MgqmRo0a1JxqJJjcWAyqF/DvE1pjfPvq6NsoBEDZ89nY24kwqWMN9j13tS/jBQkoXx+5oyzezYV7G5Fo8ImLCvBQSy+iWvsVUDyEGBjz7b1kZQ6kLIF0COvPPNXqMO7mZC/44C+LT5wu36lVJx4hYvYBRM45iGE/X8SZuFdYcui+WjoYLneT1Fd5mQVS/HntBWLmH2GjM8sKk/xYG5piO8pa19ZYdK+niNI+dvclG02qyfQLKBcU2YVSnr8QQ9VSgUdVaLPkrOzPddR15foB3XqRrSaocS9vbVe6qh8pRX9Ssovwzs8X0XfVWb3ab7+UwL7mLqBVOXQrGXe1lDNUlcmn77iu1/krmmfp+Zi9+yYu6VHvWhVDleRpKvWGDb2FCQlxmp5r5kZvIS4gIIAnpH3wwQe8gveZmZlwcxPO0k4xDBfVyEyVC7B5hA/m9KqLKqUF4uJfafc707SICfNxxSc96uDj7rUBANsuJrDCDPeBbqjPhURAMGPQJ+p0zzWlNokrwDK+REKBE49S8wQf2NroFh2o1MTpGZ2alitB8y+PYfZuzcEC3/wjXKtVW+JhU+epW3dKczJnBk2RzeYuj8zkRHvMMZdrK+1VtYpCSNMU0ZcnUSwCVL/zjPxik+VLK++vu+iA9jxhunw6GUFYF6r+lBT9eZGlELTlRHeU6J/XnrNpQQBh6wLDxNJ0QZpQvY71TXBe0Sw9dB87Lidi7OYrBu9b3hrthu7tKOB7rskKZW70FuIaNWqECxcusO+XLl3KE+LOnj2LBg10pwyg6EY1e7umC5BZGei68aru/+voFujXKASfdK8DQJkYFwCuJypU8Zo0cZq0ftybVqFUxtauVMXDWXdWG67mJyOfG3mlOKam7PaHNKi7J3dSahu5CY8HNq1msCau+ZfHkFkgxY7LiVh5LA4Rsw/otR+gEACXH3mAt3/6DyceKHz2CopLEDH7AOb9VXZfNn3QR8DV9Nwxt6fYAIH8ZdoWAx6lKQJyOfPi1Mcd0afUV+6/JxnIk5QI/uaWqonyc9ee+kM1bZDqyCI4Ue/M7zy0uTJbfZXSSGFryVJf0ZTI5Jix8zq2/veMt73vqnOImH0AEbMPYOJWpbCllr5GJsej1Dz2Pvnhzhu8zzVpebRV8mCwVDcIVRhtoq6FkpAArK8M16VugPAHBgqBEgFry2MBH3VLQG8hbt++fZg2bZrGz1u0aIGVK1capVO2Tqfa/AtR0yqEqWlqSE41AOhQyx8rhjZmUzwM5IRT55RGlnJvQtybRE6hsMCoOu1OxylWmU4qkadCZU5U4a54dnF8wZgkwRn5hmlLuOHi35eWfwIUwuK9ZIWZUh9T2loVv7bvjj00qB8zf7+B7/99hAtP0jFq42UQQvTO2VdeuL/PtYRMpOaoRwBrehiY25zarV6Q7kYcmGuOGyUd5uOKGpySXM/S89mFiqOdmBXuTRaBVsbn7N8ftAOg/hsIpUDQdr6X2eq/t5ijqWTuIbdeZOPm8yxM3nYVs/64ie+OPsSOSwk4+SAVD1JykV0grbS5yLRx+E4K9lx7gc/2KoO2ikvkuJGYxb5P5Zjw1p1+wvueFh24hy7LT+HgrRS2NBqXs49eIV1loSUpkeELLfVxGawljkGo3JsQQuMp0FML1jLSV3C7oXewVgLHMdTSU1EYp3o6gObNm6N+/frGOpxGVq9ejcjISDg7O6Np06Y4c0Z7Rv5Tp06hadOmcHZ2RvXq1bF27VqT97G86JM/C1Bmpta1WtOliq7JCXA4elfh88aNcOIKdHeTcwSjn5j7FZOTjvF1UxUwR7eL1NoXQHEzFFJdR/q5l36ueCCp+tdpgltPs0WkUnssFitNtNp8UgAgt0iKpYfu63U+fSGkbIlfy8LpB2kghODvG0nov/o83vn5olobTQ8Dc5tTDU19oZqyBlDMgQkdlBrZvKIS9joWi5U36F91XAcVSctIH7iUVmhRvU4YzVmz8Cp6HevXC8/UtENcc3NyqZBXJJVj47l4HLiVjJ1XErHyeBxm77mFkRsvo/uK02j4xRFEz/0Hr397Eu+s/w8zdl3H14fvY8uFeBy9+xK3X2TjVZ6k0kVICi1etblgfPPPA0TOOYh91xULT+b+suzIA0TP/Udwn5m/87Vzk7ddxbaLCYJtuViLUK2pwoIqQhryHD3dHN5pGcZLmM/A9bPWB1XlA6CfFckclKlXL168wLlz55Camgq5yoU8depUo3RMiJ07d2L69OlYvXo12rZti3Xr1qFnz564e/cuwsLC1No/ffoUvXr1wtixY7F161acO3cOkyZNgr+/P9566y2T9bO8qPrEaXqEMRearuS6hqwgtl9MwOL+MTytjOqkKpTK1HySmBuJv4cT4tMLkFVataGsN5hhP/+HPZPa8rZ9dfg+3mpSlX34BHg6I12H+WuIiuaPKxAY4meRYQIz2/Sd1/HXjSSjH1eIXEkJJm69ykZtPUrNw/lHr9AmSnnDY35zR3sxOtbyx5G7TBCLuQ2qfFR/U1U0lafipuFYdvQhvi6tGGEvFkMkkoMQvmuBufF1d2QffJoir5uGV8EVTqJfQHGNM/6NgZ5ObPR6nqQEnpwH0ZDmodh84Rl61g+CnViE/TeTsfL4Q7Z9h1r+CPF2QUp2IVJyJEjJLkRmgRSFUhmevMrnJQFXxcFOhEBPZwR7OSPIywXBXs6c94r//u5OaqlSrIG31pyHnUikl//ktB3XMY0TaKDqcM/lxAN+aS3VNEmaeJGp7gOXmlPEq1xgCQgtroQQsgho8/lkrnEnezHcnOwxul2kmmUlSUATrQ0nB2VfP+paC8uOPrTYwAaDhbiNGzdiwoQJcHR0hK+vr1oOMlMKccuXL8eYMWPw/vvvAwBWrFiBf/75B2vWrMGSJUvU2q9duxZhYWFYsWIFAKBu3bq4cuUKvv32W41CnEQigUSinGg5ORWfG4YxkzJokjWYSXHsXqrWmqb6qKIbVPPCzefZqFmazoG7j6qgJCjElf5n8q7Fl+ac8/NwYk20hsCkgeCSlivBs3RlKpEwHxf8r1cdDN8gnJtt3fCmaBvFX5VxvyNtEY6qPBe4UZaXihLgGFTD7mfvuYXTn3Ri3zMCd4iXM880ZG5NHACM71CdDc54t1W41rbuTvboXi8Q/2iIpAaAm8+z2AeDWAR88WY9fL7vjskCTPQ96ndDGuLOixzcS8nBhA41WCFHY4JgEVDd342XI1EmJ6xmvH6IF9JyUyEngEQqA+HM22AvF9xf2ANO9mK2qgc3XVH3ekF4pyV/cVwklSEluwjJ2UVIySlESrZCuFO8V2x/lSeBVEbwPLOwdN4Ip7wQi4AAD2cEejkj2FMp3AV5OSPI0xnBXi4I9HKCk0CVGlNyLzkHcal5WLT/Lm8eMMQ+Ex6PPuSaIHBkj0otagCYu+8O1g5vavRzlQd3LcFIXAz18Vs5tDF+PvMUk0p9n90c1c9jaDUWbunEl6WWn+ISuV61wysag4W4uXPnYu7cuZgzZw7EWrJIG5vi4mLExsZi9uzZvO3dunXD+fPnBfe5cOECunXrxtvWvXt3bNiwAVKpFA4O6mbLJUuWYMGCBcbreBnwU6nPpumaYfJnAYp0FD4azIv6lMF6t1U4PvnjJqupWn9GczRjUlYhL0AAAPuUYsw+zITV1w8CAN5sGMITbISczLkO6fZiMV6r6Q8vFwfeqjjEyxmTX49Cdw2+VK2r++LJqzw01dMUBQBLDmmPDrRGVJO7MvKLWCRC93pBuF7q72MJOXBbV/dlhThN5cG4VKviKri9V0wQDt5KQZFUjtxSrZudWASnUu23tlyDxqB1dV9ceJKu8XNvV0d81ieafZ/KeYAUFstYbSK3l0kC0YiM8CESKTT7+cUywcUc4+7A1TwwCE1dZwc7RPi58QIlVJHK5EjN5Qh3rNCneJ2SXYSXOUUokRPFtpwi3NB4NIXbBCPYKQU9vnZPW7SyIRQWy9Dr+zNlKpxuSaTmFmHzhXjU8HdXW8iaC261oeIS9Uo9DIYupFpV90Wr6kofNm8BdyRD84CG+ijvH9x9EzIKEO5rWVk4DL7yCwoKMHTo0AoV4ADg1atXkMlkCAzkh8oHBgYiJUU4KjElJUWwfUlJCV69eoXgYPXkr3PmzMGMGTPY9zk5OQgN1e2Mb0w8nOzhaCdmV9+aEh02qObNvi6vQzYjbDF1QDVFegIKUxz33IDyoVIvxBNXnmWyD0h9/QiquDpgYd/6PCFu0rarCPZyZv11ACapr+J7sS+9KWx4rxkGrlVGTp+f01nruba93xIlcgJHezGrsdHld1WWepdcPJztUVAsw/TONbHsqGEBEWXB0U6MM7M6oeXi43rvw/gxiUR804exa6eWBVfO6lpopa2KWpqeUmoHeuLgLcW1PaG0moidWMRe/6byUSQcU7U2QlWET65P6W+XEtR8SkUQoVPtALX5qnSSFyG/VHjLKCjmPZy4v6qQlqSsAS0OdmJU9XZh8/EJIZMTpOdJeMKdQuBTCH4vS7V6khI50vOLkZ5fjDsCuRQZPJztWeEuyNOJFfIYoS/Y0wWeLvY6tSi5RVKzCnBbLsRjeOuIch/nakIWa82IX9q73MczBtxrf8flBIzQMM7yKsOFfr8RbbRr77VhzxE+Tz9Mw/DWVi7EjRkzBr///ruaRqyiUJ2EutSbQu2FtjM4OTnByclJ8LOKQiQSwcfNkS3sru2+w2ihypvDhtHqMTdzTTm2AGHNHvO9MukdGE1A8wgftYLsQgxtEaYW0HHhSToCVLSSs3bfwsTSBMWM6bZZhA+GtwrHFpXwf02IxSI4lo6vbZQf/rnzEi0ifDS216csji5uze/Ovq4IIW7jqOYI9HTm+UQJUSKTsyY7riaOJ8SZX4ZDg2pe6FjbH1VcHXmCiCZ8OSk5vnorhn1dI0B5A2a+Fyd7O/ZGrY/Wujzocu7matcBwNNZOSfSOL6t3AeV0IKLEdxEIsDNUaGJK5ERpOcLm8kHNq2Gn07zte+mrKVqJxYhwNMZAZ7OaKihDSEEWQVSrabblOwi5ElKkFtUgtyiPDanoBDODmIEe7mUmmpLzbYc022Ql7PZq3Z8vu8OBjcPrXATckXArf/6UiA6nqG8QTHuHMVBrUB3rH23KS+4zVC8XBzQoZY/Tj1ME0xOb24MFuKWLFmCPn364PDhw4iJiVEzSS5fvtxonePi5+cHOzs7Na1bamqqmraNISgoSLC9vb09fH2FQ5EtBX2TjjI32ld5ErUHgKujnd6h2f6lwlKepARpuRKtq3ChvHTMtPN0UVxSl55msMENADCyTQT7elaPOvjqMD/Sk8k4/8/09ui+4jS7XdWhtUgq45hTlX384PUoJGUVYkhzw7SmjFlYW6TZZ3/e1vgZw0/Dm2LcFu11YhnaRvni3CPNJjVNuDvZ652Mlfn5dD2Toj49hI2jmuOPK89x4JYimXdydhEv2aW5U4wACo3UplEt9G7PTTXA7X/vmGB8gGu8thn5xZwAAhNp4kr/lyWnF+MPKBUQMDX9NHuuKlPzhPu64W5yDgqlMszZfUuwvZDm0pRCnD6IRCJUcXNEFTdHRId4amyXWyRlNXeMYKeq3csskKJIKsfTV/l4qis5upmp/dlhfDPQeDlXlx95gEmdogxORWVMCotl7P1F8V7zPCtv3juuBlgsEpVZgPusd138Efsc77WJwHdH4wBYVuATg8FC3OLFi/HPP/+gdm1Fln9NxdWNjaOjI5o2bYqjR4+if//+7PajR4+ib9++gvu0bt0af//9N2/bkSNH0KxZM0F/OEuid4Ng/KFHHUwmalLownd20F+Ic+as/LZciEf7Wn747VKiYNtFB+5BLBLxTDvM6ZtzNFrcmyX30pjQoTreaBgMRzsxWpSa+xhzsKsj/0bDjG9ypxpYdeIxJCVy5JcKMlw1d6CnMzaMbK7XWLkwgqA2X6gCPUzV2vxOOtfh5/1bPrgRtv33DDUC3FG/qhdcHe1AiMLc4GQvxvnH6RgvIBDun9IOXi4OeJZRoLXOpZ+7I5qEKfz99Lkdjtp4mfc+T1LCX3GaX4YzmEiOzxZ3bgjdo/7Xuy5rTr2akAW5nPByqBkTbUJ4kIZoQl2m3nohnmqmRmbtI4IyKvfE/VS2jJ3iM+UYhcyp5hbi9MXD2QEezg6ICvDQ2IYJyNBmuk3Lk1iML9zHf2iuCGMo3//7CD5ujhjZVnd6J1OhquXdevEZ5r4RLdhWRpRuHd4uDrz624ZSnt/z/deq4/3XqgNQ+vMZmqqkIjBYiFu+fDl++eUXjBw50gTd0c6MGTMwfPhwNGvWDK1bt8ZPP/2EhIQETJgwAYDCn+3FixfYvHkzAGDChAn48ccfMWPGDIwdOxYXLlzAhg0b8Ntvv1V43w2FuzLWJhwzN3Ch8Odm4VVw5O5LNp+cNrgq6O//faSz/Rf77/KEOOZB6efuhKgAdzxKzUORVM4KEdwHhkgkUnM8Z5xHueYjLlyfqNUnFUl39c07pA3mGNqcaXP1WH3ZC5RpAYBfRjZD6+p8AS/Q0xkzutXWeKzu9YIQv7Q3iqQyzP/rDnZcVgjTwd7OcLK3QxUd+fEu/q8L+wAuaxkpbpCMdTzK+XDzPOkykTrZi1EzULlav5aYZVDQiyFwTUUezvZsQmJvVwes4CSi5sJqCTn7ckukCTlyc2EEtMdpeegaHYj9N0s1IpwfVuiaMiR629LRNyBj7OYrvHJYlYX5f981qxCnGihXXLoYFwpIIRy3jj4NQvR2kxGiUah3mfflwrgr+XuY19VKCIOfgk5OTmjbtq3uhiZgyJAhWLFiBb744gs0atQIp0+fxsGDBxEernBaTE5ORkKCMjliZGQkDh48iJMnT6JRo0ZYuHAhvv/+e4vOEcfA1Uhpu5VqyxU37816+KBTFNbpEWpuJxZp9QvTRE6RFEVSmVJYEyn7JCmR6b0SYoRAL1cHDGqqno1eKL+SMTQFdqwmTvODPj1Pd444J3s7zHsjGpF+bjgwtR3uL+yB+KW98XqdQF5+MkNwdrDD0rca4PHiXoj7sqegn8zEjjVw/KMOvG3c70VVs6kv+i4iLBWuJk01Pcyqd5rw3jvaiXmLClOYTJh5wE0Vwv09V7/ThBdhx4UV4oTMqQDmvVFP43lFIuCt0vlUIiO860F14dI8gi+4mkobaak42In1ToNhLKr7V5yT/PnHuqvSmAqhRbImrTTT1k4kQpNw7zKd79iM9pjWuSY+7VO3TPurUjtIoeU1tc9sWTBYiJs2bRp++OEHU/RFLyZNmoT4+HhIJBLExsaiffv27GebNm3CyZMnee07dOiAq1evQiKR4OnTp6zWztLhap60PUMZB/RDHH8DhqreLpjZvbZaJKkmDM1InVVQjFaLj6PO54fZPDwikYgV4g7cVPZJ0xh6xyjMqv2bKOtjfjOoIc+/z04sYitBcNGk/TIE5hjaNHG6tFkb3msGABjVNhInZnZEvRAvo/qf2IlFGrWONfzdUcPfHdc+74qm4VWwfDDfTVw1MERf7Hk+cWU6hMXwuoo520HlumG+W6b6gcSE9UO5vqGTOiorSGgTlJn+/vdU6UfJXRzVCvRA12hhv2ARRKx/o1Qm512Xqn6g0cF8v7OKqiZiSegTNGNMvh/auMLO9eHO6xV2LlWEsidoWoTLOebUfo2q4qu3YnDkw/aCbTURFeCBD7vW0mjZMRSlYsLy5oTBy45Lly7h33//xf79+1GvXj0137I9e/YYrXO2jDMnb5O2FA+Fpb5LiUZIRvvB61E4rlLL08FOxPqLOdmLeRfxd0cfsj53XB8cZoWVnl8Mn9IoQU0j+PGdxiiSytW0VU3DquBRqiLSbG6faN73wVBTiw+MvjAPcE0ZwZ9oiUztFROEH99uYhaNxc5xrXA5PoMtDl/FzRG7J7ZRa1eeNBEM1irDXfu8KxIyCtBQxaTSTqUsDyOwMtfgg5Q89DByBUHG/Bnh54Yzn3RCTpEU1bxd8cV+RW1MbZVNmAVdYkYh0nIlPJMO8/OqCqZcuEEbXE2Tv0quR2fOHPT3cELbGpaRX6wiGd02EmtOPhb8zMXBTmMqp/sLeyApqxDnHqfj8726A6EYHO3FiKnqhVsV4Gv1MkdiUn9Pbfi6OaoljS/S8F0yaws7sQgikQhDmqtXY6poGK25ttRb5sJgTZy3tzcGDBiADh06wM/PD15eXrw/inHQ1wQ3vr3C8dLQjNRCNA6rgmMz+KY5bq3RgSpmzl8vqPsqiEUiTOoYBaC0pisj3Gm4b4hEIsGxzn+zHjaPboHbC7rjvTYRvAzamvpTFpjVoKboVE2RbNM618TqYU3NZnJqWd0XH7xeU+f5vx0knMDh5vxugtsZylqezJKo4uaoJsABCqHo/OzX2fdM4ABjdmXqlZqKUB9X1Avx4iXYlWiZvz3qK5NWayrbpElYF4mUQlwxJ3hndNtItd81kDPHjn3YQafvZWVEm8+TtqhJZwc7VPd3x/BW4TrnFpcQbxf8OUl98WUqRvwiXN3G1NSvqpANuFrx+FcFgm2Z79kSouIZqlZRRrzqqlVe0ZSp7BbF9ERwskLna4lo83bh52UrL25OfIGKm1hVH/OKSKQUQJ+lF5QplQagOEb7Wv7s+461/TH19Si9gi4MwaE0d5GmjN6qCrpfR7dAmxq+RgmqqAjqV/XCZ73r4q8bSXizYQgWHVBUnnDVYe7lanYs6F5qNLjmc8aE1qq6D56+yjdpLijuV6lv8IWPmyPCfV3xLL1AoyZIHyFOKpNrjVYe3joc8en5CPZyUcvZSAHebhHGFrLnIpTbL35pbxy6lYwlh+6rVUYBFOXV2tbwYzWjpz/uhPbfnDBJv7mcffQKL3OKNNYXNhWMYNanQTCuJ2YhI78YxTLha1nGCnEV1j2dcAOdpDKCCnad1IoFdYXC5TWOyUdbNQbGJy7ZwAK/mlB9GIRwcu7oI5CJoBT8uA7l5dXmiEQizOhWG34eTpi77065jsWFMaUlZxehSCrj+QwVSWW8fFsArEqAY+CGyrs62qOGv5vOoBBuYk5LqNhgbKpVccHbLcIgKZGxD2Hmt7+XbPx6yYkZ6u4OIpEI9UI88SQtn00JowkmBdCLrELeA4X5bbT9no6lZcoepeahTQ3f0nOrt3OwE+OLvka2I1shY9pFYsPZpwAUqY261wvC0bsvMblTFD7sWostf5YvKUFxiVyjH13PmGD0jAnGzedZIAToy0kLVMPfnVegPkzA59dU9Fp5BrGfd62w8wH8iNMof3dcys/QuFhiXAssKbCGe8+XlsgBCwpSpUKchcIVejT5DgCAZ6kmztCgBE1IVCYWN3FinqQEPesHafcLEPGz4nM2G4V+jati1YlHaBflr7uxHnBXpJfjM/BaTeVxt1x4xhtrVIC71QlwqqgWNNcE35xqqt6YD5FIhCUDYnjbmMXHodspagJ9eeAWTFddzPw5qS2KSmQ6HbCfZSjM+ozbhKoP3ZDmofjz2gt0rhOAp+n5eJKmaC+CCNX9lJqiM3Hmi1C0Fj7vE40JHWrA08We9YVigsOcHezgVXrPZf7rgtn38qdd8PXh+6haxQUxVc3nepSeX4wXWYVay6IZG0bzLRaLWDeCLw/cQ68Y9dKXjMHHklLc2IlFEIsUlhlLC/ix7ieSjaBNcGBuJPmSEq3O0foS4s1Xs3Mf5pISGVYPa4K/P2incX8RRAj2cuHlGSv9wCh4Ojvg/OzOWDZY2NfLUHzcHNm+qiZGXnea7+DM9Q+s7HCvuZwy5pqzNgZwIqQzC3SnldGXawmay8452ov1iqBj0o9oio5rVd0XZz7phDXvNuWnlREpfAOZXJHlLc9nK/h7OBm99JW/hxO+GdQQ07vUMruf6Yc7rlfo+eQcEynjd6jpvqKMTrUcIQ4QztdoCVAhzoIZ3746vF0deI7NqjDJEuXEOOHP9nZivNEwhH3P1WgXSeUQiUSIqeaFf1VykzEw865jLb6m7L/HZfONE8LYmeRrByqiXJnvL7tAim0Xn+GVSn64p2mWXa7HmHBTjGirRlGZqO7vDo/S+WRMvzhtdYj1hZnnquluuM+5UB9XONqLeaZwhjqlea4YLYJlPR4pADSmieFiLO3ZpfgMpJS64Dx9lW+UwDhtcM2pH5UmOtcUzMPmibMw6cRRS75Gc2JhXxOFy+yedXDt8648vzRVuA7q8enGETJ8OE7NmQVSNp/Vd0OU2q/q/u74/u3GannIGJ86J5WUIHdN4GdkLJi+MvnBGn5xBJ8K1Et1sK9c00UocpPBg+O5u/9mUgX0xjJg0mxwy1OVF3sjpGthgiCWHlLUHNamC3AUOJ8yQtWyHkAUJWuGNRHczr3HquaBLA+tlhzH4dsp6PTtSSw6cA+EEIzbfAURsw8gYvYBRP3vIC/XZ3ngauKcOQnqhYrdcwU+S4K5/6flGSeI0Fjo7UjFlLLSxYgRI8rcGQoffdTJXOfPW8+zUa2KC55nFmLtu7qrNGhibPvqbPoQqUyOT3rUwdTONdV8hN5sGII3G4Zg9u6bbGkopjeq2gBNedgsAcZpvEjHCsutjNUPLJVWkT64kZgl+Bm3HI61+wEaAmPiOfUgFcNbhRvlmFzNsayMLg8xVb2w5+oL+Lnz3RSE7hAO9upbmQcQq4mzrOcjBXxhn0t1fzc2+0CQHiUUDYFJALzpfDzebhGGI3dfsp+VyAkmb7+KpuGdNZ6XEKLXc0rGSRvCfY4sPngPn/WJ1tjWkmDyn+65+pxXI9zc6C3EjRw5Eu7u7rC3t9foeyUSiagQZwZaRPrg0tMMXtqE8kx2bgki5qavzcmbm+eNmXfd6gXyat5ZSmFpIbiauCvxGRrbjeHUiq0MqP4kw1uFI8TbBT3rB/HSX0zkVBao7Ex5PQrfHnlo1OuVK8RpSmWji9frBGDB33f18mnjJWounZCMds4SywZRtDO8VQSiAtzxep0Aows23MwHmhz2k7ML1Z4nf157jg933gAAbBzVHJ1qBwjtyiLnaNe4Jf1+PvsUtYI8MLhZKKctEwSh/zgqgvohnriakKUW/Gdu9P6a6tatC0dHR4wYMQKnTp1CZmam2l9GhuYHIMV0MKtzqUxutIePb6mzfw1/dx0tFRo5Biblgb6lviwBRhMnKZFj4NoLGts1s6DVlzFQXYyNahuBiR1rIMLPjacV4OYsrOwwrgvGNDtyc+6V1feIefDlF8tQUFyidZ5zAyWU5lTFKwtWiFMgrO13dbLDon4xeL1OoEnN4ZrcceLT8/HZ3lt4XFq9ZsWxh6wABwCjNl7WeWzCEczEYhGmvh7FfvbJHzdxJ0lZsYIxsVqaJq5/aXUcbSm/zIHeQtydO3dw4MABFBYWon379mjWrBnWrFmDnBzL9XWyFRjTZXEZV/lCfP92Y0zvUlOvqgi8MkClV5SLgOnVUpGWVmvQVG4H4Ne8rCx4qERFqppF1o9ohm8GNqjwepLmhJsY11hwH0Y5RZoTd2vDk5POYuO5eKTmluaFFHjQjSut4sJFdTiWFvlHUbBrQmve+061/fEaJ7DIlL/aB9uvCW7/cOcNbP0vARO2xOJMXBpWHIsz+NgyFcFM1XT8zT8P2NfMQsOSUowASmtUXKrmUozmwCCFZcuWLbFu3TokJydj6tSp2LVrF4KDgzFs2DBIJJbl7GdLmOLB0zbKD9O71NLop8FFKFpUtZajak4uS+LoHYUfSJ6GyhgHp76GeW9EC35mzYxWMQ+raua6RgdiEMfMYQso55JpVFZlLdnj7GDHasfTciVsOpwCgWvWVcC9oVEoLYloDdQL8UL3eooo1WAvZ2wc1YJ3D+YuCJio+ooiLjUPwzeUrWyXXEewAldDzAh8FibDsdpvmYWps8tkdXZxccGIESOwYMECtGjRAjt27EBBgXAdNIrpYTKymyv0meuDQ0q7wF3ph3g58xzlLQ1dqVmiQzwrpebCXeU3saxbk3lg55KJzFbleQCMahsBQJH8m/nt3AWSfHN9VJnzqdZBrXxXc+Vh9bCm+G5IQ2we3ULtM+7vvXdyW5P14cniXjg5s6Pe7QetPa/1c6ISrNClLj+dSmpOERIzCpCaU8S2NXYqqfLiW+q2ZKwk4MbCYCHuxYsXWLx4MWrWrImhQ4eiefPmuHPnDqpU0V42hmI6GCHqwctcs5zfz90JfRuFoG+jEMGaiwUW5kOgiouWqNO3W1RuTdSBqcrEzZYcfFJRMK4JN59na62UYgjcwunlidJmHh47LieyfQvxUk8/5OumdG9gNHaq7g0Uy8VOLEL/xtVQU0DT5ufuhMX9Y7BsUEM4O5jO818sFiHcgFJgl+M1J7QG1IMVokM8MaNrLfbz+ym5eO3rE2ix+DgyCxQR4pbmE8e4NBjrvmAs9L4Kdu3ahZ49e6JmzZq4fPkyli1bhsTERHz99deoU6eOKftI0QETcbafk9Onoi//lUMbY+XQxoKfqVZCsDRqB2k2S/RtVFXjZ5WBukGe7Gtf1SobNgg3Cu/m82wtLfVHzlHqlSc1AbdE3IkHaQCUee24cF0ZGPNtddUAJct6PlIM4J2WYXiraTW9rQNNwrz5VTz0xJjWB1WfOACY2rmmYNtlRx+otbUEmIXQ01eWlfRdbxvX0KFDERYWhg8//BCBgYGIj4/HqlWr1NpNnTrVqB2k6KZ7/SA2T5sl1XWLDvbE3eQctLPwjP/LBzdEu69OqG2vX9WTV2y8MiIWi7B/SjsUSWVqJjdbpBZH+2GsdByM7i3c1xWRfmWP9O0dE4wpv/Gdz4VcKPh1lxWfl+e8FOvmf73qYqQeEaRC1PB3w2MjVKphjqEqmL3RMAR/3+AnE08vrZRjaSlGqlVRar1vPs+ymAwMegtxYWFhEIlE2L59u8Y2IpGICnFmoFWkL/vakpLqrhzaCLuvvsA7LfQrum4uuHnxuPz9QbtK6QunSn0zFuO2RBpU88LN59lGWxAxpiR90vVoQywWoXGYN64lZLHbAjydBNu6OtqhoFjGE95CvJyRVFpqiWI7NIvw0Rkk4OPmiAhfV1zlXFsAMKBJNV7kqDYO3UpGT4GC9o9SlW4+qnFyjgKBcwXFimAdS4tO9XZVLnJf5lhOIKfeQlx8fLwJu0EpD1zziSXVdasZ6IHZPa3D1P5x99q8m9WKIY1sQoCjqMMtUZWaW4QWXx6Ho70Y1z7vWqYAHcZR2xhXk2r1DE3VNHaNb40Lj9PRr7HSHcCRk8BZRO2pNoWuIIGhzUNx8al6nlduahsA+LRXXTg72mHlsThsfb8FHqfmY/L2q4rP9t4WFOKuJyrdElTvqX4e6tp/JjLcEu+/TGJ9S7J4WW7IIEVvuBOU1kYsG6qRmsYub0OxHphFkVQmx5cH7gFQmFZPPEhFnwaG5ztk4hqM8VC6pPKgtdfwcK5f1UtNw5qWaznaA0rFou3K69+4KqZ3qYW31/+n9tmgptWw7b9nyJOU4MCU19jAtXdbKixz3ETgjAZNG6rm1JFtIrDu1BPBtpYWnQooNYeWJMTpbXW+ePEiDh06xNu2efNmREZGIiAgAOPGjaO54syESCSyyIvLmuDeMD7oFIWWkZWrOgNFfxjt1pTfrmHfdaW/Tlk9FZQ5ssrbM3X0yePI0LqG0u3CApUcFBOiuoD4sn999vWgZtXgaC8WFPScHexweHp7nJ31Oi/zAHM8bnm+Ig3lqPKKpOxrVRNpsEB0NYMFynCwL13gWVL5Or3vAPPnz8fNmzfZ97du3cKYMWPQpUsXzJ49G3///TeWLFlikk5SdEPL6pSP9jX9AQCezvaY2b22RaryKRUDk0JANeXKH7HPy3Q8xifOFJeUJk2cEJ7O6ul/KLaB6lUytHkYvn+7MSZ3qoHW1RXCfVmuT9X7JFtJhANX+Av21t/CYWnRqYDpk4GXBb3NqdevX8fChQvZ9zt27EDLli2xfv16AEBoaCjmzZuH+fPnG72TFN042IsBC0/lYcmE+bri/OzX4eVCH3S2jqacV6cfppXpeMzt3hQPJXs7A4Q4em3bLKqXnlhUWgpRoO61oQxtHspmR3ialo8AD76gxtSmDvBwgp+7cCCOEJYoxDEWL31MxxWF3pq4zMxMBAYqsyyfOnUKPXr0YN83b94ciYmJxu0dRW+yCqS89xZ4/Vs8Id4uFl1ZglIxqPpHlhfVbPXlYXF/fvk6wzRxynHR24Ntsebdprz3gpaGMl4UjcO8tX7OlDOMCjAsOtsSfeLEpX06VcYFnSnQW4gLDAzE06dPAQDFxcW4evUqWrdWFuvNzc2Fg4PpVnqZmZkYPnw4vLy84OXlheHDhyMrK0vrPnv27EH37t3h5+cHkUiE69evm6x/5qZeiKfuRhQKRSeaauiWFTnj42CEZ1LfRvzACnsDkmm1rO6ruxGlUqJPkmljiEw3nmfx3sc+y8THfyjcsFI0pLdZObSRcH8sT4Zj3ZYsyWKj9x2gR48emD17Ns6cOYM5c+bA1dUVr732Gvv5zZs3UaNGDZN0EgDeeecdXL9+HYcPH8bhw4dx/fp1DB8+XOs++fn5aNu2LZYuXWqyflkKVVxpolYKxRh4GFsTV/rfGJo4riM5YJg5lVumyRIfkBTzUtZrgpsQffHB+7zPpnKSUz/RUOnAQ6D+L2CZmjhGGNYUxGEO9L5bLVq0CAMGDECHDh3g7u6OX3/9FY6OSsHhl19+Qbdu3UzSyXv37uHw4cP477//0LJlSwDA+vXr0bp1azx48AC1a9cW3I8R8gzJcSeRSHhRtjk5OWXveAViSSsDCsWa2TymBYb89J/RItCMqIiDvZ0Y9mIRm9TbEE0c9x5BA6AoqpTVJy4qQHPZQrkeBZk1FZS3RJ84pvTW88wCM/dEid53AH9/f5w5cwaZmZnIzMxE//79eZ///vvvmDdvntE7CAAXLlyAl5cXK8ABQKtWreDl5YXz588b9VxLlixhTbZeXl4IDbWOAui5RjYBUSi2SuOwKni4qKfadgcDtF5clD5x5eoWS8faikjqJmHe8PfQ31E8wteN1W7U0VIvmGI9qPpIlgdTyEz6CHGujsK6JEsU4piE2fdTcnW0rDgMrk7m5eUFOzt1ydnHx4enmTMmKSkpCAgIUNseEBCAlJQUo55rzpw5yM7OZv+sJVija13174dCoZSdLWNaoFGoN7aOUSwepTKCf++/NPg4hM0TZ5yH0s/vNceDRT2wZ1Jbg0xO9nZi3JrfDWdndULfRlV170CxeFR//u71AoUbAqgVqD2wwFgyU9xLpYDDLU/Vs36QYHs3R6U8wdUWW6A1FbVLFz+GBBSZGrOWmJ0/fz5EIpHWvytXrgAQjqYhhBg9n5eTkxM8PT15f9aAJpU0hUIpG6/V9MfeyW1RI0CZlf6Lv+8afBxWG2HEW5WTfdnmu6ujvcZawRTrQ/Xx5+2iWZEyuVOU9mOV4wId1TaCff3VYWX5wuhg5fPz8z7Rgvu6cnxQBzerxr62RE0c43teIifKgCUzY9Z8Ch988AGGDh2qtU1ERARu3ryJly/VV8BpaWm8tCe2jKO9WeVxCqXSwi3SXZYkn6bME0eh6Isuy2Z5Ls/hrcKx8Vw8ACA+XRnA4MB5LrlpMJtyNXFcP1RLrD7EDSSSyuVwEptfeWJWIc7Pzw9+fn4627Vu3RrZ2dm4dOkSWrRoAUBRBiw7Oxtt2rQxdTetAtVC2LTANYViHLjC14usQoP3lxvZJ45CYYgO5tfH1SaI6eOfVla4z59HqXkokcnx5FU+biRmKdvYC3eOm5fxGqd9TpFUoLV5UV3QWUJaUQvogm7q1q2LHj16YOzYsVi3bh0AYNy4cejTpw8vMrVOnTpYsmQJG3SRkZGBhIQEJCUp6h8+eKBQ8wYFBSEoSNg+b62oCnEUCsU42HFW32WJAldaU6kURzEuMdW0C3HdopWWKl0uN47leIa4OPKPPezni0jL49dS1xTAwK3/y+2+psop5oT7nC2xEE2h1Tz5t23bhpiYGHTr1g3dunVDgwYNsGXLFl6bBw8eIDs7m33/119/oXHjxujduzcAYOjQoWjcuDHWrl1boX2vCMoaOUehULTDrTnaKybY4P3Z6FSrudtSrIkPSn3dZnStBa4YdOWzLljLqdTQLToQXeoG4OPuwim5Pu8TjSBPZ42+a9rwc3fCoKZKf7aLTzPwJE1pVv1uSEO9jsN1M6tpYIWHisBOLGI16sUWIsRZhSYOUES/bt26VWsboqIuHjlyJEaOHGnCXlEoFFvg4+618c0/D/DbpQS80yJMTQOiDTZPHPWJo5iAj7rVwuBmoQj1ccGne2+z21XrlNrbifHze801HifCzw0X5rxe5uv0m0EN8Xvsc8HPXPQMvKviZvlJ6+3txCgukSMhvUCtTqw5oGvDSoIlXEwUSmWlsFjGvn7jx7MG7cv4IlERjmIKRCIRwnxdjbJIMNVCw1OHG8KG95qhdXVfLO5fH6E+LgCAltV1lwozB0zwRWaBZfjsUSGukqDqk0ChUIxHv8Zlz6tm7DxxFIomzH2FRfq5CW7XVbu1c91A/DauFapVccUPbzfB8FbhmNa5lim6WG5aVWdKb8l0tKwYqBBXSdBXXU2hUAwnSsU/5/vjcbz3mnJGvcqTYGVpWxqdSqns/PxeM8HthgTeNQr1xsJ+9Q2qRlKRMM/atFyJjpYVAxXiKgnerrR2KoVSUSw/+hBJpelGTjxIRfX/HcTr355UW51/dUhZEJz6xFFMTaCned1qgsx8/oqA0aifeJBq5p4ooEJcJcHZwY6nyqbPCwrFuCwZwK9TyQhsjKD25FU+6nx+mBXuAH7eKzonKaZm7GvV0a9RCNYMa2KW8wulMfFzt/xgBUNggi90+flVFFSIq0TUDqRFrSkUU/F2izDe+9eXnULE7ANqxbDbLP0XXx64i4z8YjxKzWO3U584iqlxcbTDiqGN0bMMqXCMgZ1YxAYmMByc9ppZ+mIqGod5A6B54igmwJ7miqNQLIL1Z57i68P3edvo7KTYAmc+eZ19XS/Es9JlTmD8+8pSgs8UUCGuEsHNuE0X/RSKedlxOZH3/tg99frPFEplxq4SRvM4skIc1cRRjAwtvUWhmJb7C3uUed/49AIj9oRCsVwGllZveP+16mbuifFhLF5MvjhzYzUVGyi64ZpTaZ1GCsX46Ko/SaFQgC/718eYdpGoE1T5/LQZZcnFpxlm7okCqrqpRDjaK39OObEMez2FUtmoF+IpuP3DLpaZnJRCqWic7O1QN9izUqbVsbScrFSIq0QM5BQgzpeUmLEnFErlZd3wplg/ohmeLumFFpHKTPSv1wlAiJdmJ+7+5aj6QKFQLIMWkT4Y2LQaetQLMndXAFBzaqUiOlipIZBYiL2eQqlsVKviimpVXAEAXpxcUTHVvGAv4Jf6doswhPm4YuxrkRXWRwqFYhqcHezw7aCG5u4GC9XEVSK4qmtLqetGoVRmVI1FH3SKUmsTU9ULEzvWEBTwKBQKpTzQu0olJczX1dxdoFAqPbN61oG7kz2mdq4JABjUrJqaIEf9UykUiqkQEULvMNrIycmBl5cXsrOz4ekp7NBsSdxJykZSVhG6RgeauysUik0gkxO1fFgFxSWInvsPAOD3Ca3RPMJHaFcKhUIRRF/ZgwpxOrA2IY5CoVgGWQXFeJUnQVRA5UuzQKFQTIu+sgcNbKBQKBQT4O3qCG/XylX8m0KhWBbUJ45CoVAoFArFCqFCHIVCoVAoFIoVQs2pOmBcBnNycszcEwqFQqFQKLYAI3PoClugQpwOcnNzAQChoaFm7gmFQqFQKBRbIjc3F15eXho/p9GpOpDL5UhKSoKHh4fJ6sDl5OQgNDQUiYmJlT4C1pbGCtjOeG1lnIBtjRWwrfHSsVZOrHGshBDk5uYiJCQEYrFmzzeqidOBWCxGtWrVdDc0Ap6enlZzgZUXWxorYDvjtZVxArY1VsC2xkvHWjmxtrFq08Ax0MAGCoVCoVAoFCuECnEUCoVCoVAoVggV4iwAJycnzJs3D05OTubuismxpbECtjNeWxknYFtjBWxrvHSslZPKPFYa2EChUCgUCoVihVBNHIVCoVAoFIoVQoU4CoVCoVAoFCuECnEUCoVCoVAoVggV4igUCoVCoVCsECrEUSgUCoVCoVghVIijUCgUCoVCsUKoEFcBXLlyBUVFRebuBoVSbmwhIxGdr5TKBJ2zlRsqxJmQJ0+eoG/fvmjRogV27dpl7u6YlMTERPzxxx+4evUqpFIpgMp988jIyMCrV68AAHK53My9MR3JyckYNGgQdu7cCaByj9WW5itA52xlhc5Z24IKcSaAEIJJkyahZs2aEIlE8PLygru7u7m7ZTLmzJmDWrVqYdmyZWjTpg0mTpyIJ0+eQCQSVcqHwqeffoo6dergp59+AgCIxZV3Gm3YsAG7d+/GihUrUFBQADs7u0r3ULC1+QrQOUvnrHVji3NWE5X3SjYTe/fuhZubG2JjY3H+/Hns3bsXdevWxaFDhwBUvpXuxYsXsW/fPvzxxx84ceIEfv75Z8TFxWH48OEAAJFIZOYeGo+srCyMGTMGx44dQ1hYGP777z9cvnwZQOX7XRnOnz+PIUOGwMnJCV9//bW5u2N0bG2+AnTO0jlr3djinNUGFeKMAPeiSUtLw9atW3Hx4kW0bNkShYWFqFGjBjIyMlBQUFCpbpCAYkLJZDL07t0bzs7OePfdd7F06VLcvHkT3333HQDrnlTcvru4uCA8PBxz5szBsmXL8OLFC/z555+QSqVWr8FQ7XtJSQkAIDg4GEOGDEGbNm2wa9cu3Lt3D2KxuNKM1dbmK0DnLJ2z1oetz1ltUCGunBQWFqK4uJh9P2bMGAwYMAAAIJPJ4OLiAj8/Pzx69Aiurq5WrdZmJhJ3DAEBAXBxcUFBQQG7rVWrVpg5cyYWLlwIiURitZNK9bd1dHTEtGnT0K9fP3To0AGdOnXC6dOncfToUTP2svyojpMQAnt7ewDA5cuXUatWLfTv3x9BQUFYu3YtiouLcffuXXN1t1zY0nwF6Jylc5bO2coOFeLKwZw5c9CuXTv06dMH33//PXJzcyEWi9mLiLkRdunSBfHx8UhISLBaX4zly5dj8eLFAPj+JJ6enrC3t8fx48fZbSKRCO+99x5cXV2tdmWv+tvm5ORAJBLB09OT/X2nTp0KQgj27t2LV69eWeXKXtM45XI5Xrx4ATc3N0RERKB58+Z44403sH37djg7O+Pff//l3VitAVuarwCds3TO0jlrC9jWaI1EcXExBg0ahL/++guffPIJQkJCsG7dOrz99tsAlDdM5r9MJoOvry8SExPN1ueycvnyZXTq1AkzZ87Enj17cOHCBQBgo9kGDRqE4uJiHD58GKmpqex+wcHB6Nq1Kx4+fAiZTGY1K3tNv+0777wDQHHTYG4iYWFhGDx4MK5evYr9+/ezn1vDQ0HXOMViMTw9PeHg4ACRSIQ///wTixYtglQqRUxMDKZMmQJHR0erHmtlnK8AnbN0ztI5a1MQisHcvXuX1KxZkxw5coTddvbsWeLi4kK+/vprIpfLCSGEyGQyQggh6enpxNHRkezfv5+33RpYuHAhGThwINm4cSPp1q0bef/999nPiouLCSGErFq1itSqVYv89NNPvH3btm1LxowZU6H9LS+G/rZFRUWkV69eZPDgweTmzZtk69atZNGiRWbpuyHoGichhBw/fpwEBweT+vXrE29vb/Ltt9+SdevWkUaNGpFVq1YRQqzjWral+UoInbOE0DlL56zlj9FYUCGuDMTGxhKRSETS09MJIYS9oJYsWUKqVKlCHj58yGuflZVF2rdvTz766KMK72tZYcb07Nkzcv78eUKIYnwtW7Yku3btIoQQIpVK2fbvvPMOadSoEVm3bh3JzMwksbGxpEmTJmTHjh0V3/lyYMhvy9wo9u7dS6pXr058fX2Jo6Mj+fbbbyu+4waibZze3t7kyZMnRCqVkujoaDJu3Djy9OlTQgghSUlJZPDgwaR9+/akqKjIXN03CFuYr4TQOUvnLJ2z1jZnjQEV4srAtWvXSL169cgPP/xACFFeYMXFxSQyMpK9kJgbZklJCalZsyaZMGECuxK2Rh4/fkz69etH+vXrRzIyMgghhEgkEvazuXPnEjs7O9K0aVPi4uJCxowZY3Xj1fe3LSkpIYQQ8ujRIzJixAgiEonIxIkTSV5ennk6biDaxhkREUGmT59OCCHk5cuX7GcMd+7csZqHASG2O18JoXOWzlkFdM5WXqgQVwYyMjJIv379yJAhQ0hSUhIhRHkxLVu2jISEhLArPubGsXnzZvLgwQPzdNgIMJNow4YNpGXLlmT58uWC7W7fvk32799P7t27V5Hd0xvVm5sqhvy2hBDy8ccfk2rVqpGbN2+artNloLzjDA4OVjNJ6DqmuTDmb2oN81Xf34HOWTpn6Zyt/NDABhUSExMRGxuLpKQktc+YPDxVqlTBG2+8gfv377OlPpjwbi8vL1SpUoV1sLSzswMADB8+HLVq1aqIIeiNPmNlkMlkAICBAwciOjoa+/fvR1xcHADg6tWrABRpDOrVq4fevXujTp06Ju694aSlpfHSKnBD0Q39bZl9ly5disTERMTExFTUMHRijHH6+PioOQlboqO7MX9TwLLnK6DfeBkqw5xNTU1Fbm4u+76yzlljjNNa5qwxf1PA8uesqaFCXClSqRTjx49HkyZNMHr0aDRs2BDnzp0DoLzI7O3tUVRUhB07dmD06NFo1KgRdu7ciRMnTrDHef78Ofz9/REeHm6WceiDvmOVSqX49ddf2fdyuRyenp4YNGgQ5HI5FixYgM6dO6NZs2bIzMy02NBuqVSKcePGoW3btnjjjTcwatQotf4a+tuqRkdZAqYYp6ViS2MF9B9vZZmzJSUlGDNmDFq0aIEuXbpg2LBhSE9Pr3Rz1hTjtFRsaawVirlVgZZAbm4uefPNN0mnTp3I1atXyf3790m3bt1Ihw4deO1WrlxJfHx8SN++fQkhhNy4cYMMGzaMODo6kokTJ5Jx48YRDw8PsmbNGkKIZaqyDR3rW2+9xfrSMDx79ozUqFGDiEQiMnToUJKSklKBIzCMjIwM0qVLF9KpUydy9uxZ8tNPP5HGjRuTNm3akPv377PtrP23tZVxEmJbYyXE8PFa+5yVSqVk2LBhpFWrVuTkyZNk+fLlpH79+qRdu3bk7t27bDtr/31tZZyE2NZYKxoqxBFCLl68SGrWrEn+/fdfdtv69evJm2++yV4kP/74I4mIiCDbtm3j+R3I5XKyePFiMnbsWNKrVy9y7ty5Cu+/IRg6VtVJcvz4ceLu7k4aNWpErly5UqF9LwuHDx8m9evX5z3s7t69S8RiMZk6dSrJzMwkGzduJGFhYVb929rKOAmxrbESYvh4rX3OJiQkkJo1a5ItW7aw25KTk0nVqlXJlClTSEZGRqX4fW1lnITY1lgrGirEEULOnDlDRCIRe3GkpaWRRo0akQkTJpC1a9cSQhSh6fn5+bz9rHEVUNaxMrx69Yps3769wvpbXn799Vfi7e3N23bu3Dni4+NDatasSQ4cOEDkcrlahJq1/ba2Mk5CbGushJR9vAzWNmevXbtGXFxcSFxcHCGEsFGVP/74I6lZsyb5+++/iVwut/r7sa2MkxDbGmtFYznOARXE4sWLMW/ePOzYsYPd1q5dO3Ts2BGjRo1Cz549ERgYiKCgIDg6OuKzzz7DoEGDcPv2bbi6uvKyXVui0ygXY44VUJTh8fX1ZbNmWxpC4w0LC4OPjw+++uordtvPP/+MMWPGQC6XY9++fRCJRHBxceEdy5J/W1sZJ2BbYwWMO17A8ufswYMHAfBLfNWuXRvBwcHYunUrAKUP2+TJk+Hl5YXdu3dDIpHA1dWVdyxL/n1tZZyAbY3VIjCjAFmhXLx4kYSFhZEmTZqQnj17Eg8PD/LWW2+xJorc3FwSFxdH2rRpw0v8eP36dVK9enU2WaY1YEtjJUR4vP379yeJiYmkqKiIfPXVV0QkEpE2bdoQd3d3Ur9+fSKVSskPP/xAqlatau7u642tjJMQ2xorIbY33v3795OqVavyrAKMCa2goIB88sknpGbNmuTly5eEEEIKCwsJIYRs2bKFeHl5se8tHVsZJyG2NVZLwmaEuBkzZpDevXsTQhQX1s2bN0l4eDiZOHEiSU5OJoQQcvnyZVK7dm2SmprKqnGlUilbvsRasKWxEqJ5vBMmTCCpqamEEEL+/fdf8sMPP/DKuCxdupS0a9eOZGVlmaXfhmIr4yTEtsZKiG2N98yZM6RHjx7kgw8+ID179iTNmjVTa3Ps2DHSvHlzMmnSJEKI0qx24sQJEhAQQG7cuFGhfS4LtjJOQmxrrJZGpRfi5HI5ycrKIu3atSMzZ84khChXB6tXryZNmzYlK1asIIQQcv/+fSISiUhsbCy7/59//kmaNGlCrl69WvGdNxBbGishho1XFYlEQvr160emTJlSYf0tK7YyTkJsa6yE2NZ4mYf2w4cPyfLly8mTJ0/IlStXiKurK/n5558JIcqEroWFheS7774jbm5uZM+ePWyViUWLFpGOHTtatK+UrYyTENsaq6VSKYW42NhYtZVps2bNyPjx4wkhSqfK4uJiMmDAANKvXz/y7Nkzkp+fT4YMGUJcXV3JhAkTyIgRI4iHhweZO3euxV5gtjRWQgwfb//+/cmTJ0/Ytvfv3ycPHz4kI0aMIJGRkeTChQsV13kDsJVxEmJbYyWEjpcQZZZ9qVRKPvroI+Lv78+Om/ksJyeHfPLJJ8TDw4N06NCBDBo0iLi4uLCF3C3tPmUr4yTEtsZq6VQqIe6PP/4g1apVIzVq1CBhYWFk7ty55Pnz54QQRf4Zd3d3NvqFWQXs3r2bVKtWjS0YnZ+fTz755BMycuRIMmLECIst42FLYyWk7OMNDQ3lhaQvW7aM1KhRg7Rv316tiLIlYCvjJMS2xkoIHe/cuXNZdw65XM4+sJ88eUJCQ0PZepiqpaN27dpF5s2bRyZMmGCRpcFsZZyE2NZYrYVKI8RdvnyZ1KlTh6xYsYLcuHGDrF69mvj7+5OJEyeSrKwsNtkls9rlFsn19fVlVb8MjArYErGlsRJS/vFu2LCBfZ+cnMwzIVsStjJOQmxrrITQ8XLHm56eTghRamfkcjlZvXo1sbe3ZzWOEomEZGdnm63/+mIr4yTEtsZqTVi9EMdI/mvWrCHVqlXjXSQ//vgjadGiBVmyZAkhhJBVq1YROzs7curUKbbN48ePSY0aNcju3bsrtuNlwJbGSojtjNdWxkmIbY2VEDpe1fG2atWKLFy4UG2/9PR00qZNG9K3b18SGxtLunXrRrZs2WKx5jVbGSchtjVWa8TqhTiGTz75hLz++uu8ZIF5eXlk8uTJpFWrVuTBgwdELpeTYcOGkaCgILJgwQJy7do1Mn78eBITE0NevHhhxt4bhi2NlRDbGa+tjJMQ2xorIXS8hCjH26ZNG3L79m1CiFJzQwghGzduJCKRiIjFYtKnTx9SUFBQ4f02FFsZJyG2NVZrwuqEuCNHjpApU6aQFStWkIsXL7Lb9+3bR5ydncnjx48JIcoL6ciRI6RNmzZk+fLlbNspU6aQRo0akaioKNKkSRNy8+bNih2EntjSWAmxnfHayjgJsa2xEkLHy6BtvG3btuWNVyKRkFWrVhGxWEw6dOjACgOWhK2MkxDbGmtlwGqEuKSkJNKnTx8SEBBAhg0bRmJiYoiXlxd7kRUWFpI6deqQcePGEUL4jpSvvfYamThxIvueKSvFrUVoSdjSWAmxnfHayjgJsa2xEkLHa+h4mVxhhBCSkpJCpk2bRn799deKHYQe2Mo4CbGtsVYmrEKIy8/PJ++99x4ZMmQIL9S+efPmZOTIkYQQxapg8+bNRCwWqxXIHTZsGOnUqRP73pJt8rY0VkJsZ7y2Mk5CbGushNDxMpR1vJaKrYyTENsaa2XDKmqnurq6wsnJCSNHjkRkZCRKSkoAAH369MG9e/cAAHZ2dhg8eDD69u2L999/H6dOnQIhBCkpKYiLi8OwYcPY41lyPTZbGitgO+O1lXECtjVWgI63vOO1VGxlnIBtjbXSYTbx0UC4IffMSvXdd98lY8eO5W0rLCwkHTt2JAEBAaRbt24kJCSEtGrViiQkJFR8p8uILY2VENsZr62MkxDbGishdLyEVM7x2so4CbGtsVYmRIQQYm5Bsqy0b98eo0ePxsiRI0EIgVwuh52dHV6+fImbN2/i8uXLiIiIwDvvvGPurpYbWxorYDvjtZVxArY1VoCOt7KO11bGCdjWWK0WMwmP5ebx48ckMDCQXLlyhd3GZDmvbNjSWAmxnfHayjgJsa2xEkLHS0jlHK+tjJMQ2xqrNWMVPnFcSKni8OzZs3B3d0fTpk0BAAsWLMC0adOQmppqzu4ZFVsaK2A747WVcQK2NVaAjreyjtdWxgnY1lgrA/bm7oChME6/ly5dwltvvYWjR49i3LhxKCgowJYtWxAQEGDmHhoPWxorYDvjtZVxArY1VoCOt7KO11bGCdjWWCsFZtMBloPCwkISFRVFRCIRcXJyIkuXLjV3l0yGLY2VENsZr62MkxDbGishdLyVdby2Mk5CbGus1o7VBjZ07doVNWvWxPLly+Hs7Gzu7pgUWxorYDvjtZVxArY1VoCOt7JiK+MEbGus1ozVCnEymQx2dnbm7kaFYEtjBWxnvLYyTsC2xgrQ8VZWbGWcgG2N1ZqxWiGOQqFQKBQKxZaxuuhUCoVCoVAoFAoV4igUCoVCoVCsEirEUSgUCoVCoVghVIijUCgUCoVCsUKoEEehUCgUCoVihVAhjkKhUCgUCsUKoUIchUKxGTp27Ijp06ebuxtqWGq/KBSKZUOFOAqFQhHg5MmTEIlEyMrKMvkx9+zZg4ULFxrtPBQKxTawN3cHKBQKxdbx8fExdxcoFIoVQjVxFAqlUpKfn48RI0bA3d0dwcHBWLZsGe/zrVu3olmzZvDw8EBQUBDeeecdpKamAgDi4+PRqVMnAECVKlUgEokwcuRIAAAhBF9//TWqV68OFxcXNGzYEH/88YfO/mg7pqo5NSIiAosWLWL7Hx4ejn379iEtLQ19+/aFu7s7YmJicOXKFd45zp8/j/bt28PFxQWhoaGYOnUq8vPzy/L1USgUK4AKcRQKpVLy8ccf48SJE/jzzz9x5MgRnDx5ErGxseznxcXFWLhwIW7cuIG9e/fi6dOnrFAVGhqK3bt3AwAePHiA5ORkrFy5EgDw2WefYePGjVizZg3u3LmDDz/8EO+++y5OnTqltT/ajinEd999h7Zt2+LatWvo3bs3hg8fjhEjRuDdd9/F1atXERUVhREjRoCpnHjr1i10794dAwYMwM2bN7Fz506cPXsWH3zwQZm/QwqFYuEQCoVCqWTk5uYSR0dHsmPHDnZbeno6cXFxIdOmTRPc59KlSwQAyc3NJYQQcuLECQKAZGZmsm3y8vKIs7MzOX/+PG/fMWPGkLfffltnv4SOSQghHTp04PUrPDycvPvuu+z75ORkAoB8/vnn7LYLFy4QACQ5OZkQQsjw4cPJuHHjeMc9c+YMEYvFpLCwUGffKBSK9UF94igUSqXj8ePHKC4uRuvWrdltPj4+qF27Nvv+2rVrmD9/Pq5fv46MjAzI5XIAQEJCAqKjowWPe/fuXRQVFaFr16687cXFxWjcuLFRx9CgQQP2dWBgIAAgJiZGbVtqaiqCgoIQGxuLR48eYdu2bWwbQgjkcjmePn2KunXrGrV/FArF/FAhjkKhVDpIqYlRE/n5+ejWrRu6deuGrVu3wt/fHwkJCejevTuKi4s17scIegcOHEDVqlV5nzk5OZW/4xwcHBzY1yKRSOM2pk9yuRzjx4/H1KlT1Y4VFhZm1L5RKBTLgApxFAql0hEVFQUHBwf8999/rACTmZmJhw8fokOHDrh//z5evXqFpUuXIjQ0FADUggQcHR0BADKZjN0WHR0NJycnJCQkoEOHDgb3S+iYxqJJkya4c+cOoqKijH5sCoVimdDABgqFUulwd3fHmDFj8PHHH+P48eO4ffs2Ro4cCbFYccsLCwuDo6MjfvjhBzx58gR//fWXWp628PBwiEQi7N+/H2lpacjLy4OHhwdmzpyJDz/8EL/++iseP36Ma9euYdWqVfj111919kvomMZi1qxZuHDhAiZPnozr168jLi4Of/31F6ZMmWK0c1AoFMuCCnEUCqVS8s0336B9+/Z488030aVLF7Rr1w5NmzYFAPj7+2PTpk34/fffER0djaVLl+Lbb7/l7V+1alUsWLAAs2fPRmBgIBvluXDhQsydOxdLlixB3bp10b17d/z999+IjIzU2SdNxzQGDRo0wKlTpxAXF4fXXnsNjRs3xueff47g4GCjnYNCoVgWIqLLeYRCoVAoFAqFYnFQTRyFQqFQKBSKFUKFOAqFQjESEyZMgLu7u+DfhAkTzN09CoVSyaDmVAqFQjESqampyMnJEfzM09MTAQEBFdwjCoVSmaFCHIVCoVAoFIoVQs2pFAqFQqFQKFYIFeIoFAqFQqFQrBAqxFEoFAqFQqFYIVSIo1AoFAqFQrFCqBBHoVAoFAqFYoVQIY5CoVAoFArFCqFCHIVCoVAoFIoV8n/txCzu50jwjwAAAABJRU5ErkJggg=="
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from pytesmo.validation_framework.adapters import AnomalyClimAdapter\n",
    "ds_mask = AnomalyClimAdapter(ismn_reader, columns=['soil_moisture'], \n",
    "                             timespan=(datetime(1991,1,1), datetime(2020,12,31)))\n",
    "anom = ds_mask.read(ids[0])\n",
    "pprint(anom)\n",
    "anom.plot(title='Anomaly (wrt. 1991-2020 avg.)', ylabel='SM m3m-3', figsize=(7,2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## ColumnCombineAdapter\n",
    "\n",
    "This adapter is used to combine multiple columns of a dataset after reading the time series. It takes any function that can be applied accross multiple columns and will add a new column of the chosen name to the data frame.\n",
    "E.g. in the following example we combine the `proc_flag` and `snow_prob` and `frozen_prob` column into a new column `good` that is 'True' when all the chosen columns are 0 and otherwise 'False'. Then we apply a `SelfMaskingAdapter` that uses the newly added column to filter the dataset. This is just for demonstration, we could also just use the `AdvancedMaskingAdapter` for this example. A more common use case for the `ColumnCombineAdapter` would be to compute the average when multiple soil moisture fields are available in a dataset.\n",
    "\n",
    "This example also shows that you can stack multiple adapters together. If they depend on each other, it is important to notice that the innermost adapter will be called first!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                              sm  proc_flag  snow_prob  frozen_prob  good\n",
      "2007-06-03 03:33:38.999992  84.0          0          0            0  True\n",
      "2007-06-03 14:51:50.999997  82.0          0          0            0  True\n",
      "2007-06-06 02:31:33.000005  76.0          0          0            0  True\n",
      "2007-06-06 13:49:46.999999  75.0          0          0            0  True\n",
      "2007-06-08 03:30:13.000023  83.0          0          0            0  True\n",
      "...                          ...        ...        ...          ...   ...\n",
      "2014-06-25 13:51:37.000006  83.0          0          0            0  True\n",
      "2014-06-27 03:32:04.999979  79.0          0          0            0  True\n",
      "2014-06-27 14:50:18.000019  81.0          0          0            0  True\n",
      "2014-06-30 02:30:05.000000  77.0          0          0            0  True\n",
      "2014-06-30 13:48:17.999999  83.0          0          0            0  True\n",
      "\n",
      "[713 rows x 5 columns]\n"
     ]
    }
   ],
   "source": [
    "from pytesmo.validation_framework.adapters import ColumnCombineAdapter, SelfMaskingAdapter\n",
    "\n",
    "def select_good(x):\n",
    "    return (x['proc_flag'] == 0) and (x['snow_prob'] == 0) and (x['frozen_prob'] == 0)\n",
    "\n",
    "ds_mask = SelfMaskingAdapter(ColumnCombineAdapter(ascat_reader, func=select_good, new_name='good'), \n",
    "                             '==', True, 'good')\n",
    "pprint(ds_mask.read(1814367)[['sm', 'proc_flag', 'snow_prob', 'frozen_prob', 'good']])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:pytesmo] *",
   "language": "python",
   "name": "conda-env-pytesmo-py"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}