spectralAnalysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Spectral Ananlysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Import Libraries "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "import pandas as pd\n",
    "import sqlite3\n",
    "import csv\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "import numpy as np\n",
    "#from wordcloud import WordCloud\n",
    "import re\n",
    "import os\n",
    "from sqlalchemy import create_engine # database connection\n",
    "import datetime as dt\n",
    "from nltk.corpus import stopwords\n",
    "from nltk.tokenize import word_tokenize\n",
    "from nltk.stem.snowball import SnowballStemmer\n",
    "from sklearn.feature_extraction.text import CountVectorizer\n",
    "from sklearn.feature_extraction.text import TfidfVectorizer\n",
    "from sklearn.multiclass import OneVsRestClassifier\n",
    "from sklearn.linear_model import SGDClassifier\n",
    "from sklearn import metrics\n",
    "from sklearn.metrics import f1_score,precision_score,recall_score\n",
    "from sklearn import svm\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.naive_bayes import GaussianNB\n",
    "from datetime import datetime\n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "from sklearn.cluster import KMeans\n",
    "import matplotlib.pyplot as plt\n",
    "import pylab as pl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "ename": "ModuleNotFoundError",
     "evalue": "No module named 'spectral'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mModuleNotFoundError\u001b[0m                       Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-21-fb1dae139610>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;32mfrom\u001b[0m \u001b[0msklearn\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcluster\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mSpectralClustering\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0;32mfrom\u001b[0m \u001b[0mspectral\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mSpectralClusteringg\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;31mModuleNotFoundError\u001b[0m: No module named 'spectral'"
     ]
    }
   ],
   "source": [
    "from sklearn.cluster import SpectralClustering\n",
    "from spectral import SpectralClusteringg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Create Dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_csv('dataset/28_1800.csv')\n",
    "table=df.pivot_table(index=[\"Name\"])\n",
    "table=table.reset_index()\n",
    "# features order:\n",
    "# 0.Name, 1.AirTemp., 2.Alt, 3.DewPoint, 4.Lat, 5.Longt, 6.Pres., 7.R.Humidity, 8.Visib., 9.WindDir., 10.WindGust, 11.WindSpeed\n",
    "\n",
    "# change the order of features\n",
    "table2=table.iloc[:,[1,2,3,6,7,8,9,10,11,0,4,5]]\n",
    "#print(table2.tail())\n",
    "# specific features - new table is created\n",
    "# enter the fourth column to address feature e.g: 1. AirTemp, 7. Relative Humidity\n",
    "# e.g: table3=table.iloc[:,[0,4,5,1]]\n",
    "table3=table.iloc[:,[0,4,5,7]]\n",
    "table3=table3.dropna()\n",
    "\n",
    "# for debug print tail of table3\n",
    "#print(table3.tail())\n",
    "#table3\n",
    "title=\"RelativeHumidity\"\n",
    "\n",
    "\n",
    "# k=clustering number\n",
    "k=4\n",
    "cluster=KMeans(n_clusters=k)\n",
    "table3[\"Cluster\"]=cluster.fit_predict(table3[table3.columns[3:]])\n",
    "#table3\n",
    "\n",
    "# centroids = clusters' center points\n",
    "centroids = cluster.cluster_centers_\n",
    "weather_clusters=table3[[\"Name\",\"Lat\",\"Longt\",\"Cluster\",title]]\n",
    "#print(weather_clusters)\n",
    "weather_clusters.to_csv('ClusteredData.csv', index=False)\n",
    "\n",
    "\n",
    "# informative\n",
    "#print(\"centroids\")\n",
    "#print(table3.columns)\n",
    "#print(centroids)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Cluster0 Size: 2728 , Cluster1 Size: 1976 , Cluster2 Size: 3 , Cluster3 Size: 2537\n",
      "max_c0: 35.5  min_c0: 0.0 Color:R\n",
      "max_c1: 148.2  min_c1: 67.5 Color:G\n",
      "max_c2: 655.0  min_c2: 655.0 Color:B\n",
      "max_c3: 67.4  min_c3: 35.666666666666664 Color:Y\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi0AAAGsCAYAAAAR7ZeSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAMTQAADE0B0s6tTgAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeXRc133g+e99S9WrQi0o7OAKkBIpUiYlU4stxZaXxLbGSZy0xThz2k7sJGfGM5lJMslJp5NpO/Ecu/skc7ozmaSnz3Qn6baTdLttS+lJ7MRLEseOpdjWQomLCJISSXADCjtQ69vv/PFQEAiCJECCC8jfRwcS8areq4sCj96v7v39fldprRFCCCGEuN0Zt3oAQgghhBArIUGLEEIIIdYFCVqEEEIIsS5I0CKEEEKIdUGCFiGEEEKsCxK0CCGEEGJdkKBFCCGEEOuCBC1CiIsopT6mlDp/ndd4p1JKK6WstRrXraSU+pZS6jNXePyHlFJ60fefUko9e3NGJ8TdQ4IWIe4w8zdYPf9VV0q9opT6iRv4ep9RSn1ryeF/BPq11uEavcaygdRNDI4+CPz2Kp7/r4EPtL5RSj2rlPrUWg9KiLuNBC1C3Jl+D+gH3gT8F+DzSqkHbtaLa619rXX5Zr3ejaa1ntZa11bx/JrWevpGjkmIu5EELULcmepa67LW+rTW+v8E5oB3th5USv2iUuqUUqqhlHpBKfXOy11IKfUBpdT3lFJVpdSIUurfKaXa5h/7GPAvgHcsmt0ZWDwDopR6s1IqUkr1L7nuV5RSf3AtY7rCWC9ZllFKfVYp9WeLvh9WSv0zpdSfz7/WUaXUI0qpPUqp7yulakqpv1JKdSw656LlofnnvqiUcudfb+By41BKfRb4AeC35t+TYaXUW+fPLS057zsyIyPE5UnQIsQdTCllKKU+CJQAf/7YzwK/BPw8yUzMnwB/rZQauMxlHOBfAg8A/z3wLuC35h/7AsmszndJZnb6gXOLT9Zavwy8DuxfNK524D3z51/LmK7XrwJ/DjwIHAf+FPi/gV8D3gbsIAnGLqGUMufPPQs8NH/ep6/wWr8EPA/8G5L35xGt9feA08CHFl13kCS4+dPr+LmEuKNJ0CLEnenXlFI1wAOeAc4DX5p/7BPAL2mtv6a1PqW1/gPgWeAjy11Ia/1FrfWX55/7D8D/wXwAorVuAnXAn5/ZKWuto2Uu8wXgJxd9/+PAOPDcKsa0YX4WZOEL+Ooq35eWZ7TWf6a1PgH8DrAT+Hda629rrV8B/phFM1NLvBfYCPyc1vpVrfWXgM9d7oW01nNAANTm35+J+Yc+B/z0oqf+FPCPWuuT1/gzCXHHuyMy+4UQl/hD4P8CNpB8wv9NrfWkUioHDAJfWFztAqRJAptLKKV2k8y0PEQyY2Ox+v93/FfgE0qpzVrrcyQzDF/UWutVjGkMePuS674F+DNW7/CS6wK8uuRY92XO3Qm8rrWeWXTs+WsYw58An1FKbZ8PVD5CksArhLgMCVqEuDPNaK1fB15XSn0YeE4ptQdoBQX/lItv0gDVy1zrL4FDwIdJZkeeAP7DagajtT6qlHoV+Aml1H8Cfog3lpjaVjimaP5nWqCU2rTk+TGglhyzeePnbgkWD+8yxy43E62Wud6qaa1HlFJ/C/yUUurrwGbgi9d7XSHuZBK0CHGH01qfmC9J/gTwC0AZ2KK1/ournauU6gK2A/vnl01QSn1oydMCwFzBUP4ryRLRHHBOa/3C/PHx1YzpKiaAviXH9pAEXWvlOHCvUqpdaz07f+yRq5xzuffoP5HMYvUCf7noekKIZUhOixB3h38L/BzJDf1fAZ9WSv2MUmq7UuphpdSvK6Xevcx5M/Nf/4NSaptS6ieBjy95zhlgp1LqPqVUl1Lqcv9f+QLJzf1/m/8zAFprvcoxXcl3gG1Kqf9ZKXWvUupfsaSyZw18HRgF/kgptVsptR/46FXOOQO8VSm1cUnF0P8HdJL8bv5kjccpxB1HghYh7gJa628DJ4BfnU9y/bX5ryHgy8CjwIVlzotIloXeS7J083HgN5c87WmSnI4XSGY6tlxmDK8DB0iqg76w5LEVj+kqP+cR4JeBT86PxwD+22qusYLXiEiazQ0CLwO/AnzqKqf9a5Lg5NT8Oa1reSTvxQxJMCSEuAKVfMgRQghxKyil/pIksfdXbvVYhLjdSU6LEELcAkqpIklS85MkeTdCiKuQoEUIIW6NvwAeBj6ltT5+qwcjxHogy0NCCCGEWBckEVcIIYQQ64IELUIIIYRYF+7YnJZ0Oq27uy/XhVsIIYQQt6MLFy74Wuv0co/dsUFLd3c3588vu5WKEEIIIW5TSqmJyz0my0NCCCGEWBckaBFCCCHEuiBBixBCCCHWBQlahBBCCLEuSNAihBBCiHVBghYhhBBCrAsStAghhBBiXZCgRQghhBDrggQtQgghhFgXJGgRQgghxLogQYsQQggh1gUJWoQQQgixLkjQIoQQQoh14bYIWpRS71NKvaSUelkpdUQp9dH54z1Kqa8ppV6bP/62Wz1WcXcq18qcmDpBuVa+1UMRQoi7lnWrB6CUUsB/Ad6ltT6klBoAjiml/hz4beB7WusnlVKPAE8rpbZrrcNbN2JxN6n5NZ4++jRDk0MYGMTE7Oraxf7d+8mlcrd6eEIIcVe5LWZa5rXP/7cATAEe8CHg/wHQWr8AjAEy2yJumlbAMtg+yGBpkIHiAEOTQzwz9MytHpoQQtx1bvlMi9ZaK6U+BPy5UqoOlIAPAnnA0FpPLHr6MLDl5o9S3I3KtTLHJo8x2D6IoZL43jRMBooDHJ04SrlWpi/Xd4tHKYQQd49bPtOilLKA3wB+TGu9FfhB4HPzD+ulT7/CdX5FKXW+9VWr1W7MgMVdo+JVUKiFgKXFNEwMDCpe5RaNTAgh7k63PGgBHgQ2aK2fg4VloBFgL4BSqnvRc7cCZ5e7iNb6d7XWm1pfuZzkG4jrU0gXiImJdXzR8SiOiIkppAu3aGRCCHF3uh2ClnPAJqXUTgCl1D3AduAE8CXgf5k//gjQBzx7i8Yp7jJ9uT52de3i9OxpojgCkoBleG6Y3d27ZWlICCFustshp2VMKfVxksqgmGQJ6Oe11heUUv8c+FOl1GuAD/yUVA6Jm2n/7v2XVA/t7t7NU7ueutVDE0KIu47SemnayJ1h06ZN+vz587d6GOIOUa6VqXgVCumCzLAIIcQNpJS6oLXetNxjt3ymRYj1oC/XJ8GKEELcYrdDTosQQgghxFVJ0CKEEEKIdUGCFiGEEEKsCxK0CCGEEGJdkKBFCCGEEOuCVA+Ju86hsUOMVkfpz/ezt3fvrR6OEEKIFZKgRdw1yrUyn/jmJzg8fnihUdyenj185t2fkXJmIYRYB2R5SKxb5VqZE1MnKNfKK3p+K2DZVtrG9o7tbGvfxuHxw3zy7z95g0cqhBBiLchMi1h3an7tktb6u7p2sX/3fnKp5TfKPDR2iCPjR9hW2oapTCDZrXlb+zYOjR3i0NghWSoSQojbnMy0iHWnFbAMtg8yWBpkoDjA0OQQzww9c9lzRqujKNRCwNJiGiYGBqPV0Rs9bCGEENdJghaxrpRrZY5NHmOwfRBDJX99TcNkoDjA0Ymjl10q6s/3ExMT6eii41EcERPTn+9fszF6XplG4wSet7JlKyGEECsjy0NiXal4FRRqIWBpac2YVLzKskm1e3v3sqdnT5LT0r4N0zCJ4ohTs6fY27t3TZaGwrDGxMTTNBpDJJ8HYrLZXXR378eyll+2EkIIsXIy0yLWlUK6QExMrOOLjrdmTArpwmXP/cy7P8Oenj2cmj3FyemTCwHLp9/16TUZWytgcZxBMplBHGeARmOIycnLL1sJIYRYOZlpEetKX66PXV27GJocYqA4sDBjMjw3zO7u3VcsXe7L9fFHH/ijG9KnJVkSOobjDKLmZ4GUMnGcAer1o3hemXRayqqFEOJ6SNCyDpRrZSpehUK6IP1EgP27919SPbS7ezdP7XpqReev1XLQYlFUAdRCwNKilIlSxvzj8rsTQojrIUHLbexaSnvvBrlUjo89+LHbKpgzzQIQo3V8UeCidYTW8fzjQgghrocELbeh1s34q699lZHayEKlTBRHC6W9H33go7d6mLdcX67vlgcrLel0H9nsrvmclgGUMtE6wnWHaWvbLUtDQgixBiRouY0snllpBA1euPACOzt3sjG/EcdyLintbd2wWzkatmmTNtMAbO/Yftvc0O8W3d37F5JxlTLQOqatbTddXVdftrqdZo2EEOJ2JUHLbWRx07SJ+gTFdJGx2hgvj77MY5sfAy4u7YWkNf3B8kEmG5PMuDPYps3mwmb6cn18YOcH+MjejywsJcmN8cayrBz9/R/D88pEUQXTLFx1hkWWAIUQYuUkaLnFWje4ac+9qGmaYzugoJQpMVIbYc6do+gULyrtXdj8Txn4kY9jOQRRQMWt0Jnt5MsnvkzGzvDUrqeu68a4mpvwulIuQ6UChQL0rd3PlbxHK7ve4kBVlgCFEOLKJGi5RZY2IptuTrJRjWCyAY1DMV2kP9fPaG0UNLihSy7OLZT2jtfHOTJ+hA25DRydPIpG41gOjuUw681iGzZ+5PPChReYbc5eU27MemiWdk0BVa0GTz8NQ0NgGBDHsGsX7N8PuZv3c62ku6/MiAkhxBskaLnJWks0cfWr2NHIQl+PnCqS5RCZ4AC11OMA7Ovfx4sjL3J86jjj9XFqQW2htPe5s8+hUMQ6Jo6TRmtKqYXX8UKPlJlizp3j1YlX2de/b9U3xlbAYlkdaO2jlL3QLK2v79bOAlxXQNUKWAYHk6AlipLvn3kGPrq6n+t6ZqGutbuvEELcrSRouUkW5y44NNhiPE8xu5OHNmwkbTkUnQ6KbTuZawyBcR9YHdiGTX++n4f7H+bJe5+8KBeltZeOoQwMI7npaa0XXi9tpQniAK01aTO96huj55Wp1Q4RBJPUaq/QCgxsu5c4DiiVbm2ztMXdZ5Ok12hFAVX51CEqx/+Rwpbt9M2/b5gmDAzA0aPJktEKlorWYhZqcXffxb+flXT3FUKIu5EELTfJ00ef5oWRF+jJ9lCyNEXaKdfHOFB+mcc2JUm2D214mIMXpjlYPclkcAE3dNnTvYef3fezl+SetPbSOT31Evfk0pyO5pj2XYI4oOSUCOKAlJni8Q33MVE/g4pm0GZp4fyr3RijqEKzeRytQ8BEaw+l0gTBGGE4fUubpV1L99mFoPH4cxjmEeLmGXYF/ex39pEznCRwMYwkx2UFQcu1Bk2LXU93XyGEuBtJ0HITnJw+yZde/RJBHHBy+iRZw+XtxTIpZ4DTM6fZXNiMbdo4ps3Orh1coMT5ybM4psOF2gWePvr0JUmzYVjjV/fu4W+PP89YbZzdGZfXqx4vzuUoOAXaLIOfHOzmwc5pjjPCbO1FMumt1O2H8VTnVW+MUeTieaOE4RxRNLdw3DQLWFY7UeSu6Xu0msqma+k+u5Dw2nUPRnyWSHUwFJZ5xn2Zj2YfS5aI4jhJyr2KtWzZf73dfYUQ4m4iQctN8KWjX2KqOcW20jaiOOLM7CwvBZNsyU4w0tQMzw6ztbiZDtun7DlMmvt4qP+hKybNTkw8jR2N8E/e9PNcqI4w15zhv9MTuPQxZ7+F7anX6bRqpFL9bC+McS4YxvO+j/IOU9M7ub/rx/jgrqcum5Nhmg5BME4UVTGMdgzDII5jwnAKrQNM01mT9+ZaSn5X2332koTX/n7M0VEGOjs4Go5Q9qfpOzcDu3evaJZlLVv2X6277x1buSWEENdAgpYbrFwrc6FygUK6gFKKs3NnqfpVXo83Ug9OsNGJSOmA2WqdudQAXx+dYlNpmHwqjxu4oKAj03FR0uzST/qbCpvoSKU4VD5Hxf02F+IKDeNVitmd7MyfJw4nGeh6kGbQwPMn2JHtpacE1anL52S47nm0DjHNEuATx6AUmGYJrV1c9zzZ7I7rfn+upeR3td1nL0l4fWgfvHQAc3QUw6hQCU/St/tt8NTKZjduRMv+pd1910PllhBC3GwStNxgFa9CLpVjQ34Dw7PDzHlz5FI5RmplXq7H3JPvZGOunfFmnelgjKpX5bXZ8/zdyb/DUhbZdJZ2p51NhU2MVkfpy/Vd9Enf98epVl/i/OwQoe/Ra4UUU8eJdYrx+lmy4RQDHXsBRcZuw8KlkN3K2NgXSaf7aWvbg22XLsnJSGZT2jAMB6VKCzdorX3i2ELr4Lrfm+sp+V1N99lLEl7TDjz+ONHsNPHkSQqP/AJsW/kGijejZf9a5MwIIcSdRoKWG6x1w3yw70HmvDlOzZyiGTSZdWexlEUlspiZnWPOmyOOY7zIQ83/YxkWISFRHFH363znzHd4c/+bMc0CcewxO/sstdoB/LBCGNQomgboGCccxsAjZfXg+nPU3DGcyECFGm24VNXLNBpDaO0RBBOkUv3k8/suyslIpfqxrA7AII5rgELrEMPIYxhZUqn+635vrqfkdzXdZy+b8Kpn2H3f2+hbRcDScj0t+69mLXNmhBDiTiJByw3WumEeHDuINf92h1FIrGM0miAOMDAIogBTmej5fxzTQaPxwiSIieKI5y88z6GxQ/S09TBdP0fsHcMgBqMdkxlMfCJyhKoTR58mo4dxdEBzeho/ABUbmKGJZXagSjapVC9KpfH9MtXqyxSLjy3kZOTze8nlHqBWO4zj3ItSGq0Vvj9KPv8g+fzqb/RLrUXJ70q7z651wuu1tOxfqbXMmRFCiDuJBC03wf7d+/nGyW9wfPo4WStLPaiTs3NUvAqz7iyxjgnjkIA3llzCOAQg0hE60DSDJt84+Q2G54bpz6R5X/csOeXSaU6TNk1SNAELixqmfg3QQIxBiCYksgxUCLFloOujWEEWI2+CZWBZHfj+CEEwfVFOxuDgZzh9+hPU64cXZhPy+b0MDHx6Td6Xq5X8jtfHOVg+SH++n7291xckXS3h9VqtpmX/St2InBkhhLgTqMUNye4kmzZt0ufPn7/Vw6Dm1/jjA3/M5498Hsdy8EOfUIeUa2VGKiNERCgUmjd+D4uXTGIdYypz4Xh7pp13bNjOE501evL3kvW+hWMAhMRxiImPAjQGGo0Zh9iRhQpAm4pYRdiuRf5Umviejdj37AMUvj+B42yhvf1tl+RMVKuH8P3R+WWk5YOHpTMOK52BWK56aEN+AwdGDnBs6tjCsT09e/jMuz9z1/QuGR397GVzZiSnRQhxJ1NKXdBab1ruMZlpucGePvo0R8aP0J5up6utCy/0GJoYou7Xk6BCmSitiIkXAheNJtbxJX82lUnFrXBi5jxvzhk46QpZZdEMfQq2RRi7aB2hDDBUBBoMDcQa280Sx020jkH7OGc1YfMCfkcJ8nmiaJZc7v3L5mQkgcrywcrSKpc49uaDlTyG4XC1qpflZkA+8c1PcGzqGNtK2zBVMvtyePwwn/z7T/KHP/qHa/Frue3dyJwZIYRYryRouYFa1TH3dN7D2cpZNJrzlfM0wgZ+5JOxMgRxAAqIkwAlIgK4aOYFwDZsbNPGUAanK3McmNQ8rA9ipiMMHeLEPikzSKpjWxRgQKwiIl3F8BVGbBLmY8CgeDAkcMdpfvh+Ojp+mE2bfnHVP+PSKpfZ2WdpNI7T1raLYvHxFVe9tEp+D40d4sj4kYWABZLE3G3t2zg0dohDY4eue6loPbiROTNCCLFeGVd/irhWreqYklOiP9fP2dmznJk7Q9WrEsQBbuQS65ggDhaCFQuLlJHCxLyoisg2bbzQoxE0cEOXb5RdTtZiTHyCSBHGIbFaZhAKtAVhG4QZTewkib71LR6NHo9w+gztnztAv/nDq/75lla5hOEcQTBOJrMN30+66S6terma0eooCrUQsLS0KopGq6OrHuelYz6xorFcTrlW5sTUCcq1K19jLV4rKa/eIQGLEEIgMy031OLqmH39+zg8dph6UCdtponiCMuwMDAgWf2hx7HIGhrTKnCmWsWNk1b5URwRqxit35iJwcjxzal26v4FnugKQSUzM8kcyhIGxCaYLpheTOEU2A0Tu2HT922HtJqDk78Bv//7K+oI27K0yiWO3fnvTZKlIhcorqrqpbURZKSjiwKXVkVRf/7aSq3XolnbSrv3SmM4IYS4MSRoWWOt3Aw3dHEshw25DZyePU3JKVFIF+jJ9hBEAa7tkrWzuKGLEXi8owe2ZAM0mpQ5x0g+xTfKFvUwyXVpRI2LEnbrfp26X6cYx/xoH1iLfpPLBi6AEUPuNSgcT7rbuv0RNBqQLsKzz8Jv/Aa84x3w9rcne/EUClcMYpZWubRyWLSOgHj++9VVvfS09bClsIXjk8fZ2blzoaLo1Owp9vbuvealobVo1rbS7r3X8lqyDCSEEFcnQcsaaX0KPzh2kBOTJ5hsTtKV7WKwOIgbuYzVxnBDl7SVRilFX1sfkY6I45gf6FFsz1u4tFP3m9iY3JNXRDri7ydS1IM6rRSXnJ2jGlRJGTHv6YGf3AzdztXX+cwaOBNJ4GKECiyFCiOirAnaAMeBVAq++EX40pdg795kA8Fdu2D/fshdOkOwtDOsZRWx7d75nJb7sKziFTvFLk6+zaVyC0FBX66P41PH+f6F79OV6UIZir29e/n0u66t1HotmrWttHvval9LZmWEEGLlJGhZI60b7nRjmiiO2Na+jWl3mmpQpS/XR97OY5s2m/ObeXXyVb579rukzBQb2hx25NOEZjedqTyN4BxBFHK6FrDBieh2MjQjAx1pUOBHPgDv7FY80aXpSiX5tlfjlCE9BV4HuJ0aFYeEaTDDLLh1SKdhdBTCpD8MHR1JoDI0BM88Ax9dfoZgaZWLbXeRy6UwzTyue3rZqpflllmmGlPk03l2du7EUAa7u3fz/IXn6cx28tMP/PR1Jd9erVnbWPUkbu3KvVtW2r13tY3hpF2/EEKsnAQta6D1KbzD6eBg+SAd2Q4MZdCV6aJcK7Oraxcz7gxv6n4Tw3PDpM00mVSGulenZIUUnSKG3U25XiaIAkId0mZncawYFTcXAhV00myuZMOOnKY9BYYCrZPlnsuKwO9NEnK1A34x+bNdg/FHG/R/3cXq2gG+D11dMDMDrgvFIgwMwNGjUC4vu1R0uSqXKy13LF1mmW5O863hb/Hm7q2kdQexdsAo8ujGRxmeG6anree6fj+Xa9bmBnUOlQ/y0qkKIYUr7jC90u69q2kMJ+36hRBidaR6aA20PoW3ggtj4QakFlr0Gxg8sfUJql6V41PH2VzYTHumHYw2gihgzpulM9NFykxhKYsg9AhinxnPW3gdY/7X1WaBpcBWJNVBy4wpXvJ90AaNQWhuSL68bnBLMPbumNc/HhL2FZkfdLIs5CT5KJhmkvdy8GASuFzG0iqXy1W9LLfMEoU13laqs1W9QLr5Tdr9vyUf/COWChZmMa5HaxkrmfmZLynXEYcu/A1nmrCxfQ+DpUEGigMLOSpLtbr3np49TRQn11jcvbc1Q3O511puiWxlszJCCCFaZKZlDbQ+hafM1EIzOEMZaK2JibFNmziIMQ2TjJXh7VveThiHdGY6Ga+PMxWeIcsoM0EHPW3dBJGPEY9zfM5nKtBJ+bOhsFTy66qHLoEGNwYdg1ZgXmmAJm88QUNsgdWAqAjahNo9inJ4nk0HBmFiAjZuTGZZXBdefBGOH0/O/eY3r5jjcjmHxg4xWh2lP9+PYzmXLLP0G69h2HVmwwxdcQHTaiMVl8n6LxHTt6I9iK5m6TJWzatytqlwCu9ZCBqutsP0SvcvWmljOGnXL4QQqyNByxpYvIdOb7aX8fo4pUyJaXeavlwfM+4M20rb+OvX/prnzj1Hzs4xXh9HGYp8Ks9zdYPt6YjBtiniQLEp38fZ+kYOV6exjRqFdIFG2EDHGq00ldDkeC1iQwb6Y5I2/ldaHmpNxczPzGAnh8wmxDkItWZ20xzFCzXy9VyShAtvBCz33Qd79iQVRVfIcVm6r0+5VuYT3/wEh8cPL9zkt7dvp5QtMePO4Ic+OStgiznHqNlDwx3HNExQBn5cwvWO8abOfWvSun/pMtaZygSjF77CoJG96HlX2mF6pfsXrbQx3NJE5qXt+mVpSAghLiZByxppfQr3I5+p5hSnZk/Rne2m0+lkd/duZpoz/M2pv2GmOUPNrzHVmCKMQ7KpLF2ZLp6dbmMi6maycZZgdBLL7qESBIRxSM2vkbJS+LGPpSxCHfLdSYuMikmpmB35pE2LOR+4tHJcFkqf548vLBmZEORB+ckSU9QBRjpm5JFNFMa66X75LNap8I2A5aGH5s8zl81xWZpYW/ErtKfb+e6571Kul9nZuZO2VBtRHPHa9GtURiu0Z9ppT7fTYdcxirOkrS425jfSCOo0gyaamM3Zbn5g2xOr/l1cKVhobXBYigtU/Sqj1VGyqSzFdLI8tpIdplvde69mJZspSrt+IYRYOQla1sjST+GtPi2tm9/P/cXPYSmLvlwfE40JDGUQ6pCG36BhNbANm7P1GlXfJmWmKKWSJZSYmDAKydpZPO1RC2oAeDrFl8sxz8/Ak32wK58EIJsdcCxImfNBTFJFe0mOi2EkSbmRglQFLKsda+Me5nrGCO/vZPO5R/DSdaL7BzFjj3S0KMfFMKBSWQhaWgFLf66fFy68wIGRA5yrnKMZNUkbaWbdWTYXNy/kscx5c+zs3IkXe3iRRcUtk0pp3rf9R8mkO2mGLo5pk9Iz5J2VN5Nbaflwza/xtde/xqmZU7w48iL5dJ6N+Y080PsA5Xr5ohyVG03a9QshxMpJ0LLGlvsU/tzZ55hqTrGttA3Hcqh6VSaiCZRW+LGPF3nEOl4IZmzDppAqUE1VyaVyTDYnqfm1JFdmfvdmy7DwY5897WAZ8L1p2FuEHiCj3qgmiklijKVRizLmZ2MiiFHong5qtRfQGmrhy1Q6XyJ1/3mMwgUAsn4/3bV9WKGdJOoWkmBscWLt985/j4Plg8z6s8m2A5GHZVo0ggbn5s4RxzGNsEHWzjJYGqQ3U6A6+1f02DFpdYZj5/8zucxWdm38EVQ0s+olkpWWD7eCrPdufy+vlF9hpDLCscljTDQm+ND9H7okR83Al3UAACAASURBVOVmWMmsjBBC3O1um6BFKZUG/g3wPsAHXtZaf0QpdS/wOaALmAU+prU+eutGem2iOOL0zOmkEkaBbdoEUQA6+XM+lWfOmyOMQiq6wompE+TTeR7sfZDDE4cppUtMNidxLIfJ+iTZVJa012BrFspuMtPSYQMKmjGYRlIOvWx7XM3CJo1KgXaymLlubLsLUHjeBerqFEZXjuIZ0J0dVNNnCIw5eo9sIL370WSWpVymcv4gql6nmqlyZvZMUkGlIWNlqPv1hS6+URwx2ZzEUhZaa4rpIlH9mxSMGSb0AAXG2GjahO4pTpz/PHs6foK2ZhEvPER609V7tKy0fHhp9dLjmx9nzp2j7teZaE7wvu3vu6TcWQghxO3htglagN8mucXu0FprpVRrXeDfA/9Ba/1ZpdR+4I+Bx27VIK/F9o7teJHHRGOCnmwPylA4BYdjU8fwoyRPZcadwcQk7+RJW2lSRoogDjhfPY+pTLrauijXk6WnIA6Ybc7SmUryV/IW3JtLloPSBmTN5Higk1mYxWINUZzEMYrkX3ZbH6lUL6CIojpxXKet7UH8DVX8oI1m4zBeVKea8nAf2UV+xw66P/f/Yh09TcFoEBvP0+g6R5j3iIlRSpE205imSRAHSaBiauI4Zi6YY0fnDkopcGsXCIwuFAYX/HZK+QEc5RFOvczU61/Aj/8KdEy2sIfu934Gq3j5mYiVNnVbrklc0SlSdIp4kbdsAq4QQojbw20RtCil2oCfATZprTWA1npUKdUD7APeO//UZ4B/q5Qa0FoP35LBXqP+XD9RHFENqpjKJNLRQmVRLpXDj3yMjIFC0dvWixd5GMpgpDpCe7qd16ZeY86dQ6lk12c3cKmHyUzJjnwSqMyEEPiQciBjgheBQhHHJoowaU4XK4LIABWjtMaxLEwnQxS5GEaaIJjENHPzy1azeAMx6bAbO96GYVSxu99O4+DXmZyFvsEP0GcY7KrN8cLkYQyvgZE1iHWMG7n0ZHtoBA1qfg0dJiVM2zu280D/A/jBDMnoDGpBjWK6SMoukRr5Htr0CDOdZNiEjiMalcNM/s0n6dv/hxe9p4vzQFZaPrzSJnFCCCFuP7dF0AJsB6aATyilfghoAp8iWQ4a0VqHAPMzMGeBLcDwrRnq6lW8Cm/qeRP9+X7OzJ5Z2MG4t62XkzMn2dO7h9enXqc9086FygVmvVmCKKA/10+H00EjaDDtTpMyUygUju1Q9+vMBDDpwUPt4M3nrIQa/Dhpy+JpCAOFj6JgGRhaAxrDiAgjBYaJZWcJgnGq1RqW1Yad3sJM9SRe7QIWSUt/ZfWwoX0v6BRmIyZ1PqS+CbzZKumoyP7swwCcnn6WyI4I4oC0kaYz20k/aQzt4NidvOvep/jlt/4yTx99mkMj36YQuTTCCu1OiS3FLRjNGYx4hqaZJkOyRKMMEye7jfrcIbzzyVLR5RJuHWcQ1z19xfLhxeXpA8WBhQ0ZlzaJE0IIcfu5XYIWG9gGHNVa/7pS6gHgb4Ef4dKGr8t2JFFK/QrwK63vi8XiDRrq6hXSBSzT4q2b3sr9PffjBi6O7ZCzc1T8ClprIh1xdu4sDb+BH/gLS0CO7RDpiDa7jVwqhxu61P36wrUPzcHDHdCbhnYrCV4mPehIzc+2kCLWJp4OsJQiZcREkYOVNtFo0vm3k09ncd0LGIZJee41gsgjbShiVUTpkCCcoTz7PFt73oPVtEGZKCAyXIiK5AyHj+XextsnMvzZ1jTfrB1kvHqGe1On2NoGpXQHO7p6eWhLN45BUmV1z5N89aVJYu8kuWwPGBY6rBHaHpmonxLtuOYskdHEjDPJEk9lFNh72YRbx9m20PfkSuXDC03ihl/ECAJi22b34CO3JAFXCCHEyt0uQcsZknyW/wygtT6olDoNbAU2KaUsrXWolFLAZuDs0gtorX8X+N3W95s2bVquu/0tsfTTfTFdXPh0/+T2Jzk0fogTUydohk3SZhrLsMin83ixR5o0pUwJN3IpOSU2FDYwWZ/k4NhB3NBl2o8ZrsPQXJLXkrMhiJNcFsOwaaQe4MTsDDnbps+4QHumCKle7GiEmbDIYPEtFNIFDOMAc7XjEI6ijU48IiztYdLENlJ4oUtg9ELGQeuQUPn4ZhUzdkhHRYgitsft/NbbfpmP5+DE8B/QbB4j4wzSk+sln8rhuqeZnHyGrq6n0NWv8Za+rZyeOIvvP0+gszTjFDk/w7b4Hsr57+FZswBJtZRt0J21r5hw67qn2LTpl+noePKK5cM5Hz72CpSHNBUTCpGmb1Yn832p6/tdS+myEELcOLdF0KK1nlRK/R1J5dBfK6W2AoPAd4CXgY8AnwWeAobXWz4LXL4FfDNokrbS3Nd1H1W/SsWrMNWcwggM2qw2ppvTKBRhHDI8N0x/vp9cOkfaSoOGehRypuGzNQsvziazK44B9xZKFAs7yZoFup06db/GIXcTvZkfpM+MiL3nyOcfoD2TbEZYLD7OjK+phEeohVV67RqWEaKUxjZMYtI0/UmiYi8z28bR1SpTWQcUZN1eul/txNr9IPT1UfLKbMxEOKX3XJRf0qrkCYJZfH+EYu5NvDm/l+n6WeqNE+RzeykdGOZc5i8IrYhUXEjycKIqUUYx4n2ezdE/42oJt9nsDq5YPvz00zA0RN/gHvoM46qdfldipT1iLkeCHSGEuLrbImiZ9z8B/1Ep9TtABPyP88m4Hwc+q5T634EKcG13lVtsuRbwAL/3vd+jJ9tDW6qNre1bmahPYMwY1PwapmFStIrYho0bukw2Jjk+dZzOTCdBHNCWamNX9y5GmyMUgjG2Zj0wDLoznWzufoIZ560Y2iNVnOPVseO8MnOWRzMR1TDm4bZ+9vQ+cNEYhyvTpPHYlK5jKY3WSd9/HUekaBJUv8N0NAJ9fXQ0dmOWp9EKGs5xJvftpe/JZHnlSpU8UdSgXn+VfH7fwuMdbVsoZTfiusNYjz1JfPQvMF1FpOoorbHsEumeN1GrHcJ1z3Nd+/WUy3DsGAwOzjew4bKdfldjpT1ilrreYEcIIe4mt03QorU+BbxzmePHWWclzi3L7VGzuPnciakTKBRZO7uw0WIunUu65cYhaTNNPahTaEsqXtpSbdS8GmEckjEzeJHHkfEj5NN5DkX9VM0i77/nnZSyGzkyPcJA2sY0HLSVoz0X8OG+t/PkvU9SSBfQ1a/RaAyhrSRpdc6dxvcv0GFqTKWJUEkHOjQGGgwDU2lMs0Cp652oPgMqc6imi+PY1FMzeHaNNLmrVPK4GIZz2ZmSwJjD7ryHlO5ABy7KdjCzpeQ5roHWwfXt11OpJCVXxpJa8GU6/a7USnvELOdagx0hhLgb3TZBy51k6V48MTG7unaxf/f+ixqXtcpv8+k8/bl+ytUyHdkOsnaWC5ULnJ49Taxj5rw5gigg0hGO4ZC20mTsDLl0UiqdMlNkU1k2dz7Azz76mwBEl9mNuPX6obNk1+PGJBoD0zQhDknyn+OkNZwGG0UcN1FqUSBSKEKhmOzD6M4t9EK50kaA2ewefP8CWsdEUZU4ToIY08yhdYzjDAIxRqaAmg9W4I2ZlFSqn2LxB664X88Vl1oKhaSjbxxfHLhE0UWdfldjpT1ilrqeYEcIIe5GErTcAK2ApSPTgR/62KbN0OQQzww9w0cfeOPT8+IE3Qfml2rOzp1lqjmFF3mgIWWlCOMQU5k0gyaGbZAyUxTTRephnVKmRMbK8I6t7yDUyeaKfbm+q+5GvHTPm5Tn8p3yv2ATQyjDxARI5lvmtwAI8f0ytdprFAqPY5pv7I683NLMlTYCHBv7M6am/pI4DlDKQusQw7Dp7PwxOjrezfj4Hmq1w2TcEsoN0I5N05khn99LPp90x11uv54wrDE6+tkrL7X09cGuXUkOy8BAMsMSRTA8DLt3X9PS0Ep7xCy12mBH8l6EEHc7CVrWWLlW5tDYISYbk7xSfmVhpqO3rZcgDChvL18UQCxO0O3P9XNu7hy9bb3UvCSnxcAg0hF+lAQ/YRwSxiHKUOTsHLWgRiFdwLEcGkHjoo6uK9mNuLXnTTYLve2PElb+FtsAgyh5wkXbAMS47glGR/8j/f0/g2m2XXZpZmlQNO25VGKHsHqKZuMkvj+DUqBUsveAUvbCfkmDuV/k9D/8BPX0qyiSOZ+8t5GB9//CsmNvve/l0c9iROfpzL9p2aWWhZv+j70dzDmiUweS6ic3mwQsT11byfOVZpautGS10mBH8l6EECIhQcsaq3gVjk8dJ4xDurJdyd47WjNWH2O6OX1Jm/jFCbonp09S82sEUcDw7PBCgBLFEWEc4lgOlmERxiFaa5RSoCGMQ2zTJg6ur6Pr++/7pzz/0r8nisbfCFqMRTGLSqN1iOed40L5S5jZPaQNm57SI5f0QmkJVI7/duwLnJl+mT51gmJ8koJdJW2kMZXCtjtIpfqx7S7q9WPJksinfp/7DrdTfWgQP+eTqqXIvzSDd/h3aPzOb10009Baijs5+QKDxks0dIG+XI2H+veRthwcZ4Dq2HMEIyfwrTl0Nk29PgR7Ibt3M4bvks3fQ/c9T11XAHClmaXLWWmwI3kvQgiRkKBljbmhy1RjisHSIGq+D55Sig6ng1Ozp3BDd9nz+nLJvjgAr02/hhd5dGQ60FoTRAEz7gyWsoh1TMEpUPfrZKwMjbDB/W33M+POXHdH17QRMdj/T5ic+ga+d5JkmkUTA3GsiHQAxPhhg6b7OqfmtlOnwHZPs78Dlt7yG42TfOPgrxM3T7DPmsCIJ4lVjI4i3NDHtByUqmKaOZrN1/G8szSPfJ30kSOwbRv5igkVCK2A0SdnaPD3cMiB7q6FmYbWLNWOQg/FoJ2M6qRcG+VA+WUe634z6qUDNKPvEPk58nO9VLY20BtSYNqEaZ9C16M03OHrDgCWziytdAnnasHOcnkvVb9GI8pQn3mBUul9slQkhLhrSNCyxhzLoTPTyUxjho5sB4ZK9uKZac7QnenGsZxLzinXyoxWR/nKia/w3LnnmGkmmydONaboyHQkGyiaKbzYo6eth562HsaqY4zURuhp62FTYdNCou31SPbwydLb8xOMjv4RUVQj0g3iOGn/T7KbAlXfx42nGQ/LbO7dd0m+Tms549zon6IaQ/SaETqaJdAxKQWWgpgIFdfxwyZhWME0HbTWXPD/hOa+SbobW7HCJLNmYvAsjfYazsk0qppDbx6g0Rhi6NSnmZg8zc78DrTK0qp06sx0MVodoTLcIDM9SrAlpljbTtQNfvM1UqN96ME+fH+EKKqtaeJrOt2H5yX5Kp7HVa93tWBncd6LF7q8NHqAcm0UhUGaGQ41PsuP7/lfZWdqIcRdQYKWNVZIF9jRtYPpxjTlenkhp6Uv10dHtuOi5ZvFVUZHxo4w2ZxMYgM0A6UBhmeHmWpMYRs2KSNF0Snyts1vo5Qpsa20jS3FLfzIvT9Cf75/TfbMWbxckclsp1p9CR3rpI4ofuN5kYZmFPHqxFFm2MSj+Z0cPfYs5fyb6duWtNmvVl8kjDw8lcdhEq0j0kqjlEYBhkoqj5MqHpc4NjFNhzgdU91QoeEcoXNkI2ZgU+2awm6kiFIhVrFIFNWYnX2WWuM17lWKlPsCnrkR39hIHI5Si3P4jRqNiTpsADvOYessvlmFdBbmKijPA9Mgjl0sq3jFKp+Vup7ck1Z+TrlWplI7sZA8vTjvpRWwdGa7MLTG1BEHZ85dkuAthBB3Kgla1lhfro8Heh9gaHKIXd27CKIA27SXXb5ZqDJyOoh0xLb2bZyvnme8MY4f+WzMb2TWneWe0j28ZfNbUErx4T0fxrGcZSuC1kJrucL3x6lWD6PxUDoJMkINs0GyRUAzgiNTk/ij32Z39ThGWnPy5d9hZscAwcNT5PP9WIaNIiJOiqIxVevP8xtItbJsF44YVM3TVPbUQE0zMXiayNZoQ5MdszE3Fknlz9KYOI7vl7GsdrwgBtVOKjzPuDvKiXqGLnuYplvn9WyGfd5DmNpBE2PEFvNDQfs+ZGIMw1l5Y7qruJ7ckyuVyWezu5iYe5FybYTObDeG1lhM45sb2VB8E0cnjlKulWWzRyHEHU+Clhvgkpb9wRt9UlplyG7ocmzyGIPtg0zUJwAwDZMthS3MNGcopAt0Zjqp+3We2PoEzajJ7u7d7O3de0PH3lquyGTu5dTsWY6PfpMNmeQxL05mWVIGHJyF0aaP5Z6nPBMy1BFTaTfoKx9hw9EapYGH2eqYpMwsXjCJuRCc6IV/L45ZtA4IwybEU9CmIdYETnMhC7i2McYpFPGqL+N553GcQbRu0pbKMOc3aQaQU5M01WM8N1ehW1nUzk7TLBR4b9GmYZ/DibpI2Xk8qwxGhXRqG6aZW3ljuiu43p4rrb8vg+2DGMogiqOFZbcP37+fsVqZrHqRVGwCMb65gZr1ZkyVVJgtTfAWQog7kQQtN8ByLftzqdxFgcxkc5KRyggb8htwLGehI66hDPpz/eTTeap+lUbYYKI5waMbH72puxBnMtsZbzZ4ftpga1tMn/NGkNGM4O/Gku/rtuKrXbMMNFLsGY2x791OsXmAs5MXmLbHyVs+QQQpHZNNUlSSgCVZPeGNSl9vvh/MvCV/M2MnpOGfwoyyRFEV3z9POr2VLe2DvD79OmO1MdIpGK8OcaZewMv1ccGOOHbmr3j8zHZy99RplM5hjyn8+3N4lo0fKYLKqxdVP12pt82VLM49qXhzNEOXjOVQSL+x9FSusey1y7XyQgBrzL8hpmEyUBzg9YnnGZm9h66u93P6zDkG7B5QWSIj2cU8iiOqfpWJ+sQNm30TQojbhQQtN9DiPimffeWzF32SLjpFDo0d4sDIAR7f8vhCR9xSpoRhGDy68VHOz51nc3EzH33wozf9ZjQTwESQZXMhx7HZCsMNSBvQbsPxCoz6YGrIxRaeEfPWeANGpUJU30Bg5FGlMY42PN7a20O7lcbQEEcVIq2wDAsdx1hGON/EDogXxSyLeq0t6l4CBOj5ZGCtY7T2UbpJbyZF4KcxYh/D6mZP744kAVrNccoa5c+Lk/zzI/fiOQ2quZDvZjdxwukiRYiPZruneTJX42uvX72L8eWYZgE/9nj57LOU62MoDDQxfW297CgWeObYVzkyPbLstSteBYVaCFgAlHZpjw6QNo4xNhbQ4ZTYm/d4uTLLQD6NzRi1IOarp78PwFdOfGXVYxZCiPVGgpabYLlP0h2ZDnZ27mRoaoj7uu9jX/8+Xoxf5Pj0cbqz3cy4Mzyy6ZGLWu/fTBWvQsV8kLZMgy3+AdzIRyk4XoV/mExmWfr8FD9e38T3s1NJUq1S4PmMNNoZ7Z6gJ2MSaJPIbMc1tqNNH+0dpcdQwOQbLxZfPMmy2EJvuxgwIkBhGDm0Dmk2T+O6p4gil04jxNMG7+yY4UDTxW0qVKNGvqOfc+0lytv20Jfp5vPmYYbKxxncvQ+jvbSwDPPs2WfpzHYuuzyzkiTXdLqPI9MVKvXjdGa2YSiTOI6Yqx/ny9MOlVTbZa/d2s4h1jG2rmLgkg2OYMYzTOkC+ey9OKkc9xfrpMMDTNRdtDao18e4P5Njc/+HcVKlVY9ZCCHWGwlaboLlPkkDPLzhYaab05ycOUlXpov+Qj8PbXyIJ7Y8sWYVQdeqkC4QYNHV/ZN8a7TO8NxRKkHETNKqhZyR4tGwl76aJs5qUtqEKIBajcbmbo4FKfba57EURCqLIiLUnUxd6GCzXcPvidAKFCaaaD4z9+riuI5l9WAY4PsTJDMwFtqwmfYMOs1ZHsoc4tnKLiaVyyarj7xR5GQaTqoJXmqe435VwvB8IFmGKTklvjX8LT6464OXLM+sNMm1XCvzSrXA/W33kYrLoJP1r7q9ha8PD/PDOzuueO3dnYPMTf8lXU6IoQOcaJjpKEt/7s0U0sX5H77CRsdlW897mPF8Zi58n7d2OvjqGFUeW/WYhRBivZGgZQ1cLQ9i8SfpxYGLbdjc33P/Da8IuhatfZFeGHmBiSAild6GpWdpo0ExtCiGJicKPv1xhnuqNoeCc/S1Fejb0oN9/27unf5jtrQ5mGYnsbLQmFjlV6lXKgT5PZhhHW3VUb4mTK08aAGIozmsuAsAw8hiGDZKK4p2hVFfUTLHKaYcCpbNVqOHfwhepxI3CYk4GJynShv7bINWx5wgCjAw8EP/otdpbaOwkiTXilchJk0t9QM04zkMXGIcRhsugT6PH1352u/ogleacL4JaSKKaDrSBpvnJ9nCcI4gGMc028nZaQLSOFaWUJVIRSOY5hyRUVzVmIUQYr2RoOValcvUpkZ5evIfGAqWz1VoWbwx4kBxANMwieKI4bnhm1IRdK32797Pa1OvMePOUHJKdGQ7KMZFLGWgKlWmvDlOZ7N4aYPTTozr1GhXz/Pjc6/wWKlOWtVxoik0Dn6U5kwwDh1pvlcd497Qpw2DjB0uKnu+ighM30AbPr4xAirE1jY662AYkLVMNioXy1D84EYDbVl8f+I57LCLPfYGqpHLafccIzlQjdd5rKMbINkCgZiUlbr45eKImJVtjbA4MMUoElGcv3ZyjZR5+Wt7XpnIH+aRrR/gPr9K3Zsgqpu0OX2E4QRhOEccu/NbN4BhOGSspFdxrBRoAwOXiOKqxiyEEOuNBC2rVavB00/D0BBPW0cYYpLB3p0YDz9MZC+/mzMsUwZNvCZdbG+kXCrH++99P187+TU25jZiWzYZK0MzaOJ3+v8/e3ceHtd13nn+e+5ae2FHAQRBgJS4SZRkSmIkWYstW7Ylz1ixRXsSLy2lJ+6k3Yn9JJP96Twzk05PtmnHE7eTSTzxoyRykralpGPHiWV3FFuWLFkStVASQVGECHABCjtqX+6958wftwCCJEgCFECA5Pk8TwlE1UXhYBHwwznveQ/F0ZeZcdKkohu5pXkDVb9Kt3iVFrtA0nIxDEHZD7Coka/miUckG02HmXw3qjJDOTqOQGCb5x8LENa1BB4YBigbISSen0dUiijHxRIBEdOjLl2m/QT5hEmrP85HJNTGZkhLSXfLRkbabI7nj5Or5UjYCWaqM+zq2MV0ZZpmS2EbHp60GSos/WiEuWA6OPk8fcl2MOLUWfDc1WlSbuqMwJpJZCiXDzG38yjlpkm5aXJqnHo9G37YsooQNr4/Szy+A8tKk7Igk+giWzxBs1VnzCuAtXg/IE3TtMuFDi3L1Qgs2d4WDlYC+sVmjOwYvPgS5q23nrWmYLFt0JfCL5YtLVvojHdSC2ok3HD2yLXcsHuvYeAm0rSmujGEQasjuMq1ydcdAjlNOrqZdMSg6hdBZqEYxXTrCBEwM3ENnd2z+Ea5car04sW4Cw6YBsJuwZFRRRDzqKcCpBlgVAUCBRZYCFwnxfam25mqTPOWaZDIGExv6qdoJrgqmYLpwxycPMjhY/tpIsLOzC7uuevf8fTB3yFf+h5KGZhCcnNqF/dsff95P0fZYpbZ8ijvSIzTVD7CVHEfQoCnWrim7X7ed9Uv8O3D3z5rYF3stOdkcjf5/AtUKm9Qr49jmjESiV2YZhKlAoQwuaZ9O9nZl3hqosDLxToSya6OXbx/y/nHrGmadinSoWU5slk4eBD6+8nLibCHq2lCWxuMjEAuh5lOn7OmYOE26EtBJpHhQ9s+xDfe+AaTpcn5mQIU9KZ6idgRakENP/BxnQogQNh40iQI8ggBQhZxDA+R8pE1sLwxCMaYGa5T9iQ7OkG2guFyZm3LXBe68OxGhAJnQlHNKIQU4c4jR4HyGjuMbKp+mYmxbzAetPPW7DHyjs/rlSqTXgR/3CfjtnH9rMu/OSzZIiJk5AlGx3+H3buS1Ns/QjWoEzFtHDVDJfc46ejiO3EWdrHtEa8RY5J0fBs39b6TQHlYcpKOpgiZROaMwAowUhiZD6+nn/YshI3rdpFK3UhLy72Nc6ESpxyu+Eb2FQJrA5t7bqE7EPOdlx8ffFzvHtI07bKkQ8ty5PM0igpIqQhyriGcYYRLFtUqQTJx2dUUfPK6TxKxIjw/8jx+4GOZFttatvFS9iWeG3kOL/BQStHkSDZ2VvBUlBpxUAGBLGEoA3wP5RkYVoDRsg/LEUzlYaIEG21IeKAkYDVKXEzC/OMBCgwPhA9WBQJXYeVBOoLYqEmlC4KEhbR8am6UmhSkbY+4GzA8a5CvzvDmjEFFRfEDnzdyr7DDS3F974cpCo8BxvDyL9L6+i6it+5oVKOAUulzdrOdCyxXpVpo8QLqbA57tBgmt/bcilKdp7x9JpE5o8ngXB3Uh7d9AGDR054Xnls0d7jigbHneGL6GBuad9EcaZ4fc9pN691DmqZdtnRoWY5UKjzgT0oyZpodVhcDfpY+0YwpJYFrn1KrcLmYW9r6wFUfOGVp68G/f5DRwihSSlzLpVAXvGoV2Jauk4j0AyNIlUDkxhFS4htQrICIgOkLEhFFWoFyoVqHaCQ8JgCDk2cEGWCWw7Dix0BZ4fWp12DmHYpSxkfGQAofpUDWCkTcCJ70QI7TFfF4ZTrKGzMTOJaDIcEJAo5YZX4q/1dsttpIuGW6ExOkh5/hxg1duL2bARozHosfpLiw946lwmMYDMM8ecJ0LXdKN9y5tz9bu34hBA9eHwaSscIgJR+GKzXGBr5J3ImzZ8Memm0oVEf558EnefLEa7w2NsgbsxN0JbrY3bWbiBXRu4c0Tbus6dCyHJkM7NgBAwPQ18feyG4eLb/AwNQbGO3tSDWz7otr346FS1uD04PsG92HEAJPetRrdQxh8NSUSZvrUa8OURazGCoAV2EBQgZETBAy3Pli5CBuglmDmgluAESZDywooAYiAGmBMwuxo+BOQvQEVLrBS0HVBkuAL4CoRNYVhlvGVnk8meS1AjRZkqTt4UuLSg2mKPJk7U02qATN4zNEKJKt5Hnxh3/Hrcdvgxt3oxybYq3AcH6CZnlqDdLC3jtSrBz5ZwAAIABJREFURcLBNra023gUK8eJGuFBjEFQpVw+xHTt5HlThXqBqlclYkfm66AGpwf5wdEf8MLIC3x/6PsM54ZRKJocl/dnIryvdxumYVGvzfJj6U2M5mMk3RZGi6O8lH2JW3tu1buHNE27rOnQslx7984X4yYMg4dkF9mtN5K/505SLWvbEO5i+spLX2G6Ms3mps0EKmC2MstsZZSbm2tsjvnMVMFygACqHqSiYBphJkGArEMsCb4JsgwiBcoB6uFSkBKgIoANKgBMiA1BZBwiWai3QcsLMHMDVHqBarhCF9QhKJnYnTaWUcVWRf7NhlkcAyZqkPcEIyXFiwXI132+VnqOFkzuCQR9LYpjVsA12SGcfXX2d88yXIHsiTNb5MeMKso7xmTBw3Y6iZldOMFxDL9AixhHViUT5R9iGDGkrGMYEY7njxGvDfODI4eYrJXml4e6El20Rdv4yktfYbw0zosjLzI0O0QgA4QQ3JguERMlvjv8Gm3xDq7rvBaHGfY01fjh9Awt0RZGCiNMV6b17iFN0y5rOrQsVyIBDz0UFuXm85BKkclkuJJ+RWSLWYZzw7imi0KRq+bIlrLc1VpjcyxsunK8DGkLHLOx2iPCiRMJWBJUHmQcrATUqmBHQdTDAxQNP5xdsbIgY2AVw4kXpwjJQYgPwvi7waxD6k0obwRMMCpQTkLdKuMoE58EHc44zTYYAnqiUJeKoAV2V+BISVGZ9Xlt3GbfCTA3WrQkZzmebGKm9CSTs73EWz9Av5WZX8b5uwNf5c5WxcvHvkHaO0itWmVSNjNoXsN7WgoY/lHidhuu6eIrDynLVGpjHCo5DM0eI196FVWxyHtXs6lpE6YwGc4N8/TRp3FMB8uwODh5EIC4HafZUfRESkzWHDyZRwK5aoHmSAvbUz6jfpq3CtPkq3kGZwa5vff2y3amT9M0TYeWC5XJhLcrUL6WpznSTNyJcyJ/grqsEzdqbIwpch5sjEJPDBwjDCtNThgalALDhkCAMkFVIQhgtgjtCXBL4CfDHUICMD2InID4W+AnoPvvIXkEas2E9S4CjBpEJsCPgExAXIBhSUp1MK0JMq5CKggID3x0BdRUOMapKnQ2+yhV4alxh8JsClmZRWSySGOGhBcjUfku036CrNxGcyzD7NQ/8FyhxnjdpKflJo7lj9HmjTJR+xFv5Syubr+dq9q3YQmDQuE5LKuZ4en9TNbbaIm18sakQ3ekyv6ZYV6rF+hv6ufQ5CFmqjM4poMvfeqyjoGBK0psijm4Jvh1iSc9CvUCBycH6Eh0siUVZ0/3tWyq2QzODPLze35+3TYq1DRNWwk6tGjLlnJTxJwYV7dczVhxjEq9QpsjUQp8Bb0xiJhhnYnZaLSiCEOGL8EygSYoT4ThJF2HlrfALUO1C6QLZgmCGBj18PWm18LAAuDOQGwYypvCpSJ3EmgF6YEzBWYn2NECdUOGS1EBOAJsIxxHBEjZcGc7VAO4usXnNekzkK/xE+kYHS4cLVkcLk6Tmx2nM6oYrbzA96ot3NdZoehsoL95MwaC7mQ3hVqCDfEsNQU7O/eQcNPU62OAoBbUKdbLpJ0or04dJlctYtg1kAXyxVF+NHuQGV8gEPiBjxCCiAF3tUs2xepETY+r4hAxAt4sCKJWFGEY5KozjBsFXDvc5nx77+06sGiadtnToUVbtrnur2OFMVqjrdSCGmVfIIRiSzyc0Zg79Nk+re+KZYAXgB1AcwUis2FgCWJQb4HYcah2hDuF/DjYM5B6A9qePPV52p+EiTvD4GLnod4EpoJKR1jvgpDzm5BM82T7l4UvY2a4fJUGfnMH/I+JgO3xAjZVOlxJExWmPUFOtbElqZjwSkxVZhGeS296E0dyR8nVcgC0WQVMfN6aGmRTUwcuVaSsUFMGpvAZKQwxXjiOlBX6G5+jSgBCBAyX4ckJqEmBUIL7u6E/DicqMFJVRAzojYXbqszoJiKWC/4oB2ZNDHmMPRv26CUhTdOuCDq0aBdk78695Ko5Xsy+SDqSZrocMFkrc3MzeCr8pWxysjccC17aZrhteTwNM6/ATQI6DsPEDii1glUFuw5ODvq/DE0Hz3z/VhW6vhMuFQXxcGZm5D6o3wXuCFQ3hKHpXOcwCgF2Y4wpB368O3ylFHiYAcx40GSCqWYo1h3yNYWM+UwVRvjRkRIRO04y3oxSAbn6BNPlIk3m31DOp4g5TTRbdTxvEkMFWP4UO5NVUk0wVYeh8snPzaYY3NkGT04qPphR3NMZzgD1RCFbhSOl8LotKYt0KoFhRMnJXk6IOB+/9uO8s/edK/RV1TRNW990aNEuSMJJ8LlbPke+lufRgUe5tv1aXs1+mzvbIBMJA4tlnPl2czuZhQknJPyRD+/eADdsBDkF8QJYEWjzYNdTiweWhdwZYCYML/VOsKfDIFLt5dTEdA5zO6yNxj9SBsQtaLbDgGHUq+TrPhWZIV8s0O1UOFGcoFWVSJQ9avEZpCyzNWUQN+rY0sOrV8j7FnHLwBSCuixjES6b5b1waDEzDHd5DzbF4UM23NIKSStcyhJA2gZbCI7Vm+lSFo6MUWITY34XwpxkS8uWJX/NNE3TLnU6tGhvyy/c+gsM54Z59vizzNYNjlckG6LQ6pz9bQzCJnJ+Gfqr8PgMBJ3wnq3Q1gg6xTrs2wg9kXBW5XyCeGP3EeBFw/vEEkPL6QRhPY5lQncU2lzwJLhFxSvTBnu64myK1/BUnbIaRfg+rh0lYhoMVQUJ06I9AgZ1WpN3kigd4cR4lbJXoS8GTTbckIbEgv/7DCDWaGvrqUbtjyFIWAbXNxtsCfK02oJy7VkM9UMKVZdseRffPvztM04V1zRNu1wt8rewpi1dwknwhQ98gdt7byfnm2SrkFpCFDYFdDaB7UJNQuDDpA+5E1A+BMGrMNsLg/csbRxmKdwm7Ueg2gOIxszJ26DmnhvYEPW5sW2C5oSgWGxlcLKNQxNtvDiZxA5MLBHBMSQRQ+LYLrOeQODjBRWSkVYMI07RC4+zbnagzYGcH97KPvTFwmWilBU+3hOFTleRdhQttke341MKfGZqAVIJtsar/GT37Pyp4pqmaVcCPdOivW0JJ8Fnbv4MX33lqxzIe4suCy2mPQmYcHUK7ugFW4ILOK1g5WF2GvL94dKPO3Pu53JnIIiGDeqsWfAbsxZzQ1nsBOnzmcs8Cmhx4J6NdW72ZijUC4xWDI6XBNe4DptTkqhdIOl4NFsGUkBNCJRSmAQopXCsJGU5S8mHzY0iWwhndHYmwTXDouVmJ9waLoxwiQjk/PJVi6Voj5QQcyNTr3KT28JLE3WyxffrhnKapl32dGjRVkRHvAPHcmh1lr4iY5tw02bobYNMFOr1sABXToOZgoQD7pvh0s9c3cpc0e3pIabWDGalsRW6l7Axy1xiCc+4RF5IcmEuPIDjgBAKy/QoKbg5DhGjTtSS+EgKHsQtRUVWMKXCxUZVhpDRq5Cyii3qKKDkQ9KGiITOSNjPZqwGGyJhiFEifDm3+wkW1gcpQKHCTdI0Bc+wQVj6rCFN064IOrRoKyJfy9McbabFnV76Gwm4th1sJ5xZiJngpMOajpkitLph4SoBjL4v3N48V1wbGw63Pc/VuwTxsP1/0+sgalDeANIhnLpZIYKwv0vcgk43rDtRhFurfRRVaVKTirQVhFutVY1qeQhVPsa1sTo3J8NtzkU/fD6pwtelCj8sQ4ThZG5mZZ4RBpiFoUs0Fq8EPhnxGjFjCYU/mqZplzgdWi4h2WL2lFOW15OUmyJmxogup0pKgGOH5wUpD4QFgQGiBZpmw6Jcdxxyu042khMqnIkob4LJOyHznfCpzBLzXXKDBBgSqIK0afzGX7mP1WwkCsMMg0sAoEyiJlhITBV+cApFYAQoEdBlhuEkYUM6gOFC2D8maoRDy7hhcAnkqbuujIWfz7mPY/4+QSAhaUoSorByH6Cmado6pUPLJaBYL/LogUcZmByYP2Rv4eF960EmkWFT0yaov7qk6w3CFvythfCb0K5APd4IGS7EHbBnIf0KTN1xMrBA+DKShdKmk/Uuc11yC9sgSEJkBIr9rGqpuUG4tFMNwFcSz4OWxoHPQqgz+tMYjb4xURM2J2G8BjErDD4RMyxInmvMdzo5/5+TLyUK0zBJRFpW7WPUNE1bT/TuoUvAXGDpb+qnv7mfvnTfutw1ErNj3LGc359VsC0QLngt4cyFXQezHM6wdH4/PCQRdTKwzBGN+4L4yfvanwR3IuykKxQnu9utMseEhKlochVCnQwni5lb+nHMcBt1yW/snlJhHcvpbydlY1lokZkiwwDXTBCLbiIa1f1aNE27/OnQss5li1kOTh6kv6kfQ4RfLtMw6Uv3cWDiANlido1HeNJoYZS25dSQ2OGOH2s6XN6ZW9rBhKb9Yev+hcs+C6nGfWbp5H1WFTKPQ/INiBwHo8pFCS2CcMnIEmGX3SVRkDShLsOZmkDOldgufchSCmy3l9bWD+G662u5UNM0bTXo0LLO5Wt5BGI+sMwxDRMDg3wtv0YjO5Ov/GXt0BEynFXBArMOVjkMMPFh6PheGELmln2qmZPBRYnw9fjwmbuI3JkwtCgXZJQLai63XBfyLoSApBNOoAxXYLQabnVW6vyhRREWKxdkCif9ATo7P3EBI9A0Tbv06JqWdS7lppBIpJKnBJdABkgkKTe1hqM7VYuIU1XhAYTnMvdRxEaAKai1gxmEhydGx8DKndzmDKcejjhXiBsfPvMQxTntT8L4Haz7SG6KsKalLwbZWth51znt/8iFTX0l4ZLSdN0iUIIDZYejiXa2SlgflU3r13ouYtc0bel0aFnn5k5UHpgcoC/dh2mYBDJgKDfEzvad6+oH8CbHXXoX2ipMWlDIgGkAEXArsP0IiKYzl31OPxxxsWZzc48HNtglEBVQF2m25UIECo6Xw7OHNsZgshYW9irCglxJGGykCmdm6gGUAihKm0qgmBA3ooojPDbwGA9e/+Bafzjr0qVQxK5p2tLp0HIJ2Ltz7xk/eHe27+SBHQ+s9dBO8VpxP7e3Lu1arw5lC6IVUCWQBpTTcOhO2P3NxUPJ3OGIp/MjMHovlDeHPVqkA9VOMIJwOzWRt/NRrQ5Fow6G8BRnCDvl1mVYG+OpcKLIV43zlCQcrcBs3aQoBQV6uH7TvdiGPV/btJ4C7HqxsIjdEAaBDOaL2HXQ07RLjw4tl4CEk+ChGx5a11Pc2WKWbHGE5nMclDgvgHIApUNgNYGZDH9BWz4UXKgPLv39+hE4/Bkobgar2DituQqV7sb26dMOJVzBdi0r4tp02A34W6PhTiIFXJUIg0uzYxCxFKZS1KSg4AnyKo4T2c62jr04VpjG5mqb1tv3xFpbShG7/pxp2qVFh5ZLSCaRWbc/ZB9/83E+vimYb7x2LhUPRsowFIV0HvrGwy3AyoNSCxSjS3+/o/dCaTNER0DZYVApbwJlgDI5Y2no7ZxFtJIEYRM5m7B9/3s74PlZ+McR2NUENzRFmfAMAuWRrVocr7dSkib3bd1La/Lq+edZj7VN68VSitjX6/9PmqYtTocWbUWUxp9lc+z8ta9KwfQLkG+CuAc5F4aALTPhEpEyIDV77ufIJiDvQiwKpf7wzKFqFwSxMKjU042dRgurWBdY68CykCKsV+mLw/cnYbQGKX8b+ZxCyiLVwCAV7WPvrr24psuR2SM0yWBd1zatF5dSEbumaUujQ4u2IroNRcU4d58SpWB6Gqr/CvVbwM5AfDoMLiUHVDukj0DP8cXfvujAoztgoB0MBZFWuCEFnS2NmYtqWM8i/AVN5xakqPUUVuDUPGUK2JmCeze/i82t1+AFHnVZZ3vrdvbu3MuWli2LFpWux9qm9eJSKmLXNG1pdGjRVsSgV2PT+ZaGBEyWwUtC81MwczvYnWEL+4oLnW/CPV8LL52bTUnVIFMM75sLLP0zjSwiYSYGERea6mEoEkE4Y4N96rteb4EFTvZjmZsQandgQ/UQ//HOv1m0dulSqG1aby6VInZN05ZGhxZtRbxU8eg4XzdcBd2tMBBASwU6vwuFNpjogA8/Cze8Gc6mPHw9fC99Kx4ZbLK8K/cMdxyFg20LAgtAHqLjUO6CdBX8RFjTIuOEFzVSwXoJLIutVs0FFwEcLYOjxsjmDnDDhrvP+jzrubZpvdFBT9MuLzq0aCtie/t29h+DDbGzXxMEYNuQbAHGIBAwHsD1+8PAAvAn23vZ1/8TtHcWw5kT0cW3x+7gYORvERw9o2am9AoY26CeAlLh8tD8mUNi/QQWONlvZbG6n2oAh4rQ4gjG80dgw0Uf3mVNBz1Nuzysq56hQoj/XQihhBDXNl6/WgjxQyHEISHEc0KInWs9Rm1xP737p6kG575GNPqSJHfD0TYYaoKdE/DAQPh4NgH7Nv8EbZ0FChMJipNJChMJ2joLHOn6CfLOqSEkfg0kdoOTAtm6ILAYrMuGcoqw34o8rU9/PYB9M1CRELUdOlL9azI+TdO09W7dzLQIIXYDtwBHF9z9p8CfKaUeFkLsBf4cuHUtxqedWyaRYXeTCZw7ubge7MzBdXHIfOdkvQrAd5tvpbWzSHEiycnUIShOJGnLFBEv38qR5mfY1gwbPhb2d8EIT4dGESaadRhW5syd4lwJwu3OgYTJOgzk4XAJMhGIRLefc2lI0zTtSrYuZlqEEC7wJeAzNJb5hRAdwG7gkcZljwH9Qoi+NRiidh5D40+Qts992J9phlua0yOQbIXm04pl89bc9H2YPCKxEqn0NJFYGaVga76Nm/pgw4NhUzoAQ4Q7icKucqzr0ALhLiHXCDvfnqjC63mY9qDVhaJs5efe/bW1HqKmadq6tV5mWn4LeEQpdUSc3DO7ERhRSvkASiklhDgK9BK29tDWkXz5CCgXpcpn3fZsADiQvQ8iE9D69Knt+nfMZjlIF7g1tvQdIt0yBUAgFL4X59rpIVo2w6xLGE5OCz0s0ptlPXbBtQxwgFdy8LfHoNm2ua3vHj5/75/rugtN07RzWPOZFiHErcDNwB8v8vDpf7if9e9oIcQvCiGOz92KxeLZLtVWQSrWj2G51PzzX6uAeguM3H/q/XfPPIM6kuCaXa+QaJkkX0wzW0yjpEuzOUvrDUco9RF+F5ztO2Gdz7TMkQqenYLm2Ga++hNP8lcf/ZYOLJqmaeex5qEFuAvYDhwRQgwBPcDjwLVAjxDCAhDhFMxGTq15maeU+rxSqmfulkjoE1wvpr6Ouwm8rnOuD0kAA2QE7CmobIDCaTWnP/fqI6QoIWSEWLRGKlqFyRg3HByk1gH+3CzLEq23WZY5E1UoqTae/umnuWXjLWs9HE3TtEvCmocWpdTvKqW6lVJ9Sqk+4DjwfqXUXwAvAZ9sXPoAMKSUGlqbkWrn02v8LxRL575GQrjDJw8VA8ZO29rba47w3ree47p9x+gcyNF74Bj3DT5ORzGPskE5nLtw5hLxcg7+5Mf/m55d0TRNW4b1UtNyNj8DPCyE+A0gD+iz5Nexzs03MPLq+a+TJkxeDV4F/nUjbLwe9g5Aog5mCeoWTLceYSYBswpGRLjLqKcMdHP+mZZ12KNloWoAb1VbuLtf7xLSNE1bjnUXWhqzLXP/fgO9xfmS0dm3h+hbjXb65wkWkSTYJmxpgde6w61hD74SFua+5kO+H1pHwp1B0oBiP+Q7wtOgz2sd17XUJGRr8Gvv+n/WeiiapmmXnHUXWrRLl+tmcFQaKXJLut5wIH49bDTheeD9g+H9Lx+HawxwOpmfNUm1ATHCNjBz25sXo8KbXPOFzzMpFW51zsk2brvqk+d/A03TNO0U6/BHu3Yp6+37ImqJNSfCAncj1HfAwc3hmUOjCZAeFPfB7A8g/zyUD4Lhhr/0CTh7/zpJGFrO05n37VDqwktqhICRMlitP7eiY9I0TbtS6NCirajNOz9FfRmhwbahtxOiFvywB/5mJ8hGPUpQAG+SsHC30YNFVJmfTZkjCbvL+oTnGflGI1ysUMHuXLPdig+14MJXnxTw1WMmd1314ZUZmKZp2hVGLw9pK648ZuNu8JZ0rRAQSUC+E44qeGED9M5CyYJrJsFUUC9CTUHCA+UTftc2QoxsJAhhLGj8v5S6lyVQhCEFERbPRo0l1tSchVRwdUJyXed1KzNATdO0K4yeadFWXHPxp5Y1y6GALT50lKCpAkUHXu0ID1Q8koZDVTCHwKlC0BzWq0hOBpYLdd5v/kZBccELm++65tubZVECSsFlsF9b0zRtjejQoq24u//9n/La1NJqPyQgAwg8oHGg4FXTYCXh/Qo+XoGfPwjv+wOouOfflbRU5/vGl4RFs98fC2dI4qcfGbBMAvACeGEGnjjyxNt7Mk3TtAtQq2Uplw9Rq2XXeigXTC8PaavivcMf4VD670hY4SGBi5GNVFOswfEZKNrQHsCNW6FtG6h3QL0KU3mYPtJoLLcC5gLLufq4GIRbsu/JrExQChQ8MwUnqg5HZo5A//nfRtM0bSX4fpGJiUcplweYO5EtFttBe/teLOvS6h6vZ1q0VXHNp77E2GEo+yfDyekMEc7GPDkGA07YXO7ONti4OQwK7gTYs+AlIXcLmCsQsZfzDW+wcjM7YxX4UmNLdzqSXpkn1TRNW4K5wBKJ9BON9hOJ9FEuDzA5+dhaD23Z9EyLtjoyGQazPUQ3HT/n0oop4IEtkJ+C8Sq0d0DdhRYFdgaKMQgcCIyw2PbtMBZMsVzsbrk/mIRZH7oTbboQV9O0iyZcEjpIJNKPaPwQFcIkEumjVDpArZbFdS+d40T0TIu2ahLX38ZI9fzXCQE/tQeifeA1QYcDPR5IGzDAKoMVEO5pPgvjLLf5xxZ8p1/sxnOBgh/NgIXFB7d+kJSburgD0DTtihUEeUDMB5Y5QpgIYTQev3To0KKtmvfs/gTRJe64EQL+/bVwfQJaXMAF6YKoN5ri+mD5oM7SA0ay+OyJMf+fxnVnu3AVjVXhlTxsbt7MO3vfqQ9J1DTtojHNFCBR6tQffEoFKCUbj186dGjRVs1VHXuo+ZElbxN2XLAMUAb4cZBOY6uwDQgwA4iNQ1A9e+M4eZbpFikbgeUikwr++mj47/de9V4e2PHAxR+EpmlXLNfNEIvtoFo9gmr81adUQLU6RDy+85JaGgIdWrRVlElkaOr4KMtpTaIaVVbCD2davCT4FlACswKqBMI+R4HsIsFkLcIKhMtCw2V4dgYSdoIv3fclEs6lVamvadqlr719byO4DFGtHpkPLG1tl94fUboQV1tVH77uM3z/2b+lxfWWNOMiXbDLYE9BwQwDi1QgmsCSUFGQPE9XWsnap3FfwbEy/OMozHjwZx/8/BqPSNO0K5VlJejqeohaLUsQ5DHN1CU3wzJHhxZtVXU038LGzvcwPfltIkv4bgsEuDXIxQAbLAUiAiioifBk6PNakFrWYpbFl/D0NDw7BU9Owi0bbuHTN3364g9E0zRtgTCoXJphZc5a/0GqXQGuufr3EKKN2jl2/8wRAgIfVBOYMiy8laJxU2AvMWavQb3tvKem4IuH4TvjsLNjN4985JE1GommadrlRYcWbdUlk9exo/uX8Se6KNXPfW0QAEYYXso+qAgYNggHXHfpoWU1Ustph0sv+nhNhktCs57gt971W+z7mX1sadmysgPRNE27QunQol0UPds+w03OfVgTPWfd+aMU2M9CdAgqPrgxMI3wbCLZ2OosxNl3Dq22873bQMJMHeoKtrVs49M36iUhTdO0laRDi3ZRWFaCrnu/wN38Em/NngweC18eHoNb/j+wJ8L7THNBWCFs+x9efDFHHo5NqUZB8GnvXhHuEqoEMFYLJ3fKPnx050d1PxZN07QVpkOLdvEkErif+hx3bv9HXpptZroKJR+mq/Dim/Cu3wCzIvAC8GdBKE7tt9IIDorw5cWqWREivFmNWR7ZeP++CgOLVBCEw+N4BdriGX7ljl+5SKPTNE27cujdQ9pFt+2aD7LtmmmeeOPLZN/8PplqK3tv+1/hK5PwyCMknaeozL5JU3vYfl+pxrJQHXDC05crNbBdltZu922S6uQsjxBgEoYUQ4TvXgI2UPBhpALb2nfrfiyapmmrQIcWbc3cve3TsO20uo89e7jqsUd4df+vUOosEE+GMy7Sg5oHbj38prVtUO7JpZrVzC6LrUaZ4tR/CwEbInC8BP8yNrSKo9E0Tbty6eUhbX1JJLAe/Fnu/cR3KDy/g6kRg2oVKnVwPJOu0Y3cPPafqJsfY/+syY8mBTNVqPrgrcJ6keIc3XfnrmnUusQs2JIEU+ZWfiCapmmanmnR1qfIddfynq5fYWLw21SKxzEtSWvTzWzs/ixW9xZ2jj/BvqceR6kc2WIYGNrj4ayHsYLTLoFqJHvRCDCnPa4AJcImdkJAmwtbk2vVIUbTNO3ypkPLBcgWs+RreVJuSu8QWSUTE49SLg/Q2v8uhDDmD/gaU/8DY/QHBOUBrk+34Hs5xiswXoeIB5YJ7tsMLYEKZ21c89QAtNjTisbNaMxZRk0wg5m3NwBN0zRtUcsKLUIIE9iolBpaneGsb8V6kUcPPMrA5AAGBhLJjrYd7N25VxderqBaLUu5fJBIpJ8gKCBlFcOIEIn0MTn5D0Sj/cTju9iauJuhsb+h3a2StiSOAYYCJUFc4MKnAsoBPDcNvVHYnDi1fiVQp75+OlPANakqTxx5grv7776wQWiapmmLWnJoEULcAfwN4WaJXiHEzcBnlVKfWq3BrTdzgaW/qR9DGAQyYGBygMcGHuPB6x9c6+FdNoIgj1Ie+fyz1OujhAs0EstqwvMmSaV+DCEMYoluuoY7qFsupWAUVA2R8FF2cEHboQMVbln+wpuwPx/ed3tbnNuaK2yISFwLOhxIn+P8o0DC1Qn4waG/0KFF0zRthS3n79HfB+4CpgCUUs8Du1djUOtRtpjl4OTB+cAWtkUdAAAgAElEQVQCYBomfek+DkwcIFvMrvEILx+mmaJUGqBWG8G227DtVmy7jWr1BJ43hWnGALCae3CiG3D8Gk3RDaTcLoRtIaSFscwGdDM1+NoJk8+9YvFWOcH21u3cd/V9fPKW/8IvfOgEz6qf5IuDEWLniPkK8FS4rHRw9NsX/gnQNE3TFrWc5SFLKTUoTt1KcZ6TZC4f+VoegZgPLHNMw8TAIF/L6/qWFRVWiyilECJ8aRgmQiiUOjmPktx5P/kD/52Kd4RoqYWyK/ANiZACw1RLmnEpG310bf0DUuYz3Bs7TpPbhG3a3LzhZh7Y8QAJJ8FfP/DX/MzDAttYvCAXwvsNAbaAsWJ+ZT4NmqZp2rzlhJaqECJBo22FEOIaoLoqo1qHUm4KiUQqeUpwCWSARJJyU2s4ustLEOSJxbbj+1OnLA9FIpsxjCiVymEsazdCmIhIAnfn7SRr99Ikr8Wb/UNy5WdRgQ/m+adbakYPd/zYM6SjGW7bsvecRdY94QTPWXvCGIRdcwMBhUB3E9A0TVtpywkt/wl4HOgWQjwMfAD45GoMaj3KJDLsaNvBwOQAfek+TMMkkAFDuSF2tu/UsywryDRTGIZFKnXLKYW4ppmgXH6TSGQT1epQY1eRJB7fSdumB5iYeBS71IpbjqOCMlWrds4FUMOAu/Z8j1j05Ncuk8ic9WuZl00oZs/ZyM6TUKjDja16Y56madpKW/JPVqXUd4QQbxKGFQH8tlLq8KqNbB3au3PvGbuHdrbv5IEdD6z10C4rrpshFttBuTxAJNKHZaXxvGkKhRdJJK6lp+dz1GpZgiCPaaZw3cz8jqN4fCde+TkqxTIkOXtoMaC9/SFisS1LHtcn7niM8YPvwTpHajEE1CTERZ7jM/vpab5uWR+7pmmadnZCqWVWLF4ienp61PHjx1fluXWfltXn+0UmJh6lWHyFSuUQnjeJbbcTjW4lkbiO9va9WNbJbeblo08zOvYXEMlw6K0vYliz2NGzP780BO+6PX/KcyzFX35T0BM/++OBhBkPDpdg0+bf5RO7f3VZz69pmnalE0KcUEr1LPbYeRfehRATQojxs91WfrjrXyaRYWvrVh1YVpFlJejqeohIpA/H6aat7SO0tX2IWOxqyuUBJicfCy8sFuHhhzG//NfI117i0Jt/hLTLmOeZQ/Sl4Pj0c8seV2ty6zkfNwSYBrQ68N3D31r282uapmlnt5RqwZuAm4E/Bb4O3AO8D/hb4L+s3tC0K12tlqVeHyWZvBHbbgYgCIoYRpR8/nlqtSw8+igMDOBu2EW+yYd4BVGzUfLcbXENFLP5l5c9pg1td5zzcSUgakBdQr5WWfbza5qmaWd33poWpdQwgBDiTqXUXQse+qwQ4kng91ZrcNqVLQjygEAIAymrFAovzu8m8v0Zsm9+kQ1vVrD6t1KzC/gKPN/CdU2QgsXPZwYMkNIkHmld9pi2bPosk+N/fta0rxT4CsZrcE1mz7KfX9M0TTu75Wxx6BZCtCmlJgGEEG1A1+oMS9PCXURBkKdWG6VUOkAQzGLbbSilUCqgVh4i2zNGi9tM3SzgqhhVL4qSUUzPIRqZPPOEZqPR5t/K0Nv+/mWPKZm8jlnjKlrkmTXogQqbyxV9eD0Pf/zB37zAj1zTNE1bzHJCyxeAl4UQcwv19wH/18oPSdPCQtzp6W9TqRwhl3sGz5vBdTONIFPAcdrxZZ2xzAHKygdAODmitSQVM48T2JTFFmIMMjfnIggDi2+0sqX3s7ju8muSssUsXz3Rzbtjh9mWDOtX5sesoBrAy7Ngp/5nXfOkaZq2wpaz5flLQogfELbyF8B/VUq9umoj065oc6c8t7S8j+np71Cvj1OrnSAIiiSTN1Gu56jJCUwnhXMiwE72ULEn6XYcRnIRqqkqBjk8YliGxMICM04yfRsbOu6nre3Ctqnna3mOF8v82lvwsR64tRViFkgJeR+emoDvTrg8/7N/uMKfEU3TNG1ZHbCUUvuB/as0Fk0DTp7y7DhdFIsvI2UFwzAAi0DWGJgZI6gO4BHFcqucSNZ4x9QELWNNTLcPs83dzEynpC6niEZvJdP2XoKgTKUySDy+i0zmwg+3TLkpmmPNVKXgkaMG3xgN6IqAEDBSCbc7f2zH/WxpWXr/F03TNG1plnPK8xEWqWxUSm1e0RFpV7y5Atxi8WXq9SyO04WUNXx/llxlHD84QMyOERUCjy5OtMRR1QK3jjsk8x0ki+1Id4rojo9jx8MlGsOIYFlpqtUharXsBS0NQbjd/eO7Ps7TR5+m6lcpeBZ5TxEQAJByUvz5j//5Sn0qNE27TOj+XitjOTMt/9OCf0eATwEzKzscDeD4zH7ylVFS0a4rsqPqXAFuvT6KZbUihEEk0ku+VMOTHknLQ8kCdbOHurWJ1okhRivj5DZej+PGINEJs4OIV1+HWzKQz0GliohGEI7RCEUX/kNj7869HJ46zOef/TxVPzx+SyBojbbyzZ/8JglneQ3rNE27fBXrxTM6qe9o28HenXv1z4oL8LY64goh/kUp9Z4VHM+KWc2OuKslV8nynf3/kcncPqRSGELQlr6R913326SjV1YyHx7+XSYn/zvR6BYaJbTkykc5kp8k6rQTiARCSPxqFA4eYjKmuDYeI2IoXFyC6ChmUeKanSTHkhjKQSmfao9Nz11fxG1++8s32WKWL+/7Msdzx7mp5yY+vfvTb/s5NU27vDz88sMMTA7Q39SPIYxTzqx78PoLX6q+nJ2rI+4Fn+omhGgG+i94VNoZvvXyr3J08gdMey5CGMwUj9I6/To/OvoE//dH31rr4V1U7e0fZXb2Ser1CQzDRClJxN3AjLJxhUXevo2oPIyqD1CNzjJrVXhTuliBS7EaYxMm3aU8NTeAro2k6p1UzUniw+B+4yl48O2Hlkwiw2/epbc1a5q2uGwxy8HJg/OBBcA0TPrSfRyYOEC2mNVLRcu0nJqW5zlZ02ISBhbdEXeFHJ/Zz7GJp5jxI/heGVsOsz0OpoCIeYQHviz48NV/wCff9UtrPdSLIhbbQkfHxygUXsC22zDNOKaZYEPhmwxXQBpxnhyD2RmD8Wwe0zH5aJegxW+l1TA5OjsN1Tobm1soO2OYyiVV66fNvQ4OHIBsFjJL+GGRzUI+D6nU0q7XNE1ryNfyCMR8YJljGiYGBvlaXoeWZVrOTMvC35Y+cEQpNbLC47liDU29TtmvEHdaqdYH2JwIA4sCYia8pwO+eeiXuergLLd88tcgcfmvhba37wWgXB5AqRqeN8kNG+8nN6n4bwe+QbaYxTRMCkGc22tRSjJPTU3Sb6dpy/sMlyNsrm0jqvJ0FPfQVL0qPLjCMMIgcpYQki1myU+Pkvruk2QOjYTXSwk7dsDevVfE517TtLcv5aaQSKSSpwSXQAZIJCk3tYajuzQtJ7RsUUp9ZeEdQoh/e/p92gUy4ggBY/nDvCMZBpaqhJQFURPaXbi/B37/2H/G+uNv8kefefyyT+hzhybWalmCII9ppnDdDEF8kL989evErBie9KjHIxyVdW5xyyijiCULRDdaqIpPWRVpkVFMGaFsj2F6Nq6U4czJaU4pmHv1NeTkJN0dm7gzspUulSYzMACPPQYP6nVoTdPOL5PIsKNtBwOTA/Sl+zAN85Salsv9Z/hqWMqBiXN+bon3aRfgqo49lFQbrWaOmAllCWkLmuywDNUQkLbhngxUk/v5zW/9b2s95IvGdTPEYlvntyl//cDXKdVL9Db1kklmqAQ1kk2z5C1J3YRyMkCkBcmmKrXkPnLuINnk04wmn+K4+feM3jyF33bmbMlcYOkXLXRN1hhNCv66+jy/mf8H/rD6rzzcOUrxwMvhkpGmadoS7N25lx1tOxjKDXFk5sh8YHlgx4U1uLzSnXemRQhxE/BjQJsQ4jMLHkoDzkoMQggRITw1eidQBrLAzyqlhoQQHcBfAluAWuP+p1bi/a4nmUSGqzb+W15949exjSoJMwwpUoU31wBLQLsDH9oAX37rX9g/tp/rOq+sLdHZYpYT+ROk3BRCCCZKEzTZkt4Y7M/b3NYSEDWqVD2PRCSJ8DyCUg6/Mk56thvVtZ3yNUkmJx87pclctphl38g+2mJtFHJTvG6MMaYkm802plWZZhFjQI7zmDXNg+dYWtI0TVso4SR46IaHdJ+WFbKU5aENwE1AHLh5wf154KEVHMufAf+slFJCiJ9rvP4+4HeBZ5VSHxBC3Aw8KoTYolTjwJnLyCeu/2n+a2mE4ek/IG2DL0ECETOcaVEKhAEdEbi+eYzXx16/4kJLvpYn4SToTnYzNDtErppjU6odQ1QwqBIxFZ4yMVG4loUUdUy7h3qzg/+OH8Nq7iGiAkqlA/NN5or1Ig+//DDPjzxPc6SZSjnHpMiyk15Mw8CQ4BHQJ5o5wFtkrerb6PKiadqVKJPI6LCyAs67PKSU+gel1E8Bn1RK/dSC2+eUUi+sxCCUUlWl1D+pk01jngXmOu1+DPhS47rngTHg9pV4v+tNwknwa+/6fYbrm6gF4WF8UTP8IkkFlQBqQbhctCsNj+z7wloP+aKbK2y7IXMDLdEW8rU8b81mqfs1tsR9HENS8hVKSWrlSQKvQKn6BoXZN8kdfBa/PI7nTSJludFkLlwWOpY/FrbojzaTTLRSdhTHy6PIQCEVRJSFOTWN0dZOPh1Z48+Cpmnalem8oUUI8c6T/xT3nX5bpXF9FvimEKIVMJRSEwseGwJ6FxnnLwohjs/disXiKg1t9f36vf/CieBqqkE4w+I1Tg+uBuAYUJHhy2Ozz7F/7Mo6CmqusG04N4wpTEpeifGaz3hN0emCwiJuSSw8fBlQU1Ay6hQiRQ4bT/Pagf+X6Znvkc8/x9TUPzOSG+Tg5EF6073E7TijhVEswyKSbmPGUYyWR+kuG6SnSwRdGeS2rbriX9M0bY0sZXnoIeBp4JcXeUwB/7SSAxJC/AZwNfCzQJQzzzsSi72dUurzwOfnXu/p6bnwVr9rrDu9hd954BC/+FWTuzskNRnOtFhGuFw0UQuvCyT8076vcd19V9YS0d6de3nq6FMMzgxiYKBQvFF02NMCEVOScgI8X6GEQUCAawqUBCvuMzVTgNwI/Z17qNdHmBj9Kq+Ov4ovfaSUTJQnmChPIAyTqYhkQ6afLa27CdIZhtTMRan4n55+grHcy3iilY629+spZU3TtIbzhhal1KcbL9+92oMRQvwS8BHgvUqpMlAWQiCEaF8w27IJOLraY1kP7OaPM1p9hKQFngzTW94Lt0J7EkarUC3n1nqYF12xXqQt1sbd/XfzzLFnyNVyFGsjnKgqRqs+e5oEEeGjhEJKkwAwMBAioCnhM1EOaKMLK4gwmP0W8SBPIr6NOgk64h0czR0lX8vjWi4zssQ/Tj9LS6WF+7ffv6oV/+XyIK+9/pPkS68TBAFKCN4cbKXc9PN8eNd/0OeUaJp2xVtOR9znlFJ7znffhRJC/CLwk4SBZXbBQ18H/gPwfzQKcTPAZbd7aDG/d99fcf+ffpU72hS2Gc6seI3dRE9NQaVucmf/XWs9zIturstkZ6KThJugM9HJ6xMeeVmm2y1QVgoVFAmUgWMKhBfDNsNTmP26waBUZIefpN2cQfjHuDEVJ+AlakYXw94m6kGd8fI4d/Tewc3dN1Oql5ioTBCxIqsWHHy/yCuvvI9KbQghJbZpAhKLUczZP+CxgYw+p0TTtCvecvq0nBJwhBAmsCI/wYUQPYRHAjQB/yqEeFkI8aPGw78K3CaEeBN4GPjU5bhz6Gz+3Z1f47sTMFiE41U4UYHvTcD3x+AGt4+7b9y71kO86OaKcZNOkq5EFyWvRHOkmedno7yRh7JfxbEUgfApeQLfiyKEJKgazAiDGhU2ReukbA8fl0KQYqJuEVFjpIP9jJfGaYu2sbVlK+lImu5UN9e2Xzt/VshqGBr6LWq1YaSUKMMmXAhUGEhS5hSRyf+T49NXVv2Spmna6ZbSp+WXgV8B0kKI8QUPxYCvrsQglFLHOXutyhjh1ucr0gd37GV74hXuf+ROSiJH1Qffi3BjpJevfPLraz28NbGwy+T1ndcD4AUeg9ODDBUURbWBj3ZNgawSSElZ5XFKDnXbpmxU6HRaiYkyvjIoygiO04TnlWlObKIlOsZgpRnPSJKOpuff52qeFVKrZcnlfghIFBYmPuFioGCupCslxjlx9Lfpafnair5vTdO0S8lSlof+jHCJ5k8Ii2Pn5JVSM6syKu0UWzZex2u/PssT+x7lyNgb9HduuyJnWBba+/+3d/fBdd13ncffX0nWg21ZSfxQy1USK88xidOU0m2XJgF2+gfLdNohZrow7SQs7LRDswuTYfiHDjC0w7DbnSzDlllaCDUP3U5psjt0SztlWCgNNJ3Qpo3dxq5rJ07ixKJxnFi2ZUmW9OWPe2RfKbIsP+ge/aT3a+aO7nm4937Pd46tj8753XO27ThzBdv+tf1EBFM5xVvf+Fb6VnWyOv6GtXmYY2MjnBgf5fSqK+loH6NjAjb1rGJq6iin4w0MT63h6OirtEcb7W2rWN3eS+ZRNq+9kb6us6FlMe8VMjk5TMT0P8XGEZbZGX6MHnJ8P8eP76K3d2UNvJakaQsZiHsMOAb85OKXo/ms9KDSbPZVJl8++TJf2PcFBq8cpPf012ifXMfpqTZ6Vx8n2ofZ1NNJR9ft/K/9r9E9PMK2Hjhw4kUmq/c7dfoU7dHGDX29DF51G6s6r2JyarIl9wppb19HT88gw8OP09423sgtbTB9lCWnYFXnFjo7uhgfPwwYWiStTBcyEPd64PeAO4AzV9fKzE2LUJe0INNXmZwe5xKTr9I5OcREbGKiYzNMjTAcP+C2Tffw9NDXOXLqBMdOt7Nqqo31q8Y5Mt5BW1sHN155Peu7JujuvonfvfsTZ2+cSBtTTC3qvUK6ujZzxRX3cOzY1xgdfYaptglyKs+cHJpsu4qBdTcycfoFOjv7F6UGSSrBhdzl+Y+BP6RxpdqfAv4zjQu9SbU7c9G5I//Eld0JbW1MZfLK6Ek2rN7KPw/t5bsvfZNjp1Zz4PgY3z9ykndu7uL63gSmuH7daqLzOr47chUnxk+0/F4hGzfuYGLiGIcPP8zIyPeg7TTQQfuqfvp672R8/AV6e7d7akjSinYh3x7qy8zPAlOZuRv4APDOxSlLunA7tu3guqtu49j4axwZOcIrp47Q37uFJPnOv+wiSHp73sAb170R2nv4yitdfH34ap6bvJVDvJ2TXT9K0s3wWOPy/pvXbuam9Te15OJuHR1rufrqX+aOO/6Gm276JH19d7N69c2s6dnC6dOH6e3dztatH1n0OiRpKbuQIy2nq5/HI+JaGvcAuvbylyRdnLWda3n/m3+ZvQfhxInvsGb19US088V9X+DKjlH2j67h+MQq2tugu72bJDk0MsIYaxg4DZw6umiDbReqq2szW7bcx5Yt93H8+C7Gxw/T2dnvERZJ4sJCyz9ExFXAx4FvAGPAlxalKukS3DDwC7z88iOMjOzh6KnX6OI1Xpq8goPjqzk5foI1nWvoWdXDa6OvceTkEYLgqaGneG3sNW7fdPuSufJsI6gYViRp2oJDS2b+WvX0f0fEYzQuBPfXi1KVdAk6OtbS338/Y2NDrDp+gGeeG2f3yIts6r2C8TzEsbFjdLV3MT4xTltbG29Y8wYIuGXDLfR29fLonke9+qwkLUEXcqTljMx8AXghIua8IJy0FHR1beaars3cvuX77Hn1L3l19FUGrxhk5PQILw6/yPo163nrG9/K9k3b6V7VTV9XH5NTk2eufOuNCiVpabmo0NKk2Dspa+XYsW0HoxOjfP57n+fAqwcA6FnVwy0bbuGua+6iu+PMN/gX9cq3kqRLs5DL+G+7lNdLdVvbuZYPvuWDvOeW93DgaCO09Hb18undn6azvXPGuot55VtJ0qVZSOiYb9zK6OUqRFps0xeim/bk4SfZc2QPW/u2tuTKt5KkS7OQy/gPtqIQqdWa71/UiivfSpIuTWQuz2EpAwMDeejQobrLUAFaeeVbSdL8IuLFzByYa5ljUrTizT5tJElami7kMv6SJEm1MbRIkqQiGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopgaJEkSUUwtEiSpCIYWiRJUhEMLZIkqQiGFkmSVARDiyRJKoKhRZIkFcHQIkmSimBokSRJRTC0SJKkIhhaJElSEQwtkiSpCIYWSZJUBEOLJEkqgqFFkiQVwdAiSZKKYGiRJElFMLRIkqQiGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopQRGiJiBsj4msRsS8inoiIbXXXpBoNDcG+fY2f0sXatQu+/OXGT0lF6Ki7gAX6BPDJzNwZETuAh4G311yTWu3ECXjkEdizB9raYGoKbr0VduyAtWvrrk6lGBqCD38Ydu8+ux/dfjt89KOweXPd1UmaR2Rm3TXMKyI2AfuADZk5EREBHAbelpkHz/W6gYGBPHToUIuqVEvs3NkILIODjV82k5Nw8CBs2wb33Vd3dSrFL/5iI7Bcdx20tzf2o2eege3b4Y/+qO7qpBUvIl7MzIG5lpVweuhq4KXMnADIRsp6Hrim1qrUWkNDsHfv2cACjV84W7fC0097qkgLs2sXfOc7ZwMLNH5ed11jmaeKpCWthNACMPtwUMxeISIejIhD048TJ060qDS1xPAwRJwNLNPa2xvzhofrqUtlOXy4sR9NB5Zp0/vR4cP11CVpQUoILS8AAxHRAVCdHrqaxtGWMzLzocwcmH6sdYzD8rJuXWPswdTUzPmTk41569bVU5fK0t/f2F8mJ2fOn96P+vvrqUvSgiz50JKZPwC+BbyvmnUvcHC+8SxahjZvbgy6ffbZs79wmse0OIBSC7F9e2PQ7TPPzNyPpse0bN9eb32S5rXkB+ICRMTNwE5gPTAM3JeZ353vNQ7EXYbm+vbQtm1w771+e0gLN9e3h7Zvh498xPArLQHzDcQtIrRcDEPLMjY01BjDsm6dv2R08Xbtaoxh6e/3CIu0hBhaJElSEUr/yrMkSZKhRZIklcHQIkmSimBokSRJRTC0SJKkIhhaJElSEQwtkiSpCIYWSZJUBEOLJEkqgqFFkiQVwdAiSZKKYGiRJElFMLRIkqQiGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopgaJEkSUUwtEiSpCIYWiRJUhEMLZIkqQiGFkmSVARDiyRJKoKhRZIkFcHQIkmSimBokSRJRTC0SJKkIhhaJElSEQwtkiSpCIYWSZJUBEOLJEkqgqFFkiQVwdAiSZKKYGiRJElFMLRIkqQiGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopgaJEkSUUwtEiSpCIYWiRJUhEMLZIkqQi1h5aI+J2I2BMRT0XEExHxE03L2iLif0bEgYjYHxG/VGetkiSpPh11FwA8BnwkM09FxB3AVyKiPzNHgfcB24CbgD7gyYj4u8zcW2O9kiSpBrUfacnML2XmqWpyN9AObKim3wv8YWZOZuZR4C+B/1BDmZIkqWa1h5ZZfh44kJmHqulrgOealh+s5kmSpBVm0U8PRcRjwK3nWHxnZr5QrffvgN8E3jlrnWx+u3k+50Hgwenpvr6+i6pXkiQtTYseWjLzrvOtExH3AJ8C3pWZ32ta9DywFfjnavraat5cn/MQ8ND09MDAQM61niRJKlPtp4ci4m7gz4F3Z+ZTsxZ/DvhARLRHxFU0xrh8ttU1SpKk+i2Fbw89DHQBn4o4c/bn/Zm5m0aY+RFgXzX/Y5m5p/UlSpKkutUeWjLzxnmWTQIfamE5kiRpiar99JAkSdJCGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopgaJEkSUUwtEiSpCIYWiRJUhEMLZIkqQiGFkmSVARDiyRJKoKhRZIkFcHQIkmSimBokSRJRTC0SJKkIhhaJElSEQwtkiSpCIYWSZJUBEOLJEkqgqFFkiQVwdAiSZKKYGiRJElFMLRIkqQiGFokSVIRDC2SJKkIhhZJklQEQ4skSSqCoUWSJBXB0CJJkopgaJEkSUUwtEiSpCJ01F2AJEla2sbGhpicHKa9fR1dXZtrq8PQIkmS5jQxcYKXX36EkZE9NE7OTLF69a1s3LiDjo61La/H00OSJGlO04Glu3uQnp5Buru3MjKyhyNHHq2lHkOLJEl6nbGxIUZG9tLdPUhEIy5EtNPdvZWTJ59mbGyo5TUZWiRJ0utMTg4DcSawTItoJ6KtWt5ahhZJkvQ67e3rgCkyp2bMz5wkc6pa3lqGFkmS9DpdXZtZvfpWRkefJXMSaASW0dGDrFmzrZZvEfntIUmSNKeNG3ecGYwb0UbmFGvWbGPDhntrqcfQIkmS5tTRsZb+/vu9ToskSSpDI6jUF1amOaZFkiQVwdAiSZKKYGiRJElFMLRIkqQiLJnQEhE/FhGTEfFA07zVEfGZiNgfEfsi4qfrrFGSJNVnSXx7KCJ6gf8KfGnWol8FxjLzhogYBB6PiL/PzFdbXqQkSarVUjnS8hDwMeDIrPnvBf4AIDOfBb4KvLu1pUmSpKWg9tASET8JXJGZj8yx+Brguabpg9U8SZK0wiz66aGIeAy49RyL7wR+F3jnPG+RzW83z+c8CDw4Pd3X13cBVUqSpKVu0UNLZt51rmUR8Q6gH3giIgA2AO+KiI2Z+ZvA88BW4OXqJdcCXzzH5zxE4zQTAAMDAznXepIkqUy1DsTNzH8ENk1PR8RO4BuZ+fFq1ueADwH3VwNx7wE+2Oo6JUlS/Wof03IeHwN6ImI/8GXgQ5l5tOaaJElSDZbEV56nZeb9s6ZP0vgGkSRJWuEic3kO/YiIMc6OhVkp1gIn6i5iibEnM9mPmezH6yde+BoAAAYYSURBVNmTmezHTK3ox8bM7JprwbINLStRRBzKzIG661hK7MlM9mMm+/F69mQm+zFT3f1Y6mNaJEmSAEOLJEkqhKFleXno/KusOPZkJvsxk/14PXsyk/2YqdZ+OKZFkiQVwSMtkiSpCIYWSZJUBENLgSLiP0bE7oiYiIgHZi3bGRGHIuLb1eNjTctWR8RnImJ/ROyLiJ9uffWL4zw9Oed2L+eeTFup+8R8IuLGiPhatc1PRMS2umtqtYg4GBF7m/aL91bzV0xvIuL3qz5kRNzWNP+cPVjO/ZmnH3PuK9Wy1vYjM30U9gDuoHHn7D8DHpi1bOfseU3LfgPYWT0fBIaAK+venhb05JzbvZx7stL3ifP05O+A+6vnO4DH666phh4cBG5byb0B7gYGZvdivh4s5/7M048595U6+uGRlgJl5lOZuQeYusCXvhf4g+o9ngW+Crz7MpdXi/P0ZL7tXrY9WaAVt/0RsQl4M/AX1axHgcGI2FpXTUvFSutNZn41Mw81z5uvB8u9P3P1Yz519MPQsjw9GBG7IuILEfGmpvnXAM81TR+s5i138233SumJ+8RZVwMvZeYEQDb+RHye5b/dc/l0dVr1jyNiI/YG5u/BSu7P7H0FauiHoWUJiojHIuLIOR5Xn+flvw7ckJnbgYeBL0XE2qblzd9xj8td+2K5xJ7A/NtdZE+mLaA3y3KfuESzr/WwUra72d2ZeQeNv5RfAf60mm9v5u/BSuzPufYVaHE/DC1LUGbelZkbzvF44TyvfTEzp6rn/xcYBm6uFj8PbG1a/dpq3pJ3KT1h/u0utifTzteb5bpPXIIXgIGI6ACIiKDxF+Ny3+4ZMvP56udp4PeAu7A3MH8PVmR/zrGvQA39MLQsMxEx0PT8bcB6YH8163PAh6plg8A9wOdbXWMN5tvuZd8T94mZMvMHwLeA91Wz7gUOZubB2opqsYhYExFXNM36WeBb9mb+/WMl9udc+wrU9G9pMUf5+licR7WDHAJOAq9Wz++slv0tsBv4NvA48ONNr1sDfJbGL6x9wI66t6VFPTnndi/nnjRt44rcJ87Tk5urXuwDvgH8UN01tXj7r6Pxy2ZXtW/8FbB1pfWGxiD0Q8AEjW/O7T9fD5Zzf+bqx3z7Sh398DL+kiSpCJ4ekiRJRTC0SJKkIhhaJElSEQwtkiSpCIYWSZJUBEOLpMuiuhPsbedf85I/59sR0VM9/5Xq/icX8z73R8Qjl7c6SYvJ0CKpKJn5psw8VU3+CnBRoUVSeQwtkhZNRLwlIh6vbtb4RET8aDV/a3VvpN+OiG9GxP6I+PdNr7s3IvZGxLci4sMRkdP3S5p+HhG/AWwBHqmOvrwpIn4rIv570/s8EBE7q+edEfGJiNgXEX8P/JtZtf5qVeOTEfHFBd7TSlILGVokLYqI6AT+D/Bb2bhZ44M0AsaaapX1wDcz84eBB4D/Ub1uE/BJ4F2ZeSdwYq73z8zfBl6icRXfN2Xmt89T0geAQeCHgJ8CfqSp1p8DbgLenplvBj4DfPzCt1rSYuqouwBJy9bNwHhmfhkgM/8xIn4AbAcOAycz86+qdR8Hrq+evw14MjO/X01/iirQXKIfB/40Gzd9Ox0RfwG8o1r2HuAtwDcb93yjHZi8DJ8p6TIytEhaLMHrb1tP07zRpnmTNILCfK9biImm9wHonlXPuQTw0cz8k4v8XEkt4OkhSYtlL9AVET8BEBH/lsag2d3ned3XgR+OiBuq6fvmWXcY6GuaPgC8JSLaImI1jbvOTvv/wPsjoqP69tHPNS37PPBLEXFVVeuqiLjzPHVKajGPtEi6nP42Iiaapn8W+G/VOJZR4Gcy82REbDzXG2Tmv0TEB4G/johXgP8HnAZG5lj994FPRcQIcD/wKLADeBo4SOPO1j3Vup+kcWrqaRp3sn0MuLb6zD+PiPXAVyIiafzf+DCNu9tKWiK8y7OkJSciejPzePX854FfyMx3nOdlkpY5j7RIWor+S0T8DI3/o44C/6nmeiQtAR5pkSRJRXAgriRJKoKhRZIkFcHQIkmSimBokSRJRTC0SJKkIhhaJElSEQwtkiSpCP8KC8+m0A5KW5cAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plotting clusters\n",
    "plt.figure(num=None, figsize=(8, 6), dpi=80)\n",
    "if k==2:\n",
    "    x_0=weather_clusters[weather_clusters.Cluster==0][\"Longt\"]\n",
    "    y_0=weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1=pl.scatter(x_0,y_0,c='r',marker='o',alpha=0.4)\n",
    "    x_1=weather_clusters[weather_clusters.Cluster==1][\"Longt\"]\n",
    "    y_1=weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    c2=pl.scatter(x_1,y_1,c='g',marker='o',alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\",len(x_0), \", Cluster1 Size:\",len(x_1))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "\n",
    "elif k==3:\n",
    "    x_0 = weather_clusters[weather_clusters.Cluster == 0][\"Longt\"]\n",
    "    y_0 = weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1 = pl.scatter(x_0, y_0, c='r', marker='o', alpha=0.4)\n",
    "    x_1 = weather_clusters[weather_clusters.Cluster == 1][\"Longt\"]\n",
    "    y_1 = weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    c2 = pl.scatter(x_1, y_1, c='g', marker='o', alpha=0.4)\n",
    "    # for sensor fault visibility in figure\n",
    "    # c2 = pl.scatter(x_1, y_1, c='b', marker='x', alpha=1, s=300, linewidths=4, zorder=10)\n",
    "    x_2=weather_clusters[weather_clusters.Cluster==2][\"Longt\"]\n",
    "    y_2=weather_clusters[weather_clusters.Cluster == 2][\"Lat\"]\n",
    "    #c3 = pl.scatter(x_2, y_2, c='b', marker='x', alpha=1, s=300, linewidths=4, zorder=10)\n",
    "    c3=pl.scatter(x_2,y_2,c='b',marker='o', alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\", len(x_0), \", Cluster1 Size:\", len(x_1), \", Cluster2 Size:\", len(x_2))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "    max_c2 = max(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    min_c2 = min(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    print(\"max_c2:\", max_c2, \" min_c2:\", min_c2, \"Color:B\")\n",
    "elif k==4:\n",
    "    x_0 = weather_clusters[weather_clusters.Cluster == 0][\"Longt\"]\n",
    "    y_0 = weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1 = pl.scatter(x_0, y_0, c='r', marker='o', alpha=0.4)\n",
    "    x_1 = weather_clusters[weather_clusters.Cluster == 1][\"Longt\"]\n",
    "    y_1 = weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    #c2 = pl.scatter(x_1, y_1, c='g', marker='x', alpha=0.8, s=169, linewidths=3, zorder=10)\n",
    "    c2 = pl.scatter(x_1, y_1, c='g', marker='o', alpha=0.4)\n",
    "    x_2 = weather_clusters[weather_clusters.Cluster == 2][\"Longt\"]\n",
    "    y_2 = weather_clusters[weather_clusters.Cluster == 2][\"Lat\"]\n",
    "    c3 = pl.scatter(x_2, y_2, c='b', marker='o', alpha=0.4)\n",
    "    x_3=weather_clusters[weather_clusters.Cluster==3][\"Longt\"]\n",
    "    y_3=weather_clusters[weather_clusters.Cluster == 3][\"Lat\"]\n",
    "    c3=pl.scatter(x_3,y_3,c='y',marker='o', alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\", len(x_0), \", Cluster1 Size:\", len(x_1), \", Cluster2 Size:\", len(x_2), \", Cluster3 Size:\", len(x_3))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "    max_c2 = max(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    min_c2 = min(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    print(\"max_c2:\", max_c2, \" min_c2:\", min_c2, \"Color:B\")\n",
    "    max_c3 = max(weather_clusters[weather_clusters.Cluster == 3][title])\n",
    "    min_c3 = min(weather_clusters[weather_clusters.Cluster == 3][title])\n",
    "    print(\"max_c3:\", max_c3, \" min_c3:\", min_c3, \"Color:Y\")\n",
    "\n",
    "pl.xlabel('Longitude')\n",
    "pl.ylabel('Latitude')\n",
    "pl.title(title)\n",
    "pl.savefig(\"plot_output.png\")\n",
    "pl.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Graph Laplacian"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "ename": "AttributeError",
     "evalue": "'SpectralClusteringg' object has no attribute '_validate_data'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-27-6e7da53b3de3>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m     26\u001b[0m \u001b[0;31m#table3[\"Cluster\"]=cluster.fit_predict(table3[table3.columns[3:]])\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     27\u001b[0m \u001b[0;31m#table3\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 28\u001b[0;31m \u001b[0mcluster\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfit\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtable3\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mtable3\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcolumns\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     29\u001b[0m \u001b[0mtable3\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m\"Cluster\"\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mcluster\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_labels\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     30\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m<ipython-input-24-d4c77ceb39b1>\u001b[0m in \u001b[0;36mfit\u001b[0;34m(self, X, y)\u001b[0m\n\u001b[1;32m    171\u001b[0m         \u001b[0mself\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    172\u001b[0m         \"\"\"\n\u001b[0;32m--> 173\u001b[0;31m         X = self._validate_data(X, accept_sparse=['csr', 'csc', 'coo'],\n\u001b[0m\u001b[1;32m    174\u001b[0m                                 dtype=np.float64, ensure_min_samples=2)\n\u001b[1;32m    175\u001b[0m         allow_squared = self.affinity in [\"precomputed\",\n",
      "\u001b[0;31mAttributeError\u001b[0m: 'SpectralClusteringg' object has no attribute '_validate_data'"
     ]
    }
   ],
   "source": [
    "df = pd.read_csv('dataset/28_1800.csv')\n",
    "table=df.pivot_table(index=[\"Name\"])\n",
    "table=table.reset_index()\n",
    "# features order:\n",
    "# 0.Name, 1.AirTemp., 2.Alt, 3.DewPoint, 4.Lat, 5.Longt, 6.Pres., 7.R.Humidity, 8.Visib., 9.WindDir., 10.WindGust, 11.WindSpeed\n",
    "\n",
    "# change the order of features\n",
    "table2=table.iloc[:,[1,2,3,6,7,8,9,10,11,0,4,5]]\n",
    "#print(table2.tail())\n",
    "# specific features - new table is created\n",
    "# enter the fourth column to address feature e.g: 1. AirTemp, 7. Relative Humidity\n",
    "# e.g: table3=table.iloc[:,[0,4,5,1]]\n",
    "table3=table.iloc[:,[0,4,5,7]]\n",
    "table3=table3.dropna()\n",
    "\n",
    "# for debug print tail of table3\n",
    "#print(table3.tail())\n",
    "#table3\n",
    "title=\"RelativeHumidity\"\n",
    "\n",
    "\n",
    "# k=clustering number\n",
    "k=4\n",
    "#cluster=KMeans(n_clusters=k)\n",
    "cluster=SpectralClusteringg(n_clusters=k,assign_labels=\"discretize\",random_state=0)\n",
    "#table3[\"Cluster\"]=cluster.fit_predict(table3[table3.columns[3:]])\n",
    "#table3\n",
    "cluster.fit(table3[table3.columns[3:]])\n",
    "table3[\"Cluster\"]=cluster._labels\n",
    "\n",
    "\n",
    "# centroids = clusters' center points\n",
    "##centroids = cluster.cluster_centers_\n",
    "weather_clusters=table3[[\"Name\",\"Lat\",\"Longt\",\"Cluster\",title]]\n",
    "#print(weather_clusters)\n",
    "weather_clusters.to_csv('ClusteredData.csv', index=False)\n",
    "\n",
    "\n",
    "# informative\n",
    "#print(\"centroids\")\n",
    "#print(table3.columns)\n",
    "#print(centroids)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Cluster0 Size: 2678 , Cluster1 Size: 2846 , Cluster2 Size: 3 , Cluster3 Size: 1717\n",
      "max_c0: 34.5  min_c0: 0.0 Color:R\n",
      "max_c1: 148.2  min_c1: 34.666666666666664 Color:G\n",
      "max_c2: 655.0  min_c2: 655.0 Color:B\n",
      "max_c3: 100.0  min_c3: 70.75 Color:Y\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi0AAAGsCAYAAAAR7ZeSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAMTQAADE0B0s6tTgAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeXSc13nn+e99l6q3duwocBEBUhIFiot227Es20m8TDJjJxaTzMRJ5CRnxj2Z7mzTk0lmnG73sacnPSeT7umkc2Y6cUdJnDi2pcRx0ondcbzKtiyJlERSBEmRBEiCQBE7an2r3uXOH28VCIAACZDg/nzOwRFR60WRR+8P9z73uUprjRBCCCHErc642QMQQgghhFgLCS1CCCGEuC1IaBFCCCHEbUFCixBCCCFuCxJahBBCCHFbkNAihBBCiNuChBYhhBBC3BYktAghllBKfUQpNXqNr/EupZRWSlkbNa6bSSn1daXUJy9z/w8qpfSi7z+ulHrhxoxOiLuHhBYh7jDNC6xuflWUUq8ppX7sOr7fJ5VSX19283eAPq21v0HvsWKQuoHh6EPAb63j8b8NfKD1jVLqBaXUxzd6UELcbSS0CHFn+ndAH7Ab+HPgM0qpfTfqzbXWDa114Ua93/WmtZ7RWpfX8fiy1nrmeo5JiLuRhBYh7kwVrXVBaz2stf6/gHngXa07lVK/qJQ6rZSqKqVeVkq9a7UXUkp9QCn1olKqpJQaU0r9vlIq1bzvI8D/Drxz0exO/+IZEKXUw0qpQCnVt+x1/1Yp9btXM6bLjPWSZRml1LNKqU8v+n5EKfW/KKX+svleR5VSjyul9iilvqeUKiul/rNSqmPRc5YsDzUf+4pSym2+X/9q41BKPQu8HfiXzc9kRCn11uZz25c971syIyPE6iS0CHEHU0oZSqkPAe1Ao3nbzwG/BPwC0UzMnwB/p5TqX+VlHOD/APYB/y3wbuBfNu/7LNGszneJZnb6gHOLn6y1fhU4CexfNK424D3N51/NmK7VPwf+EngIOA78KfD/AL8GPAncTxTGLqGUMpvPPQs82nzeJy7zXr8EvAT830Sfz+Na6xeBYeDHF73uAFG4+dNr+LmEuKNJaBHizvRrSqkyUAeeB0aBzzfv+xjwS1rrL2mtT2utfxd4AfiplV5Ia/05rfXfNB/7TeBf0QwgWusaUAEazZmdgtY6WOFlPgv8xKLvfwSYAL69jjFtas6CLHwBf7/Oz6Xlea31p7XWJ4B/A+wEfl9r/Q2t9WvAp1g0M7XMe4HNwM9rrd/QWn8e+OPV3khrPQ94QLn5+Uw27/pj4GcWPfSnge9orU9d5c8kxB3vjqjsF0Jc4g+AfwtsIvoN/19oraeUUmlgAPjs4t0uQJwo2FxCKbWLaKblUaIZG4v1/7/jL4CPKaW2aq3PEc0wfE5rrdcxpgvAO5a97luAT7N+h5e9LsAby27rXuW5O4GTWuvZRbe9dBVj+BPgk0qpHc2g8lNEBbxCiFVIaBHizjSrtT4JnFRKfRj4tlJqD9AKBT/J0os0QGmV1/oicAj4MNHsyFPAf1zPYLTWR5VSbwA/ppT6I+AHubjElFrjmILmz7RAKbVl2eNDQC27zebiz93iLR7eKretNhOtVni9ddNajymlvgL8tFLqy8BW4HPX+rpC3MkktAhxh9Nan2huSf4Y8M+AAnCP1vqvr/RcpVQXsAPY31w2QSn148se5gHmGobyF0RLRPPAOa31y83bJ9YzpiuYBPLLbttDFLo2ynHgPqVUm9Z6rnnb41d4zmqf0R8RzWL1Al9c9HpCiBVITYsQd4ffA36e6IL+r4FPKKV+Vim1Qyn1mFLq15VS37/C82abX/+9Umq7UuongI8ue8wZYKdS6gGlVJdSarX/r3yW6OL+y80/A6C11usc0+V8C9iulPoflVL3KaX+Nct29myALwPjwB8qpXYppfYDz1zhOWeAtyqlNi/bMfQFoJPo7+ZPNnicQtxxJLQIcRfQWn8DOAH882aR6681v4aAvwGeAM6v8LyAaFnovURLNx8F/sWyhz1HVNPxMtFMxz2rjOEkcJBod9Bnl9235jFd4ec8AvwK8JvN8RjAX63nNdbwHgFRs7kB4FXgV4GPX+Fpv00UTk43n9N6rTrRZzFLFIaEEJehol9yhBBC3AxKqS8SFfb+6s0eixC3OqlpEUKIm0AplSMqan4/Ud2NEOIKJLQIIcTN8dfAY8DHtdbHb/ZghLgdyPKQEEIIIW4LUogrhBBCiNuChBYhhBBC3Bbu2JqWeDyuu7tX68IthBBCiFvR+fPnG1rr+Er33bGhpbu7m9HRFY9SEUIIIcQtSik1udp9sjwkhBBCiNuChBYhhBBC3BYktAghhBDitiChRQghhBC3BQktQgghhLgtSGgRQgghxG1BQosQQgghbgsSWoQQQghxW5DQIoQQQojbgoQWIYQQQtwWJLQIIYQQ4rYgoUUIIYQQtwUJLUIIIYS4LUhoEUIIIcRt4ZYILUqp9ymlDiilXlVKHVFKPdO8vUcp9SWl1JvN25+82WMVd6d6vUC1eoJ6vXCzhyKEEHct62YPQCmlgD8H3q21PqSU6geOKaX+Evgt4EWt9fuVUo8Dzymldmit/Zs3YnE38f0yk5PPUa0OEWX8kGRykO7u/VhW+mYPTwgh7iq3xExLU1vzv1lgGqgDPw78BwCt9cvABUBmW8QN0wosjjNAIjGA4/RTrQ4xNfX8zR6aEELcdW76TIvWWiulfhz4S6VUBWgHPgRkAENrPbno4SPAPTd+lOJuFC0JHcNxBlAqyvdKmThOP5XKUer1AvF4/iaPUggh7h43faZFKWUBvwF8UGu9DfgB4I+bd+vlD7/M6/yqUmq09VUul6/PgMVdIwiKgFoILC1KmShlNO8XQghxo9z00AI8BGzSWn8bFpaBxoC9AEqp7kWP3QacXelFtNa/o7Xe0vpKp6XeQFwb08wCIVqHS27XOkDrsHm/EEKIG+VWCC3ngC1KqZ0ASql7gR3ACeDzwP/UvP1xIA+8cJPGKe4y8XieZHIQ1x1G6wCIAovrjpBK7ZKlISGEuMFuhZqWC0qpjxLtDAqJloB+QWt9Xin1vwJ/qpR6E2gAPy07h8SN1N29f6EYVykDrUNSqV10dT19s4cmhBB3HaX18rKRO8OWLVv06OjozR6GuEPU6wWCoIhpZmWGRQghriOl1Hmt9ZaV7rvpMy1C3A6ioCJhRQghbqZboaZFCCGEEOKKJLQIIYQQ4rYgoUUIIYQQtwUJLUIIIYS4LUhoEUIIIcRtQXYPibtSqXSIRmOcWKyPTGbvzR6OEEKINZDQIu4q9XqB4eGPUakcJppoDEml9jAw8EnpvyKEELc4WR4St7XoJOYT1OuFNT1+ePhjlMuHcZztJBI7cJztlMuHGRn5zes8UiGEENdKZlrEbcn3ywvt9VszJsnkIN3d+7GslQ/LLJUOUakcIZHYjlImEJ3YnEhsp1w+RKl0SJaKhBDiFiYzLeK21AosjjNAIjGA4/RTrQ4xNfX8qs9pNMYBtRBYWpQyUcpo3i+EEOJWJaFF3HaiJaFjOM4ASkX/hJUycZx+KpWjqy4VxWJ9QLhwYnOL1gFah837r12hXODE9AkK5bUtWQkhhFgbWR4St50gKBLNmCzN3K0Zk+j+S4tqM5m9pFJ7KJcPLywRaR1Qq50mk9l7zUtD5UaZ544+x9DUEAYGISGDXYPs37WfdGzlJSshhBBrJzMt4rZjmlmiGZNwye2tGZPo/pUNDHySdHoPrnsa1z2F60aBpb//E9c8rlZgGWgbYKB9gP5cP0NTQzw/tPqSlRBCiLWTmRZx24nH8ySTg82alv6FGRPXHSGV2nXZrcvxeJ4HHvjDDe/TUigXODZ1jIG2AYzmDJBpmPTn+jk6eZRCuUA+LVuqhRDiWkhouQ3U6wWCoIhpZqWXSFN39/6FYlylDLQOSaV20dX19JqeHwWVjdspVKwXUaiFwNJiGiYGBsV6UUKLEEJcIwktt7Cr2dZ7t7CsNH19H7llAl02niUkJNThkuAShAEhIdn46ktWQggh1kZCyy2odSGenv57Go2xhV0yWgcL23rz+Wdu9jBvCVFQufkzGPl0nsGuQYamhujP9WMaJkEYMDI/wq7uXVecZSmUCxTrRbLxrMzICCHEKiS03EIWz6wEQZVS6WUSiZ3E45tRyrlkW29rZqFVn6GUjWHEAUgkdlwy83CrzErcqfbv2n/J7qFd3bt4enD1JSvZcSSEEGsnoeUWsrhhmudNYpo5PO8CpdKr5HJvA5Zu663XW23pX8fzpvC8WQzDJh7fSiyWp6vrA/T2/tSS15ZlpusnHUvzkYc+sq5Zk8U7jgxlEITBwo6jZ/bJbJoQQiwmoeUma13gkoZLfVHDNMNwUAosq51GYwzfn8eycku29bbO0QEDrRsYhoPWHkFQROtOpqf/BtNMoLVeCENXu8x0xy5fFApQLEI2C/mN+bny6fyaPiPZcSSEEOsjoeUmWb4s4DDFvtQY+zZvIm45WFaOWKyPRmMcrSEMXbROL2zrbTQmqFSOEIttolY7itYaw3AAB9+fQymbMGwwO/t1TDNFKrV71e6xl1squtWXL646TJXL8NxzMDQEhgFhCIODsH8/pG/MzyU7joQQYn0ktNxgrbqSvzr29wzNjF38LdvPMV05xMHxA7xt69sByGQeoVh8hVrtOJ43QRCUF7b1zs9/m6gr7MUma0qphffRuo5hxAjDGkrF1t09tuW5o8/xyvlX6E51k4wlSdvpW2L54prDVCuwDAxEoSUIou+ffx6eWfvPdS0zULLjSAgh1kdCyw3i+2XGxz9FpXIE1w/w5o6yL7WTCpvROGB14Dg7mSsfY94dJOd0oJRNPN5HJvMYnZ3vX1JAe/EcHWMhkGitF95PqThh6GHbCQzDQutwSXBZS/fYUzOn+Nwbn8MPfU7NniIkpC/dx77efTd9+eJqa0Hq9QLBhVOYp18hPrA7CiwApgn9/XD0aLRkdIWloo2YgbrWHUdCCHG3kdByA/h+mZMnf5lK5TCm2UbVrxLXZeLheZRvULKjIttq7DHc2jTl6ilM/zxauySTe+jr+7lLCmZb5+jMFQ/iaQf0PEq7aO1h2+1o7WEYMYzEw1S9Km7pCJ2ZB9fVPfbzRz/PTG2GjkQHdb9O3IozXh5HKUU+lb9pyxdXUwuypOfNzCz0HyKZKNNdfgRLO9GDTDMKMcXiFUPLRhXQXs2OIyGEuFtJaLkBxsc/Rbl8ZOGQvkBViKsRrMBAaZ+62oo2bLzQZlzfT9LpIPTPYhgOjcZ5Jiefu2SnT7lR5qXqbii9RFZNYOkaCbNBzExjmllCbXCy2s3BiWEMDLLBEJudwzzQNUibk7ti99hCucCJqROcmz/Hm9NvLiw95eI5al6NlJ3a0OWL9SyzXE0tyOKdWSrTha4OU82NM5V+lXwpCo0EQVTbkr38z7WRBbRXs+NICCHuVhJarrN6vUC1+gaW1YZSJmHogV8gYwb4/hkc4yxOMExNbWXO93hb0iGmHyGZefSyO32eO/ocQzPjDLT9AmEwhhHOMlaZYCDZxzu3voV/HD3JoXKZvnQfrxVe43zJQPslegtDvP++D/GjW5/GstKrXiyL9SKvjL9CPahjmzaO5aC1Zs6do1Qv8cP3/fCGXFyvZpllvbUg0d/BxZ1ZZHOo/Gac8fNU8uepm/PEG2kYGYFdu644y3I9CmhX23EkYUYIIS6S0HKdBUERpeIopQFNvX4W3y+SS/RSrh4nDH0gxKZKb2yALnMG1x3BNDOEoQuAbXcs2emz/Df9hrWFUi1GzjzHzNw3GHFK1ItvsC+1k68WRhkvT9GT6kHrbiYrU7xYOEFV/xlxM75qWBgtjjJVnWJrdivz9XlcvzkWw8b1XXb37N6Qz+dqllnWWwsSFRurpcXIjz6COgBq7hjB+Ekot0WB5ekrL8vciALaW33XlhBC3AwSWq4z08ximklsu5d6/SxBMI9ppmk0ClgEmLFuMLPooIJlFPAaFeYao8zM/CNKKUwzg223EY9vodEYJx7PL/lNf6IywYGxA9xjDZGP1RmtBcz5xzANm63WWXLBFLXkXlTNRXk+ZuDTnejmc0c+R1+6jz29e2hPtF8SFrzAI2ElaIQNupPdeKGHH/q4vks6TGMZ1/5P51qWWdZTCxIVG4dLi5HjDvptb0FPpzHt/w7ad6y5T8uNKKCVpnNCCHEpCS3XWTyeJ5kcJAwb+H4R1z1NEFTxvGkMw8EyY2g9jx8W8UNFGFYABSi0VmjdaDaMqzA39y0ymYfJxrPU/TovnH2Bg+MHscIiO3vKvDlvEGjNSDBCI6gzm+4lyyznRl7FnWlgKpMGLq9OzTBkXqAe1JmsTdKX7uORvkeWhIW+TB8diQ4MZVD2yigUGk1HooMgDOjL9F3zZ3MtyyzrqQVp/R1ENS39S4uRex4nnn/7usd+PQtopemcEEKsTELLDdDdvZ8wdCmXDxGGDbTWzd4pmjCsEXW0DYla7Ovml9OsaYnCDmiKxZcolQ7RHutB++cYmR4i1CGb0m00wlmKDY9ULEVnspORuRFOzQ7TZdZ5ozHFpB0nhkFSm3TOu8QTAb3pXuJmnEKpwKu8ytu2vm0hLOzt3cu+/D4OTxzmvo77mmNWjJfHeajvIfb27r3mz2UjllnW2n22u3v/QjFu9LmGVyxGvpzrWUArTeeEEGJlElpuAMtKE4Z14vE+DMPG94t43iT1+hhB4DYf5QHhomf5RG1XQqCB788wM/MlqtUTaDPDg84sm/INXp2e5o15Rb9Zw8CkXC/zZuNNNBrCgNAMmQo1pSDADqCGgfYt0iULp+ajMgk6kh2MlceYqc0sCQuf/P5P8rGvfozDE4cXZhP29u7lE+/+xIZ8LldaZpmoTPB64XX6Mn3XHJIsK01f30c2/NDItYam9ZCmc0IIsTIJLddZq6ncxMRnMM0c4KOUQRjW0boOBLQOMbxIEc22tP4MoAmCIq57AsxukkaMeGYHD9QncAyY9JJsSviMVjyCMMRQil5HM1wCt2aSCaBhKho6oGaE7JlSTJ0ZouvBx6MLo4ZTs6d48p4nFy7C+XSeP/zAH3LowiHGS+OXDQ8rzTisZRZipWWWTZlNfPX0V/n9l39/4bY9PXv45Pd/8poDQhRUbu1ZCmk6J4QQK1OLu6jeSbZs2aJHR0dv9jAYH3+W+flv47rniMW6CMM61eoQ9fo4QTBHq34l+vK4GFgUF8NM63sLsFBGkorv4hk9jJYmcAyPb09ZlNw6PU6Ih0YBI1X4zgQoz2Srn6Tm1QhDHwX8yEmb2Y4E4299ENIZ5upz/OSen+TnH/75de1OWWmXy0DbAADDc8Nr3vmyOOC0Zne2t2/HVNEF+/Tcafb27uUP/ps/uLq/iNvMSp9rq2ZGdg8JIe5kSqnzWustK90nMy3XUas/SCJxL/X6WaItz6MEQRWtGyiVQGuv+ehWKPGb32suBhaAGIZhN3u91IgbDYpelUCl8ajRF6syVgn4x/OQsqDiwWwjeqapAmb8Eh0NAwOTQiLEdSy+b9hnhglOPf0gP9zzw/ziW35x3T/jSrtcvnj8iwB8YOcH1rWNOZ/Oc+jCIY5MHFkILBDVcmxv286hC4c4dOHQhtTT3Oqk6ZwQQlzKuPJDxNVq9Qex7XZisT5c9yyue4YgKKG1h9YuUTDxgAbRX0es+WVycWko1jy12W0+t4FJSLvtkjVKJJWmUYf7m6UOo5WLgQUgMGDKgbIZUrRDlNZ8r8vlr7fVuDA1wpNfOMjPbf7hdf98K+1yKXtl/NDHD31K9RJw6c6XyxkvjaNQC4GlpVWEOl4aX/c4l4/5xPSJK47jSqJAeoJ6ffXX2Yj3yqfz3N95vwQWIYRAZlo23OJCz8X9QTKZRyiXDxMEFQwjjtYBYDV3srSeHV34DStH6M8DF4t0o8cEXJyBsUhZGTYnJ2n4HpszIbYF7+uBL4yBGywdl29A1YKequYdkwa7pw1OdRoMVE2e+eY8/M+/Af/+36+5VwmsvMvF9dyFwOH6LjlywNp3vvRl+ggJCXSwJLi0ilCvdqv1RjVrW3KGUXP5LpkcXHLMwtW8l8yoCCHElUlo2SDLL2ZBUCQe34pldeK6w1hWO5aVJRbrIQw9TNPFMJLN2ZMGIRCG0fRI0JjFMGxiKt5cPop2EC2tc/Hx/fPYGIQqoCcBhTp0x+GpLvgvF5YNUEVx5+ECPHE+xFEGuQaMpBsU7Dr5F16A3/gNeOc7Yf9+KJejgwOz2VWDzEq7XBzbwdc+CoVjOQuPXevOl55UD/dk7+H41HF2du5cKEJt1bRc7dLQRjVrW3KG0SrHLKznvaTzrRBCrJ2Elg3SupjFYn2Uy69Rr49RKh3AsnLE41tpNC4Qhi6GEQcUhpFH6wCtQwLACy200Y5JBYOAMAwJjACDOFBvvotCqThaVwEI0XhBgNaQs8EP4esV2JaC9tjSJSI0dLhgh+BEkzyYXoChNcU45B0H0ml4/XV44QXo7IxOPA5DGByMgkx66UV0pV0uaTuNbdgACxfd1Xa+LJ6V8lR64eKdT+c5Pn2c753/Hl2JLpShrmmr9UY1a7vkDCNAKRPH6V84ZmHWY13vJZ1vhRBi7SS0bIDFF7Ni8UUajQK23UUs1k2jMYFhONh2HqVs4vGtVCpvUCx+F8OIgZGgEcYIjR60ymDreTRgGCZeGJCwM+jQbzafA62jQt2AZqDQYKpoOcgxYXMSfA0pE2YXjdHQsGMWxtIwmgE79LEDCK0Y2bkahAra26NDA48fhw99KPo+CGBoCJ5/Hp659CK60pblDz7wQbTWS3YPLe4Wu9ISy8GpaU6UMgy07cRQBru6d/HS+ZfoTHbyM/t+5pqKb9fSrK3d5or9W1Y8w4gouCgVza4V66y5MZx0vhVCiPWR0LIBWhezIChFZwpZHQsXNsOwsO0uwrBEKrUb1x3BMOKYZgLfrxCEATXaMM0eYsF5IEQBAQ4+ijAoEdWytPq4RMUqOgwXynS9EGwTTAP2tcFsHYJlO9mtAI53AhqmkxDzoRKH3ZMN0pPzcN8WyGRgYiIKK43W1iMT+vvh6FEoFC5ZKrrcLpfV6jSWL7HMuzMUy99ge+wexisdOJZDzsnxxOYnGJkfoSfVc01/P5dr1qZwCUt/z+jcGKvVqLSseIYRLMyYmWaWbJw1N4aTzrdCCLE+sntoA7QuZkEQLdtcvKDp5sUshVIGbW1PEQQlarXjxONbse02LCuJqRuY4Sy+0UmIg8bEpI6Jx9LAohbes9Xw31RgNFu7GER/1sCe3KIBhtEzz+TgTDsc7YRDvTCRgKEO+OUf8Clv6gLXBa2jL+diPQrlMszNwalTq34GK+1yWem2lZZYivUyp4plJude5rXzX+Urw1/hO+e+gxd6Cxfva9FaxhqeGyYIo9DXWrJ6JFvCDsZwnAESiQEcp3+hRmW51hlGrjvcLKTm4hlGqV3E4/nLvtfy5bHFYWox6XwrhBArk9CyAVoXM8+bQOuwuZSj8bwp4vFNmGa6eZuJYSTI5d5BMrmbZPJBks492FYWIxjHDIs06KahuvFDTcxoEAUWE7CbX7GF922dUhQzIGFC3IRN8aieZWcm+i8ABtRtCI1mCa8BRgieBXMOvL5J8anc6Si0zM1Fsym5XPT9d74D//APUa3LZz4Dzz4bhZgVXG4bcKl0iOnpL1MsvsTyJZYT029S8irErQS9yQydiU4KpQIHxg5s2MV7/679DHYNMjI/wvDsMCPzI+zu3MTu9uxla1SW6+7e3wwuI7ju8EJgWXyG0UrvtdJhiusJOEIIIWR5aMN0d+8HoFYbxnVPY5oZHGcLqdReXHcEx9nOzMzfMT//bUwzTaMxgVIK08zQFjeouQH1cBqNhUuSlL0ZS88RhkWUAqVizUCkAQsMHxU2Z1maFJC04N40dMXg7wrLinGB9jhkFAQNUPNQdMBA8+3UNO+eGWLvnj3RMlEQwMGDMDYWPXFwEB58cMX6lsttAw6CMsPDH6NSOQy0ji/w6e39MEoZFN1JKvUJ8qke6v4EjdDCUAbtiXaOTR/jkU2PbMjFe6VlrKxRZHz8jy5bo7K85f9azjBaT2O463latBBC3Gmkjf8Gq1ZPMTn5eer1USwrs3CasOfNMj39N1Srp5uHJk4Thn60PGR24XnjKLMbzxvHMGJYhkOjMUnUeC5OLNaO78836yfqgE1AHRWuPA4NvDoLv3kU6n40C/NUV7SzSBNFi7EivHgB3BAG/DRv2fVe3v7Ae9h/FNIHDsOBA9GW582b4eGHoyWjIIiKdX/lVxbqW8bHn11So+J5M5TLh4nFuqhWT+N550kkdqKUidYB8/MvAh7J5CBV3+V8cZSYleBIKcc3Z7MLF2/bsPnEuz/Bw30Pr/vvYS2BoVQ6xNjY7+I492Lb7Rc/u+aSz5Ytv7Ihhypu1HiFEOJuIG38b6Bkcgfbtv36kt/EAY4d+3nAIhbL4/uTaG2gtY/vV1GqilI2SpcxlYFtZXGcfoLAJQjmMQyTICgRhj4XG875mADNJR+WhRcF7MrCf5WHvxprBpYkFGrNk4w0bMrAowa8NgK5ZBu99wzy7YkDFHft4Rd3/iSFYJ7iwCay2W7yZrPGxTSjrdDFIuTzS2pUtG4wN/c9yuVXaTQmCIIKEGCaHShlN0ONhW134ron8f0ytpFAh3WCMGRrup8fbH83rudimzaz7uy6m8mtpe/J4pkh1x2jVDpEMrmTbPYxlLKX1KjcKNfjtGghhLjTSGi5ThafJjw39218fxrH2Y5pOs1W/JOAQusGWtfROsT3JwED35/H98tYVgLLSuN5E83AsvhEaJPWTqLlgaXFNuCJdnhtLpphaQUWAK1g3I1uPx2HelcHL4+9DBpeLbzKi8mvke8cJW5MElZg0Opjv/MIaW1HW62zURhbvA24WDxIpfI6vj+LUjammWgWJ7vU6+dQyiQW60PrGpbVQcwZ4OzcMWq+i1uvkLJeIKxMkW37EGPV2auq61hL35PFu5fi8c0Ui69Qqx3H92dIpR68pEZFCCHEreGWKcRVSsWVUr+nlHpTKfWGUurTzdvvU6BnUBMAACAASURBVEp9Ryl1Qin1klJq180e69WIlhyGqdVOAhqlbKLts6CUjW13NB/n4/tFqtVjaB2QTA5i250kk3uIxfqIx7djmh1YVtsV39MAehzojEdhJWFFxblJqzWm6PZMNkk2101XsovuVDeu5/JG6TTTbTEGZqFfdfBy4wz/ofw1CiNHYNeuaGmoUMA8MwmVEp43i+ueWejqa5pxlDJQyiAMA7QO8LwpfL9MGEYlxGfnjlJtzJJJDGBZXbiBDY0Rzl/4DJt1hr21HIXTh9b8Ga+l78ny3UuG4dDW9iRdXT9KLNZHT8+HyeefuWS7sxBCiJvvVppp+S2iOYP7tdZaKdVaF/j/gP+otX5WKbUf+BTwtps1yKuRSOwgDOv4/iSxWA9KKQzDoVo9BjRQysLzZgET08xgmnEMIwkY1GpvorWBbedoNM4TBCWUshe23K4yyQKAUmCpqLHc9mS0qgNgaphpwMlSc19SoodsPItCUWlUqHgVHso/RCFRYsJOcXLqMGNGhQO6zrmeQR7rv5/9/+n/JX18mLhhkOw7zXz/YcK+1m6naAbJsnJE/WvKhGG0o0rrgCCYxbDy1GpzxOxuDKXoTnYSBha+2YM/dZQzh/6CcS9DqEMG83vY/+FPku64/KzLWvqeZA1YqUGcbXcQBPOYpoMQQohb0y0RWpRSKeBngS26WRmstR5XSvUAjwDvbT70eeD3lFL9WuuRmzLYqxSP9zUv2KWFgtRYLI/nzTa3RDeaF1KFbfeidZ1YbHNzR47Cdc82l11iGIZJGLpXektCok65D3fYKMPDACq+wsKgOx6Qi8E3Jw1GGx6lqeNkYhkMDNKxNI2gwXR9jq+2hcQy3XQZ2zHDEl3b38HQwS/z/CQ8M/ABMAy6wx780b+lFBsmzEIY1rGsHJbVgWlmqdWG0drF80oYxgyJxP00VBbKcxhKgdYYVAitTmZGxmnoOr3tXbSFWwiCgKHCYZ7/89/kmX/6B0t+vuXFq5drItfaOm2a0SdzuQZxQgghbk23RGgBdgDTwMeUUj8I1ICPA3PAmG72rm/OwJwF7gFGbs5Q1y8IiqRSu4nF+nDdM0S1KCaxWC+12ilSqT3UaiexrDbq9fP4/hxae5hmlkRiO7FYP3Nz/4hpZptLLjHCcB64WOGyEq0hayve0gFnKiY6DMnY0SnRIWArk2mdRxkKjeZC5QJtThvTtWnGSmN4oQdAd7Kbfb37CP0YqXpI75TP0XYo6BJ5clhGiq16P5x6jslHq3jmDNEcDjT8CsrqxrY66Op4B729P0MisZ3jI/8Wa/YARlgCpfBVjvlalmLtONlMHJsUAKZp0t+xnaPjhyicPkR++97LFtsuPwtppb4nyeRgs6alfyFA3oziWyGEEOtzq4QWG9gOHNVa/7pSah/wFeC/5mLtaIta/mQApdSvAr/a+j6Xy630sJsiChsW2exbSaUebB6c6GCaaXy/SGvZxHXPEobVZnGuh+dNY1lZgmAWrV1Ms50wrDSLX4PLvmdIdAbRrG8TYtEe85n34xyvKLpjaRoYdCQc3rLlIQp1m7HyGJ2JTl4ZewXXj2Zxck4OL/CYq8/x0vmXeM+97yEX2qBMDAVF7ZIn1/oh6Tv/EEY8y6T6DtXqCYq1MeohVGnngs6Rte/l6Xu2Y1lpHrz3N3mtcJBKZYh4vB8n3oVbP4eTqBPU86TCdubDGpNhGTTUlE9xZpz89r2XLbZdS9+T7u79UTHu5CuouoeO26R6Ht+w4lvZviyEENfHrRJazhBdZ/8MQGv9ulJqGNgGbFFKWVprXymlgK3A2eUvoLX+HeB3Wt9v2bLllmlA0+qY2/rt3rJyC7/dd3S8n0rlELXaCcKwBsQxDAvTzBCGdaLlok7AwLKyxOMP4HlTVKtDhKGD1pe2uNdEZw9VfCD2IJvjFerESFZHaRhdZJN9nC+eIxNP0db2CDtUloPjBzkxfYI5d452px0/9HE9l5pXI2bGqOkaPckeMB2C0KcUNJgMS2QDh7yZgyDA8iw2932UrvaP8vyh32WsOkRXZjuB0YtFmqPTwzD0PE8PPs1zR5/jjfJmVO0UieohEmaCjOFQrSTorg/yDe9NXvfOUaWBDjXKDNgy9QL7Z+69bLHt+3a874qN3SwX+r4M9VOaIAZmQxPfoWE/cA31t2vZbr0aCTpCCHFlt0Ro0VpPKaX+EXgf8HdKqW3AAPAt4FXgp4BngaeBkdutngUW/XZfHUIpY6HpXBDUMIw4icQDBEGJICji+9MEgYFppmg0LqCUBUS7j+LxPgwj1tx1ZAFptFFHN5dyFNAIoRw4zKs+AqubcugRo8zm7BYe7f9ByoHJG0aZXGYvZbsHC/i+rd9HGIYcunCIGXeGht8g0AEaTdWvEjNivDT2Eu073sc3cxNQKvG3yiE0YNDoZf+FTtK7HoJ8ntlygUPzAf1t76XRChawECzmanOMlcfY0bUbQ+2lVD3D+fk3yXbuZdfwCJ+b/SpnYjVQ0KlTVLwyyknxzZmDMJRe0yGDl+178txzMDREfGBPVJ18hZOs12ot262Xu5agI4QQd5tbIrQ0/RPgPyml/g3R2sf/0CzG/SjwrFLqfwOKwNVfVW6ildq/A4yO/jtsuwfLSpFIbKPRmMR1DXy/jGGYQBzb7kJrn2r1OPPz38M044Rh1LjNMDIk7M2UG7MEQQ2NSZ0ETnqQR7Z8EC+EamMeXT8O3lkytk8uHuNC7hFeLWXYagULtR/Dc8PU/TpKqSWH+Gk09bDOickTTFen2ZLv473OLpIXZggMGNLHeX7rXp55OlpeudwunqpX5Y3JN3ik75GF+zPJbdznbGFkfoQffd9PUn/+O1TdaZLapopHW7KdLf17mfVLHJs6RspOrekU5RUVCnDsGAwMLNpOdfmTrNdiLdutVwpRVxN0hBDibnXLhBat9WngXSvcfpzbbItzy0rn0yxuOletngAUppmkdSK0ZaWbMzE+EEfrCr4/g2W1YZrpRUtISaKNNwbQIB1rA9VDLPd+kvFu4pRxrDQJZZKJpXEtj1jsHXR2vh/TzNK7NU1l0W/4pUYJANuwqQU1AAwMFIqQENUsJQp1yDvv/X6STjvMz2O6Lv1xm6N6lgJl8qQvu4vH9V0c01l1pmTe8Lh39zsIxw+TUQnsmEMiE7XYt4IaMSPGluwWhueGL1tsu6piMdoLbixrUbSs0+96rWW79fKxXW3QEUKIu9UtE1ruJJc7QHBx07JotiXENDPEYn00GgUsqwOlkgTBeer1YbQO8bz5ha65UXDxsO0uVLMkOQxd4vE+LKud+/p/iVisd8WlqK6upxfePw5Laj8mK5NMViY5FjuGW3PR6OYeo4iBQblRJh1P0/CbpzDmcpDLYQLG7PySpZnVdvHs6d7D+fJ5Qh1SapRwPRfHdkjbaUJC+jJ9WKaFnUiTTnYtXMy11vjaxzZtfmzXj/Gts99atdj2svUh2WzU0TcMlwaXIFjS6Xe91rLdermrCTpCCHE3k9ByHbQCg2V1NPuv2FSrQ0xNPU8+f3HKf3GBbjq9j3KZZj+WaaDe7BwbI9QNFCZa1wnDGrFYH7bdRhCUMYwEppkkk3lrczkpWNNJxC2tkJGNZ0nGkjSCBgYGurlpqxVcAgJKXglv3uPg+EHePfBuHCtqxLbShXm1XTxPbn2STx/6NJ9/4/MYhoGlrCiMGDYffOCD7O3dy+ObHmd4dpipC8N0BQ4qmWAi5hEzYzy++XF2dOxgR8eOJeGk3YZS7U3+6sg3OTIztnp9SD4Pg4MUhl6muKWbrJUiTzo6BLLV6fcqXC6orTYDtN6gI8W6Qoi7nYSWDVavFyiXD+F5U5TLr9GaabHtXsLQo729sCRALC7Qjcf7cN1zWFYX1cY0fggaAyP0MA2NhdU8g8hDa72wZdo0M9h2B74/u6Q52uKlqCvJp/M81PsQXzz2RSCqY9HLdptrNH7o88LZF0hYCd418K5VL8zpWHrJTA7+DC+Nfp1nD/4OL42f5OTMSUxl0pnsJGklseP2wnP3d74D93t/yBf9IU7HfAA64+184H2/uGTrcj6dp8tJMzn5HKPVIQ5PHKFRm2Jf6gEq9qP42r6kPuTUzCk+v/UM56aOkb1wgFDBYNDB/gc+SPrpa9vyvJbt1ss/87UEHSnWFUKIiISWDRYERWq142jtY9tdRPt5NJ53Ad+fafZYuXhxXzwrUqudwvfLDE0OYYUWlhnDjPrqQejiGXHihgXECIIKhpEkDOuYZg7fn73m5mhPbXuKh/oe4pXzr1BsLN1K3W5D1jbxQpM5r8rrY19jW8rAsnPs6n5i1Qtzl5NmbuqzHDjzl8zWZuhF0RWcJ/Bh2qsz486QT+d5OPswx6aORXUc/+r/5J8c9viR+97FKacKQciOk9Mw9xpj+8aWzDS0Al9DdTBaDehMbCemL6D8VynZb1uoD3lx9EW+PvJ1/uHUP1DxKmR7s2zatI196R0MhUWeH3B4Jn1tAWB5UFvLjMhago4U6wohRERCywYLAhfPm8ZxBrjYB09hWR247mmCYOX2+/F4niAoUm1UmHWn6TMtQtJoFQA+phFSC+PETbCsHGFYxfOmUMogHu/dkJOJ+zJ9PLbpMXZ37+bzRz/PTHUGRZ2nuuGBNLTFFAnTw0BTDEPe2l7knlwb3W0aZ4WjN0/NnOJrR36dSvkor00XqPl1vKBOZ8xndxq+PhUn1CGTlUkOjB/gnso9vP6dv4ITB8lv30kek7ybpqw8nttVYqj2NYwvO4Q9XQx2DfLBe9+xcPhhsTrZrA8x8XUHsWAM05zHDeO8MfEGL46+SKleYrY2Sz6dp81po+DOo8Jxntj8xE0rfL1S0FmpWDdGmfszCU5OvkSh/D5ZKhJC3DUktGww03SwrE58f7ZZVBsVwvr+LLbdveKBfPV6gUZjnKmpv6VSPkCKEooQS0/h0wEqjsYmhothP0BPz36CoILrniGVeoDe3mc2pP18a7nilfOvcF/HfZwIT/BYtsA7uqDPgZztEzeiKFYOXGzmaU/vuKRep7Wc8Z+H/pQH40Ocq3qUGhVsZeOHIQUXtqUgNVtn3ouWRMr1MoVygT8rab46UGAwGWd/+R7S2uZT2VMcic1xb9Gi3U8T5PoZmhrCDi+ww5qkVKphGmazdFhHF3dtYODyythhzpfOk4vl2JLZQt2vU6wXGS2O0t/Wz1h5jHKjvCGFr9eyjLNaX5nFxbpKu6T9g8SCcUAxYMxRGH+WroF/KqdSCyHuChJaNphpZkkm78fzZvC8Aq2allgsv3CAYMviXUbl8hE8bwrDiIEKcY1tOOEZLKbR2kZjUdYZ8umH8P0ZtA5pb3/Xkh1BG2H/rv0AnJ47jaOqPNkF/SlIm2AaYEbnG5K2waofplTaQi73JJWJl6nP3ku8fQfPFb7EK+dfwaJOOpYmrMwQhgG+wUKlTKghacGMF6DD6Na4jlOOQXvd4B86xnnDnierLf42NUZ7YHG2w6UvGOGRcAftTju//8qzvLtzmmKQQhNt1y7WS3Qm2smYLjOuy4npE/Tn+qn5NWzLBgVJO8l8fZ66X8fAoNKorK3HyxVsxDLO8uLpxcW6Of8gsbCArzoJUVR1gBGMXlLgLYQQdyoJLRssHs+TTu+jWh0ilRpEaw+l7BVrThbvMoKARGI7ljdJ2p6n5BVR5iZs5qkY2zlW38aW7vcweN/7rrgj6Fq0lise6HyA3/7az7ItUcQxwIhKcxZKc5UGX5cpTX+X9MuzqOoFLkzNca6iOdBdoHvPWyjoGJoQQylMw8ILPDQaRdQqpdIs12kV/FqGxWu1YQ7cO0fDcynHwAkVmTBGZ8XGybQxbru8WniVN6ffZLg0z9u78+xIBMz6MYbnzzJdnWAgneC1RpKjtZA2p41HNz3Kt85+i7gZJxfPMV+fBw2NoEEQBkzWJnli8xPXNMtyrT1XLrdNfrBrkFNTL9OVGF8ILNO1KfoyW+jMPEilcpR6vSCHPQoh7ngSWq6D1Vr2d3U9vVC7kDRc6s16DM+bBECp6OTnngy4tRTj1RniJBjRD9Df/SRPDz5NPJZmrTuCrkVHsoO3b7qfNt4EAiyD5tEB0SyJoaDR8JhrnMU5W2Y4E3KgUmYi28Hrc2fY9nqJud4slTDJtsQkZ0LFdBiFnbwDZyow64EiOmHaxKTiVSg1Shh2SKhCfAVVBbO2z0ynRXd7O6nSBSYrk0xVpsg5Oc56A6TsCVRwlr44BFrTkXmAtHU/lOY4Wy7TleyizWnj3Pw5upJdhDqk5p6n4Z4il+jjic1LC4nXslV8uWvtudL69+I4A81/M8HCstv+Xfv5wuEC83MHcDHQhPRlNvFI/mGUMlHKuKTAWwgh7kQSWq6DlfqkeCrNp49crHdwmGJfaox9mzdhGg6tjrhKGdiGzcN9j1ELbcrVU/zA5n/Glva9N/RnyMazZNQEShkQRidKa6J6ltasCxo8Ql7ZVqF8IcfmMZNs1WM4n6NcKfCANUWCgGxckY81mGmEjLkGx4sh35yi+ZrRLEtAQBAGtNuQMqBiw2wdQjN6XMn08WuTOH6Jht9AKx01orMyHKpYjM3N0OO00WONcy9FNsfO0N8Z8l1vjC8c/Qx1bTDnzjFRPsMTbTV2diXY3t7Dzq572dQRFRKvtSngap9XaxnH1iUMXEIcGqQXlp5WC0P1emGhoFg1Q49SJo7Tz8Tsy8xyL+/a/kPMFEbxzW6SdpJsPDpdW+sArUOCwKVaPXHdZuCEEOJWIKHlOlrcJ+Uzrz27pN4BP8d05RAHxw/wtq1vX9QRt70ZXmxiepatXU+Sv8GBBaItzt1OkqqbIk4Ri5DWHEIrsyiA0OBIWbNptItEOkFursymnhSd3ZO0x2LUVB+T9QTl2jhZs8x0I8Y3p8EPfcBfeD/HgKe6YVvi4mufqcA3p8FtNrCtBy7JWDI6VkAr6l6dUr3EZHWSsZrLYNqnzQrQVh7faAcd0mW/yf2OFy2vZbdwf6xBV0yzueMx3nXv+5fMaGitV53tuFLNSD6dZ1fnAPMzX6TL8VHNGZEx12J35/vQpS8xukoYimZJ1EJgAaj7LgfGDzJXPsb4eY8ybTyWcdnTPofWcS5ULuCYNlZwgSAoMTHxZyu+thBC3EkktNwAK9U7YHXgODuZKx9j3h0kk3mEYvEVarXj2Hb3Qg3MtW5jvlpBUOT+7oc4OVGlUn0dgxq2cTGw+AEU3SSVEw4vTsEPhZADUIonjCx+Ek42LHoSNrl4jh3t24GQTe5Z3HgvXzlzcKEXjKlMnuoO2OZAoXYxtGxLwlPAf4lWz/B1QM2r0Z3qJgxDTsye4OTsSTztkTV8YgGMe0kebHbqrfgu467mrb157jN2EwRF8qqItu5h2p2hWJ8nG8/hOP0Uiy8BBqnUg5fMdqy1ZuSdXfBaDUZrF5e97kkqdsRepVrtXTUMtY5z0M2jDWq+y5ELR5ivz9DnZFGxe2mQ5nszZU5MHyTvuChtoFVId9zh3s7dVxW0hBDidiOh5QZYrd6hGnsMtzZNuXqKOF3E431ks4+Syz1FLNZ3U6f5TTNLzLR4aNtPcPa8S6U6hK/rC0tDgW9gN7qZHmsjNIaJhc2fLQhI1CskUj2cNUM0mpSdJNABnckulO3TmJyi6lWj91EmOStgWwIK9YuFvhoWtka3z8Js87gj13eZrEzSk+rB8z0CAkxlkoiFUYmvMnn9wiGe2PQEk9VJcrEUWWOGhHoJLEU8HKOhNSYpar5LNp5DKbN5QKVeMtsBrLlmpF4vEDRGeHzbB3igGTwSlkPCCJia+gLZ9OBlw5AZG+DlM19ktObR8H1G5ofZkUlTcfYRGNH5ThPVEt+advmJne+h3YmjA5fJ6lfw5su8LXt1QUsIIW4nElo2wJUKN1c7Y8bXNqN6N32bP0xH3Lml6hFa5yKVSi9jGw2yyXuYcUu4fhnTSqOCAFuXKW3ezOMzHcxQIJhX2PEcufw2xuzvEgYOm7ObMVDU/Qbni6OkmWBTqp2BTI7RShXLtEiZpVaJzBKt21ImzDZvMw2TMAyZrc0SEtLhdGCaJgkzxDLm0V7AhfIFhqaG2JrdSt6axgjL+Na9oAzscAozmCWtiiSaMzJaByhlAWqhrmhhDM2akcVb1VeyeIknG88t1Jw0Ghear9tY8vjlYegbUzBfg82OouIF1B3F2ariSEPx+BaYr88zUZmgzWmjEsZJG73YXCATa2O8XFiYNVrptYUQ4k4hoeVqFQr4c+NM8k2q9hiXqye40hkzN7rIdq26u/dTrb6J581i2+20J9qYb6SY8WxIQEKNsa99jKAz5LSuUbdDzmLyYul1dhdrPJCDl859F42i00nT70xR82o8lEmQ3eQxXDL56qRH1Y8WU1S4NLio5lclaH4fgklIaCjKjTIA9bCOjU1JOVyox9mU8JmsRz1S/GCOHFNMel3E7CSGUniqDc8fp9NOkzSjUOK6I2SzTyyqaelvzr5E963leITFSzyLQ49SNhCiVGzJ4xeHoUK5wND0CP1tH2BOl6joSV4uG8RjfcyUJrm/Po/ruWitQbFwUGWIg1LRUlRr1mj5awshxJ1EQst6lcvw3HMwNMTkliNUk1M4uZ2oRx9Dx1Y+zRnWf5jercCy0nR2/hAzM18iHt+MYdjkMgl6/FrU46RWYybdwZligv70o9QtTap4BFeVOFzqIes0aLeraA1bY1VMFeO1eZstSuMbXWxLT/AU8DfnFWcqmm3JaEmoVdOysDW6OUlh/P/svXl0XNd95/m5b6lXewGFwg6SAEGJiyRK1r5Z3te44zimnT1SPO2kO504y3RmMn0mp8+k+3QnPmfSWScdp+V2j5N0piM5TtLxHseWxGjfaJmgKO4kgMJeqP2td/64VQC4gqBAEqTuRwenlnfrvVslEPjit3x/EmQQElghMSOGR0Dge3jCo+E3+KYb8lAX3JA1uSGXJmeFNL0UU16OZGMWgYEgyabEJgoxl0bjEJbVcVrt0KFTjzK/8BKOFSdlJy+6rqgdmZpZfIHAKJC0U2Ri6VZt0i0EwTyWlT2nGCpXDy6lD0ORI57IkU1OU6wUAWj6TWzTpuSW2Nm1k1xciZPQyNEUvaT4HlZUIgjimGb6ooWWRqPRXGto0bJWWoLF3Zanng+Je1sRxSl46WXEffedt57gUobpbQQSiVFisd7WYEYVPXJMBxHOUgPGwwxd3SoFRNBgzre5qSPFC/MlvjWVZii9mZxZxTamaYhOIjFLzasx2nkDL4+X2BR36TRVC/RDBVXDclr30OzyXkIJARE9FUEz5uNZIYEQ2AgCQ9AIBV+fEhzzctyy7YOcai6QNZ7DMEzu6r8Ji5B4LEdHop9a7TV6rA+QcLtwUqNUo7aj7QRx4hg02Zrfxsc2re44XKwWmaxM8tSxaWT9CGnxApGEQiLPbZs+ytDQrzA//7Vz+vbAudOHt/ffzgvRC7w+/zrT9WmSdpJbem4h42QIo1ClycIa840pRp0Qr7qXmSpYVp7u7o9etQJujUajuZxo0bIWikU4cABGRgitliGcYUJXASYnoLyIyOYuWE9wvhkzGxXH6aNQ+EFmZ/+2NaDRbKUfQNibCBsJ7KiEQZNmGAIGthGjI+ZzpFblWCWgw6yyyfZY8KawENQai0xVX2FzEOK40BvAGz58Ywo6Y6qGpRYuR1gAkGCioi1ddcmMEWGaAgMD34yIwojIhISZYMEX/PErf0fGyXBnep7h2AQT0xNsyvTjNiOmSxY9x4bpeH4vGAZEEY8NzTG2NcNI7/YlC/7X5o8RXcCCf+WsodemX2O2PsuOrh3c038/jgh4sTLD4nychzf1neXbA+B5E4Rh9pzpQ9uw6c/2c8fgHXxo24fIOlnSsfRp0bp+8SqbkxY7Bj6FJXzCsNYaBRHX7c4ajea6RIuWtVAuK0tYw8CMWoZwRAjDAGFAo4nMpK+7eoLe3p/EMOJUKs8jZYAQFonEdmZLz9NX+TaZoIlEkCRCxEMWZQ/pRA+yehIRLpC2BfHIx2xECBHRmwnIx1WEJCMi7m5CeRxmJSxGUAbCVkGL7YPZcuK1lCZizoHdU4KjeQPfkMwnDDwT4oFgMDPAbFBlojzBUG6ItJ3GNgxEOMdCZZ6knUDUGiw2Zujd+mlCy2cqLHJk8kVG2I3Rd/EW/G0BkU/kCaOQrZ1bmapO8XzR4L5N9zGQ6z3t9Y7TRxCkz2lg97HtHwQ4Z/rQllUldiRL0brZ0nP4pZPkM7dg20lATf+OxXp155BGo7lu0aJlLWSzyuUsinDIkfT6qceKxL1OhIyQcfu6rCdoO/zm8x88rUtqZuYd5Iwy9dDBMhMIDAr2FHgTDKa20WecpNsOCV2P/ljIYFdEw5MkbfBrkmQMFmpQiMP7OuFLCxAaINsVuKj7gQGJAEYXUPOIAtg+FTGeivAtiGJgB+o14/PHadrQneolCmYh9CDWgWPW8KIGiSiO6UUs9hQ50PwiiaDAfNRgODFLYvppqqf6kUNbgQtb8K/03pmpzSAQmMIkn8wzUZ1gsblILp476/Xns+sXQiwJkoXqYVIWFJJ9TI3/J2q172GaGXwSBGY3cSuB5b2B6x6mFM4Qi/WTydyOYcR155BGo7mu0aJlLfT1wc6dMDYGw8N0V29nJvkC9ebriE3dyNjVNYS73Kx0+C2VnqHZPIFjpoAqYVgDJJFh0BULeBf/SCPu8/0SzLgREw24pyDZ3AlTVbASsFiRFEtQEdCbg0IFTrb7nFuiJTLACSHpg2fDXRNQSsD3e6BhgW+CGYEdKoHjGiFmZFKvL1K3qzhGhn6ngomPlB4Nv4RthQgzopR4g8R8F8nxRQKrgmOX4NCXqJy4H+64ndC2iYhIGmdb5K/03onbcSKUJ40hDPzA51T5FKEMlyz822LEX3yRrsxNhGGFKGpiGPGlOqh6/TCy8iSx+hjNsM4rB7+MHywi1poe/gAAIABJREFUSeJGPtUojsQilCZW8g6G4yksK4/nTVKtvkw2e5/uHNJoNNc1WrSslT17lopxLcOgP+rHvfkOwgcewsxcXUO4K8ns7F8hpUsisZkEAY3mOGFYxZcSpE+vI5EChuJQbcLJJkw0IGXD4Ro0mmrKc8qAtA8+YDkQq4FnQcJTIsUKYaACsQiIIO3BveNwMA/3nIJnB+GfNoNv0LLrlSSwiQWSsudh0SAKFnBlhBDgS4G0JKYBUgRMJvcSGxIMRQY1GZB2FmnsP4z7ouTkzjx3ZSu483/O5Bkt7Vkni+sVKVWq2FYn/el+TpZPUnWrTNemiWTEM+PPsLOwky8f+DJHS0fpZJxB8SodzrMMpLNYrXboWKwf2y4wOfl5fH+aRGIbc/P/QMObxYsMmuEi5cAgbzeJWRa+GOBUfYaYDNlsLGBZeVx3At+fP+c0cY1Go7le0KJlraTT8Mgjqii3XIZsFqfvrfULwnWL+P4CQjhIGeL7c8ioRihdzChcTu0AlgFOEjoSsC2jOoBEBHYNzBhUHCVaAgucCnQEUGldp1BWUZZ8E1I+1G2VInr7CZjMwHAZxroh4atoSzxU6aWK6VE3IBN0MNmosysZ4QNBBLaQiNbE6gBJ3WhiOjHyVYMwTBNm6zB6kok3quxghIF4Fk/kycU7l1I5JyY+z4uT++gJvkttvoljOgybfRzxAk6UT1GIF3Bshy2dW5isTPLVg1/iR0b6iYenSAWTBJ7PTLmD/tzNGEaCZvM4i4t7MYwYhpGjUnmRcu0AQWQgDBszapC04jQij1hUxzAb5OMZTjbqDKZzEMwThmWazcPkcg9et5E+jUaj0aLlUunrU19vQcKwTCzWSzw+SL1+gCiSRJGPaE2DbiNW3DcFpG3wIrizF/aWVMrHj+AuF741D2UPkiE0bMg0lWCxI3jHcbV2JgkP71PniwQsOCpV1FMDz4SmDU0TzCACE+a8Es/NhLy7S3UeIZSI8qX6xjcMaBBhCA+ZC+kzcjStgHxHk9xokSCY5bWZYZg5RHcyz009u3GsTl47/l+oeSaDnXdQrByl7s1j+m+wOSbo2/QubizcuOSlMlmZZJszgfB9AjNPhEHM8JHBBKVyhZiVx/fniCIXjBiGmFQmclETyzBxI6Ww4sLFMCIMAoScIx5NIMgjEjfTEbNpNA4zMPCLZDIb06hQo9Fo1gMtWjRrpu3+msu9k2bzKGFYYuXE5gthC8jEoLsXrBJIA+In4bAFt81DZMK4q4RLwoe5pKpXadhw9wT0KSNcds7A3k0Qa0ViSnEVbclImE9K6haEhMx6qnU6ZoAlgNatoe6STkAYRRhGREmOk/UGKXo+M5HLFrtBzjqMBKLG9zg6/hK5xCDSP0EqdicxcYrOTIMgsmn6CbrsMjLXQzLhEAET9SZZM6Bg16lFCXo5hM0iBj4SSRgu4IYlIiRhJJBRQIiNYUgMJEQhZlRFCAhlRIiNxMDHRsg6GTzipk0QLNDR8aAWLBqN5rpHixbNmmm7vy4u7iUevwEhTuB688DiyhracyIExARMNqB2CqJZ+MA3Ifsg1GMqDdRTh6mUEiymhOkUPHgSPj62fJ49Y7DowL5eFZUpO2BLmE/AYrzVMh3BZBNONqDHgW5HRVxCljdpAKYJQoI0JHPWKaZdgchY+BFYQR3D7EAYHXhBlaY7hYVHyt9PFKWIjDRGLI0lHbpj01jR39FobCWMArpkHseoIaUkKccxo1MIAsJIqhEFEUStwQWmIQEfiwAvkK0ZBhKrlcoShBgiohI4ROZmZFCiEDMxg5OksnfrlJBGo3lLoEWL5pLo7t5DECxSqTwLSOJON9VmBRFFq77WMWFLDuwbYLsD/QUYKcHLvXCwqyUgBHTX4MY5+MXnYPf06edIe/BLzwECXuuG9x6Dp4fgxX4lYqqOShkt+PDkLPzwoHpdM2pFXAREEgzR0i+t27gFW3OSUPoEUgksL1xAmDYGAQTTiNAjJqfxGkmiegXiKWI5F1/AbG2OI00TacTJGifJRU0838NrLuBaTSwjwmhFelq1vYpIPZaRVJO01ZYII3UoQlAPYNLPgtXLQHo7Q9kUPT0/TkfHA2/q/6VGo9FcK2jRorkkLCvNpk2/RBCUmZl5jHh8K01vjjBaWPW1oYSCA30J1eo8dSPsnoexg7B9BnxL+a4sJJU/y5mCZSX/y8vw2E54fgCO51QNTE8Zvl9oLYjgpQV4Tw8kTSVUOlqzC43zhIQEqgYnkEo0mASIsIgQNobvgCsIHRfTkUjDp1ifIu3WWfAFtTBivl5Emh2EToYOu0pdGsQNb2lq9dI4xTP13YqBkRKoeDDvq+hULUoQEGdLJk1nroeRzlF8f5ZEYnTVz1uj0WiuF7Ro0bwpNm36FVz3OIuLz0DUvGBqqI0B5GMw14SBeYjiUOiG9zvwn/e3OpeBgTLEAiiml2tZziTtwSOvwg3zKkV0KqtqYQ50La8xDTjSen3KVOfPx5SJ8fn2K1B1MM0QGiFk7YiGb5J3bVKijwWO0ZABgeUSxUO8QFAN46RtE8uK0QxDHDNOwspwrNlJzJ4lxZQqCF7lsxFCXdfHYNazqAQhnXYTVwY03TpHi3/N7HyS4fxuMvNfO2uquEaj0VyvGKsv0WjOj2Wl2bbtd8nlHkQIC+MivqMMQ0VbQkPVpfgGTHvQn4f7I7hhDm6bhIwPxzuUGFmN0QUYqEJfBfZ3K9M50Qpb1HxwpSrITVtQClRdy2oCS6LEgylUXYnEpxqFNEWTtN9DR9hLw4sxU03iujaOCCj5JpIEhmFQ9haRMmQh7MA1Ckgjvvobae3LEtAIIsbKPmPlkBk3YijukzZCQgkHF5t8c2Jhaaq4RqPRvBXQkRbNm8ay0gwO/jyzs39HEFRWf4GEeHsoYhysGMx60GXCbD+ciKCzqSItJ3OqjXk1+qqqo+iVnqVLYERKGC14cLwKHTbUQyjEVlfrEpUaciNwDCW0ZOCz6JSAElGYwIgcIuFRkwZlfLwwYt+iRIpAGQ86kqKXIZQhgTQJhEVs9beCFynhkjRhwZdI4HgdtqYkr9R89ld86iGY4lUcK80PjPp0dupZQxqN5vpHixbNupDJ7May0ngX0fkcCSUKTMBNqsiHK1X780mUiZwAyjEYXlSmcW2KaRV5ybpnp4zefgL+ahcMVWA6Ca6lhEtkwhOzyo13W0qJAVDRHvMCdS3NUImFXke1amdtlbpBQkSDIGpwoARfPmXjYbE9HzGUCIAIA9i/EHEi6CZjF3Hrk8QTFyHoUF4ybeGyPQ0lH3I2zHnw6qLymVH7D3lu/HkcSvxYz09p0aLRaK57tGi5hnDd4mkDCzcasdgm6s2Dq65zA0hZqlak7MF8EwYTcLQCb9TBTMPmEsSAwbISKNWYKrgd6wZDKuGzc0a1Pqc9dd7QUEW7dgjzMTiSB9dWx5oh/PUE3JuHoSSUA0ibYJjnTxM5JmzPtKIe1gqB0+o0MgwYTUMx8gE4UI5xuGpQDTxqASz4VWzxOruzIVuGo4vKxXqt4lwDVYvzYgneqCrR9N5eVRy8EilDipUJ/sf+v+fn73vbRVxBo9Forl20aLkGCILq0nRgzpiBs5EKMG275+yOmHMwVla1JWarELY7Accq8BdvKJ+Wpg1dddX2vHVBRVS+cKsSLCML6hMIhXr8OPDwq+q8WVeJmYEKpAPl8bLSOKbTVmLFbbU9G2K5W6etR9rb9yJVKDzrKeFinmOtRKWc/tVWdd566IFQ0Zm9MwYCSdry+edbzt+pdCaWoU7sRiqF5odwqqGOzbpwQxrGGyq1ZgL9CZj1U0zPHqdYLZ41jVqj0WiuJ7RouQZoC5Z4fAQhjKUZOLOzj9PX9/DV3t4SicTIRa2b9+Df7od+B/pMmAlgclFFSjCVGKjGYFNZRVKKaThQWBYsoATJcEkV3ba7i9p1Lc8PQMZTxbMrwygpCxZ9VePS0YrAiBVLIrksTCyhuosSVsvX5fRTLbVFS+CmLEx5Bou+5GBFsiUJVk9EIGF3Fkazqxf9tjFaJ0+YKn11eye8VoZ39cC2NHTHYGsKFgOYrAtONkwWjF5u6YxTdstatGg0musaLVo2OK5bpF4/sCRYAIQwiceHqdX247obpwAzDGurrpEo0dAMVTrohFD2+zFUq3LTAiOEXTPw4cOqK+hgVysdg+o2aloQDyDnqpqVsrNc39IWOU9saVn1t1JJoKZKu5G67bTbAxQ5bcBjm3bExWxfmNPXLGWKBHgSMlZEtwNbknCoCqMpFRFJmeevm1ntc5IS7uiE4ZRqv64Gal9TTWWSN9mUPLOY4/3bdpC0k2Sd7NovpNFoNNcQWrRscMKwDIglwdJGCBMhjNbxjSFaGo1Dq66RUhWWgtICEjU/yJTgtqIsqQDuP7Vs25911bGnhmAqvWwk21uFQkMdb9P2bfmnIVWM65kQtb7LF3yVYrGEEhT5mBIUsVZdSzuF40VKIFR8laKh5dC7MirTJmxFc2JCnTdhwc05df9EHfriLRv+NQqX9rW6HchYsK+sBJcAchYQwGjGIJEZIZ/Is6t7l46yaDSa6x4tWjY47eGEUkanCRcpQ6SMWsc3CG5z1SWhhOdbprl2CDkPZamPEi+bF+G2KfjQ4eUC276q6iR6vQBb51vzg4R67Eyd3UXUV4XdU/CNrRCc0S69bxHu7GxFP1q+LYSqwNUU0IigHqjHjVDVjgwk1GNLnC5aotYaAxUJaXf1NEPojcO7u9Ve1ypY2kihIlCGCf1xNUMJVFHwpiR40ma7vUg2V+F9N37g0i6i0Wg01xDaXG6D0x5O2GweRUrV+ytlSLN5jFRq14ZJDQHEotUF1Jy3HGnZOatmDlkRWFJZ8DuBSv2sjJ4U0+rxjlmYT6pBivMJ9TjtqeNnknVbkZzTVAbMu3CsBi/PQ8kDN1SCyWyHNiTEDbANlYIZr0PJhUag1vqR+nJD9T4WfZUC8iP1WoGqnRFSnSM8e2sXjSHU/iWqTbvHUe3XaRNCaZGJd3P38Md4W1cXjcWvv4krvTVQqdaDuG7xam9Fo9FcIjrScg3Q3b1nqRhXFeJGpFK7Ntxk38iSFzzut+YApUzwmuBEyvF2cHE5enIkr0TMyuhJ2VFrHzh1dk3L0dxyTUvbw6VpwqEu1UV0PHd6tGXBg+M11RHUl1CRFF8qHxSkiph4Er41AU/NKZHzRhXe29MSI6hoUdJS6SZDKvM521CRmLA1P8hvRWbsN/kvTMpWgbCEDgsQKhLU4dgMF+4mn9qMlOGGq2/aSFwr3XcajWZ1tGi5BrCsNP39j2x4n5Zq+PoFjy96qjW4FkJfTXUIbS5BNc5Sa/LmEtRisK9neVBiu5U5QgmVXCsKEwqIDFXI+4Vb4YUBVcNSt+FETombWHB2iuiJWeULk7OhEqiIymJdFehGEnIxeKEEh2sQb0VfGhEkUfUlkuVUTbejIiL1UEVjKgFsSqlzdtgX3zV0LiKp9mQbqp6lI2YgkcRNk57cLgr5dwAbs75pI3GtdN9pNJrV0aLlGkIJlY35S8l1izTdmQuuqYQqygFwgw0NE4amoGGp+UN2pKY6H+6EycyyaGm3Mo91qzZnUyrBcqxDdRl9cyv87XYIhIrYlG21tmkpEXMmzRCenYe78zDnqqhI29RNGJBqRTY6bfhAjxImT8/CrR1qOnUo1evGysphd1vLx6UUqG6kZqDqYN6sYFkI0gRRE0eGSBHDlwamACexk80DP43RmmW0IeubNgjXUvedRqNZHS1aNOtCo3GYxfoCmQsMN5xpql/qnxyAu1DREFJw/DAEbVdboaIq/Wc43u8Zg8dQYiSeACMBN0/A/SfgFz+oBEvag9CEZlK9JhDnrymZbKpUUaLViYOpxEvWUo699+ShK6bajWuBispkbdVuLIC8o157qAZNCWVfFd96kao5uZQ25zYS8EiSsGOUfZNyJJBWgWlRYCg3wvZ8L0Ioo5mNWt+0UbiWuu80Gs3qaNGiWRcOH30Bwzh/2amU0B2HBwpw4iTcuQ9mtkFjADYDRw60alo6VedPO8rSJu3BT74Oh7ZAZVgV7KZyML4JFnMQejCehcBQhbqZpuo4OrPUvDOmampCVGTkltTyknZLsSHgblsJkLgB8yEMJtX9BW/Zq8UxoI4SO9+dUWmvlAkDcdgqL94F90y6Cx/BC0Ncb5YeEWBYeWKZ99HTs4eeVO81Ud+0Ubimuu80Gs2qaNGiWRdmFxeIX+C7SQpl5pbz4T5TRVN+6AB8GSVcphah7CrB8u/+8dznmHkI7F4YOKq6c6QAboSRXnjpqHLBbVpghTCTaDnstoib8FABfnwTdMTUXo7XYWwRCnFVNDuUVOKqGaqCXAMVbcnYqrU5b8PmJEy0OrvbIkfQmljtATElYtrmdGttz/NCSCRuIG06BEEFx9lEd/ceksnRpTXXQn3TRqHdfadqWoYRwtTRKY3mGkaLFs26kEvOUK2e/7iUKrrh1SGRhlgMsjX46X1wogoPvAJ9h5YjLN/uvI9TqT6GakXevfA0bifUt0C8qAQLqFtnEgZugDfiIFsppkVHTXheWVTyma3w/v7TPVM6HRivwddnYFcG/GA5nRRJZT5ntLxSSiE0DBWpsQ04WFGFt31xVaez0Lp23lbzg84w0r1ohAA7vJ+ewd0XFCQbub5po3GtdN9pNJrV0aJFsy4U8lmqlWXn2HNhAOU4jFQg1XpOCuhuwNvGwFmAN+Kb+cOHfhSxpYoEXqKfvzn+dn7u8Bdx5OSSYGkTC6GjrtIykw7MJVrOuis2MZw8W7C0GUzBnR7k4+BYqubGalnzGygRE0ZqFlDQKth1TBVl6UsowfLE7PL5HuyCW3OXXtNiGXB86h8Y3rHn0k6gOYtrpftOo9GsjhYtmnUhk7kDLiBYTAENH1IxiEdguEpYNPsgdVwJFoDfe+hHESMVFuYygEAiyY9UeDT1w/y8+CMlilrC5UAe3igAAvLzMJZRpnPyjE38x1su7Ep7S6dKy1iGWidbBnNhyzAOlCfLggv7y9ATh2fmVEt0O8ICKgrzAwOq1uWSkSCthTdxAs350NEpjebaZ0M54goh/q0QQgohbm49vkEI8U9CiINCiOeEELuu9h4156a3dw+RYZ/3uERNJg6BGQfqncuCpfCEWvPtzvsQw1UW5jJYkcCKwI4E83MZvB7J8dlRmn1qptBn74M/uQv23QRP2vD4IBQzqph3pXIaSakalgshUPOH2lOb/Zahm90yjTMN1fl0sKIEieBswRI34ZEtkLtA99TFIAX0Dd735k6i0Wg01ykbRrQIIW4H7gVOrHj6T4DPSSlvBD4LPHo19qZZHdct0qBwwTXZGHx3Cp79Lrx4GIYeg75vgNUqbB3r6CMSYLXGMkdCFdOaUhACr4/tInkcvvx28PtUi/FEBZ6aZrmI5IyISlcMgovw0m+/zGi54oIqpg2kqnuxBNzaCfcXYCQN7+xuDVNs8aE+JZB4E11DAKGRZqT/Ry79BBqNRnMdsyHSQ0IIB/gj4MeBf2w91wPcDry/texx4A+FEMNSymNXY5+a8xOGZVwvTfwMj5L2hOMggq8W4VszcK8NngPv8U8P1meDIoh+JBLPFEvdPxJJyoDIn+U/BfDns2DFlAvtnH92Omglcx4sBJA8fxCodY2WcBHLT4RSRVMKcbBXvK9GqETKjjT88RF4Vze8p1f9Y3oz/iyR4TCy+X/X9RYajUZzHjaEaAF+E/gzKeVRsVx8sAmYkFIGAFJKKYQ4gbL1OHZVdqk5L6aZhSiNiypabdeGtP93Vnx4cQEQ8Owm+F4vdNXh3+xdnub8voWn+fvi28n3VqiXMpiRKmDJd1Y4NZNh8ZanmcxAMYKoPVDxTJEgT3/uaE0NSOyPXzgCcuahdjmLQI0fmPMhY6r31htXEZi35eE3E8qoLmEoA7o1axYDDAR2fBcjm3+Rnp6fWOsZNBqN5i3DVRctQoj7UAapv36Ow2dO4Dvv7wQhxK8Cv9p+nMvl1mV/movDcfpIJm6gUnsZUEWt7eyJH8J4Y3mtlMpP5bFdsH0eHn5VPd9Xhe2n/pJX7B9lsKtK1Eq1TM5kwP1LjvQq/5XoQsrgjGNxU0VkLqQmVuqcdklMzFD7zMeUGNnUPt6qewkj1UG0uVUzk221Oq+VTPrtjG79DVKpW3SERaPRaFZhI9S0vAPYARwVQhwDhoCvAzcDQ0IIC0CoEMwmTq95WUJK+TtSyqH2Vzqtp7deaYa3/hiTNZh1VToIUGkWVHTiB/qW60CEhIoD39mipjO3+aFTJwiszzI28TIHT00yNvEyHcFn2V07Qd2GUnu44kXyoT4YTSvDuPO9TLAcWTnteaHSPXFTFeo6pirMNVDpopwNyVbkJWGuvWvIsLdx151PkM+/TwsWjUajuQiueqRFSvlbwG+1H7eEy0eklK8JIX4G+EngC8DHgWO6nmXj0pHeipiKcbTbwzHVdGWH1i/8GLy7R9WI/H0RumzoFyBzUHZUlAWgvwo3z0C+/jSyA5ISUmWY8dXww4rNRedgOmPKNK4/roTHuV7WFirijMcSJU7aKaWVkRj7DHEStUI1a/0LYGj0P67xFRqNRvPW5qqLllX4OeALQoh/A5QBPUd+AzPUuZuqazGY8EgaKmWyUihsTsI/61UipuCAGcLWPESLELS6iPqqsGsRFh+EgfSyXX+tCptL8HrX6vsQIUhTGc7d2qEiJOcSLBFKaMiW6JCo64kL+M2ci0vpFmrKONsGtIGcRqPRrIUNJ1qklMMr7r8OaNOKa4iF2DY6wn3kzuGNIgRsy8FxF15bVCIhsuHgPWDHYPQxte4d3fBKEk41QURKgPRloLcXqLGqomh3E92QgVzs/MuFVL4usuXL0gzVnKErQbrn167MhTQajeY6YsOJFs21Tc/gXTjGvvMeFwJGU3CgAjekoDsGiyl4o19Z8nc9A+Eg3LUPdsSgYUEigP+xC7IFKHgw63N+JSKV10oE3JdvRVJYTvecuZcwUm64pqHqU64EVV+wreuBK3MxjUajuY7YCIW4muuIR277DIlVfvlvScFtHdCbUB4qM1IFUKZuheIHAamiIFkXemtQjsHJrBIjHZKz5g+Bql8ZiqshiAbQaatCWQlLhSrRGa8LJRSb8GoZyt5Z3dKXBQmUza1kE/2X+UoajUZz/aEjLZp1pW/rbqafhcHk+dcYQg0VnPdVx40Vh7iAhgczgxArQWwezPbU5nirKBaIVdVgw2ZLbsdNeKgAW5It0SGgWFdtyLmYirjYQgkUibq1Wo/rATw7D8MpqIXQwZURLsn0boY6d1/mq2g0Gs31h460aNadzuIdq66Jm9AI1PDBoRREJpgOhElwCzB313JtStaF7gRMlmE8hEQShhIqutIWLMVm66sB9xVgdweUfGUsF0rVutz+ClHDD49WVUrIFqojyI8uv2BBwp2bPnq5r6LRaDTXJTrSoll3PvWez/P1qVuxL5AmEkDKWq4jcRvg1JW/i1UGrwAzN4LvQkJC4gSMhXDPqEovSZQ3Sn8cXphbblVOmKo+JQSqgTK4O1yDTQnVtRREymzuxQX44nElcHZkVcRntbTWeiAEGLJ++S+k0Wg056BYLVJ2y2SdLH3pa88fSosWzfqzezf+7zvYN7vnXTLVAC9S0Y1AQs6FMA52Haw5mMnBUwFMN8CswY7nYecnoSOpoikS5VZ7Yxq2ZWB/WZ3XMZYFzHhdzTgKUamgLkdNav7SKVUIDPDX45A04Ec2v7lBhxeLYYDvz13+C2k0Gs0Kql6Vx/Y/xtjsGAYGERE7CzvZs2sP6di1Y8aq00Oay8Kd3X9AFJ3/eG8SPF+lZhwgCVh1cKbgUCdMpuGZTvhWCr7WC3/xAPTmVAqoLUqakapF6XJWRGxWpHhqoRIzL5Xg1UUVXfnjI8uCpc2LpfObz10OTPMCBT8ajUZzGWgLlpGOEUY6RxjODTM2O8bjY49f7a2tCS1aNJeFvh/7NAefNs5rnW8J2NGpBIgpWwW2U9A0YaofFn1w5yHjQtIFNw3+iigKqOjJvKtM5NoW+o3W2IBIQqPlHtcIlcHc6xVY8E/fR9yEH9v85qYzXzQGgInjbL4CF9NoNBpFsVrkwOwBRjpGMIT6YWkaJsO5YfbP7KdYLV7lHV48Oj2kuWwkt32QYv0r9J8nsCCEqiOZD6GyGfqT4FahIeDVCThkQSOuxgF0SpCGioasFC6HqtBhq24hx1Ln3DsDCBhILE+aPl6DJ2ZbLzKg04K8DT8zArd3XN7P4fT3nCWT0Z1DGo3mylF2ywjEkmBpYxomBgZlt3zN1Ldo0aK5bGy+/YcoHfnKknA4F0Ko2pSZBhzKQWoWvjsNr06pyAoCnADKERyvq1EA7RRRyoLBBHy/DC+XlKCZcJejKZ22WlMLlp+LG/BQtzK4u7NTFeieb2/rjQFkMrdgmtkrc0GNRqMBsk6WiIhIRqcJlzAKiYjIOtfOzySdHtJcNt69/dOUfHNVUSCEamk+7sLv1+BrUzBnqAGJZqsuJhTwxDScrMKmFNzXpez+d2ThgQL8xGZ4qAfuyithAkqonGqcLmJ+aAC2p6ErpiI0V0qwqH9pNtnsvXqis0ajuaL0pfvYWdjJ0dJRwigElGA5tniMXd27rpkoC2jRornMJAd+8ywn2nMhBNyWg8MSYoESLE0LQkNFXCIBogkvHIPJOsy5MNlQU6NP1VVNStZSni0PFU4/d9yA9/fCT25SQmdnBt7WoepgBFfmH4EBxOMjDA//xhW4mkaj0ZzOnl172FnYybHFYxxdOLokWD6+8+NXe2trQqeHNJeV20f28Jf/8H/xQMFbtTsnY0MyBo0aFGowk4ZMDNIWGDUoN6HXhEJcueluTqkoigTKvjKbO1ZXPi6d9nKE5aFuGE1AIaaEjWOqtJG5xmnOl4phAOTYvfsrWNa101qo0WhJythPAAAgAElEQVSuH9KxNI/c9oj2adFoLkTWyXLS+QSz7p/T7Vx4rRBwUwc8UYItLty3BQYyYIXKQ2W6BCdqMJxUoqPbUWmeQMJUU50jZqji25SlREunDVsSkLWVKPIjFXkJI7DMK2Db3wrj3HHHEySTo5fzShqNRrMqfem+a1KstNHpIc1lpS/dx2jXzXxxfCuN4MJrBXBPXs0L2j4MvVlVdDsDlCUMdsF7tqpBiLVARWDyDvTF4eYcjCSVA64Q6jgo8eIYapBiNVDdSknryrU4G0Bn54d1x5BGo9GsA1q0aC47n9j1CaLEdr47N0rzAsJFolI32RhkuqDUhJ1ZuC0P2zthIAejXTDrwmi6NUcogmaoBE/MgPsLqr25nRqq+UqodMVgR0b5sphCRW4Ey1GW9kBFeRH1NxeD0RIshcIn2bHj0fU5qUaj0bzF0aJFc9kZzY/yyZs+SdO6l2NzsfOuE8D3FiFnqcLbG1OQScBCAHOhanu2LBVhsVumcRI1U0iiHHATJpxYMdpnwVNGdh22uoApWvN/aIkUln1fogu0Zl8K27b9Z26++f/T3UIajUazTuiaFs0VYc+uPQCMfXuekcJXcc4xnNAN4S+OwkeLkBuCVApKrciMgYqGBJGqT/FCWAxaXTmmitBMN8E2Ib3i3J0xNe15wYNNSRVhCaLlrqG2gDHWuyjX6Ka7W09z1mg0mvVEixbNFWGpcr3/7Tz5t/8Khr+OYy471jZD+M1X4GdfVNGQ+ToMdC2LFoESJrOuiqa0hxsWHPU4kmpwYiTBk6jCGEPVvAwllmcStb9CqQp45z11jnZ6ab3o6PqYjrBoNBrNOiPkeiXxNxhDQ0Py1KlTV3sbmvNQPXmYzz33YQL3ICcX4Y2X4vzYTB8PHmzwuZ11vj9Q4cGHlIBpd/jMeXCsBrd2Qp8DvYnl/KYQqghXAuMN+JVXoeTBx4bgPT2qMPfGtCrMFWK5dqUZqXVdMRWxWQ98w+Y9D87r9maNRqO5BIQQ41LKoXMd05EWzVUhvWmUX930+rJnwE+1PAOqVXb+Pz/LxAt/xd5iwNYOqHgqEuP70JWA70yD34RfukkJDUOowIobQcWHnjj80jb4/UNK3ABsTamITLtmRbTSRDEB/Yn1TQ2lkndrwaLRaDSXAV2Iq7mq9KX7uLHrxmXfgHSaPb/8Od5z00cYP2ZzaLHVsmzBoAPTFdh3CmwXZptKqCz6UPbUbYSqjRlJwS1Z5Xybt0FIMFvVtxJ1G8llg7n1FC3V+vF1PJtGo9Fo2uhIi2bDkY6l2fOLf8yJ//bzfHHim+xzqsRtcH1Brp7h1/PvYWh7PzL2KJYIaUYBIcvCQ6DqWu7Jq+nPM03lzZIIwTKWO4esVjpovQ3mavWFdTybRqPRaNpo0aLZkDx24mvUbrmB/1Xcw8TCcU6VJylFNd655SE+9Y5f5tGXHkXMJLjBLgFgtmpUYhZUAhVFSVoq+hI3VSpopaFcO9oiLoOVv0dtnc+o0Wg0GlhjekgIYQohhi/PVq4ditUiB+cOUqwWr/ZWrkuK1SIHZg8w0jFCJWFgd/ezc8f9vO/WjzEuKjz60qNMVCd4YfEuFprKjC5hqa9qAPUQxuuqcHfBU0W4ttFqbxZK0ASR6iASLPu0rBdJE779+p+u81k1Go1Gc9GRFiHE24H/jiob2CyEuAv4jJTypy7X5jYaVa/KY/sfY2x2DAODiIidhZ3s2bWHdEwXXq4XZbeMH/o8c+oZJquTS591f7qfuBVnobHA7f23Yy9U+IOXX+ITOyoMJHwaoaQJHKnAY+Pwvj41DTrvQOAqzxZTtOYTgaprYXXlLtdoOpe2YHr+eeDTl/oRaDQajeYcrCU99FngHcBjAFLK54UQt1+WXW1Q2oJlpGMEQxiEUcjY7BiPjz3Ow7c+fLW3d92QdbKMzY4RRAHdqW4EAikl45VxvMBjd+9uDGGQyxTIehmeOjRER7xKJZxkv1/neE1QNyTHa2Cl1Dk3J5Xw8FsdRpNNcCVsTkBPq8Oo7f1irFAxZuxGnpqxuDF+gJwdXVRoUgig+Z11/EQ0Go1GA2sTLZaU8rA4/U9Ob533s2FZmbIwhPrVZRomw7lh9s/sp1gtXtOTMzcaovWflBIhBBKJQGCbNhERkYzIdQ/Rnx1ksjyOQxcz0mTOn8A3a5iG4IlZlfjZkoKgpiIuzRCO1CBApYYmWh1IjplgIJXEoAkY2HYPAwM/x+Dgv+RtkRKsE8d+nZuyUxe1/87wjcv22Wg0Gs1blbWIlqYQIk0rsi6EuAloXpZdbUDKbhmBWBIsbUzDxMCg7Ja1aFknym6ZHYUdzDXmTksPDWYGySfydKe6OVo6ynBumNvv+SgvPPNlXq8coZsUN8lOXrQDKrg0Q/jGlEoLpUw42YDbOpTBXBs/AjI/zIfu+CwiKOL7c9h2F4nE6JKjbRp45LZH+OL4Zy76PdgGVCr79HRnjUajWUfWIlr+HfB1YEAI8QXgg8BPXo5NbUSyTnbpL/yVwiWMQiIisk72Ku7u+iLrZLFMi3uH7qXiVmgGTeJWnHQszbHFY3xi1yd48sSTS7VF/bc/yJ3yQ9yRGOWr889y6PBfE1SmqcsmUqhi3AXgbyeU++32rGp9DiIwnV382gP/rWUGN3rhjYm1WeYePfUou3f+3iV/DhqNRqM5nYsWLVLKbwgh3kCJFQH8eynlocu2sw1GX7qPnYWdjM2OMZwbxjRMwijk2OIxdnXv0lGWdeTMzzoXz532WY/mRxnNjy676TrKTfcLr3yBqYUam3KboFTCD5t4K77D25GX5xcgZauuot/+4C9ftHtttvuTUPncRb+Phje31reu0Wg0mgugZw+tgXN1D+3q3sXHd35cdw+tM2d+1mWvzObcZvbs3MNo/uyISPHIPn73hT8kmezgLw59ieL8CQzHJ2VCLVTRltMwoDPWyfz/Mb+mff3PbwiSFyn1E4Wf5b7df7Km82s0Gs1bnTc1e0gIMcMFrCyklD1vYm/XFEuTis/4C1+z/rQ/68Pzh/mr/X9FxatQapT43EufO73NvFqFxx6jfGAvwnqNI1RwnBLv7IvoSS273R6vwROzKtrSbgH6vR9Ye+rmuPk+dvDNizKkM6LFNZ9fo9FoNOfnYv5mvLN1+8+BPPA51O+BTwHjl2lfG5q+dJ8WK1eIJ088yUJzgZt7bl5qM39h/AWmqlM8fNvD9D32NRgbI7t5G/PVA7zgzXFnZ41kMuSUCzJS36xbkvBQQaWH2oyX1/7te/8Nv8Ab3/8mBWf1tQuLz635/BqNRqM5P6uKFinlcQAhxENSynesOPQZIcQTwG9frs1p3tqc2WbeDJq8NPkS45VxXph8gROTB7jr4DR7tryPPivJvKwRWFU2Z20OVT2kEUHLSK7oqtbnTgcWfHX+pr/25rddvXfz/L4keeqrerZEUX3N59doNBrN+VmLjf+AEKLQftC637/+W9JoFO0284pbYao6xd4Te5msTtKT6qEz3km3kWbMnOPR5tPs9Q4RJ8aAncAXEZZhq+/uSH3J1nTndruzJSweGn5ozXtynD62DXySMFp9bVycMyWr0Wg0mktkLS3Pvwu8IoT4+9bjDwP/Yf23pNEoTGHy6tSrBFGAH/kcWzhGX7qPbEy1n1vxJJOywncaz/G8d4xD4TQdRpIEHkMyTlWE+EYAqBSREKpjCODOgTt598i7L2lfD+3+A/7xu19Ydd1N/T93SefXaDQazbm56EiLlPKPUELlNeD7wIellH98uTam0Tx54sml+2k7TdJOUnJLHJg9wEB6gEP+FNMZQWcTBowsWSOO7ztMlwWprKQz0dHy1YW+uCrGXfDh9t7b+bMf/rNL3tfU1H/HNsUFBy0afpzMpndcYIVGo9Fo1spaIi1IKfcB+y7TXjSaJdr1LO8ffT9Pn3yaY6VjVLwKKTsFAoaSSQ7N72NgZITxI0fpLvkMmBYTVDhS6yPdZ7PFmGdLehADwYQbY6Ge4H974MP89nvfXBlWvb4fMFvDFoOzOomCCLJdH8I0teGgRqPRrCdrmfJ8lHO0Pkspt67rjjQalic9v1J8hZJbIhlLErfixETIu3sEQ+IZEol5hDjFyO4hOuLv5fapWSi/zvf8cb42YzCSGWVzpovtPbfTY/dwS8ug7s3OiUomdwEQITAME4iQkVwSLyWxiS2Fdy6NAdBoNBrN+rCWSMtHVtyPAz+FckfXaNad9qTnMArpSnZhCIPuZDdpby9WVGPe38WkC33pXrYnLcKTX0OOpbjHAJ8Akevg1u3vI5ddthFarzlRg4Of5ujR3yAKZvAjE9MwEEYASMLIwsh8gNGhT63Dp6DRaK4nXLdIGJYxzaz+o+YSWYuN//fPeOpFIcQ/rPN+NMAr499munyUnuwItw1eWrHo9UJ7ujNA1gwYydiU/D4e3Pwgx0rHWXRLhKdKxIJp6L2VY36TndEox0vHifbtgwffq060uEjYqBH5lXWZE3XzzV/me6/9IARzyFAiBHgyTinzc3zktn9/0aMBNBrN9U8QVJmZeYx6fYx2W2MyuZPu7j36Z8UaWVNNy0qEEJ3AyDru5S3PxOJh/uveT+G5byAk5Kwmz79W4B03/QY3Dv3U1d7eFaXsltlZ2HnapOcOq0I2mSGT6CftZHhw8wO8eGQvxcobJHOCE2KKmhERi5ksCpcvzT7F9u/73FlKYhenOWZW2BXm6St9DfbsgfSl/7Do6LiXtz84zfj4nzKz+AKROcTA4Ke16aBGozmLtmCJx0cQwkDKkHp9jNnZx+nre/hqb++aYi01Lc+zXNNiogTL/305NvVW5b/u/RRe8yDJWIHtzjHiRhVDLnDk9Z+mOv1f2LHj8ySTq0wivk44bdKzV6HpN0lbPv08Q7EJCSuOY8W5O7OFResodSci1tnEtUp0iAS3YDA2J/jW8eeY9+Lc1DnKLnsnH7d3w9gYPP44PPzmf1gMDn6awcFPr8M71mg01yOuW6ReP7AkWACEMInHh6nV9uO6RZ0qWgNribT86xX3A+ColHJinffzluWV8W/jN9/AjA1yo/MqGVN1pUjAErCw8ATPPXc/73zn1Gqnui44a9KzoyY9T5QtNicFKcticfGfcIOjBIWThA7Yho9V3UwlAl/63OjAYDrBs14fP5G4l92xltnb8DDs3w/FIvSd/cPirLxzsQjlMmSz51yv0Wg05yMMy4BYEixthDARwmgd1z9XLpa1iJZRKeXnVz4hhPjUmc9pLo3p8lEAOowaWWO5jVaI5dsomuYrf9fDvXd9i3zf7quz0SvInl17zpqqfVPho9zWJZmf/1uCYB5hO9ixbpr1IsKKKDunOFQxCTwXacKmziY9CwPEDXv5xKYJhqGEyAoRclbe2XNJvlYm930DTIHpWTijd73p1JJGo3nroKwPIqSMThMuUoZIGWlrhDWyFtHyC8CZAuVcz2kugZ7sCFJAF8eWhMq5iKdmeP6V+3jnuw5f9yHF803Vdt0izdpLCBEjDBcJe9KECz7NWhVMQRTGqdgGODHGw0VeDl/nn4X3cRDIijh9pCGKVORkBWfmncNX/pEZ+SzTdyRJe0MgI5ITR+n+kov109rtVqPRrI7j9JFM7mz9bBlGCBMpQ5rNY6RSu677n+PrzaqiRQhxJ3APUBBC/PyKQzkgth6bEELEgb8EdgF1oAj8CynlMSFED/D/AqOA23r+qfW47kbitsF38/fODTjGk0g4y7BsJXaszqsv/yh33/udK7S7q8uZU7XDsIzrjgMhhtXNqdoRfFNiOGDLkGoiRJhpLGFgYXGiucj/6X6Zu5wRIhmxc8Fiz46Pkl4RZWnnnS0rj+/PYNR8ao3XiDriQIgdpTAjh/rALLPjf0Nf8aM6VaTRaC6K7u49S38UqULciFRqF4XCx6/21q45LibSMgjcCaSAu1Y8XwYeWce9fA74qpRSCiF+ofX4/cBvAc9IKT8ohLgLeEwIMSqlDNbx2huCn3ng8/zNk/eQj82turbefIpS6Rk6Ou69AjvbWIRhkyCYJR7fypHSccpeg5SZwTabWFFIEDVZ8HwGEyaHXItm5ND0a+QbHunIZKxb8PhOWFmG63mT1GrfQ31bGUTVefxskVS0hcBoEAkfiwTxsEAtcwR34TCOFi0ajeYisKw0/f2PaJ+WdWBV0SKl/Bvgb4QQH5JSfvVybEJK2QS+suKpZ4Bfbt3/JK3Wainl80KIKeBB4DuXYy9Xk4HcKP/yI7N8/R8E9oVCLUAYhhw98Vu8rePLV2ZzGwjTjGPbBWpukbJbIm07yDAGMsAQkpGMQ96DcS/D3tkGUcwg1jVAs/8WSOZIOAbPze/jAyuccUulJ/B9JYSEMAliAU3TpynmsaIEhlQ1MSISiAjC1NX8BDQazbWIEiparLwZLiY99ICUcq+6Kz585nEp5VfO8bI3y2eAvxNCdAGGlHJmxbFjwObLcM0NQ//I55g6+rOYFxAuQsDc/P98S7bLmWaWRGI7tcobZNlHInSRBMSsECmhYAeUPUG9vkjD80mZdRrlMV53m8x09SBMi1KzxBde+QK/cPcvYMsqnjdJIrED35/CsvIYiSzCSuAGc8TDbVhRAiKJnJtGburC7H1rtJ5rNBrNRuJi0kOPAHuBXzvHMcnpEZI3jRDi3wA3AP8CSHD2vKNz/ioXQvwq8Kvtx7lcbj23dUXZvfXTnOq8h/0v3Yp1IeEShRw6+afctO03rtzmNgCO00c6vZuFyquAJCRBhIckwJUSSwQYEjYl4We3QSUAZJ0kBxgJqpxI30koQ04tnuLxscf5xA33EYY14vFNRFGTIJgDDMhmoOpiTQGNOaQMaA7FSN36g5dVKH776Lc5unCUkc4R3j3y1nZE1mg0mpVcTHro063bd13uzQgh/jXww8B7pZR1oC6EQAjRvSLasgU4cY59/g7wO+3HQ0NDZw13vJYY6tzNK2Y/VjR5wXXVmVdg2xXa1AYil3s7xeIXiRsCNwqIGT5NaePKGJao0OtIyqagx4GnZiEILUpE3BArUgi+RzK3m7cVChyd2cvR1AmalecwzU6EkBhGCil97FgBPxNS3uHhigaJ5HYyvfdftuK5/7+9ew+S6yzvPP59untmdJmbrh7ZY0sj2xjJxthgBy9gE5KwVYHisrFSLLfg4NrCFXuXlIvd/WPZJAVsajdsOalAqhJXvJg1hCLIm8K7hkotAYIhphwwRvZKQljWWJKtsTS6jUbS3Lrf/aN7xIw0F8nSdM878/1QXep+z+nu5zx1kH865/R7dh/Zzce+8TF29O+ABBUq9HT28Pl3fp7buhfftUuSdLbC7KtURcRT5zP2atWOlHwAeEdK6diERV8H7q2tcyvVE4IL7tdDU3nt+vtnXB7A8uiuTzHzTplly65mRdv1FJq6GK60cLrSwnB5lGIUaS0U6GpuYqRcpLPQQgISBfaPjHFjaz+/3vkynWM/YnPh2xzq/yYjhTWMVIYplVZy+vQuTp/eQ7HYRkfHHay94kM0r7iWpas209X10Tm7V8jv/N3v8HTf04yURzg1dooTwyd4pu8Z3vPV9/DwMw8zODI4J98rSbm4kHlaJq0bEUXgkvztHRHdVG8J8ALw3ahOVDKcUnoT8B+BRyLiF8AI8JGF+MuhqVzT80lefvG/MFY5NuXyymiRa7vvrnNV80Ox2E5EM6XiEq7s3MSJwQodS5s4NnyUwdOHqJAYS2Wao4k1hU5OxSAJuH55mc6WJkaLHZyuBCPDA5Qro1RODpAIVjYfoKN0mkgFmpo6aWu7mUJhCaVSB0NDvXN2DdHjux7n6QNPU05lKqlCOZXPLOs/1c+nv/dpjg8f5xNv+sQl/25JysWsR1oi4t9HxCHghog4OP4AjgNPXIoiUkr7U0qRUro6pXRT7fGm2rJXUkr/MqV0bUrp+pTSP16K78zFjTd9i0Kl6ZzxymjQfeKdtHQv/Jlxp9LS0kVb2y0UCiXK5RO0NK+gKcZooUIUShwcDpojMTBcoFxupkSRVsp0LA1SqZUyzew5+nMSsLRlHR3N0LTkNRwYaWdgrIVly3pYvvwGCoUlwNlTbl96j2x7hJHKCBFBJVXOjAdBInFs6BgPPf0QfYN9c/L9kpSD8zk99CDV+Vn+ofbn+OPKlJLTgs6xzs7beOute+kafDNNp9oonmqj9ZX13HDibq55x4ONLq+h1qzZwqpV7yGiiZTKnBg+wljlNMta1pOa1vLiqRIHh0c4MHKIK8aWsnrJcpYuLTFYbmH/4FHKlSJrl3ZSSkNEZYwiY7Qu6WZwdIThscqZwAJzO+V232Af/Sf7SSlVHySi9r9Uuw59efNy+k/189T+S3ZGVpKycz4X4h6nelTlN+e+HE2l1NHF5vf8kOH92ygPHKDYvm7RHmGZqFRq5Yor7mH16vfx0tGn+L/P/S2XtW1gpAK9sZ1/PvYcm1peoat1hMKyItdfdh3rWqCp9c1cOdzPgf6XaI5BqBwlgHKljaFCDykVKZMoFqtnP+d6yu2B4QHWta2judDMaGX0zPh4YAmC1pZWjg8d5+ToyUv+/ZKUi/O+piUirgb+DHg9cOafoCmltXNQl6ZQDSqGlbO1tHQx1vRaBriSVaWr+NG+f+LAyZPc0HUrQ6ND7Dy9j42Fy1i/7De4cm0Lhw8/RunUbtriICk1127F2EIpHWFpeYheruZ1K3+DoaHeuky53d7STseSDq7suJJ9x/cxXBk+E1gAlpaWMloeZUlpCdevvX5OapCkHFzIhbh/DfwlsBF4F/BvqU70JjVce0s7FSocHTpK32Afq5etJghaSi0MllZw1eq3sOvwLm5sqTA89AqRTlEqLmOkPExTcSkjhcsYiss5NfIynSvfxaaNn6jblNtdrV3ccvkt/OLwL2guNrPn6B6GykMEQVOhibbmNlJK3LH+Dm68zNAqafG6kNDSkVL6WkR8KqX0bER8HPgO8MdzVJt03rpau9i0ehM/3PtDUqpeE1JJFY6cOsJVbWvYUHyWUuEHHD98mmKhmUrlNJ1L1tE/2sThkSGay8c5nDq4fNkVvGXjHUB9p9zesnkLQ2NDPPbzxyhGkf0n9jNaHqVjSQerlq3iV674FT7z9s/UpRZJmq8uJLSMn2w/ERHrgVeoTvQmzQtbNm/h+NBxth3cdmbe5MvbLuf2NSMsH/0ZI4xRauqkVFjO8PBexsaOsG7ZBta1b2Ro5BCvbbuFZaUibUvW1b321uZW7rnlHt732vex+8huAIbLw4yWR1nXts4jLJLEhYWWf4yIlcAXgB8Dw8Cc3EBRejVam1v5xG3VeUyeO/QcV6+4mjUtRZYPPc7psUGWN6+jVDhNodBEqdRBuXyckZFDLC120hJjlCqHWb78Vxt6L6eu1q4zN3GUJE123qElpfQfak//JiKeADqBx+ekKuki3P2Gu9m6fSs7+ncwOnSMUhxjTUs73Z3XMDq8l3L5OKXSCsrl04yO9lMuH6dYbGN0tDo2NjY4Z7PeSpJevQs50nJGSmkfsC9qU9dK80lrcyt33XQXfYN9HB3czfDhERh9iVKhSHHJVQwN7WVs7BgpjQJFli9/I52db6WpaSVDQ3vo73+Urq6PNnozJElnOe97D00j65sSamHrau1iU9dbuGzl2ygUSoyMHCSiyJIlGyiVOigWl7JixdtZu/Zf0dy85syykye3MzzszLOSNN/MeqQlIjZfzPulRluzZguVyhD9/Y9x+nT1IteIFpYv30xn5+2T1p08Xb/XlkjSfHI+oWOm61aGLlUh0lyZOHPueGgpFts4ePArRDRPWncup+uXJF2c85nGv6cehUhzraWla9IvgwYHn+bUqR0sWbKBiOKcT9cvSbo4nt7RorVmzRYOHdrKqVM76jJdvyTp4hhatGiVSq2sW3dX3abrlyRdHEOLFr16TtcvSXr1LvYnz5IkSXVhaJEkSVkwtEiSpCwYWiRJUhYMLZIkKQuGFkmSlAVDiyRJyoKhRZIkZcHQIkmSsmBokSRJWTC0SJKkLBhaJElSFgwtkiQpC4YWSZKUBUOLJEnKgqFFkiRlwdAiSZKyYGiRJElZMLRIkqQsGFokSVIWDC2SJCkLhhZJkpQFQ4skScqCoUWSJGXB0CJJkrJgaJEkSVkwtEiSpCwYWiRJUhYMLZIkKQuGFkmSlAVDiyRJykKp0QVIF6yvDwYGoL0duroaXY0y1vfCNgaOHKB95Tq6Nt7Y6HIkzSKL0BIR1wJfAlYDx4C7UkrbG1uV6m5wELZuhR07oFCASgU2bYItW6C1tdHVKSODR/rY+pVPsaPvWQpRoJIqbOp6HVs+9FlaVxqEpfkql9NDfwU8mFJ6DfAnwEMNrkeNMB5Yenqqjw0bqq8ffbTRlSkz44GlZ+VGelZfzYaVG9nR9yyP/s1/bnRpkmYw70NLRKwF3gB8uTb0KNATERsaVZMaoK8Pdu6shpVCbbctFqvBZfv26nLpPPS9sI2drzxHz8qNFIpFAIrFIhtWbmT7gW30vbCtwRVKms68Dy3AlcDLKaUxgJRSAvYCV01cKSLuj4j944/BwcEGlKo5MzAAEb8MLOOKxerYwEBj6lJ2Bo4cIIgzgWVcsVikEAUGjhxoUGWSZpNDaAFIZ72Oc1ZI6YGUUvf4o9VrHBaW9vbqNSyVyuTxcrk61t7emLqUnfaV66ikCpVyedJ4uVymkiq0r1zXoMokzSaH0LIP6I6IEkBEBNWjL3sbWpXqq6uretHtnj3VoALVP3t7YfNmf0Wk89a18UY2db2OPUdeoFzbl8rlMr1HXmDzuhv9FZE0j8370JJSOgj8FPhwbehOoDel1NuwotQYW7ZUg0tvbzW8jAeWO+9sdGXKzJYPfZZNXa+j98gL7OnffSaw3PnBzzS6NEkziOolIvNbRFwHPAysAgaAj6aU/t9M7+nu7k779++vQ3WqO+dp0SXiPC3S/BMRL6WUuqdclli1H8UAAAvGSURBVENoeTUMLZIk5Wem0DLvTw9JkiSBoUWSJGXC0CJJkrJgaJEkSVkwtEiSpCwYWiRJUhYMLZIkKQuGFkmSlAVDiyRJyoKhRZIkZcHQIkmSsmBokSRJWTC0SJKkLBhaJElSFgwtkiQpC4YWSZKUBUOLJEnKgqFFkiRlwdAiSZKyYGiRJElZMLRIkqQsGFokSVIWDC2SJCkLhhZJkpQFQ4skScqCoUWSJGXB0CJJkrJgaJEkSVkwtEiSpCwYWiRJUhYMLZIkKQuGFkmSlAVDiyRJyoKhRZIkZcHQIkmSsmBokSRJWTC0SJKkLBhaJElSFgwtkiQpC4YWSZKUBUOLJEnKgqFFkiRlwdAiSZKyYGiRJElZMLRIkqQsGFokSVIWDC2SJCkLDQ8tEfHHEbEjIn4WEU9FxK9NWFaIiM9HxO6IeD4ifq+RtUqSpMYpNboA4AngMyml0xHxeuB7EbEupTQEfBjYDLwG6ACejojvpJR2NrBeSZLUAA0/0pJS+lZK6XTt5bNAEVhde/1+4C9TSuWU0hHgb4F/3YAyJUlSgzU8tJzld4HdKaX9tddXAS9OWN5bG5MkSYvMnJ8eiogngE3TLL45pbSvtt6vA38IvOOsddLEj5vhe+4H7h9/3dHR8arqlSRJ89Och5aU0u2zrRMRbwO+CLw7pfTzCYv2AhuAf669Xl8bm+p7HgAeGH/d3d2dplpPkiTlqeGnhyLiDuAR4L0ppZ+dtfjrwMcjohgRK6le4/K1etcoSZIabz78eughoAX4YsSZsz8fSSk9SzXM3Arsqo1/LqW0o/4lSpKkRmt4aEkpXTvDsjJwbx3LkSRJ81TDTw9JkiSdD0OLJEnKgqFFkiRlwdAiSZKyYGiRJElZMLRIkqQsGFokSVIWDC2SJCkLhhZJkpQFQ4skScqCoUWSJGXB0CJJkrJgaJEkSVkwtEiSpCwYWiRJUhYMLZIkKQuGFkmSlAVDiyRJyoKhRZIkZcHQIkmSsmBokSRJWTC0SJKkLBhaJElSFgwtkiQpC4YWSZKUBUOLJEnKgqFFkiRlwdAiSZKyYGiRJElZMLRIkqQsGFokSVIWDC2SJCkLhhZJkpQFQ4skScpCqdEFSJKk+a9vsI+B4QHaW9rpau1qSA2GFkmSNK3BkUG2bt/Kjv4dFChQocKm1ZvYsnkLrc2tda3F00OSJGla44Glp7OHnhU9bOjYwI7+HTy649G612JokSRJU+ob7GNn/056OnsoRDUyFAtFNnRsYPuh7fQN9tW1HkOLJEma0sDwAEGcCSzjioUiBQoMDA/UtR5DiyRJmlJ7SzsVKlRSZdJ4uVKmQoX2lva61mNokSRJU+pq7WLT6k3sObaHcqUMVANL7/FeNq/ZXPdfEfnrIUmSNK0tm7ec8+uhzWs2c+emO+teS6SU6v6l9dDd3Z3279/f6DIkSVoQ6jVPS0S8lFLqnmqZR1okSdKsulq7Gjap3DivaZEkSVkwtEiSpCwYWiRJUhYMLZIkKQvzJrRExK9GRDki7pswtiwivhoRz0fEroj4rUbWKEmSGmde/HooItqA/wZ866xFnwSGU0rXREQP8GREfDeldLTuRUqSpIaaL0daHgA+B/SfNf5+4C8AUkp7gO8D761vaZIkaT5oeGiJiN8EOlNKW6dYfBXw4oTXvbWxqT7n/ojYP/4YHBy89MVKkqSGmfPTQxHxBLBpmsU3A/8VeMcMHzFxyt6YdqWUHqB6xAaozoh7AWVKkqR5bs5DS0rp9umWRcRbgXXAUxEBsBp4d0SsSSn9IbAX2AAcqr1lPfDNOS1YkiTNSw29EDel9ANg7fjriHgY+HFK6Qu1oa8D9wJ31S7EfRtwT73rlCRJjdfwa1pm8TlgaUQ8D/w9cG9K6UiDa5IkSQ2wYO/yHBHD/PK00mLRCngF8mT2ZDL7cS57Mpn9OJc9mWyu+7EmpdQy1YIFG1oWo4jYP93tvBcrezKZ/TiXPZnMfpzLnkzWyH7M99NDkiRJgKFFkiRlwtCysDww+yqLjj2ZzH6cy55MZj/OZU8ma1g/vKZFkiRlwSMtkiQpC4YWSZKUBUNLhiLiYxHxbESMRcR9Zy17uHbTyGdqj89NWLYsIr4aEc9HxK6I+K36Vz83ZunJtNu9kHsy0WLdL2YSEddGxD/VtvmpiNjc6JrqLSJ6I2LnhP3i/bXxRdObiPjzWh9SRNwwYXzaHizk/szQjyn3ldqy+vUjpeQjswfweqo3ofyfwH1nLXv47LEJy/4AeLj2vAfoA1Y0envq0JNpt3sh98T9YtaefAe4q/Z8C/Bko2tqQA96gRsWc2+AO4Dus3sxUw8Wcn9m6MeU+0q9++GRlgyllH6WUtoBVC7wre8H/qL2GXuA7wPvvcTlNcQsPZlpuxdsTy7AoutBRKwF3gB8uTb0KNATERsaVdN8sdh6k1L6fkpp/8SxmXqw0PszVT9mUu9+GFoWpvsjYltE/J+IuGnC+FXAixNe99bGFrqZtnsx9cT94peuBF5OKY0BpOo/Efey8Ld7Kl+pnVr964hYg72BmXuwmPtz9r4Cde6HoWUeiognIqJ/mseVs7z9PwHXpJRuBB4CvhURrROWT/yNe1zq2ufKRfYEZt7uLHsy0Xn0Z0HuFxfp7PkeFst2T3RHSun1VP+lfBj4Um3c3szcg8XYn+n2FahjPwwt81BK6faU0uppHvtmee9LKaVK7fnfAQPAdbXFe4ENE1ZfXxub9y6mJ8y83dn2ZKLZ+rNQ94uLsA/ojogSQEQE1X8xLvTtniSltLf25yjwZ8Dt2BuYuQeLsj/T7CtQ534YWhaYiOie8Pw2YBXwfG3o68C9tWU9wNuAx+pdYwPMtN2LoifuF5OllA4CPwU+XBu6E+hNKfU2rKg6i4jlEdE5YegDwE/tzcz7x2Lsz3T7CjTg/0tzeRWyj7l51HaO/cBJ4Gjt+c21Zd8GngWeAZ4E3j7hfcuBr1H9j9UuYEujt6VOPZl2uxdyT87qz6LcL2bpyXW1XuwCfgxc3+ia6rz9G6n+x2Zbbd/4BrBhsfWG6kXo+4Exqr+ce362Hizk/kzVj5n2lXr3w2n8JUlSFjw9JEmSsmBokSRJWTC0SJKkLBhaJElSFgwtkiQpC4YWSZdE7S6wN8y+5kV/zzMRsbT2/Pdr9z55NZ9zV0RsvbTVSZpLhhZJWUkp3ZRSOl17+fvAqwotkvJjaJE0ZyLiloh4snajxqci4i218Q21+yJ9OiJ+EhHPR8Q7J7zvzojYGRE/jYhPRUQav1fS+POI+APgcmBr7ejLTRHxRxHx3yd8zn0R8XDteXNE/FVE7IqI7wJvOqvWT9ZqfDoivnme97SSVEeGFklzIiKagf8F/FGq3qjxfqoBY3ltlVXAT1JKbwTuA/609r61wIPAu1NKNwODU31+SunTwMtUZ/C9KaX0zCwlfRzoAa4H3gXcOqHWDwKvAf5FSukNwFeBL1z4VkuaS6VGFyBpwboOGEkp/T1ASukHEXEQuBE4AJxMKX2jtu6TwNW157cBT6eUflF7/UVqgeYivR34Uqre8G00Ir4MvLW27H3ALcBPqvd7owiUL8F3SrqEDC2S5kpw7i3rmTA2NGGsTDUozPS+8zE24XMAlpxVz3QC+GxK6X+8yu+VVAeeHpI0V3YCLRHxawAR8WaqF80+O8v7fgS8MSKuqb3+6AzrDgAdE17vBm6JiEJELKN6x9lx/wB8JCJKtV8ffXDCsseA34uIlbVamyLi5lnqlFRnHmmRdCl9OyLGJrz+APAntetYhoDfTimdjIg1031ASumViLgHeDwiDgP/GxgFTk2x+p8DX4yIU8BdwKPAFmA70Ev1rtZLa+s+SPXU1Haqd7F9Alhf+85HImIV8L2ISFT/bnyI6p1tJc0T3uVZ0rwTEW0ppRO1578L3J1Seussb5O0wHmkRdJ89O8i4rep/h11BPg3Da5H0jzgkRZJkpQFL8SVJElZMLRIkqQsGFokSVIWDC2SJCkLhhZJkpQFQ4skScqCoUWSJGXh/wNEKbK6x/UBmgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plotting clusters\n",
    "plt.figure(num=None, figsize=(8, 6), dpi=80)\n",
    "if k==2:\n",
    "    x_0=weather_clusters[weather_clusters.Cluster==0][\"Longt\"]\n",
    "    y_0=weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1=pl.scatter(x_0,y_0,c='r',marker='o',alpha=0.4)\n",
    "    x_1=weather_clusters[weather_clusters.Cluster==1][\"Longt\"]\n",
    "    y_1=weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    c2=pl.scatter(x_1,y_1,c='g',marker='o',alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\",len(x_0), \", Cluster1 Size:\",len(x_1))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "\n",
    "elif k==3:\n",
    "    x_0 = weather_clusters[weather_clusters.Cluster == 0][\"Longt\"]\n",
    "    y_0 = weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1 = pl.scatter(x_0, y_0, c='r', marker='o', alpha=0.4)\n",
    "    x_1 = weather_clusters[weather_clusters.Cluster == 1][\"Longt\"]\n",
    "    y_1 = weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    c2 = pl.scatter(x_1, y_1, c='g', marker='o', alpha=0.4)\n",
    "    # for sensor fault visibility in figure\n",
    "    # c2 = pl.scatter(x_1, y_1, c='b', marker='x', alpha=1, s=300, linewidths=4, zorder=10)\n",
    "    x_2=weather_clusters[weather_clusters.Cluster==2][\"Longt\"]\n",
    "    y_2=weather_clusters[weather_clusters.Cluster == 2][\"Lat\"]\n",
    "    #c3 = pl.scatter(x_2, y_2, c='b', marker='x', alpha=1, s=300, linewidths=4, zorder=10)\n",
    "    c3=pl.scatter(x_2,y_2,c='b',marker='o', alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\", len(x_0), \", Cluster1 Size:\", len(x_1), \", Cluster2 Size:\", len(x_2))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "    max_c2 = max(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    min_c2 = min(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    print(\"max_c2:\", max_c2, \" min_c2:\", min_c2, \"Color:B\")\n",
    "elif k==4:\n",
    "    x_0 = weather_clusters[weather_clusters.Cluster == 0][\"Longt\"]\n",
    "    y_0 = weather_clusters[weather_clusters.Cluster == 0][\"Lat\"]\n",
    "    c1 = pl.scatter(x_0, y_0, c='r', marker='o', alpha=0.4)\n",
    "    x_1 = weather_clusters[weather_clusters.Cluster == 1][\"Longt\"]\n",
    "    y_1 = weather_clusters[weather_clusters.Cluster == 1][\"Lat\"]\n",
    "    #c2 = pl.scatter(x_1, y_1, c='g', marker='x', alpha=0.8, s=169, linewidths=3, zorder=10)\n",
    "    c2 = pl.scatter(x_1, y_1, c='g', marker='o', alpha=0.4)\n",
    "    x_2 = weather_clusters[weather_clusters.Cluster == 2][\"Longt\"]\n",
    "    y_2 = weather_clusters[weather_clusters.Cluster == 2][\"Lat\"]\n",
    "    c3 = pl.scatter(x_2, y_2, c='b', marker='o', alpha=0.4)\n",
    "    x_3=weather_clusters[weather_clusters.Cluster==3][\"Longt\"]\n",
    "    y_3=weather_clusters[weather_clusters.Cluster == 3][\"Lat\"]\n",
    "    c3=pl.scatter(x_3,y_3,c='y',marker='o', alpha=0.4)\n",
    "    # Numbers of Elements in Clusters\n",
    "    print(\"Cluster0 Size:\", len(x_0), \", Cluster1 Size:\", len(x_1), \", Cluster2 Size:\", len(x_2), \", Cluster3 Size:\", len(x_3))\n",
    "\n",
    "    # Print Cluster Max, Min Points to determine Cluster Seperation Point\n",
    "    max_c0 = max(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    min_c0 = min(weather_clusters[weather_clusters.Cluster == 0][title])\n",
    "    print(\"max_c0:\", max_c0, \" min_c0:\", min_c0, \"Color:R\")\n",
    "    max_c1 = max(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    min_c1 = min(weather_clusters[weather_clusters.Cluster == 1][title])\n",
    "    print(\"max_c1:\", max_c1, \" min_c1:\", min_c1, \"Color:G\")\n",
    "    max_c2 = max(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    min_c2 = min(weather_clusters[weather_clusters.Cluster == 2][title])\n",
    "    print(\"max_c2:\", max_c2, \" min_c2:\", min_c2, \"Color:B\")\n",
    "    max_c3 = max(weather_clusters[weather_clusters.Cluster == 3][title])\n",
    "    min_c3 = min(weather_clusters[weather_clusters.Cluster == 3][title])\n",
    "    print(\"max_c3:\", max_c3, \" min_c3:\", min_c3, \"Color:Y\")\n",
    "\n",
    "pl.xlabel('Longitude')\n",
    "pl.ylabel('Latitude')\n",
    "pl.title(title)\n",
    "pl.savefig(\"plot_output.png\")\n",
    "pl.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Clustering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "class SpectralClusteringg():#ClusterMixin, BaseEstimator):\n",
    "    \"\"\"Apply clustering to a projection of the normalized Laplacian.\n",
    "    In practice Spectral Clustering is very useful when the structure of\n",
    "    the individual clusters is highly non-convex or more generally when\n",
    "    a measure of the center and spread of the cluster is not a suitable\n",
    "    description of the complete cluster. For instance when clusters are\n",
    "    nested circles on the 2D plane.\n",
    "    If affinity is the adjacency matrix of a graph, this method can be\n",
    "    used to find normalized graph cuts.\n",
    "    When calling ``fit``, an affinity matrix is constructed using either\n",
    "    kernel function such the Gaussian (aka RBF) kernel of the euclidean\n",
    "    distanced ``d(X, X)``::\n",
    "            np.exp(-gamma * d(X,X) ** 2)\n",
    "    or a k-nearest neighbors connectivity matrix.\n",
    "    Alternatively, using ``precomputed``, a user-provided affinity\n",
    "    matrix can be used.\n",
    "    Read more in the :ref:`User Guide <spectral_clustering>`.\n",
    "    Parameters\n",
    "    ----------\n",
    "    n_clusters : int, default=8\n",
    "        The dimension of the projection subspace.\n",
    "    eigen_solver : {'arpack', 'lobpcg', 'amg'}, default=None\n",
    "        The eigenvalue decomposition strategy to use. AMG requires pyamg\n",
    "        to be installed. It can be faster on very large, sparse problems,\n",
    "        but may also lead to instabilities. If None, then ``'arpack'`` is\n",
    "        used.\n",
    "    n_components : int, default=n_clusters\n",
    "        Number of eigen vectors to use for the spectral embedding\n",
    "    random_state : int, RandomState instance, default=None\n",
    "        A pseudo random number generator used for the initialization of the\n",
    "        lobpcg eigen vectors decomposition when ``eigen_solver='amg'`` and by\n",
    "        the K-Means initialization. Use an int to make the randomness\n",
    "        deterministic.\n",
    "        See :term:`Glossary <random_state>`.\n",
    "    n_init : int, default=10\n",
    "        Number of time the k-means algorithm will be run with different\n",
    "        centroid seeds. The final results will be the best output of\n",
    "        n_init consecutive runs in terms of inertia.\n",
    "    gamma : float, default=1.0\n",
    "        Kernel coefficient for rbf, poly, sigmoid, laplacian and chi2 kernels.\n",
    "        Ignored for ``affinity='nearest_neighbors'``.\n",
    "    affinity : str or callable, default='rbf'\n",
    "        How to construct the affinity matrix.\n",
    "         - 'nearest_neighbors' : construct the affinity matrix by computing a\n",
    "           graph of nearest neighbors.\n",
    "         - 'rbf' : construct the affinity matrix using a radial basis function\n",
    "           (RBF) kernel.\n",
    "         - 'precomputed' : interpret ``X`` as a precomputed affinity matrix.\n",
    "         - 'precomputed_nearest_neighbors' : interpret ``X`` as a sparse graph\n",
    "           of precomputed nearest neighbors, and constructs the affinity matrix\n",
    "           by selecting the ``n_neighbors`` nearest neighbors.\n",
    "         - one of the kernels supported by\n",
    "           :func:`~sklearn.metrics.pairwise_kernels`.\n",
    "        Only kernels that produce similarity scores (non-negative values that\n",
    "        increase with similarity) should be used. This property is not checked\n",
    "        by the clustering algorithm.\n",
    "    n_neighbors : int, default=10\n",
    "        Number of neighbors to use when constructing the affinity matrix using\n",
    "        the nearest neighbors method. Ignored for ``affinity='rbf'``.\n",
    "    eigen_tol : float, default=0.0\n",
    "        Stopping criterion for eigendecomposition of the Laplacian matrix\n",
    "        when ``eigen_solver='arpack'``.\n",
    "    assign_labels : {'kmeans', 'discretize'}, default='kmeans'\n",
    "        The strategy to use to assign labels in the embedding\n",
    "        space. There are two ways to assign labels after the laplacian\n",
    "        embedding. k-means can be applied and is a popular choice. But it can\n",
    "        also be sensitive to initialization. Discretization is another approach\n",
    "        which is less sensitive to random initialization.\n",
    "    degree : float, default=3\n",
    "        Degree of the polynomial kernel. Ignored by other kernels.\n",
    "    coef0 : float, default=1\n",
    "        Zero coefficient for polynomial and sigmoid kernels.\n",
    "        Ignored by other kernels.\n",
    "    kernel_params : dict of str to any, default=None\n",
    "        Parameters (keyword arguments) and values for kernel passed as\n",
    "        callable object. Ignored by other kernels.\n",
    "    n_jobs : int, default=None\n",
    "        The number of parallel jobs to run when `affinity='nearest_neighbors'`\n",
    "        or `affinity='precomputed_nearest_neighbors'`. The neighbors search\n",
    "        will be done in parallel.\n",
    "        ``None`` means 1 unless in a :obj:`joblib.parallel_backend` context.\n",
    "        ``-1`` means using all processors. See :term:`Glossary <n_jobs>`\n",
    "        for more details.\n",
    "    verbose : bool, default=False\n",
    "        Verbosity mode.\n",
    "        .. versionadded:: 0.24\n",
    "    Attributes\n",
    "    ----------\n",
    "    affinity_matrix_ : array-like of shape (n_samples, n_samples)\n",
    "        Affinity matrix used for clustering. Available only if after calling\n",
    "        ``fit``.\n",
    "    labels_ : ndarray of shape (n_samples,)\n",
    "        Labels of each point\n",
    "    Examples\n",
    "    --------\n",
    "    >>> from sklearn.cluster import SpectralClustering\n",
    "    >>> import numpy as np\n",
    "    >>> X = np.array([[1, 1], [2, 1], [1, 0],\n",
    "    ...               [4, 7], [3, 5], [3, 6]])\n",
    "    >>> clustering = SpectralClustering(n_clusters=2,\n",
    "    ...         assign_labels=\"discretize\",\n",
    "    ...         random_state=0).fit(X)\n",
    "    >>> clustering.labels_\n",
    "    array([1, 1, 1, 0, 0, 0])\n",
    "    >>> clustering\n",
    "    SpectralClustering(assign_labels='discretize', n_clusters=2,\n",
    "        random_state=0)\n",
    "    Notes\n",
    "    -----\n",
    "    If you have an affinity matrix, such as a distance matrix,\n",
    "    for which 0 means identical elements, and high values means\n",
    "    very dissimilar elements, it can be transformed in a\n",
    "    similarity matrix that is well suited for the algorithm by\n",
    "    applying the Gaussian (RBF, heat) kernel::\n",
    "        np.exp(- dist_matrix ** 2 / (2. * delta ** 2))\n",
    "    Where ``delta`` is a free parameter representing the width of the Gaussian\n",
    "    kernel.\n",
    "    Another alternative is to take a symmetric version of the k\n",
    "    nearest neighbors connectivity matrix of the points.\n",
    "    If the pyamg package is installed, it is used: this greatly\n",
    "    speeds up computation.\n",
    "    References\n",
    "    ----------\n",
    "    - Normalized cuts and image segmentation, 2000\n",
    "      Jianbo Shi, Jitendra Malik\n",
    "      http://citeseer.ist.psu.edu/viewdoc/summary?doi=10.1.1.160.2324\n",
    "    - A Tutorial on Spectral Clustering, 2007\n",
    "      Ulrike von Luxburg\n",
    "      http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.165.9323\n",
    "    - Multiclass spectral clustering, 2003\n",
    "      Stella X. Yu, Jianbo Shi\n",
    "      https://www1.icsi.berkeley.edu/~stellayu/publication/doc/2003kwayICCV.pdf\n",
    "    \"\"\"\n",
    "    #@_deprecate_positional_args\n",
    "    def __init__(self, n_clusters=8, *, eigen_solver=None, n_components=None,\n",
    "                 random_state=None, n_init=10, gamma=1., affinity='rbf',\n",
    "                 n_neighbors=10, eigen_tol=0.0, assign_labels='kmeans',\n",
    "                 degree=3, coef0=1, kernel_params=None, n_jobs=None,\n",
    "                 verbose=False):\n",
    "        self.n_clusters = n_clusters\n",
    "        self.eigen_solver = eigen_solver\n",
    "        self.n_components = n_components\n",
    "        self.random_state = random_state\n",
    "        self.n_init = n_init\n",
    "        self.gamma = gamma\n",
    "        self.affinity = affinity\n",
    "        self.n_neighbors = n_neighbors\n",
    "        self.eigen_tol = eigen_tol\n",
    "        self.assign_labels = assign_labels\n",
    "        self.degree = degree\n",
    "        self.coef0 = coef0\n",
    "        self.kernel_params = kernel_params\n",
    "        self.n_jobs = n_jobs\n",
    "        self.verbose = verbose\n",
    "\n",
    "    def fit(self, X, y=None):\n",
    "        \"\"\"Perform spectral clustering from features, or affinity matrix.\n",
    "        Parameters\n",
    "        ----------\n",
    "        X : {array-like, sparse matrix} of shape (n_samples, n_features), or \\\n",
    "                array-like of shape (n_samples, n_samples)\n",
    "            Training instances to cluster, or similarities / affinities between\n",
    "            instances if ``affinity='precomputed'``. If a sparse matrix is\n",
    "            provided in a format other than ``csr_matrix``, ``csc_matrix``,\n",
    "            or ``coo_matrix``, it will be converted into a sparse\n",
    "            ``csr_matrix``.\n",
    "        y : Ignored\n",
    "            Not used, present here for API consistency by convention.\n",
    "        Returns\n",
    "        -------\n",
    "        self\n",
    "        \"\"\"\n",
    "        X = self._validate_data(X, accept_sparse=['csr', 'csc', 'coo'],\n",
    "                                dtype=np.float64, ensure_min_samples=2)\n",
    "        allow_squared = self.affinity in [\"precomputed\",\n",
    "                                          \"precomputed_nearest_neighbors\"]\n",
    "        if X.shape[0] == X.shape[1] and not allow_squared:\n",
    "            warnings.warn(\"The spectral clustering API has changed. ``fit``\"\n",
    "                          \"now constructs an affinity matrix from data. To use\"\n",
    "                          \" a custom affinity matrix, \"\n",
    "                          \"set ``affinity=precomputed``.\")\n",
    "\n",
    "        if self.affinity == 'nearest_neighbors':\n",
    "            connectivity = kneighbors_graph(X, n_neighbors=self.n_neighbors,\n",
    "                                            include_self=True,\n",
    "                                            n_jobs=self.n_jobs)\n",
    "            self.affinity_matrix_ = 0.5 * (connectivity + connectivity.T)\n",
    "        elif self.affinity == 'precomputed_nearest_neighbors':\n",
    "            estimator = NearestNeighbors(n_neighbors=self.n_neighbors,\n",
    "                                         n_jobs=self.n_jobs,\n",
    "                                         metric=\"precomputed\").fit(X)\n",
    "            connectivity = estimator.kneighbors_graph(X=X, mode='connectivity')\n",
    "            self.affinity_matrix_ = 0.5 * (connectivity + connectivity.T)\n",
    "        elif self.affinity == 'precomputed':\n",
    "            self.affinity_matrix_ = X\n",
    "        else:\n",
    "            params = self.kernel_params\n",
    "            if params is None:\n",
    "                params = {}\n",
    "            if not callable(self.affinity):\n",
    "                params['gamma'] = self.gamma\n",
    "                params['degree'] = self.degree\n",
    "                params['coef0'] = self.coef0\n",
    "            self.affinity_matrix_ = pairwise_kernels(X, metric=self.affinity,\n",
    "                                                     filter_params=True,\n",
    "                                                     **params)\n",
    "\n",
    "        random_state = check_random_state(self.random_state)\n",
    "        self.labels_ = spectral_clustering(self.affinity_matrix_,\n",
    "                                           n_clusters=self.n_clusters,\n",
    "                                           n_components=self.n_components,\n",
    "                                           eigen_solver=self.eigen_solver,\n",
    "                                           random_state=random_state,\n",
    "                                           n_init=self.n_init,\n",
    "                                           eigen_tol=self.eigen_tol,\n",
    "                                           assign_labels=self.assign_labels,\n",
    "                                           verbose=self.verbose)\n",
    "        return self\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "#@_deprecate_positional_args\n",
    "def discretize(vectors, *, copy=True, max_svd_restarts=30, n_iter_max=20,\n",
    "               random_state=None):\n",
    "    \"\"\"Search for a partition matrix (clustering) which is closest to the\n",
    "    eigenvector embedding.\n",
    "    Parameters\n",
    "    ----------\n",
    "    vectors : array-like of shape (n_samples, n_clusters)\n",
    "        The embedding space of the samples.\n",
    "    copy : bool, default=True\n",
    "        Whether to copy vectors, or perform in-place normalization.\n",
    "    max_svd_restarts : int, default=30\n",
    "        Maximum number of attempts to restart SVD if convergence fails\n",
    "    n_iter_max : int, default=30\n",
    "        Maximum number of iterations to attempt in rotation and partition\n",
    "        matrix search if machine precision convergence is not reached\n",
    "    random_state : int, RandomState instance, default=None\n",
    "        Determines random number generation for rotation matrix initialization.\n",
    "        Use an int to make the randomness deterministic.\n",
    "        See :term:`Glossary <random_state>`.\n",
    "    Returns\n",
    "    -------\n",
    "    labels : array of integers, shape: n_samples\n",
    "        The labels of the clusters.\n",
    "    References\n",
    "    ----------\n",
    "    - Multiclass spectral clustering, 2003\n",
    "      Stella X. Yu, Jianbo Shi\n",
    "      https://www1.icsi.berkeley.edu/~stellayu/publication/doc/2003kwayICCV.pdf\n",
    "    Notes\n",
    "    -----\n",
    "    The eigenvector embedding is used to iteratively search for the\n",
    "    closest discrete partition.  First, the eigenvector embedding is\n",
    "    normalized to the space of partition matrices. An optimal discrete\n",
    "    partition matrix closest to this normalized embedding multiplied by\n",
    "    an initial rotation is calculated.  Fixing this discrete partition\n",
    "    matrix, an optimal rotation matrix is calculated.  These two\n",
    "    calculations are performed until convergence.  The discrete partition\n",
    "    matrix is returned as the clustering solution.  Used in spectral\n",
    "    clustering, this method tends to be faster and more robust to random\n",
    "    initialization than k-means.\n",
    "    \"\"\"\n",
    "\n",
    "    from scipy.sparse import csc_matrix\n",
    "    from scipy.linalg import LinAlgError\n",
    "\n",
    "    random_state = check_random_state(random_state)\n",
    "\n",
    "    vectors = as_float_array(vectors, copy=copy)\n",
    "\n",
    "    eps = np.finfo(float).eps\n",
    "    n_samples, n_components = vectors.shape\n",
    "\n",
    "    # Normalize the eigenvectors to an equal length of a vector of ones.\n",
    "    # Reorient the eigenvectors to point in the negative direction with respect\n",
    "    # to the first element.  This may have to do with constraining the\n",
    "    # eigenvectors to lie in a specific quadrant to make the discretization\n",
    "    # search easier.\n",
    "    norm_ones = np.sqrt(n_samples)\n",
    "    for i in range(vectors.shape[1]):\n",
    "        vectors[:, i] = (vectors[:, i] / np.linalg.norm(vectors[:, i])) \\\n",
    "            * norm_ones\n",
    "        if vectors[0, i] != 0:\n",
    "            vectors[:, i] = -1 * vectors[:, i] * np.sign(vectors[0, i])\n",
    "\n",
    "    # Normalize the rows of the eigenvectors.  Samples should lie on the unit\n",
    "    # hypersphere centered at the origin.  This transforms the samples in the\n",
    "    # embedding space to the space of partition matrices.\n",
    "    vectors = vectors / np.sqrt((vectors ** 2).sum(axis=1))[:, np.newaxis]\n",
    "\n",
    "    svd_restarts = 0\n",
    "    has_converged = False\n",
    "\n",
    "    # If there is an exception we try to randomize and rerun SVD again\n",
    "    # do this max_svd_restarts times.\n",
    "    while (svd_restarts < max_svd_restarts) and not has_converged:\n",
    "\n",
    "        # Initialize first column of rotation matrix with a row of the\n",
    "        # eigenvectors\n",
    "        rotation = np.zeros((n_components, n_components))\n",
    "        rotation[:, 0] = vectors[random_state.randint(n_samples), :].T\n",
    "\n",
    "        # To initialize the rest of the rotation matrix, find the rows\n",
    "        # of the eigenvectors that are as orthogonal to each other as\n",
    "        # possible\n",
    "        c = np.zeros(n_samples)\n",
    "        for j in range(1, n_components):\n",
    "            # Accumulate c to ensure row is as orthogonal as possible to\n",
    "            # previous picks as well as current one\n",
    "            c += np.abs(np.dot(vectors, rotation[:, j - 1]))\n",
    "            rotation[:, j] = vectors[c.argmin(), :].T\n",
    "\n",
    "        last_objective_value = 0.0\n",
    "        n_iter = 0\n",
    "\n",
    "        while not has_converged:\n",
    "            n_iter += 1\n",
    "\n",
    "            t_discrete = np.dot(vectors, rotation)\n",
    "\n",
    "            labels = t_discrete.argmax(axis=1)\n",
    "            vectors_discrete = csc_matrix(\n",
    "                (np.ones(len(labels)), (np.arange(0, n_samples), labels)),\n",
    "                shape=(n_samples, n_components))\n",
    "\n",
    "            t_svd = vectors_discrete.T * vectors\n",
    "\n",
    "            try:\n",
    "                U, S, Vh = np.linalg.svd(t_svd)\n",
    "                svd_restarts += 1\n",
    "            except LinAlgError:\n",
    "                print(\"SVD did not converge, randomizing and trying again\")\n",
    "                break\n",
    "\n",
    "            ncut_value = 2.0 * (n_samples - S.sum())\n",
    "            if ((abs(ncut_value - last_objective_value) < eps) or\n",
    "                    (n_iter > n_iter_max)):\n",
    "                has_converged = True\n",
    "            else:\n",
    "                # otherwise calculate rotation and continue\n",
    "                last_objective_value = ncut_value\n",
    "                rotation = np.dot(Vh.T, U.T)\n",
    "\n",
    "    if not has_converged:\n",
    "        raise LinAlgError('SVD did not converge')\n",
    "    return labels\n",
    "\n",
    "\n",
    "#@_deprecate_positional_args\n",
    "def spectral_clustering(affinity, *, n_clusters=8, n_components=None,\n",
    "                        eigen_solver=None, random_state=None, n_init=10,\n",
    "                        eigen_tol=0.0, assign_labels='kmeans',\n",
    "                        verbose=False):\n",
    "    \"\"\"Apply clustering to a projection of the normalized Laplacian.\n",
    "    In practice Spectral Clustering is very useful when the structure of\n",
    "    the individual clusters is highly non-convex or more generally when\n",
    "    a measure of the center and spread of the cluster is not a suitable\n",
    "    description of the complete cluster. For instance, when clusters are\n",
    "    nested circles on the 2D plane.\n",
    "    If affinity is the adjacency matrix of a graph, this method can be\n",
    "    used to find normalized graph cuts.\n",
    "    Read more in the :ref:`User Guide <spectral_clustering>`.\n",
    "    Parameters\n",
    "    ----------\n",
    "    affinity : {array-like, sparse matrix} of shape (n_samples, n_samples)\n",
    "        The affinity matrix describing the relationship of the samples to\n",
    "        embed. **Must be symmetric**.\n",
    "        Possible examples:\n",
    "          - adjacency matrix of a graph,\n",
    "          - heat kernel of the pairwise distance matrix of the samples,\n",
    "          - symmetric k-nearest neighbours connectivity matrix of the samples.\n",
    "    n_clusters : int, default=None\n",
    "        Number of clusters to extract.\n",
    "    n_components : int, default=n_clusters\n",
    "        Number of eigen vectors to use for the spectral embedding\n",
    "    eigen_solver : {None, 'arpack', 'lobpcg', or 'amg'}\n",
    "        The eigenvalue decomposition strategy to use. AMG requires pyamg\n",
    "        to be installed. It can be faster on very large, sparse problems,\n",
    "        but may also lead to instabilities. If None, then ``'arpack'`` is\n",
    "        used.\n",
    "    random_state : int, RandomState instance, default=None\n",
    "        A pseudo random number generator used for the initialization of the\n",
    "        lobpcg eigen vectors decomposition when eigen_solver == 'amg' and by\n",
    "        the K-Means initialization. Use an int to make the randomness\n",
    "        deterministic.\n",
    "        See :term:`Glossary <random_state>`.\n",
    "    n_init : int, default=10\n",
    "        Number of time the k-means algorithm will be run with different\n",
    "        centroid seeds. The final results will be the best output of\n",
    "        n_init consecutive runs in terms of inertia.\n",
    "    eigen_tol : float, default=0.0\n",
    "        Stopping criterion for eigendecomposition of the Laplacian matrix\n",
    "        when using arpack eigen_solver.\n",
    "    assign_labels : {'kmeans', 'discretize'}, default='kmeans'\n",
    "        The strategy to use to assign labels in the embedding\n",
    "        space.  There are two ways to assign labels after the laplacian\n",
    "        embedding.  k-means can be applied and is a popular choice. But it can\n",
    "        also be sensitive to initialization. Discretization is another\n",
    "        approach which is less sensitive to random initialization. See\n",
    "        the 'Multiclass spectral clustering' paper referenced below for\n",
    "        more details on the discretization approach.\n",
    "    verbose : bool, default=False\n",
    "        Verbosity mode.\n",
    "        .. versionadded:: 0.24\n",
    "    Returns\n",
    "    -------\n",
    "    labels : array of integers, shape: n_samples\n",
    "        The labels of the clusters.\n",
    "    References\n",
    "    ----------\n",
    "    - Normalized cuts and image segmentation, 2000\n",
    "      Jianbo Shi, Jitendra Malik\n",
    "      http://citeseer.ist.psu.edu/viewdoc/summary?doi=10.1.1.160.2324\n",
    "    - A Tutorial on Spectral Clustering, 2007\n",
    "      Ulrike von Luxburg\n",
    "      http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.165.9323\n",
    "    - Multiclass spectral clustering, 2003\n",
    "      Stella X. Yu, Jianbo Shi\n",
    "      https://www1.icsi.berkeley.edu/~stellayu/publication/doc/2003kwayICCV.pdf\n",
    "    Notes\n",
    "    -----\n",
    "    The graph should contain only one connect component, elsewhere\n",
    "    the results make little sense.\n",
    "    This algorithm solves the normalized cut for k=2: it is a\n",
    "    normalized spectral clustering.\n",
    "    \"\"\"\n",
    "    if assign_labels not in ('kmeans', 'discretize'):\n",
    "        raise ValueError(\"The 'assign_labels' parameter should be \"\n",
    "                         \"'kmeans' or 'discretize', but '%s' was given\"\n",
    "                         % assign_labels)\n",
    "\n",
    "    random_state = check_random_state(random_state)\n",
    "    n_components = n_clusters if n_components is None else n_components\n",
    "\n",
    "    # The first eigen vector is constant only for fully connected graphs\n",
    "    # and should be kept for spectral clustering (drop_first = False)\n",
    "    # See spectral_embedding documentation.\n",
    "    maps = spectral_embedding(affinity, n_components=n_components,\n",
    "                              eigen_solver=eigen_solver,\n",
    "                              random_state=random_state,\n",
    "                              eigen_tol=eigen_tol, drop_first=False)\n",
    "    if verbose:\n",
    "        print(f'Computing label assignment using {assign_labels}')\n",
    "\n",
    "    if assign_labels == 'kmeans':\n",
    "        _, labels, _ = k_means(maps, n_clusters, random_state=random_state,\n",
    "                               n_init=n_init, verbose=verbose)\n",
    "    else:\n",
    "        labels = discretize(maps, random_state=random_state)\n",
    "\n",
    "    return labels\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Validation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}