Skip to content
Snippets Groups Projects
dev_samir_xarray.ipynb 738 KiB
Newer Older
{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Build a validation ndvi and weather dataset"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\"\n",
    "Notebook for development and testing of the SAMIR Pixel model.\n",
    "\"\"\"\n",
    "\n",
    "import xarray as xr\n",
    "import os\n",
    "\n",
    "data_path = '/mnt/e/DATA/DEV_inputs_test'\n",
    "\n",
    "# Original sets\n",
    "ndvi_path = data_path + os.sep + 'ndvi_' + str(size) + '.nc'\n",
    "\n",
    "# Validation sets\n",
Jeremy Auclair's avatar
Jeremy Auclair committed
    "val_ndvi_path = data_path + os.sep + 'xls_NDVI_' + str(size) + '.nc'\n",
    "val_weather_path = data_path + os.sep + 'xls_weather_' + str(size) + '.nc'\n",
Jeremy Auclair's avatar
Jeremy Auclair committed
    "val_outputs = data_path + os.sep + 'xls_outputs_' + str(size) + '.nc'\n",
    "xls_file_path = '/home/auclairj/GIT/modspa_pixel/SAMIR_xls/SAMIRxls_Reference.FR-Lam2019.xls'\n",
    "additional_outputs = {'Zr': 10, 'Dei': 100, 'Dep': 100, 'Dr': 100, 'Dd': 100, 'Kei': 1e4, 'Kep': 1e4, 'Ks': 1e4, 'W': 1e4, 'Kcb': 1e4, 'fewi': 1e4, 'fewp': 1e4, 'TDW': 100, 'TAW': 100, 'FCov': 1e4, 'Tei': 1000, 'Tep': 1000, 'Diff_rei': 1e4, 'Diff_rep': 1e4, 'Diff_dr': 1e4}\n",
    "# additional_outputs = {'Dr': 100, 'Dd': 100}\n",
    "additional_outputs = {}\n",
    "\n",
    "normal_outputs = {'E': 1000, 'Tr': 1000, 'SWCe': 1000, 'SWCr': 1000, 'DP': 100, 'Irr': 100, 'ET0': 1000}\n",
    "\n",
    "additional_outputs.update(normal_outputs)"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Writing datasets: 100%|██████████| 365/365 [01:13<00:00,  4.95days/s]\n"
   "source": [
    "modspa_excel = pd.read_excel(xls_file_path, sheet_name=1, header=10, index_col=0)\n",
Jeremy Auclair's avatar
Jeremy Auclair committed
    "modspa_excel = modspa_excel.loc[:, ~modspa_excel.columns.str.contains('^Unnamed')]\n",
Jeremy Auclair's avatar
Jeremy Auclair committed
    "dates = modspa_excel.index\n",
    "\n",
    "# Open empty dataset to get structure and reindex with correct dates\n",
    "empty_dataset = xr.open_dataset(ndvi_path)\n",
    "empty_dataset = empty_dataset.reindex(time = dates)\n",
    "\n",
    "# Transpose dimensions\n",
    "empty_dataset = empty_dataset.transpose('time', 'y', 'x')\n",
    "\n",
    "# Get the numpy array for 'ndvi'\n",
    "zero_values = empty_dataset['NDVI'].values\n",
    "\n",
    "# Transpose the numpy array for 'ndvi'\n",
    "zero_values = zero_values.transpose([0,2,1])\n",
    "empty_dataset['NDVI'] = empty_dataset.NDVI.transpose('time', 'y', 'x')\n",
    "\n",
    "# Assign the transposed numpy array back to 'ndvi'\n",
    "empty_dataset.NDVI.values = zero_values\n",
    "\n",
    "# Drop ndvi to get empty dataset\n",
    "empty_dataset = empty_dataset.drop_vars('NDVI')\n",
    "\n",
    "# Datasets\n",
    "ndvi_val = empty_dataset.copy(deep = True)\n",
    "weather_val = empty_dataset.copy(deep = True)\n",
    "outputs_val = empty_dataset.copy(deep = True)\n",
    "\n",
    "# Inputs\n",
    "ndvi_val['NDVI'] = (empty_dataset.dims, np.ones(tuple(empty_dataset.dims[d] for d in list(empty_dataset.dims)), dtype = 'uint8'))\n",
    "weather_val['Rain'] = (empty_dataset.dims, np.ones(tuple(empty_dataset.dims[d] for d in list(empty_dataset.dims)), dtype = 'uint16'))\n",
    "weather_val['ET0'] = (empty_dataset.dims, np.ones(tuple(empty_dataset.dims[d] for d in list(empty_dataset.dims)), dtype = 'uint16'))\n",
    "# variables\n",
    "variables = list(set(list(additional_outputs.keys())).intersection(set(modspa_excel.columns[0:-2])))\n",
    "variables.sort()\n",
    "\n",
    "    outputs_val[var] = (empty_dataset.dims, np.ones(tuple(empty_dataset.dims[d] for d in list(empty_dataset.dims)), dtype = 'int16'))\n",
    "\n",
    "# Progress bar\n",
    "pbar = tqdm(total = len(ndvi_val.coords['time'].values), desc = 'Writing datasets', unit = 'days')\n",
    "\n",
    "# Loop in time values\n",
    "i = 0\n",
    "for t in ndvi_val.coords['time'].values:\n",
    "    # Inputs\n",
    "    ndvi_val['NDVI'].loc[{'time' : t}] = ndvi_val['NDVI'].loc[{'time' : t}] * np.uint8(np.round(modspa_excel['NDVI'].values[i] * 255))\n",
    "    weather_val['Rain'].loc[{'time' : t}] = weather_val['Rain'].loc[{'time' : t}] * np.uint16(np.round(modspa_excel['Rain'].values[i] * 100))\n",
    "    weather_val['ET0'].loc[{'time' : t}] = weather_val['ET0'].loc[{'time' : t}] * np.uint16(np.round(modspa_excel['ET0'].values[i] * 1000))\n",
    "\n",
    "    # Outputs\n",
    "    for var in list(outputs_val.keys()):\n",
    "        if var == 'NDVI':\n",
    "            NDVI = np.round(modspa_excel[var].values[i] * 255)\n",
    "            NDVI = NDVI / 255\n",
    "            outputs_val[var].loc[{'time' : t}] = outputs_val[var].loc[{'time' : t}] * np.uint8(np.round(NDVI * additional_outputs[var]))\n",
    "        outputs_val[var].loc[{'time' : t}] = outputs_val[var].loc[{'time' : t}] * np.int16(np.round(modspa_excel[var].values[i] * additional_outputs[var]))\n",
    "    \n",
    "    i+=1\n",
    "    \n",
    "    # Update progress bar\n",
    "    pbar.update()\n",
    "# Reorganize dimension order\n",
    "ndvi_val = ndvi_val.transpose('time', 'y', 'x')\n",
    "weather_val = weather_val.transpose('time', 'y', 'x')\n",
    "ndvi_val = ndvi_val[['time', 'y', 'x', 'NDVI']]\n",
    "weather_val = weather_val[['time', 'y', 'x', 'Rain', 'ET0']]\n",
    "\n",
    "# Save datasets\n",
    "# Inputs\n",
Jeremy Auclair's avatar
Jeremy Auclair committed
    "ndvi_val.to_netcdf(val_ndvi_path, encoding = {\"NDVI\": {\"dtype\": \"u1\", \"_FillValue\": 0}})\n",
    "weather_val.to_netcdf(val_weather_path, encoding = {\"Rain\": {\"dtype\": \"u2\"}, \"ET0\": {\"dtype\": \"u2\"}})\n",
    "ndvi_val.close()\n",
    "weather_val.close()\n",
    "\n",
    "# Create encoding dictionnary\n",
    "for variable in list(outputs_val.keys()):\n",
    "    # Write encoding dict\n",
    "    encoding_dict = {}\n",
    "    encod = {}\n",
    "    encod['dtype'] = 'i2'\n",
    "    encoding_dict[variable] = encod\n",
    "\n",
    "# Outputs\n",
    "outputs_val.to_netcdf(val_outputs, encoding = encoding_dict)\n",
    "outputs_val.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Compare `modspa-pixel` and `modspa-excel` outputs"
   ]
  },
  {
   "cell_type": "code",
Jeremy Auclair's avatar
Jeremy Auclair committed
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABhEAAAHSCAYAAADi0+oXAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABnS0lEQVR4nO3dd3gU1f7H8c9ueg+EDgGEKC0gohEivQiIgBfFi1cRsAJXFL2I7QqCitiwITaKCpZ7VRQBqaEXQ7kQekSCYggQejohZX5/7C8Ly05CElI25P16njxPMufMme85s3uSzHfnjMUwDEMAAAAAAAAAAACXsJZ3AAAAAAAAAAAAwDWRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAACQj9zcXP3+++9asWJFeYdS5ipz3wEAwAUkEQAAAAAAMBETE6Phw4fruuuu08yZM8s7nDJVmfsOAAAckUQAAAAAAFRIy5YtU3BwsH755ZdSab9169Z66KGHJEndu3cvlWO4qsrcdwAA4IgkAgAAAACgQkpMTFRSUpIOHz5casfIW8qnR48epXYMV1WZ+w4AAC6wGIZhlHcQAAAAAAAUx19//aX69euXWvvdunVTfHy8fv/991I7hquqzH0HAAAXcCcCAAAAAKDCKs0EQkZGhjZu3FgpP4lfmfsOAAAckUQAAAAAALiE+fPna+DAgWrcuLGmTZumH374Qd26dVNgYKDatWun2NhYSdLWrVv10UcfafTo0Xrttdfs+//666+aMGGCatasqbvvvluSlJWVpUceeUQhISFavXq1vW5ubq5mzZqlESNGaOLEiQoLC9PixYsd4lm3bp0yMzPL5JkABcVTlH5J0tmzZzVu3Dg9++yzeuutt9S3b1/9/PPPhT6eVLZ9BwAArs29vAMAAAAAAECS+vfvr9jYWM2dO1cJCQlq0KCBVq5cqTVr1qhLly7617/+pUWLFqlFixbauHGjPvjgA02dOtW+f2RkpCIjI/XXX39p9uzZSkpK0jvvvKN7771XERERCg8PlySdPn1ad955p2rWrKk5c+bI09NTv/zyi+bOnavbbrvN3t6KFStktVrVrVu3Uu335eIpbL8kaffu3erXr5/ef/999e/fX5L03XffacyYMbrjjjsK3f+y6jsAAHB93IkAAAAAAHAZW7ZskdVqVffu3dW3b19JUocOHSRJ+/btkyT5+PgoPj5ektSvXz+nNvr06aOcnByNHj1at912m7p27apHH31U1apVU0ZGhnr37q3jx4/riy++kKenp6Kjo7V9+3b9/e9/d2gnKipKN9xwg6pWrVpq/S1KPAX1S5ISEhLUtWtXjRgxwp5AkKTRo0frjTfeKNLxyqLvAACgYiCJAAAAgAopKytLR44cKe8wAJSgnJwcRUVFqUOHDg7L6CQlJUmSqlevbt+2cOFChYeHq0GDBk7t5CUdcnJy1K5dO4ey8ePHa8uWLXrllVfk7e2t7777ToMHD9bs2bPVs2dPe71Tp04pJiamwGcCDBs2TBaLpdBfJ0+edGqjsPFcrl+S9Nhjj8lqterJJ5902D548GDdddddhT5eYfoOAAAqD5YzAgAAKEOGYSgmJkY///yzrrvuOt17771lctxNmzZpwIABql27tqKiolSlSpUyOW5JO336tJYvX66FCxfql19+0RNPPKEJEyY41ClKXzdu3KgpU6YoNjZWcXFxCg0N1YMPPqjnnntOFoulyPUAOLNYLHJzc1N2dvZl627evFlnz57V/fff77A9OjpakuwXuQ8ePKjY2Fg999xzpu1Uq1ZNVapU0ZkzZxy2p6Wladq0afLy8tL+/fs1atQoNW3aVFu3blVwcLBD3VWrVik3N7fAC+kPP/ywunTpctl+5QkICCh2PAX1S5L+/PNPzZ8/X0OHDpWXl5fp8Qt7vML0vSBFOecAAMD1kUQAAABXva+//lpfffWV1q5dq/T0dKdyi8Uif39/hYSEqFWrVurWrZsGDx6skJCQEoth165deu2117Ry5UodP35ckvTSSy+VWPuXM2fOHB09elRHjx7VypUr7Z9IrWhmzZql6OhozZ07N986he3r999/r6FDh2rLli1q3ry5/QGkL7zwgmrXrq1hw4YVqR6AK7d06VJZrVb72v15pk+froCAAD3xxBOSbHchSLIvd3SpSZMmqXnz5lq3bp1ycnLk5uYmyZZkzMjI0MiRI/X8888XGEtUVJS8vLzUvn37fOt06NDBfndAcRQlHin/fklSTEyMDMNQWFjYFR+vMH0HAACVB8sZAQCAq959992nxYsX2y865Zk2bZr279+vHTt2aM6cOerXr59WrlypJ598UvXr19d7771XYjE0bdpUs2fPdngAaFkaNGiQatSoodatW6tz587lEkNJePrpp/X9998rKCgo3zqF6evp06f18MMPq3379mrRooUsFosmT56sHj16yN3dXXXr1i1SPQAlY+nSpYqIiHBYtmjp0qX6+eefNWvWLNWoUcO+LSQkRJGRkdq1a5dDG8uXL1fDhg01YsQIJScna9u2bZJsydyMjAxJUmBgoNOx8563kCcqKkrt27eXj4+PFi1aJMMwSrSvkooUT0H9kiQPDw9JMv30/++//16k45VF3wEAQMVBEgEAAFQabdu2dfh5+PDhuvbaa9WyZUvdcccd+uCDD7R37141btxY6enpeuqppzRp0qQSObaHh4c8PDycPl1bUjIyMgpMUHTs2FGJiYnavn27/QGcFZXFYjFd5iNPYfr6+eefKzk5WbVr17Zvc3Nz07Jly3T69GndeuutRaoH4MqdOXNGW7ZscVjrf/v27Ro2bJhmz56tgQMH2rdv2bJFnTt31pIlSxQbG6uYmBhNmTJFK1eu1KpVqzR06FB169ZNku0uhqlTpyo7O1sRERHy9/fX119/bb+wHhcXp1GjRikrK8ve/okTJxQXF6eIiAitX79e3t7epbJ02eXiKWy/JNvvuKCgIH3//fdKTU2VZFu+aPz48fYkQGH6X1Z9BwAAFQdJBAAAUGn4+vo6/HzxMhB5QkND9dFHH9l/njRpkpKTk0sshvzWqb5So0eP1v/+979SaftqtGrVKkmS1er457DFYnFYs7yw9QBcuaioKOXk5Oivv/7S9OnT9eqrr+rTTz/VmjVrNHjwYIe6TZo00fr163Xo0CHdfffd+t///qcJEyZo1qxZmjhxoiSpTp06GjhwoObNm6fQ0FDdcMMNqlGjhn766SeFhITo+uuv1w033KCPPvpIEydOVKtWrezt+/n5qX79+vr555+VmJhov3Bf0i4XT2H7JdmelzBv3jz5+fmpVatWuuuuuzR+/HiNGjVK1113XaGOV5Z9BwAAFQfPRAAAALhE9+7d5ePjo4yMDGVkZCg6Otr+ME9XNGPGDE2fPl1Dhw4t71AqjPj4+BKtB+DKLV26VO7u7pozZ478/PwKrLtu3TqHnx966CE99NBDTvW+//57p209evRQTExMge37+vrq0KFDlw+6BBQUT1H6JUldunTRpk2bin08qWz7DgAAKgbuRAAAALiEm5ubw3I5J0+eLL9gLmP27NkaMWJEoeunpKTo6NGjpRiR6yior2lpaYVqo7D1AFy5pUuX6qabbrpsAgEAAABliyQCAADAJbKysnTq1Cn7zxc/4PNS6enpmjx5sm655RY1btxYAQEBuvHGGzVlyhSdP3++WMffsGGD/vGPf6hly5Zq0KCBqlatqt69e2v16tUO9X744Qe98sor9rWqf/zxR4WFhSksLExPPvmkvd758+f1888/6+9//7tq1qyppUuX2su+//57ubu7y2Kx2L98fHz0/vvvOxzrb3/7m73cw8PD6VOwpTEOkvT111+rY8eOatiwoRo2bKi77rrL6SGqFyuor0uWLLGPT96nbC8es7CwMK1bt67Q9a60/3v27NGYMWP0xBNPSJIWL16sli1byt/fXy+++GKx2zcMQ+vWrdPDDz+sPn36SJKOHTumESNGqFatWvLz89Ptt9+uP//8s8CxX7Zsmfr376+mTZvqmmuuUdOmTfX444/rjz/+MK1fWq+BefPm6W9/+5vCw8MVGhqq6tWrq2fPnvr6669NH/ZaUv3Pzx9//KEJEyYoLCxMki3RNHHiRDVr1kw+Pj5q1qyZ5s6da6+fnZ2tadOm2S+O16tXT6+88kq+D6otzjgWds4oqzEqjl27dunw4cPq0qVLmR0TAAAAhWQAAABUIpLsX/mZP3++vU5wcLCRmppqWu/QoUNGs2bNjI8++sjIzs42DMMwYmNjjaZNmxqSjPbt2xspKSn5xvDSSy85lY0bN86QZDz44INGenq6YRiGsXnzZiMoKMhwc3MzFi5c6LRPgwYNDEnG0KFDncq2b99uDBo0yAgNDbUf9/PPP3eoc+DAAaN69er28i1btpj2d9y4cYbFYjHWr19fIuNQkPPnzxsDBw40goODjXnz5hm5ubmGYRjG8uXLjYYNGxre3t5OY1iYvuYpaMyKWq+o/f/kk0+Mtm3b2mMcNGiQMX/+fHufJBm+vr7Fan/x4sVG+/bt7e20bdvW+N///mfUqlXLqFmzphEQEGAva9KkiZGZmWk69vfff7/RrFkzY8OGDfbtL730kj225cuXX9EYFEZSUpJx2223GcHBwcbChQvtr4Ho6GijWbNmhiSjW7duxunTp+37lET/87N9+3aja9euhsVisbdx9OhRo2XLloaHh4dRv359w83NzZBkuLm5GRs3bjQSEhKMm2++2QgICDBCQ0MNq9Vq3/fdd991OkZxxrGoc0ZpjlF+8sYkP5s3bzY6d+5sSDL69etnrFu37oqPifJ1uXMOAAAqFpIIAACgUrlcEuHgwYPGNddcc9mL0OfOnTPCw8ONsWPHOpVdnIQYN25cvjFcmkTYu3evvWz16tUOZc8++6whyYiIiHBqrzAXulevXl1gn6ZOnWovnzt3rmkbY8aMMe644w6HbVcyDgV5+OGHDUnG999/71S2cuXKAhMxl+urYZRcEqE4/U9JSTFSUlKMatWq2S/i9u/f3zh16pSxePFio0GDBsbAgQOL3X5ubq5x4403GpKM+vXrG23atLFf9M/JyTFGjBhh3+/HH390aveee+4xqlevbhw7dsxh+6lTp+wXwVu2bHlFY3A5OTk5Rs+ePQ1Jxtdff+1UnpiYaNSoUcOQZHTt2tV+wb0k+p+ftLQ0Izk52Rg4cKB9/w4dOhiffPKJkZycbBiGYcTFxdnPa9euXY22bdsa//nPf+zxHTp0yGjcuLEhyWjQoIFD+8UZx+LOGaU1RvnhgnLlwzkHAODqQhIBAABUKmZJhIyMDCM6OtoYM2aM4e/vb7/AV9DFs7feesuQZOzdu9epLCMjw/Dw8LBfoMsvhksvgC9dujTfC/kzZswwJBl+fn5O7RXmgnhcXFyBF9YzMjLsdyP06dPHqTw7O9uoW7eusWzZMoftVzIO+clLEoSFhRk5OTmmdYKDg/NNIlyur4ZRckmEK+l/3qfB/f39jePHj5do+/fdd58h2e6kOXz4sEPZqVOn7J+mf/bZZx3KfvzxR0OS8eabb5rGM2DAAHvioyTGID9z5swxJBmhoaEOCYKLffbZZ/bzPGvWLIey4va/MD755BP7cVeuXOlU/sQTT9jLFy9e7FT+zjvv2MsTEhLs24szjsWdMwyjdMfoUlxQrnw45wAAXF3cBQAAUEk1adLE6eG7Tz31lG677TZ17dpV7u75/6k0Y8YMSdKAAQNMy4ODg5Wbm6vMzEzl5OTIzc3tsvF0795dzz33nLKzs3X77bc7lOXtn56eftl2zFitBT8Ky9vbWyNGjNArr7yiJUuW6M8//1TDhg3t5UuWLJGPj4969OjhsF9pjMOUKVMkSb1798437qCgIJ09e9a07HJ9LUlX0v+811fr1q3zfe5GcdvPazsoKEh169Z12Kdq1aqqVq2aTpw44fTQ8HfeeUeS1LdvX9PjffPNN1q3bp3atWt3xTEW5LPPPpMktWzZMt/69957r0aPHq2MjAxNnz5dDzzwgL2suP0vDC8vL/v311xzjVN53rMSJKlp06YFlh8+fFh16tSRVLxxvJI5ozTHCAAAAFcXkggAAKDS+u2335SRkaH27dtr+/btkqTExETdeuutBe53/Phx/fbbb5Kk3bt3F5hsKAo3NzdNnjzZYduBAwc0Y8YMLVmyRJLyfRhrSRgxYoRef/11ZWVl6bPPPtNrr71mL5sxY4aGDx9uf4izVDrjkJmZqRUrVkiSQkNDr7i90lRS/c/vInlpvc4kydfXV5Ltob95kpOT9euvv0oyvzgu2ZJNF78/SiPG7OxsRUdHS5LTxe2L+fn56frrr1d0dLS2bNmi7OzsQh/frP8lJSAgoNDlGRkZkoo/jqU5Z5TmGAEAAKBiKbuPaQEAALggHx8f/fTTT6pWrZok2yetL754biY+Pt7+fVpaWqnEtWbNGvXv319vvfWWhg4dqieffLJUjnOxOnXq6K677pIkzZw5U1lZWZJsiZWoqCiHT3pLpTMOBw4c0Llz5yTZPg3tykr7dVAWr7OL/fXXX8rJyZEkh2RRQUojxlOnTtlfe+fPny+wbqNGjSTZLnRXxE/M513gL4lxLI85AwAAAJUDSQQAAFDpNWjQQP/5z3/snwh/8cUXNW/evHzrZ2Zm2r//66+/SjSW48eP6/bbb9f999+vl19+WZ9++qmaNWtWoscoyBNPPGGP48cff5Qkffnll7rjjjsUEhLiULc0xuHiJYqSkpJKpM3SUpqvg7Jo/1IXj3dCQkKh9imNGC9ejupyiYGLP9Xv6elZIscvD1cyjuU9ZwAAAODqRxIBAABAtucRvP7665Jsnw4ePHiwdu7caVq3SpUq9u9XrVp12bbzlma5nFOnTikyMlKLFi3St99+q9atWxdqv5IUGRmpm266SZL08ccfS7LdlTBixAinuqUxDn5+fvbvd+/efdn65am0Xgdl1f6VHq84+xQmxpCQEPtSOvm9B/PkXXwPCAhwiKWiKe44usKcUZIsFgtfFeALAABUPiQRAAAA/t/TTz+tQYMGSbItKdK/f38dP37cqV5YWJh8fHwkSR988IHDp4gv9dFHH2nr1q2FOv6rr76qgwcPqkaNGmrfvn2R48/NzS3yPmYef/xxSbblUT799FN5enqqQ4cOTvVKYxzCwsLsd4T88ssvl13OpjyV1uugrNq/1LXXXmu/eD916lT70kZmpk6dqtjY2FKJ0Wq1qkuXLpJsy/zs378/37p//vmnJKlbt24V+uJmccfxSucMV2MYBl8V4AsAAFQ+JBEAAEClUdBF0TwzZ85Uy5YtJUmHDh3SnXfe6XQh28PDQ7fffrskKS4uTk899ZRpW8uXL9cbb7zh9CyB/OQ9UNjMwYMH7d/nrRefJ+/Cb0FLvxTlws+gQYNUo0YNSbbljczuQpBKZxz8/f3VrVs3SdKJEyf09ttvO9XJycmxrxufnp7uVF6YvubVuVzdguqV1uugrNo3O17eMzF27NihMWPGmNabP3++1q1bp6ZNm5ZajHnLaknSrFmzTOukp6dr06ZNTvUrouKOY3HnDAAAAKAoSCIAAIBK49KL7GZr7vv5+emnn36yLy+yYcMGPfroo071/v3vf8vDw0OSbdmf22+/XVFRUTp06JCio6P11FNPqW/fvpoyZYrDEj1nzpyxf3/pBfC8ZMDx48f1ySef2LcvWLBAX331lf3n7du3KzU11d6fOnXqSJLWr1+vI0eOyDAMvf/++9q3b599n5SUFNPvzXh5edn77OHhofvvvz/fusUdh4KMHz/evi7+uHHj9NZbbykjI0OS7ZwNHTpUycnJkqQff/xR33//vXbt2lWkvuad+4ufwVCcelfS/4v7lJ/itl/QJ9klKTU1VZJzYm38+PH2dt5//31169ZNP/zwg3bs2KGVK1fqscce0/Dhw/Xee++VyBjkp1evXho8eLAk2yfz8+44uNi0adOUkZGhBx980J54utL+F8bFDz42u/sn78Hgkvm5zTu25HhxvzjjWNw5QyrdMbra5Obm6vfffy8waXO1qsx9BwAAFzEAAAAqgZycHOPVV181JNm/PvzwQyM3N9e0/uLFiw2r1WqvO3z4cCMlJcWhzsyZMx3qXPzl5uZmvPvuu04xvPzyy/Y6N9xwg5GUlGQvf/PNNx3aqFmzplG1alWjU6dOxsaNG+3b/f39jZCQEOO3334zDMMw3n//fXuZt7e3UbNmTeOhhx5yOPbEiRPtdXr27Gmkp6cXOF4JCQmGm5ub8eijj152bIs6DoXx/vvvGxaLxd6Oj4+PUb9+fcNqtRoTJ040GjRoYN/+wAMPGDt37ix0XxcvXmwvr1KlihEXF2caQ2HrFaf/0dHRho+PjyHJ8PDwMDZt2pTvWBSn/datWxuSDF9fXyMtLc2h7NChQ/axbdeundO+v/zyi+Ht7W16vKpVqxq//vpricR4OefOnTPuueceQ5LRrFkzY9u2bYZh2N5Hc+bMMTw9PY0hQ4YY58+fL9H+FyQjI8O47bbb7H37/PPPHcpzc3ONoUOH2sufeeYZhzkmMzPToXzs2LEO5UUdx+LOGaU5Rmby4q+Itm/fbjz88MOGJOMf//hHeYdTpq6k7xX5nAMAAGckEQAAwFVvypQpRlBQkOmFucDAQKcL7nkmTZrkUNfLy8sYPHiwQ53NmzcbAwcONGrUqGG4u7sbtWrVMgYNGmRs2bLFoV5UVJTh7+/vdHxPT09j4MCBhmEYRnZ2tjFhwgSjfv36ho+Pj3H99dcbH3zwgZGdnW0YhmH885//NHx9fY3WrVsbq1atsredmZlpPPTQQ0ZwcLBRv35947XXXjNycnIMwzCMDz/80AgODnY6rre3t/Haa68VOG5Dhw51uOhYkMKOQ1GsXbvW6NevnxESEmJ4e3sbbdq0Mb7++mvDMAyjV69exrvvvmucOXPGXv9yfZ03b55RvXp109dBrVq1jG+++cYwDKPQ9Yrb/5tuusm07Zo1axorV668ovFds2aN0bt3b4d227RpY0ydOtUwDMOYMWOGce211zqU33777caGDRsc2vntt9+M+++/36hdu7bh4eFh1K9f3xg1apRx5MiRfM9XabwGDMMwfv75Z6NPnz5GtWrVjPr16xvNmzc3Bg4caCxbtsypbkn138xHH31kT/xcet7i4uKMr7/+2nSeCQgIMBYtWmRs2bLF8PT0NJ2DFi1aVKxxLM6cUZpjlJ/SvKC8dOlSIygoyFi4cGGptG8YhvHrr78akowZM2aU2jFcVXH7ThIBAICri8UweDISAAAAAKB0WCwWubm5KTs7u8TbnjNnjoYMGaJPPvlEw4cPL/H2JWnSpEl68cUX9eeff6pBgwaXrf+f//xHK1eu1GeffVYq8ZSlovY9T2mecwAAUPbcyzsAAAAAAACK4/7771fnzp1Vv379UjvGihUrFBYWVuiL6LGxsVq2bFmpxVOWitp3AABwdeLBygAAAACACqs0EwgZGRnauHGjevToUWrHcFWVue8AAMARSQQAAAAAgEuYP3++Bg4cqMaNG2vatGn64Ycf1K1bNwUGBqpdu3aKjY2VJG3dulUfffSRRo8erddee82+/6+//qoJEyaoZs2auvvuuyVJWVlZeuSRRxQSEqLVq1fb6+bm5mrWrFkaMWKEJk6cqLCwMC1evNghnnXr1ikzM1Pdu3cv9b4XFE9R+iVJZ8+e1bhx4/Tss8/qrbfeUt++ffXzzz8X+nhS2fYdAAC4NpYzAgAAAAC4hP79+ys2NlZz585VQkKCGjRooJUrV2rNmjXq0qWL/vWvf2nRokVq0aKFNm7cqA8++EBTp0617x8ZGanIyEj99ddfmj17tpKSkvTOO+/o3nvvVUREhMLDwyVJp0+f1p133qmaNWtqzpw58vT01C+//KK5c+fqtttus7e3YsUKWa1WdevWrVT7fbl4CtsvSdq9e7f69eun999/X/3795ckfffddxozZozuuOOOQve/rPoOAABcH0kEAAAAAECp+fzzz2W1Fv4m+C1btshqtap79+72T8F36NBBkrRv3z5Jko+Pj+Lj4yVJ/fr1c2qjT58++vzzzzV69GiNGDFC7dq1U9euXSXZlunp3bu3UlNTtXjxYnl6eio6Olrbt2/Xq6++6tBOVFSUbrjhBlWtWtU01j///NNp29mzZ5WdnW1aVq1aNfn7+ztsK0o8BfVLkhISEtS1a1c9/fTT9gSCJI0ePVo+Pj5FOt7l+l6Qop5zAADg2kgiAAAAAABKzbBhwwpdNycnR1FRUerQoYPDMjpJSUmSpOrVq9u3LVy4UOHh4aYP/c1LOuTk5Khdu3YOZePHj9eWLVv0ww8/yNvbW999951eeOEFzZ49Wz179rTXO3XqlGJiYjR27Nh8473mmmuKVPb55587jUdh47lcvyTpsccek9Vq1ZNPPumwffDgwUU6XmH6XpCinHMAAOD6SCIAAAAAAFzC5s2bdfbsWd1///0O26OjoyXJfpH74MGDio2N1XPPPWfaTrVq1VSlShWdOXPGYXtaWpqmTZsmLy8v7d+/X6NGjVLTpk21detWBQcHO9RdtWqVcnNzC3yw8E8//eS07T//+Y9Wrlypzz77zKmsTZs2xY6noH5Jtrsi5s+fr6FDh8rLy8s03sIerzB9BwAAlQdJBAAAAACAS1i6dKmsVqt97f4806dPV0BAgJ544glJtrsQJKlv376m7UyaNEnNmzfXunXrlJOTIzc3N0nSpk2blJGRoZEjR+r5558vMJaoqCh5eXmpffv2+db529/+5rQtJiZG0dHRpmWXKko8Uv79yjuuYRgKCwu74uMVpu8AAKDyYJFCAAAAAIBLWLp0qSIiIhyWLVq6dKl+/vlnzZo1SzVq1LBvCwkJUWRkpHbt2uXQxvLly9WwYUONGDFCycnJ2rZtmyRp165dysjIkCQFBgY6HTvveQt5oqKi1L59e/n4+GjRokUyDKNE+yqpSPEU1C9J8vDwkCRlZ2c7tfX7778X6Xhl0XcAAFBxkEQAAAAAAJS7M2fOaMuWLQ5r/W/fvl3Dhg3T7NmzNXDgQPv2LVu2qHPnzlqyZIliY2MVExOjKVOmaOXKlVq1apWGDh2qbt26SbLdxTB16lRlZ2crIiJC/v7++vrrr+0X1uPi4jRq1ChlZWXZ2z9x4oTi4uIUERGh9evXy9vbWxaLpcT7fLl4CtsvSWrbtq2CgoL0/fffKzU1VZJt+aLx48fbkwCF6X9Z9R0AAFQcJBEAAAAAAOUuKipKOTk5+uuvvzR9+nS9+uqr+vTTT7VmzRqHBwNLUpMmTbR+/XodOnRId999t/73v/9pwoQJmjVrliZOnChJqlOnjgYOHKh58+YpNDRUN9xwg2rUqKGffvpJISEhuv7663XDDTfoo48+0sSJE9WqVSt7+35+fqpfv75+/vlnJSYm2i/cl7TLxVPYfkm25yXMmzdPfn5+atWqle666y6NHz9eo0aN0nXXXVeo45Vl3wEAQMVhMbgvEQAAAABQzh5++GF9+eWXOnv2rPz8/Mo7nGL75JNPNH/+fC1atKi8QwEAACgRJBEAAAAAAOUuNDRU9erV06+//lreoQAAAOAiLGcEAAAAAChXu3bt0uHDh9WlS5fyDgUAAACXqLBJhFWrVmn+/Pkl2mZycrKWL19eom0CAAAAAPK3ZcsWPf7445KkPXv2aP369eUcEQAAAC5W4ZII27ZtU69evdStWzdt27atWG3ExsbKarXKYrE4fAUFBclisZRwxAAAAACA/ERERGj16tUyDEPz589Xhw4dyjskAAAAXMS9vAMorLNnz2r69OlKTExUdHT0FbU1ZcoUmT0KolmzZurevfsVtQ0AAAAAAAAAwNWiwiQRgoKCNHbsWEnS4cOH9d///rdY7SQmJmru3Lk6ePCgqlSp4lDm7+9f5DsRcnNzdeTIEQUEBHAXAwAAAAAAAACgQjAMQykpKapTp46s1vwXLaowSYSLL9B7e3sXu52pU6dq8ODBuuaaa0oiLB05ckShoaEl0hYAAAAAAAAAAGUpPj5e9erVy7e8wiQRSkJaWpo+/vhjNWjQQM8884y6du2qbt26ycvLq9htBgQESLINdGBgYEmFCgAAAAAAAABAqUlOTlZoaKj9Gnd+KlUSYebMmTp9+rROnz6t7du366233lKVKlX0zDPPaOzYsXJzcytym3l3SAQGBpJEAAAAAAAAAABUKJdbpj//hY6uQpGRkZo+fbqef/55tWvXTpJ05swZPf/88xowYIBycnLKOUIAAAAAAAAAAFxHpboTISIiQhEREfafY2JiNGrUKG3YsEELFizQ+PHjNWnSpALbyMzMVGZmpv3n5ORkSVJKSopDxsbd3V0+Pj7Kzc1VWlqaUzt5t4ikp6c7JS+8vb3l4eGh8+fPOxxLktzc3OTr6yvDMJSamurUrp+fn6xWqzIyMpSdne1Q5uXlJU9PT2VlZencuXMOZVarVX5+fva+XMrX11dubm46d+6csrKyHMo8PT3l5eWl7OxsZWRkOJRZLBb5+/tLklJTU2UYhkO5j4+P3N3dlZmZqfPnzzuUeXh4yNvbWzk5OUpPT3eKKW8M09LSlJub61BW0Bhe7tzkPWDb7NwUNIZ550YyH8OCzk1BY3jxuTEbw7xzU9QxvPjcmI1h3rkp6THMOzcFjeGVvL6vZAwLen1fyRianRvmCBvmiAuYI2yYI2yYI2yYIy5gjrBhjrBhjrBhjriAOcKGOcKGOcKGOeIC5ggb5ggb5ggb5giZ/pyfSpVEuFTr1q21du1aPfDAA5o9e7Y++OADjR07VsHBwfnuM3nyZE2cONFp+44dO+wvaEmqUaOGmjVrpszMTG3bts2pfufOnSVJsbGxTieradOmqlmzpk6cOKEDBw44lFWpUkWtWrVSTk6OabuRkZHy9PRUXFycTp065VDWqFEjhYaG6uzZs9q7d69Dmb+/v2688UZJ0vbt253eODfddJP8/Px06NAhHTt2zKEsNDRUjRo1Umpqqnbs2OFQ5unpqcjISEnSrl27nN50119/vYKDg5WQkKD4+HiHslq1aqlJkyY6d+6cU18tFos6deokyTaGl05OzZs3V/Xq1ZWYmKiDBw86lIWEhCg8PFzZ2dmmY9i+fXu5u7vrwIEDOnPmjENZWFiY6tatq9OnTys2NtahLCAgQG3atJEk03Zvvvlm+fj46M8//9Tx48cdyho0aKCGDRsqOTlZu3btcijz9vZW27ZtJUk7d+50mhBbt26toKAgxcfHKyEhwaGsTp06uvbaa5Wenu4Uk5ubmzp06CBJ2rt3r9Ok16JFC1WrVk3Hjh3TH3/84VBWrVo1tWjRQllZWaZ97dixoywWi/bv36+kpCSHsuuuu061a9fWyZMntX//foeyoKAgtW7dWoZhmLbbrl07eXl56eDBgzp58qRD2TXXXKP69evr7Nmz2rNnj0OZr6+vPYEYExPjNNG2adNGAQEBio+P15EjRxzK6tatq7CwMKWmpiomJsahzMPDQ7fccoskaffu3U6TdMuWLVW1alUdPXpUhw4dcihjjrBhjriAOcKGOcKGOcKGOeIC5ggb5ggb5ggb5ogLmCNsmCNsmCNsmCMuYI6wYY6wYY6wYY6wMUtWmLEYl458BTBs2DB9+eWXeumllzRhwoQrbi8jI0NNmjRRfHy81qxZYz/RZszuRAgNDdXhw4cdnolAVs+GrN4FZP5tyPzbMEfYMEdcwBxhwxxhwxxhwxxxAXOEDXOEDXOEDXPEBcwRNswRNswRNswRFzBH2DBH2DBH2DBH2CQnJ6tevXpKSkoq8Hm/JBH+39NPP60pU6bo22+/1T333FPo/ZKTkxUUFHTZgQYAAAAAAAAAwFUU9tp2pXqwckHCwsIkiUQAAAAAAAAAAAD/jyTC/8u79aRp06blHAkAAAAAAAAAAK6BJML/W758uSIiItSoUaPyDgUAAAAAAAAAAJdQIZMIeQ/CuPTBFHmioqIUHh6uKVOm2LcdOXJEr776qhYuXGhaf9WqVZo2bVrpBAwAAAAAAAAAQAVU4ZIIp06d0vr16yVJGzZscHqatyS999572rNnj8NDl3/44QeNGzdO/fr1U8+ePbVv3z7l5uZq7ty5evjhh7VgwQJFRESUVTcAAAAAAAAAAHB5FSaJkJOTo5tuukkNGjTQkSNHJEkrV65UnTp1NGDAAIe6gwYNUkBAgIYOHWrf9sADD2jkyJGqV6+eVq9erbZt26pdu3basWOHNm/erFtvvbVM+wMAAAAAAAAAgKuzGIZhlHcQFVlycrKCgoKUlJSkwMDA8g4HAAAAAAAAAIDLKuy1bfcyjAkAAAAAAABABbGgcePyDgEoFf3i4so7hAqlwixnBAAAAAAAAAAAyhZJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMVdgkwqpVqzR//vzyDgMAAAAAAAAAgKtWhUsibNu2Tb169VK3bt20bdu2YrWxatUq9ezZU23btlWbNm3UvXt3rV69umQDBQAAAAAAAACggqswSYSzZ8/qrbfe0jfffKPo6OhitzNz5kx1795dAwYM0KZNm7Rt2zb1799fPXr00OzZs0swYgAAAAAAAAAAKjb38g6gsIKCgjR27FhJ0uHDh/Xf//63yG1s2rRJI0eOVPfu3TVy5Ej79tGjR2vevHkaPny4br75ZjVt2rTE4gYAAAAAAAAAoKKqMHciWCwW+/fe3t7FamPChAnKysrSgw8+6FT20EMP6dy5c3rppZeKHSMAAAAAAAAAAFeTCpNEuFIJCQlasmSJJKlLly5O5Z07d5YkzZs3T0lJSWUZGgAAAAAAAAAALqnSJBE2btwoSQoODlbt2rWdykNDQxUYGKjz58/zkGUAAAAAAAAAAFSJkgh79uyRJNWtWzffOnXq1JEk7d69u0xiAgAAAAAAAADAlVWYBytfqdOnT0uSAgIC8q0TGBgoSTp+/Hi+dTIzM5WZmWn/OTk5WZKUkpLi8NwGd3d3+fj4KDc3V2lpaU7t5MWRnp6unJwchzJvb295eHjo/PnzDseSJDc3N/n6+sowDKWmpjq16+fnJ6vVqoyMDGVnZzuUeXl5ydPTU1lZWTp37pxDmdVqlZ+fn70vl/L19ZWbm5vOnTunrKwshzJPT095eXkpOztbGRkZDmUWi0X+/v6SpNTUVBmG4VDu4+Mjd3d3ZWZm6vz58w5lHh4e8vb2Vk5OjtLT051iyhvDtLQ05ebmOpQVNIaXOzf+/v6yWCym56agMcw7N5L5GBZ0bgoaw4vPjdkY5p2boo7hxefGbAzzzk1Jj2HeuSloDK/k9X0lY1jQ6/tKxtDs3DBH2DBHXMAcYcMcYcMcYcMccQFzhA1zhA1zhA1zxAXMETbMETbMETbMERdcyRwBXK1SU1OZI0x+zk+lSSLknVBPT89863h5eUmS6UnKM3nyZE2cONFp+44dO+y/9CSpRo0aatasmTIzM7Vt2zan+nnPYIiNjXU6WU2bNlXNmjV14sQJHThwwKGsSpUqatWqlXJyckzbjYyMlKenp+Li4nTq1CmHskaNGik0NFRnz57V3r17Hcr8/f114403SpK2b9/u9Ma56aab5Ofnp0OHDunYsWMOZaGhoWrUqJFSU1O1Y8cOhzJPT09FRkZKknbt2uX0prv++usVHByshIQExcfHO5TVqlVLTZo00blz55z6arFY1KlTJ0m2Mbz0D5jmzZurevXqSkxM1MGDBx3KQkJCFB4eruzsbNMxbN++vdzd3XXgwAGdOXPGoSwsLEx169bV6dOnFRsb61AWEBCgNm3aSJJpuzfffLN8fHz0559/OiWqGjRooIYNGyo5OVm7du1yKPP29lbbtm0lSTt37nT6o6l169YKCgpSfHy8EhISHMrq1Kmja6+9Vunp6U4xubm5qUOHDpKkvXv3Ok16LVq0ULVq1XTs2DH98ccfDmXVqlVTixYtlJWVZdrXjh07ymKxaP/+/U7PGLnuuutUu3ZtnTx5Uvv373coCwoKUuvWrWUYhmm77dq1k5eXlw4ePKiTJ086lF1zzTWqX7++zp49a7/zKI+vr68iIiIkSTExMU4TbZs2bRQQEKD4+HgdOXLEoaxu3boKCwtTamqqYmJiHMo8PDx0yy23SLLdwXTpJN2yZUtVrVpVR48e1aFDhxzKmCNsmCMuYI6wYY6wYY6wYY64gDnChjnChjnChjniAuYIG+YIG+YIG+aIC65kjgCuVrt27WKOUMHXwR36ZVw6O1cAw4YN05dffqmXXnpJEyZMKNQ+jz/+uD788ENFRkban49wqcjISEVHR+vxxx/XBx98YFrH7E6E0NBQHT582H4ng0TmPw+Z/wv4dJANnw6yYY6wYY64gDnChjnChjnChjniAuYIG+YIG+YIG+aIC5gjbJgjbJgjbJgjLriSOWJlq1ZO7QFXg647djBHyHZtu169ekpKSnK4tn2pSpNEmDRpkl588UW1bNlSO3fuNK3TsmVL7d69W6+++qr+/e9/F6rd5ORkBQUFXXagAQAAAAAAgIpkQePG5R0CUCr6xcWVdwguobDXtivN4mbNmzeXJB09ejTfOnllLVq0KJOYAAAAAAAAAABwZZUmidCpUydZrVadPHlSJ06ccCpPTEzUqVOn5ObmZl+3CgAAAAAAAACAyqzSJBFCQkJ06623SpLWrVvnVL527VpJUs+ePVW1atUyjQ0AAAAAAAAAAFdUIZMIeQ96ufTBFHmioqIUHh6uKVOmOGx/4YUXZLFY9PnnnzvtM3PmTFmtVr3wwgslHzAAAAAAAAAAABVQhUsinDp1SuvXr5ckbdiwwemp8pL03nvvac+ePU4PXe7UqZMmTJighQsXasaMGfbtn376qZYuXapXXnlFHTp0KNX4AQAAAAAAAACoKNzLO4DCysnJUdu2bRUbG6u0tDRJ0sqVK1WnTh21b99eP/30k73uoEGDtHbtWg0ZMsSpnfHjx6t169Z65513NHPmTBmGIW9vby1YsEB9+/Yts/4AAAAAAAAAAODqLIZhGOUdREWWnJysoKAgJSUlKTAwsLzDAQAAAAAAAErEgsaNyzsEoFT0i4sr7xBcQmGvbVe45YwAAAAAAAAAAEDZIIkAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgqtSRCWlqajh8/rtzc3NI6BAAAAAAAAAAAKEXuJd3ggQMHNGbMGF1zzTXy8fHRsWPHlJWVpaFDh+rWW28t6cMBAAAAAAAAAIBSUugkwu7duxUeHn7Zes8995w++eQT1a5d277t/Pnz+uCDD/TTTz/pww8/lNXKKkoAAAAAAAAAALi6Ql/N79ixo0aMGKETJ04UWC8+Pl7VqlVz2Obp6amnn35ad9xxh8aNG1e8SAEAAAAAAAAAQJkqdBLhyy+/1PTp0xUZGanXX39d58+fN63XpEkT3XPPPUpNTXUq69WrlzZt2lT8aAEAAAAAAAAAQJkpdBKhf//+Gjt2rO666y65ubkpIiJC3333nVO9f//731qyZImaNWumd955R0ePHrWXnTx5UgkJCSUTOQAAAAAAAAAAKFVFejjBc889pw0bNmjs2LFaunSplixZog4dOmjz5s32Ok2aNNG8efOUkpKisWPHql69err22mvVvn17XXvttRo0aFCJdwIAAAAAAAAAAJS8IiURgoOD5ebmJkmqVauWZs2apXfffVf/+te/dO+99+rw4cOSpFtvvVW7du3So48+qrp16yohIUHJycl64403NGHChBLvBAAAAAAAAAAAKHnuRd0hNDTU4eeIiAitX79eX331lbp06aJBgwbphRdeUGhoqD7++OMSCxQAAAAAAAAAAJStIt2JIEmnTp3SyJEj9dxzz+mbb77RmTNnJEmDBw/Wjh07JEnh4eGaOXOmDMMo2WgBAAAAAAAAAECZsRhFvNLfpEkTbdy4UVarVX/88Yd++eUXGYahMWPGyM/PT5J08OBBPf300zp48KDeeecddevWrVSCdwXJyckKCgpSUlKSAgMDyzscAAAAAAAAoEQsaNy4vEMASkW/uLjyDsElFPbadpGXM3r//fcVEhIiSapSpYratGmjs2fP6plnntFLL72kGjVqqFGjRvrxxx+1YsUKPfnkk2rQoIGmTJmi6667rvg9AgAAAAAAAAAAZarIyxkdP37caVtwcLDGjx+vf//73w7bu3fvrpiYGPXq1UudOnXS6NGjdfr06eJHCwAAAAAAAAAAykyRkwi//vqrJk6cqMzMTIftNWvW1G+//eZ8AKtVo0aN0r59+5Sdna2mTZsWP1oAAAAAAAAAAFBmivxMhCNHjqhNmzZyc3PTfffdp7Zt2yowMFDz5s3TyZMn9d///rfA/Xfv3q3w8PArCtqV8EwEAAAAAAAAXI14JgKuVjwTwabUnolQp04drVixQgMHDtTbb78ti8UiSapVq5bWrVt32f2vpgQCAAAAAAAAAABXsyInESSpRYsW2r17t5YuXardu3erZs2a+tvf/qagoKCSjg8AAAAAAAAAAJSTYiURJMnNzU19+vRRnz59SjIeAAAAAAAAAADgIor8YGUAAAAAAAAAAFA5kEQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkES6SnJys5cuXl3cYAAAAAAAAAAC4hAqVRIiLi9OQIUN0ww036JZbblGbNm00ffr0IrcTGxsrq9Uqi8Xi8BUUFCSLxVIKkQMAAAAAAAAAUPG4l3cAhRUTE6POnTurf//+2rx5szw8PLRhwwb16tVLW7Zs0WeffVbotqZMmSLDMJy2N2vWTN27dy/JsAEAAAAAAAAAqLAqxJ0ISUlJ6t+/v9zd3fXRRx/Jw8NDktS+fXuNHTtW06dP1xdffFGothITEzV37lwdPHhQZ86ccfjauXMndyIAAAAAAAAAAPD/KkQS4ZNPPlF8fLwGDhyogIAAh7IHH3xQkvTiiy8qKyvrsm1NnTpVgwcP1jXXXKPg4GCHL3f3CnNjBgAAAAAAAAAApa5CJBFmzJghSerSpYtTWWhoqBo1aqSEhARFRUUV2E5aWpo+/vhjrV+/Xs8884wWL16szMzM0ggZAAAAAAAAAIAKz+WTCCdOnNCBAwckSc2bNzetEx4eLklatmxZgW3NnDlTp0+f1vbt2/XWW2+pT58+ql27tl5//XXl5OSUbOAAAAAAAAAAAFRwLp9E2LNnj/37unXrmtapU6eOJGn37t0FthUZGanp06fr+eefV7t27SRJZ86c0fPPP68BAwaQSAAAAAAAAAAA4CIu/xCA06dP27+/9HkIeQIDAyVJx48fL7CtiIgIRURE2H+OiYnRqFGjtGHDBi1YsEDjx4/XpEmTCmwjMzPTYQmk5ORkSVJKSorDQ5nd3d3l4+Oj3NxcpaWlObWT15f09HSn5IW3t7c8PDx0/vx5p+WW3Nzc5OvrK8MwlJqa6tSun5+frFarMjIylJ2d7VDm5eUlT09PZWVl6dy5cw5lVqtVfn5+9r5cytfXV25ubjp37pzTsyc8PT3l5eWl7OxsZWRkOJRZLBb5+/tLklJTU2UYhkO5j4+P3N3dlZmZqfPnzzuUeXh4yNvbWzk5OUpPT3eKKW8M09LSlJub61BW0Bhe7tz4+/vLYrGYnpuCxjDv3EjmY1jQuSloDC8+N2ZjmHduijqGF58bszHMOzclPYZ556agMbyS1/eVjGFBr+8rGUOzc8McYcMccQFzhA1zhA1zhA1zxAXMETbMETbMETbMERcwR9gwR9gwR9gwR1xwJXMEcLVKTU1ljjD5OT8un0S4+ER4enqa1vHy8pIk08EtSOvWrbV27Vo98MADmj17tj744AONHTtWwcHB+e4zefJkTZw40Wn7jh077L/0JKlGjRpq1qyZMjMztW3bNqf6nTt3liTFxsY6naymTZuqZs2aDks55alSpYpatWqlnJwc03YjIyPl6empuLg4nTp1yqGsUaNGCg0N1dmzZ7V3716HMn9/f914442SpO3btzu9cW666Sb5+fnp0KFDOnbsmENZ3nMpUlNTtWPHDocyT09PRUZGSpJ27drl9Ka7/vrrFRwcrISEBMXHxzuU1apVS02aNNG5c+ec+mqxWNSpUydJtjG89A+Y5s2bq3r16kpMTNTBgwcdykJCQhQeHq7s7GzTMWzfvr3c3d114MABnTlzxqEsLCxMdevW1enTpxUbG+tQFhAQoDZt2kiSabs333yzfHx89OeffzolvBo0aKCGDRsqOTlZu3btcijz9vZW27ZtJUk7d+50+qOpdevWCgoKUnx8vBISEhzK6tSpo2uvvVbp6elOMbm5ualDhw6SpL179zpNei1atFC1atV07Ngx/fHHHw5l1apVU4sWLZSVlWXa144dO8pisWj//v1KSkpyKLvuuutUu3ZtnTx5Uvv373coCwoKUuvWrWUYhmm77dq1k5eXlw4ePKiTJ086lF1zzTWqX7++zp4963AHk2T7ozMvgRgTE+M00bZp00YBAQGKj4/XkSNHHMrq1q2rsLAwpaamKiYmxqHMw8NDt9xyiyTbnVCXTtItW7ZU1apVdfToUR06dMihjDnChjniAuYIG+YIG+YIG+aIC5gjbJgjbJgjbJgjLmCOsGGOsGGOsGGOuOBK5gjgarVr1y7mCBX+errFuHR2djFz587VwIEDJdnuAjBLJDz//PN6/fXX1bJlS+3cubPIx8jIyFCTJk0UHx+vNWvW2E+0GbM7EUJDQ3X48GH7HRESmf88ZP4v4NNBNnw6yIY5woY54gLmCBvmCBvmCBvmiAuYI2yYI2yYI2yYIy5gjrBhjrBhjrBhjrjgSuaIla1aObUHXA267tjBHCHbte169eopKSnJ4dr2pVw+ibBhwwb7JxhOnTqlqlWrOtV5/PHH9eGHH6pHjx5avnx5sY7z9NNPa8qUKfr22291zz33FHq/5ORkBQUFXXagAQAAAAAAgIpkQePG5R0CUCr6xcWVdwguobDXtl1+cbNmzZrZnzVw9OhR0zp521u0aFHs44SFhUkSiQAAAAAAAAAAAP6fyycRqlatqlb/f+vUpevm5clbk7Br167FPk7erSdNmzYtdhsAAAAAAAAAAFxNXD6JIMm+vNC6deucyhITE7V//35VrVpVPXv2LPYxli9froiICDVq1KjYbQAAAAAAAAAAcDWpEEmERx99VNWrV9d///tfp4dVfP7558rNzdWYMWPk4+MjSYqKilJ4eLimTJlir3fkyBG9+uqrWrhwoVP7UVFRWrVqlaZNm1a6HQEAAAAAAAAAoAKpEEmEqlWras6cOUpJSdGoUaOUlZUlSdqyZYsmT56sPn366Nlnn7XXf++997Rnzx5NmDDBvu2HH37QuHHj1K9fP/Xs2VP79u1Tbm6u5s6dq4cfflgLFixQREREWXcNAAAAAAAAAACX5V7eARRWr169FB0drUmTJqldu3by9fVVWlqaXn75ZY0aNUpubm72uoMGDdLatWs1ZMgQ+7YHHnhAsbGxWrBggVavXq22bduqadOm6t27tzZv3qwaNWqUR7cAAAAAAAAAAHBZFsMwjPIOoiJLTk5WUFCQkpKSFBgYWN7hAAAAAAAAACViQePG5R0CUCr6xcWVdwguobDXtivEckYAAAAAAAAAAKDskUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgyr28AwAAlI4FjRuXdwhAqekXF1feIQAAAAAAUClwJwIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgiiQCAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGCKJAIAAAAAAAAAADBFEgEAAAAAAAAAAJgiiQAAAAAAAAAAAEyRRAAAAAAAAAAAAKZIIgAAAAAAAAAAAFMkEQAAAAAAAAAAgCmSCAAAAAAAAAAAwBRJBAAAAAAAAAAAYIokAgAAAAAAAAAAMEUSAQAAAAAAAAAAmCKJAAAAAAAAAAAATJFEAAAAAAAAAAAApkgiAAAAAAAAAAAAUyQRAAAAAAAAAACAKZIIAAAAAAAAAADAFEkEAAAAAAAAAABgqkIlEeLi4jRkyBDdcMMNuuWWW9SmTRtNnz69yO2sWrVKPXv2VNu2bdWmTRt1795dq1evLvmAAQAAAAAAAACowCpMEiEmJkZt2rSRYRjavHmzNm7cqKlTp+qpp57So48+Wuh2Zs6cqe7du2vAgAHatGmTtm3bpv79+6tHjx6aPXt2KfYAAAAAAAAAAICKxWIYhlHeQVxOUlKSWrZsqbS0NP35558KCAiwl02cOFETJkzQ559/rmHDhhXYzqZNm9SxY0d17txZy5cvdyjr2rWroqOjtX37djVt2rTQsSUnJysoKEhJSUkKDAwsUr8AoDQtaNy4vEMASk2/uLjyDgHAVYDflbhaVdTfk7wncTXjfQm4lor6nixphb22XSHuRPjkk08UHx+vgQMHOiQQJOnBBx+UJL344ovKysoqsJ0JEyYoKyvLvs/FHnroIZ07d04vvfRSyQUOAAAAAAAAAEAFViGSCDNmzJAkdenSxaksNDRUjRo1UkJCgqKiovJtIyEhQUuWLMm3nc6dO0uS5s2bp6SkpCsPGgAAAAAAAACACs7lkwgnTpzQgQMHJEnNmzc3rRMeHi5JWrZsWb7tbNy4UZIUHBys2rVrO5WHhoYqMDBQ58+f5yHLAAAAAAAAAACoAiQR9uzZY/++bt26pnXq1KkjSdq9e/dl28mvjcK2AwAAAAAAAABAZeFe3gFczunTp+3fX/o8hDx5D304fvz4ZdvJr43CtpOZmanMzEz7z3lLHx05ckQpKSn27e7u7vLx8VFubq7S0tKc2smLIz09XTk5OQ5l3t7e8vDw0Pnz5x2OJUlubm7y9fWVYRhKTU11atfPz09Wq1UZGRnKzs52KPPy8pKnp6eysrJ07tw5hzKr1So/Pz9JcuhHHl9fX7m5uencuXNOz57w9PSUl5eXsrOzlZGR4VBmsVjk7+8vSUpNTdWlz/H28fGRu7u7MjMzdf78eYcyDw8PeXt7KycnR+np6U4x5Y1hWlqacnNzHcoKGsPLnRt/f39ZLBbTc1PQGOadG8l8DAs6NwWN4cXnxmwM885NUcfw4nNjNoZ556akxzDv3BQ0hlfy+r6SMSzo9X0lY2h2bspijki/JB7gapKQkOC0rSLMEes7dChehwEX12H9eoefK8rfEfyuxNXq4t+TFel/Dd6TuJolJCRUyOsRvC9xtTpy5AjXLC/6+dL+X8rlkwgXnwhPT0/TOl5eXpJkOriXtpNfG4VtZ/LkyZo4caLT9mbNmuW7DwAAKGH16pV3BAAuxnsScC28JwHXw/sScC0FrFZTGaWkpCgoKCjfcpdPIvj4+Ni/z8rKMk0C5H0qMC+jUlA7l36CsKjtPP/88/rXv/5l/zk3N1enT59WSEiILBZLvvsBJSk5OVmhoaGKj4+330EDoPzwngRcC+9JwPXwvgRcC+9JwLXwnkR5MQxDKSkp9mX+8+PySYRatWrZv09NTVXVqlWd6uTdrlyzZs3LtmN2a3NR2vHy8rLfsZAnODg43/pAaQoMDOSXC+BCeE8CroX3JOB6eF8CroX3JOBaeE+iPBR0B0Iel3+wcrNmzeyf8D969KhpnbztLVq0yLed5s2bF9hGYdsBAAAAAAAAAKCycPkkQtWqVdWqVStJ0t69e03r7NmzR5LUtWvXfNvp1KmTrFarTp48qRMnTjiVJyYm6tSpU3Jzc1OnTp1KIHIAAAAAAAAAACo2l08iSNI999wjSVq3bp1TWWJiovbv36+qVauqZ8+e+bYREhKiW2+9Nd921q5dK0nq2bOn6ZJJgCvx8vLSSy+95LS0FoDywXsScC28JwHXw/sScC28JwHXwnsSrs5iGIZR3kFczunTp9W0aVNZLBb98ccfDg8+fv311/X8889r0qRJeuGFFyRJUVFRevLJJ/XAAw9ozJgx9rpr165Vly5ddPvtt2vBggUOx+jdu7eWL1+uNWvWqEOHDmXTMQAAAAAAAAAAXFiFuBOhatWqmjNnjlJSUjRq1ChlZWVJkrZs2aLJkyerT58+evbZZ+3133vvPe3Zs0cTJkxwaKdTp06aMGGCFi5cqBkzZti3f/rpp1q6dKleeeUVEggAAAAAAAAAAPy/CnEnQp6dO3dq0qRJOnDggHx9fZWWlqahQ4dq1KhRcnNzs9ebM2eOHnvsMQ0ZMkQffvihUzvz58/XO++8o8zMTBmGIW9vbz399NPq27dvWXYHAAAAAAAAAACXVqGSCAAAAAAAAAAAoOxUiOWMAAAAAAAAAABA2SOJALiYxMREzZgxQ2+++abmz5+v7Ozs8g4JQCGsWbNGt99+u95++23FxMSUdzgAAAAAABe0YsUKPfbYY1qyZEl5hwIUGssZAS4kKipKd955p9LS0uzbWrRooUWLFqlevXrlGBlQ+Tz44IP279u2bavhw4dfdp/ff/9dTzzxhJYtW6Zq1aopMTGxNEMEIGnEiBGaM2eOhg0bpmnTppV3OAAAuJwVK1boxx9/VL9+/dS7d+/yDgeo9KpWraqkpCR17NhRq1evLu9wgEIhiQC4iOTkZDVp0sT0ouONN96oX3/9Ve7u7uUQGVA5Wa1W1alTR9OnT9ett94qd3d3rV271rRup06dHH6+55579P333ysnJ6csQgUqtcDAQKWlpenhhx/Wp59+Wt7hALhEbGyszpw5o8jISEnSrl279NVXX+npp59W9erVyzk6oHLggiXgWtq0aaP9+/crOjpa4eHhBdZdunSpevXqVUaRAfljOSPARcyZM0eJiYm644479NVXX2nZsmX67LPP1KpVK23btk3ffPNNeYcIVDrPP/+8brvtNnsC79dff9WwYcPUtWtXDRgwQN98841+++03p/1ee+21sg4VqLTuvvtuhYSEaOrUqZetO2rUqDKICIAkLVq0SE2bNlWLFi00ePBg+/aWLVuqR48eatOmjaKiosoxQqDyaNiwoXx8fPThhx9etu7SpUvLICKgcvvyyy9Vt25d1ahR47J1x4wZUwYRAZfHnQiAi7j99tt1zTXXOP1hd/78eXXq1El169bV3Llzyyk6oPKxWq1avXq1010GCQkJql+/vnbt2qXmzZvnu3/dunWVkJBQ2mEClV52drYee+wx9e7dWwMGDCiwXr169XTs2LEyjA6onFavXq1bb73Vfkdew4YNdfDgQYc6X375pf75z39q69atatasWXmECVQau3bt0sCBA7Vu3brLXrQMDw/X7t27yygyoHKKjY3VsWPHNG3aNA0cOFC1a9d2qpOdna3o6GiNGzeOO9zhEkgiAC6iQYMG2rx5s2rWrOlUtmbNGv3rX//S//73v3KIDKicrFartm3bptatWzuVtWrVSjt37ixw/7Zt22rTpk2lFB2APJMnT9aZM2e0fPly1apVK99/wrZt26Z9+/bxTxhQBrp27ao//vhDw4cPV2RkpMaNG6d169Y51MnNzVVgYKD69+/PHbdAKeOCJeBa+vXrp0WLFhW6Pu9JuAIWWAdcREBAgGkCQZLat2+v9PT0AvdfvXq1unTpUgqRAZWX1Wq+6l+1atUuu6+vr29JhwPAxNmzZzVlyhQZhqEdO3YUWNdisZRRVEDltm/fPm3fvt1+odLDw8OpjtVqVY0aNbRixYqyDg+odMaOHWu/YPnjjz+WczQARo8erV9++UV+fn75/m+Zm5urxMREZWVllXF0gDmSCICLCAoKyrfM3d39shctx4wZw50KAIBK55///KfeffddDRkyRKGhoaZ1cnNzFRMTo4ULF5ZxdEDlFB4ebvpJ54ulpaXp6NGj4sZ4oPRxwRJwLT169FCHDh20fPlyeXl55VvvyJEjuvbaa8swMiB/JBEAF5GcnKzs7Gz7A1wvld8nojMzM7V69WrWrQTKUGEueKSkpJRBJAAaNGigYcOG6bPPPiuwXm5ururXr19GUQGVW3BwsJKTkxUYGCjJ/Pfmq6++qszMTC6OAGWAC5ZA+enWrZv9+/r16+uzzz6Tp6ennnnmmQLfj5JUp04dHqwMl2F+VRJAmduzZ4+8vLzk5uZm+rV+/XrT7b6+vurTp4+ys7PLuwvAVef48eOm2y+3JMrevXu1b9++0ggJwCVWrFihnJwcLV26tMB6VqtV8+bNK5uggEpu+PDhuv322+0fcrn49+aJEyf02GOP6c0335TFYtGQIUPKK0ygUuGCJVA+Vq9eLYvFonfffVdffPGFPD09JUl9+/Yt1P4vv/xyaYYHFBoPVgZcRH53GhSWxWLhYTtACbJarVe8fjrvSaD0Va1aVUlJSerYsaNWr15d3uEA+H/vvfeexo4dqzp16iglJUUtWrTQqVOnFBcXp+zsbBmGoV69emn+/Pmmz0wAUHZyc3P1119/qWHDhuUdCnDV8fLy0rFjx1SlSpXyDgW4IiQRABfh5eWle+65Rw0bNixSQiEjI0Pbt29XVFQUFyyBEkRiD6gY2rRpo/379ys6Olrh4eEF1l26dKl69epVRpEB2Lp1q958801FRUXp7Nmzkmy/X6+//no98sgjevTRR6/49y0AR2FhYTp9+rT95+DgYPXs2VOffPJJvvukpKSob9++mjFjBssZASWsVatW2rlzZ7H337Rpk9q2bVuCEQHFQxIBcBEjR47Uxx9/XOz9GzdurLi4uBKMCKjcrFarOnXqpMjISHl7exd6v9zcXB07dkzfffedzpw5U4oRApCkXbt2aeDAgVq3bp1q1KhRYN3w8HCeIQSUk1OnTun8+fOqWrXqZZdUAVB8b7/9tp555hn5+vpqwoQJGj16dKHu9jlw4ID++c9/atmyZWUQJVB5dO7cWWvWrCn2/m3bttWmTZtKMCKgeHiwMuAiBg0alG9ZRkaGkpOTVa1aNbm5uZnWeeyxx0orNKBSqlWrllauXFnsT0iGhYWVcEQAzHh4eOjTTz/VY489poEDB6p27dpOdbKzsxUdHc2zSoByFBISYv/eMAwdOnSIpVOAUnD69Gn5+PhoyZIl6tChQ6H3CwsLU/fu3bVo0SL16dOnFCMEKpf8ruEUxsGDB7Vnz54SjAYoPu5EAFzU6tWr9emnnyoqKsp+O6rFYlHNmjXVp08fDRs2TO3bty/nKIGr17///W9NmjSp2Pv/+uuvioyMLMGIAJjp16+fFi1aVOj6LDMGlCyWTgFcy2233aaePXvqqaeeKvK+J0+e1Pjx4/XRRx+VQmRA5eTj42P6IZfLOX/+vI4dOybDMPj7FS6BJALgYs6fP6+hQ4fqu+++k2T7pNal8h72+vDDD+v9998v0lIrAABcTaKiotSzZ0/5+fmpWrVqpnVyc3OVmJiorKws/gkDShhLpwCupVWrVtq8eXOx/0e888479eOPP5ZwVEDlxbP2cLVgOSPAhRiGob59+2rFihX25EGdOnXUokULVatWTampqTp27JhiYmKUlZWlGTNm6MSJE/yRBwCotHr06KEOHTpo+fLlBa6zfuTIET7xDJQClk4BXEvVqlWv6ENmXKwESlbNmjXVu3fvIu+Xlpambdu26Y8//iiFqICiI4kAuJC3335bUVFRkqTevXvr5Zdf1k033eRULyUlRf/5z380btw4/fzzz3rjjTf07LPPlnW4AACUqWHDhumLL75w2v7MM89c9kGtderU0V133VVKkQGV1/bt2/Xqq68WKYGQ56GHHtL48eNJIgAlKCMj44r2T0hIKKFIAEhS06ZN9fnnnxdr38zMTIWGhpZwREDxXNk9NQBKTEZGhiZPniyLxaI33nhDixYtMk0gSFJAQIAeeeQR7dixQ40bN9bkyZOVlJRUxhEDAFC21qxZY7q9b9++hdp/9erVJRgNAMl2wXHkyJHF2rdatWo6duxYCUcEVG5nzpxRZmZmsfZNTU1VYmJiCUcEVG5XcnePl5eXmjVrVoLRAMXHnQiAi/juu+909uxZjR07VmPHji3UPjVr1tSCBQvUqlUrzZ07Vw8++GApRwkAQPk5dOiQhg4dqo4dO8rT07PQ+6WmpioqKopPVwKlgKVTANcSGhqqb7/9VsOGDSvyvl9++aUaNGhQ8kEBldjx48evaP/8PkQDlDUerAy4iAcffFBLly7VoUOH5O5etPzeI488ouTkZP33v/8tpegAACh/VqtVFoulWPsahsGD6YBS0LZtW23atKnY+990003aunVrCUYEVG7vvvuuJk+erI0bNyosLKzQ+/32229q3769Hn/8cb300kulGCFQubi5uWnfvn267rrryjsU4IqwnBHgIrZu3ar77ruvyAkESbr77ru1Y8eOUogKAADXYhhGsb4AlA6WTgFcy8MPPyxJ6tixo+bPn1+ofRYsWKDOnTsrPT1dw4cPL83wgErHMAx169ZN77//vn7//ffyDgcoNpYzAlzE8ePH1bZt22Lt26pVK/4BAwBc9dzc3DRnzhzdeOONRVrOKDMzU+vXr+fCCFAKWDoFcC0BAQH6+OOPNWjQIA0YMECtWrXSoEGDFBERodq1a8vb21spKSlKSEjQ1q1bNXfuXO3evVuS9Nprr6lWrVrl3APg6pJ3Z8/Zs2f17bff6rnnnivS37GAq2A5I8BFeHl5afXq1YqMjCzyvjk5OfL29lZWVlYpRAYAgGto3bq1YmJiir1/q1attHPnzpILCABLpwAu6rPPPtPjjz+u7OzsAuvlXRIaPny4Pv7447IIDQBQAbGcEeAisrKy5OfnV6x93dzcZLXydgYAXN2u9E6CkSNHllAkAPKwdArgmh599FGtXbtWERERBS73V7NmTX366ackEAAABeJOBMBFWK1WTZo0Sffdd1+R9929e7f69evHwyIBAABQ5ubOnatBgwbJMIwiL53y7LPPlnP0wNVv06ZNWrp0qfbt26fTp0/L399f9erVU5cuXdS7d2/5+PiUd4gAABdHEgFwEVarVRaL5YraIIkAAACA8sDSKYBr+Ouvv1S/fv3yDgMAcJVh/RPAhRR0m+nlvgAAAIDywtIpgGu45ZZb+P8QAFDiuBMBcBFWq1WPPPKIbrrpJnl6ehZ6v+zsbP3+++/64IMPlJ6eXooRAgAAAJfH0ilA+bFarXrwwQf18ssvq06dOuUdDgDgKkESAXARTZo00W+//Vbs/e+77z59/fXXJRgRAAAAUDCWTgFci9VqW3DCzc1N3bp10/33368777xTvr6+5RwZAKAiYzkjwEWMHDnyivYfOHBgCUUCAAAAFA5LpwCuJSAgQH/99Zf+/PNP/eMf/9AXX3yhOnXqaOjQoYqKiuL9CgAoFu5EAAAAAAAUC0unAK7lscce07Rp0xy2xcfHa86cOZo9e7ZSU1N177336h//+IduuOGGcooSAFDRkEQAAAAAABQLS6cAFcvGjRs1e/Zsfffdd6pWrZruuece3X333WrZsmV5hwYAcGEkEQAAAAAAxRIUFKQ9e/bIYrFo+fLl+uqrr7R161bdcccduv/++9W9e3dZLJbyDhPAJTIzMzVu3Di9/fbbslgsatKkiT2h0KxZs/IODwDgYkgiAAAAAACKhaVTgIolNTVVX3/9taZPn67t27dLkv05CbVr19bgwYP1xhtvlGeIAAAXRBIBAAAAAFAqWDoFKFtRUVHq0aOH0/Zff/1V06dP1/fff6/09HR74sDLy0t33HGHHnjgAd166632JcoAALgYSQQAAAAAQKli6RSgbLRq1Uo7d+6UJJ05c0azZ8/WjBkztHfvXkkX7jq44YYb9OCDD+q+++5TcHBweYULAKggSCIAAAAAAEoFS6cAZcvNzU3PPvus9u/fr0WLFikzM1OS7X1XpUoVDR48WA899JBatWpVzpECACoSkggAAAAAgGJh6RTAtVitVvvDzC++3FOvXj19+umnuu2228orNABABUYSAQAAAABQLCydAriWvMRclSpVdO+99+qBBx5Qbm6uvvzyS3377bfy9/fXfffdpyFDhqhJkyblHC0AoKIgiQAAAAAAKBaWTgFci9Vq1SuvvKIxY8bI29vboSwrK0sLFizQ7NmztXjxYrVu3VpDhgzRPffco5CQkHKKGABQEZBEAAAAAAAUC0unAK6ladOmio2NvWy9kydP6vnnn9esWbPk4eGh2267TUOGDFHfvn3l4eFRBpECACoSkggAAAAAgGJh6RTAtWzbtk1t2rTJtzwlJUXffvutZs2apS1btki6kAC0WCzq27evfv755zKJFQBQcZBEAAAAAAAUC0unAK7l1KlTpu+vtWvXaubMmZo7d64yMjIkXUgetGjRQoMHD9a9996r0NDQMo0XAFAxkEQAAAAAABQLS6cAruXmm2/W5s2bJUlHjx7VF198oc8//1xxcXGSLiQOatSooXvvvVdDhgxR69atyytcAEAFQRIBAAAAAFAsLJ0CuBZPT0898cQT2rVrl1auXKnc3Fz7e87b21v9+vXT0KFD1atXL7m5uZVztACAioIkAgAAAACgWFg6BXAtlz7s3GKxqGPHjrr//vt19913KzAwsJwjBABURCQRAAAAAADFwtIpgGvJe9h548aN9dBDD+m+++4jWQcAuGIkEQAAAAAAxcLSKYBrsVqteuONNzRmzBh7QgEAgCtFEgEAAAAAUCwsnQK4lmbNmmnfvn3lHQYA4CpDEgEAAAAAUCwsnQJUXEePHtXhw4fVsGFDVa9evbzDAQC4MPfyDgAAAAAAUHGxdArgOrZv367ExEQdO3ZMiYmJSkxM1MmTJzV79mx7nQMHDmjkyJFauXKlJMlisei2227T9OnTVatWrfIKHQDgwrgTAQAAAABQLCydAriW7t27a/Xq1ZKk9u3ba9SoUerVq5eCgoIkSfv371eHDh106tQpXXo56LrrrtPWrVvl7+9f1mEDAFwcHxUBAAAAABRLURIIR48e1ZYtW3TixIlSjAio3AYMGCCLxaIPPvhAa9eu1d///nd7AuH8+fO66667dPLkSUlSzZo1NWXKFC1cuFATJ05UQkKC3n777fIMHwDgorgTAQAAAABQLCydAriWwYMHKzw8XM8995xT2YQJE/Tyyy/LYrGodu3a2rhxo+rXr28vX7ZsmZ566int2bOnLEMGAFQAJBEAAAAAAMXC0imAa7n55pu1YcMGeXh4OGz/448/1KJFC507d05Wq1UrVqxQ586dnfZv3ry59u7dW1bhAgAqCJYzAgAAAAAUC0unAK7F19fXKYEgSY8//rjOnTsni8WiRx991DSBIEkNGjQo7RABABWQe3kHAAAAAAComKKjo/Xqq6/qsccecyp77bXXtGfPHtOlU/r06aO2bdvqqaee0oQJE8o4auDqdfToUeXm5spqvfCZ0S+++EKLFi2SxWJRcHCwJk2alO/+WVlZZREmAKCC4U4EAAAAAECx7N+/X2PGjHHa/scff+jNN9+UZHsGwtdff+2w9rok9ezZ02mJIwBXplGjRhozZoyys7Ml2RIIFyf5XnnlFVWpUsV030WLFik8PLxM4gQAVCwkEQAAAAAAxcLSKYBrGTdunKZNm6bAwECFhITooYce0rlz5yTZEnf//Oc/TfdLSEjQyJEjdcstt5RluACACoIkAgAAAACgWPKWTrkYS6cA5eeWW27Rjz/+qOrVq+vMmTMyDENubm569NFHNW/ePNN9li1bpnbt2ik+Pl6LFy8u24ABABUCz0QAAAAAABRL3tIpb731ltzd3Vk6BXABffv2Vd++fbVv3z6lpKTo2muvzfd9OGrUKO3du1e33367mjRpoiZNmpRxtACAisBisAglAAAAAKAYNm7cqC5dusjd3V0+Pj46e/asJMkwDPXq1SvfTzUnJCTolltu0VtvvaW///3vZRgxAAAAiorljAAAAAAAxcLSKQAAAFc/7kQAAAAAAFyxoiydct1119mXTunTp08ZRwoAAICiIIkAAAAAAAAAAABMsZwRAAAAAAAAAAAwRRIBAAAAAAAAAACYIokAAAAAAAAAAABMkUQAAAAAAAAAAACmSCIAAAAAAAAAAABTJBEAAAAAAAAAAIApkggAAAAAAAAAAMAUSQQAAAAAAAAAAGDq/wBrz4RtDYv/SgAAAABJRU5ErkJggg==",
Jeremy Auclair's avatar
Jeremy Auclair committed
      "text/plain": [
Jeremy Auclair's avatar
Jeremy Auclair committed
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import os\n",
    "import xarray as xr\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import rcParams\n",
    "import matplotlib.dates as mdates\n",
    "\n",
    "# Settings for plots\n",
    "plt.style.use('default')\n",
Loading
Loading full blame...