Created on Cognitive Class Labs

ML0101EN-Clus-DBSCN-weather-py-v1.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a href=\"https://www.bigdatauniversity.com\"><img src=\"https://ibm.box.com/shared/static/cw2c7r3o20w9zn8gkecaeyjhgw3xdgbj.png\" width=\"400\" align=\"center\"></a>\n",
    "\n",
    "<h1><center>Density-Based Clustering</center></h1>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Most of the traditional clustering techniques, such as k-means, hierarchical and fuzzy clustering, can be used to group data without supervision. \n",
    "\n",
    "However, when applied to tasks with arbitrary shape clusters, or clusters within cluster, the traditional techniques might be unable to achieve good results. That is, elements in the same cluster might not share enough similarity or the performance may be poor.\n",
    "Additionally, Density-based Clustering locates regions of high density that are separated from one another by regions of low density. Density, in this context, is defined as the number of points within a specified radius.\n",
    "\n",
    "\n",
    "\n",
    "In this section, the main focus will be manipulating the data and properties of DBSCAN and observing the resulting clustering."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h1>Table of contents</h1>\n",
    "\n",
    "<div class=\"alert alert-block alert-info\" style=\"margin-top: 20px\">\n",
    "    <ol>\n",
    "        <li>Clustering with Randomly Generated Data</li>\n",
    "            <ol>\n",
    "                <li><a href=\"#data_generation\">Data generation</a></li>\n",
    "                <li><a href=\"#modeling\">Modeling</a></li>\n",
    "                <li><a href=\"#distinguishing_outliers\">Distinguishing Outliers</a></li>\n",
    "                <li><a href=\"#data_visualization\">Data Visualization</a></li>\n",
    "            </ol>\n",
    "        <li><a href=\"#weather_station_clustering\">Weather Station Clustering with DBSCAN & scikit-learn</a></li>   \n",
    "            <ol>\n",
    "                <li><a href=\"#download_data\">Loading data</a></li>\n",
    "                <li><a href=\"#load_dataset\">Overview data</a></li>\n",
    "                <li><a href=\"#cleaning\">Data cleaning</a></li>\n",
    "                <li><a href=\"#visualization\">Data selection</a></li>\n",
    "                <li><a href=\"#clustering\">Clustering</a></li>\n",
    "                <li><a href=\"#visualize_cluster\">Visualization of clusters based on location</a></li>\n",
    "                <li><a href=\"#clustering_location_mean_max_min_temperature\">Clustering of stations based on their location, mean, max, and min Temperature</a></li>\n",
    "                <li><a href=\"#visualization_location_temperature\">Visualization of clusters based on location and Temperature</a></li>\n",
    "            </ol>\n",
    "    </ol>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Import the following libraries:\n",
    "<ul>\n",
    "    <li> <b>numpy as np</b> </li>\n",
    "    <li> <b>DBSCAN</b> from <b>sklearn.cluster</b> </li>\n",
    "    <li> <b>make_blobs</b> from <b>sklearn.datasets.samples_generator</b> </li>\n",
    "    <li> <b>StandardScaler</b> from <b>sklearn.preprocessing</b> </li>\n",
    "    <li> <b>matplotlib.pyplot as plt</b> </li>\n",
    "</ul> <br>\n",
    "Remember <b> %matplotlib inline </b> to display plots"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Notice: For visualization of map, you need basemap package.\n",
    "# if you dont have basemap install on your machine, you can use the following line to install it\n",
    "# !conda install -c conda-forge  basemap==1.1.0  matplotlib==2.2.2  -y\n",
    "# Notice: you maight have to refresh your page and re-run the notebook after installation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np \n",
    "from sklearn.cluster import DBSCAN \n",
    "from sklearn.datasets.samples_generator import make_blobs \n",
    "from sklearn.preprocessing import StandardScaler \n",
    "import matplotlib.pyplot as plt \n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2 id=\"data_generation\">Data generation</h2>\n",
    "The function below will generate the data points and requires these inputs:\n",
    "<ul>\n",
    "    <li> <b>centroidLocation</b>: Coordinates of the centroids that will generate the random data. </li>\n",
    "    <ul> <li> Example: input: [[4,3], [2,-1], [-1,4]] </li> </ul>\n",
    "    <li> <b>numSamples</b>: The number of data points we want generated, split over the number of centroids (# of centroids defined in centroidLocation) </li>\n",
    "    <ul> <li> Example: 1500 </li> </ul>\n",
    "    <li> <b>clusterDeviation</b>: The standard deviation between the clusters. The larger the number, the further the spacing. </li>\n",
    "    <ul> <li> Example: 0.5 </li> </ul>\n",
    "</ul>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def createDataPoints(centroidLocation, numSamples, clusterDeviation):\n",
    "    # Create random data and store in feature matrix X and response vector y.\n",
    "    X, y = make_blobs(n_samples=numSamples, centers=centroidLocation, \n",
    "                                cluster_std=clusterDeviation)\n",
    "    \n",
    "    # Standardize features by removing the mean and scaling to unit variance\n",
    "    X = StandardScaler().fit_transform(X)\n",
    "    return X, y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use <b>createDataPoints</b> with the <b>3 inputs</b> and store the output into variables <b>X</b> and <b>y</b>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "X, y = createDataPoints([[4,3], [2,-1], [-1,4]] , 1500, 0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2 id=\"modeling\">Modeling</h2>\n",
    "DBSCAN stands for Density-Based Spatial Clustering of Applications with Noise. This technique is one of the most common clustering algorithms  which works based on density of object.\n",
    "The whole idea is that if a particular point belongs to a cluster, it should be near to lots of other points in that cluster.\n",
    "\n",
    "It works based on two parameters: Epsilon and Minimum Points  \n",
    "__Epsilon__ determine a specified radius that if includes enough number of points within, we call it dense area  \n",
    "__minimumSamples__ determine the minimum number of data points we want in a neighborhood to define a cluster.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0, 1, 2, ..., 1, 2, 2])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "epsilon = 0.3\n",
    "minimumSamples = 7\n",
    "db = DBSCAN(eps=epsilon, min_samples=minimumSamples).fit(X)\n",
    "labels = db.labels_\n",
    "labels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2 id=\"distinguishing_outliers\">Distinguishing Outliers</h2>\n",
    "Lets Replace all elements with 'True' in core_samples_mask that are in the cluster, 'False' if the points are outliers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ True,  True,  True, ...,  True,  True,  True])"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# First, create an array of booleans using the labels from db.\n",
    "core_samples_mask = np.zeros_like(db.labels_, dtype=bool)\n",
    "core_samples_mask[db.core_sample_indices_] = True\n",
    "core_samples_mask"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Number of clusters in labels, ignoring noise if present.\n",
    "n_clusters_ = len(set(labels)) - (1 if -1 in labels else 0)\n",
    "n_clusters_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{0, 1, 2}"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Remove repetition in labels by turning it into a set.\n",
    "unique_labels = set(labels)\n",
    "unique_labels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2 id=\"data_visualization\">Data visualization</h2>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.61960784, 0.00392157, 0.25882353, 1.        ],\n",
       "       [0.99807766, 0.99923106, 0.74602076, 1.        ],\n",
       "       [0.36862745, 0.30980392, 0.63529412, 1.        ]])"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create colors for the clusters.\n",
    "colors = plt.cm.Spectral(np.linspace(0, 1, len(unique_labels)))\n",
    "colors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD8CAYAAACVZ8iyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsvWmUpNlZ5/e77xZ7RO5L7VVd1dVdvarVqFsSI5pBkoWQ0cAwCBgbMTYjzMDBxmZ8wB8GD7aPzeEcy+IIBIjBwDFiACFoCZoGLUgtoVard3Xta1ZV7hmZGXu8273XH25E5FJZS1dlVVZVv79z6nRXrG9kZdzn3uf5P/9HaK1JSEhISHhrY231BSQkJCQkbD1JMEhISEhISIJBQkJCQkISDBISEhISSIJBQkJCQgJJMEhISEhIIAkGCQkJCQkkwSAhISEhgSQYJCQkJCQAzlZfwJUYGhrSe/bs2erLSEhISLhjePnll8ta6+E3+7zbOhjs2bOHl156aasvIyEhIeGOQQhx/nqel6SJEhISEhKSYJCQkJCQkASDhISEhASSYJCQkJCQQBIMEhISEhJIgkFCQkJCAkkwSEhISEjgNu8zuNuIIsnUhQqNWkC+mGL7rj5c197qy0pISEhIgsGtYnGhybNPH6XZDBGABnI5jw98+BCDw7mtvryEhIS3OJuSJhJC/IEQYl4Icfgy9z8lhKgKIV7r/PkPm/G+dwpRJHn26aNorRkeyTM0kmd4JI/WmmefPkocya2+xISEhLc4m1Uz+EPgA1d5zNe11o92/vzaJr3vHcHUhQrNZkgun1pzey6fotkMmbxQue7XjiLJxJlFDr86zcSZRaIksCQkJFwHm5Im0lo/J4TYsxmvdTfSqAWIy9wngEY9uK7X7aae6vWAdjMkjCTFYpof/vFHGN1WvO7rTUhIeOtxK9VE7xRCvC6E+DshxAO38H23nHwxhb7MfRrIF1KXudew0e6/m3pqtULmpussllvUqz7nzy7xyd94jrnp2qZ/joSEhLuXW1VAfgXYrbVuCCE+CPw1cGCjBwohPgZ8DGDXrl236PIMN0vts31XH7mcR7MRrEkVNRsBuZzHjl19l33/WsXfsPD84KPj1OsBc9N1NIpszu29bq3S5nN/+hr/9hfejZOolRISEq4BofXl9qxv8oVMmuhvtNYPXsNjJ4DHtdblKz3u8ccf17fKwnqz1T7rF/Z0xuXLz5y47Otv9P6ZjEsQxKTTzpogUq+2mZqsEvgxrWZEqS+NsFYSUa1mSLGU4Ud/8m3suWfwBn8yCQkJdxJCiJe11o+/2efdkpOBEGIMmNNaayHEOzDpqcVb8d7Xwnq1T5dmI+DZp4/ykY8+9qZ22JcLLO/9gftot0Ia9YB0xkEgmJmssrzU4vmvnkVYYs37z0xVuTCxzHc9uXJCatR8jh9ZoFJpImONsARBEDMwmMX1Vq7Rda3rrkUkJCS89diUYCCE+FPgKWBICDEJ/CrgAmitfwf4EeBnhRAx0AZ+TG/WkWQT6Kp9Vi/EYNQ+C/MNTp9cYKncYnGhyeBwjke/awfZrLfha10psHzpb4/zkY8+RnVd6qdRD5ibrfPQo9vWvJZtW8hYUav69A1kqdd9XvzmBaRSKAVKaVSscB2LpcUWI6N5olAipaJaaXPk9RnSGZfd+wZwXTtpektISLgsm6Um+vGr3P9J4JOb8V43gyupfepVn9//zW+CWHnE3/31UT76M0+w/75LJ8tdLbBMnF3kha+fXxMslNKgNadPlHnwkTEs29T1vZSNBsLOAn/8jTmUVKQzDlEo6R/IsDDXoNUK8TyHpcUWvh8TxxKhoVEPefXFSfYfHOap9x3ghW9MJE1vCQkJG5J4E3F5tU8USc6cWMASgoHBbO+PVpo/+t0X8FvhJc+5moz0/JmlS3oOvJSN09m516p+7/ZCKY3jWEipqFd9wjDGsgVSaizbIl9IsW1nCdu2UUrRaoaoWOK6NuM7SgwMZclkXM6dXuQPP/UtpFRJ01tCQsKGJMGAtWqf1Zw/u4TUmlJ/Zs3t2bxHux3xyouTl7zW1WSkGi4JFoVSGte1iWNJGK4szH4rYv/BYXI5j/MTy7RbkblfawaGssSRorrsY1kghEBpjUYwNJzv1Q9czyYKJY1GgIzVmvdNZ1zmZut89Yun1khWN7uJLWmMS0i4/Um8iQDXtfnAhw/x7NNHWZhv9NIovh+RzXk9pY5WmiCIkbEmDiULc43ea3Tz8ZWlFjJW1KttCqWVINKVke69Z5Czp0ztXKrOjj+QjIwVOH/WnBrK841eGuep9x3g+efOYVsCNEipCENJuxVSq5pTiGVbpNMOQgj8dsTykqkfdK9bKYVSek2gaTVDTp8o06j7xJHk7KlFrN7j9drC9wcP0m5H11VrSDyZEhLuDJJg0GFwOMdHPvoYkxcqNOoB+UKKi+eXefrP3wAgCiVLiy2UNPv7IIg5dWyOxYUmwJoFL4wkE99ZpFBKk8t5ZHIehUKKD3z4EJmci4wlp48vsLzcQgiBEBBHEsexeep9+/FSDumMedwX/vIwrUZArdZGKo1tCdqtiFYzQgiwbIFtWRw4OMzMVBXLEihpri+dMb0HlmVhWQqvc1qQUnH6RBmNIpV2GB7NUyilee3bkyDgbY/v6NUtygsNPvkbz7FjVx+WJd7UYr7ZKq2EhISbRxIMVuG49hpd/sh4gX/4wnGa9YBGIwQ0rmcRhpJUymVoJM8zf3UEMPXl4ZE8rWbI9GQNIWB5sYXr2qQzivd+8CAAf/WZ1wlDxdnTZeJYYVmCfCFFobPjPvLaDIce3caX/+4Ei+UmM5NV4tjs1D3PxvdXdvdam/TQ4HCOaqVNJuvRbsUoqZDSJKuiUOJ6Np5nYztmga9XfaJI4roWrmtTLKVNrUIAiJ56SUrF9MUaURjjpWz6+rPAtS/mVyumT16oJH0QCQm3CUkwuALZrMdHf+YJPv2b/0S7FWLbFmB28A8/tg2t4ezpMo5tcfDB0TU77mJfhlYzZGQsj+va/MPfHAdM0PA8GyEErmsjY0mzYU4iaHjt5Sm+/fx5hBDUq4FJ2QhAQLsdA3TSORoNpNMOjXpAFEnGtheRsWJhoUGzERAEMY5tMba9yP0PjnH21ALtVkS10sb3I1w3xf77hrEsizDoBpmVdFI3aDiuver+a1/Mb5YnU0JCwuaTBIPL0K0B+O2I733/fl58/iKeZ5PJuZT6s0ycXuypf7SGOFaMjBWIIrnGGiIMJH39WSbOLKKBXXv6+c4rUyilcBwbJ+2itKbVDHnhn86D1mitieOVMrQ2ytMVBFhCoDE1AoQgCiWzokY67XLw0CgPPjKO3464eL6C7VhcOLeERlCv+RRLaeq1gKGRPGEQk0o7eCm79+LddJIJAOaNV+7vXcJVF/Mb9WRKSEi4dSTBYAPWFz0b9YBmM2TfgW2kMg5HXpsBIJvz1uyYJ84s4rhrBVrdRbSroFmYa1CrBkSRQsYKrU3B1m/HyFiCEJ26xOVRUoMFWmuEBShNjCaOJY26xHMtbFdw5jtlUimHoeEcgR/3CsanT5bRSlOeb1Dqz5DJuOw7MNSROmmKpXTv2uNIkcm6vdu6aIwaaeLM4mULy9fqyZSQkLD1JMFgHdVqm//3t5+n3Qop9WUY3VakfyjL4kKTY4dn2bt/sLP79wjDuHMKEGjMjj4KjXwzCmUvHy+lIgolYSSZnaojhMayhNH3C4FSGik7QUVcW2O26qhEdedpGk29EuC4Fuery1z4Ty8jLMjn01ycqGDbAtsx9Q7bEpQGs9SrPvVqgNaaI9+Z4Z4DQwhLsLjYQgBSaVzPZtvOEpa1EuSajQDLEjz/tXO029FlVUKXU2l1H5cUjxMSbh82zajuZnCzjerW2zO02xF/+NvfYnmp1akPgOPYPPL4djJZlzdenSaT8ajXfdJpF9e12X9wCIDTJ8rUa22k1KRSNpmsZzqUNRw7PEscGUloo+b36gCx1Gil2ax/Ate1kEpjCXPa0Bpsx8JxLeJIMTKWp1EzhfBiX4Zs1qVa8RkdzyMsix/8Vw+xZ9/AGkVVJuPypXUGe5cz0Gs2AoQQfOSjj6Gh97PNZF005vSTL6TYsasvCQQJCTeJ29qo7nZkfSooihVHX5/BdW1SKafXtBWGMa+/NMV7vm8f+/YP0jeQ4eypRbJZF9d18NsRGjq6fti2o8SFiSW0humLVcoLTRzHYve+AcoLDVrNgCiWaAW2I5CazuJ9458pisyLKFYVmbVGSYWSiqVyC9ezEYBtC4Ql8FI22ZyH34459sYsjmNdku5ZL7mNY8VXnj25JhB0Tz/z8w2+/OxJ5mfqVzw1JGwdMghZPjJBsFgjNVik/4E92KmNvbZuJ6523Xfq57pdeEsGg4307xcnllZsGVZJYDzPod2KmJmqkc563HPvMG+8OkNlqQ0C6lVTRPXSNu1myMxklXTGRWvN4kITIcwufOpihTCQCCHwXIcgiLEti25fslKbe0Lrvp6MNXQK0FEo0RoyWZdU2vzTx7Hi4kSVMDJKpZmp2iUL93rJ7eFXp9eohLoNbFEkabdDnn36KLmcx/0Pj/UM/ZLegtuDxoU5jnz8s4TLDSNt0xqvP88Dv/gj5HeNbvXlXZarXfed+rluJ96SdhQbzSRuNyPTsdvZysp1RdxqpU0m43L0O7Ps2ttPOu1Qr/od9Y/pJLYtC8exUMoEmUzWodmIWF40zWWOa2HbFpYtsGxBKuNgWeKS99p0Oiu3UhoZK/IFDyEEgR/TrIcIoSkUUuza23+JZ9FGVhKrVUKr5bTZnIttGcms69mcPr6A6hx5NmPec8KNIYOQIx//LFpq8nvGyO8eJb9nDC01Rz7+WWQYbfUlbsjVrjusN+/Iz3W78ZY8GVSW2jTrAVpqvJRNoZQm05GDCiBXMPWDqKO3l7JrGa2YOLPE9t19bN9VotWKcD2bOIpZXGwTxZK04yKlIghi07GrNX47YnGhiW2LjmxUoaQmimRPjWQ7IOOb83llrEmlTW+DbVvUKj6NeoBWRoHUbtts29nXW+C7fQRHXp/h8Gszl85l+ODBnkooCmVPThuFshf0XM+m1QxNA1unWS3pLdhalo9MEC43yO8ZW3N7aqBAY2KW5cPnGHrs3i26ustzteu+8IXn78jPdbvxlgsGiwtNnn/uLFOTVdJpBzA72d37+nAcYxZnAkSKwI+pVwNT+G34vPzCRfx2xPxsnUzWxXFtXMeiPNcmCiRCmMAhhOikhFYKupYtVobPdHfq3ROBuHmBoEsUSVIpcxLJFVI06wGtVjenH3L0OzPMTlW57yGT2tFa89UvnmJwKHfpXIZnTvDeDx7kS8+cYH6+ge9HgMbtpJOmLlZ7j18tvU16C7aWYLG2xop9DUKY+2+Qm5G3v9p1Ny/M3/TP9VbgLRUMurWCbN54BYFx9QzDmPNnKxy4f5iTR+dpNUNarQi0JgxjiqU0ubxHJqNZVmZhbzZDs7sPZS/NozXEkUYITa3qUyimsCxhCsSxwpcapTWOLRAYt1KljALJzDS4eZ9dSZCxMaxbKsdYloXn2YShsaWQsWKx3OLk0Tkefmw7rUaIVObzzc/Ueyeo7qmh3Y74yEcf49vfPM83vnKW4dEcxVKaWCounl+mutxGa/Ozqiy3qFV8vJTN6Hjh5n3IhCuSGixyWema1ub+G+Bm5e2vdt2Z8QHm/+kIKo5xsmnSw31Yjt27/0Y/11uFt1QwWO2Vs/++YU4fX6DVNDMJfN9YR/zK//FfcP7sIksLTaTSnDg6T63i47o2fhybKWMa0EYquRFaQxjEhIGNJQReykJYFlLKzphKSRRKWg0TUDqGpOhuLuYm0fVCUhq0kkbNFGuadWO1obRibqbOhXPLxLGivNA0hfKO8XZXSttN9ziuzXe9azfnTi2adJgfd+oEJhhKqXnlhYtksh6plElFfe4zryeqoi2i/4E9eP15gqU6qYGVoBws1fH68/Q/uPe6X3t9Xn/1ax/5+Gd5/Nc/hu25V3iF67tuYVuUv32c+tlp0GB5DnbaY+Sdh0wa+E1+rsudbN4KSqXNGnv5B8CHgHmt9YMb3C+ATwAfBFrAT2mtX9mM934zrPbKyWY9Hnx0nFrHQrrVDHnX9+xlZCxPFMYUi2nmZmooqYiimOXJlukQxpww4ujKq3YcK4JAUupLEYQKtMa2LWOBLU3nse9HKAnCuvzGZzNRChzXQnU6n7ufwZxMTKCKY8XkhQrjO0pEocSxLWxHkEo5RLHk9Ikyo9sKa7qPH3h0nNdfmuTY63OA6c0YGcvjt2JiqUilHN72ju04jpOoirYQO+XxwC/+CEc+/lkaE7OX7N6vd7GGm1uPuNx1O4UMGrBsh+3vf5z5bx1FtkPCWovJv3uRgUfvYeeHnmT58LlrWrwvd7LZ+6NPce7Pv3rXK5U262Twh5ixln98mfu/HzjQ+fME8KnOf28p671yLMvqFTcX5s1sgj//o1fW2FDMTtWo1wKEJXppSRPbrrx6u67Nu753L9UlnwvnlnBci/nZRud9jeWE6nYPb0KPwbWyOoe/Gq1BS2ORDXDh7JIphPuRscC2zbS3divEb0U8/7VzNJsh7WZIGEkEUOzLMDiYxfNslNKcO7NEKefSaoY06iF9/U7iWLrF5HeN8vivf4zlw+dWdrkP7r2hQAA3vx6x0XWrKObEp77QOy1se+/b8ecr+OUqi6+dJmq0mPnyq8x86ZWrLt6XO9m0Fyp88+c+wci7Htj0E8/txmbNQH5OCLHnCg/5MPDH2rQ7f0sI0SeEGNdaz2zG+18rV/LKyWRcvvPKdM+KGqDYn+bE0TniWJFK273hL8FlFtTVuK7Fk9+9l/mZOmdOlVmYa9x8CekmIKVisdzEsS0yuc5OSmviSDE/WyeX96hV20RhzPJS2wRGoalVfIIgptS3ncFSjsX5JqsD5uoglKiKthbbczddXXOz6xFw6XVPffGlNQHIsm3Sw30svX4GO+WRGiiS320W/6st3pc72ehYEtda6HXT+e5GpdKt6jPYDlxc9ffJzm23lK5XjhCChfkG5fmG8cwRgoce20a7Ha0JEq16SDbrYdsmtx5FiihSiMsaM68wPJpnz74BBoZz1KvtOyIQgPkZCYyPkSl+C/oGsvQNZEllHBqNkMmLFU4eX2BxsUmt2qZWCZCdYvobrxg5qukvWPk5rXY9TVRFdxcyCJFhRFhrUjlxASVXFs7NqEdcjo0CkL9QQfohtufgZFZ+x1IDBcLlBsuHz234Wpc72cQtHy0EcXuDzctdplS6VQXkjVbPDVdHIcTHgI8B7Nq1a9MvZKOJZqPjBZ778mlmJqudcZVpbMuiPNdAo/FSDvlCaiVNZAnmZupXzBQ16gGnT5R5/ZVJwuAW5oFuECkVwjIqoG5ntBDmpDMz2eoUz81jhYDQj7Edi0IhRRhaaK0Jgoi52TqOY9FqhLie3XM9TRxLt57NLIauzrML22b59bMsf+csfQ/swc2mN6UecTk2KizHLR8VxrilLOmRdb9jV1i8L3eycbJphNZrAkuPu0ypdKuCwSSwc9XfdwDTGz1Qa/17wO+BMaq7GRez2l5hcaHJn/ynlzj2xiz1mo8lBFFsRkRmcy7tdkwUSaRU5pTgmALSlaoGtmPspT/+v38F4xR05xDHxmivfyBLsxESx5JmI6RZD5CKngy1+wNQCoiVSaWlHKTUhIEkjhSDw3maUjEwnGOp3EocS28DNlP+uT7PrqTEK+WonZnCL1c59B9/iOEn7l8TCDYKRMB1BaeNCsthpQkCRp48hGWv+x27wuJ9OcWScGycYhax7ve1Xa4gw5j27BLlV07eFeqiWxUMPg/8vBDiP2MKx9VbXS/YiCiSPPO5w0ycWTQGbJGi1TTHwbDTQWy+L0Yq2TVe8/2o083LmiE0YL5ftm2xWG5vwSfaHKLISF9zORfLgmzOZXmp1Yt+Wpuid9fVuut7lMq4lPo8in1pWs2QPff08/DbtjE9VUMAu+8ZYM++wSQQbBGbKf+UQci5v3yO5aPnKewZw1+qUX7xOLIdghDE9RbHf+fzxhqiE2Q2CkTYltlXSAVCoGOJUpId3/8EfYd2byjtdAoZBGaXH1QbjD71CHHTJ1XK4w0WOPuZr6DXpWWDpTpOIYOKYqa++NIlQedKSqt3/dZ/z7k//2rv9qjRonZmmuK+bUz8xXN3jbpoUyyshRB/CjwFDAFzwK8CLoDW+nc60tJPAh/ASEv/jdb6qt7UN9PCOookL37zPF/8mxM0GwGl/jS1is/cbN08QAOWmSiWSpmYmc44aAX1uo9WGtdzCPz4lshCbyVdR1PXsUz/g9Y9Oazn2Wb3H8o1R6NMxqFvMEuzHlLsS9FuxSa95Nls39lHKu0kzqVbTPmVkxz7zb+6pEgK0JiY5f5f+KFrKoZ2F/XlwxPUJ2Zxsyna8xUyY/14JSO+CCsN8nvGKew1KiC05qVf/jQyitGRJG4H2J7LwovHEZbFzg89SdxoM//8UaJaCwQMvv1e0oPFNdLOqOlTOTKBimPA2P0KzyG3fQi3kOXAT3+Q0oEdHP/tp9cEHWGb32U6QedyC7gMow2VVt3b27NLnPmTL5HqK5IeLvWeZ3oexG2hLtpSC2ut9Y9f5X4N/NxmvNdm0LWvvnB+mbnZKlGoaDYCpFI4timcGhsJ4yWEgFTKZmSsgEDQPBUQRhpCubIe3uSGsVuJ0hrXsnBcm0zWpVb1jRxWmRqCbQtSnk0YGStuYRnZbrMe0tefxnIswsB0YEulmJms8ejj2whDmfQYbCFvRv55pear7umisG8cf6GCcGx0LPHLVdxiriO9hsxo35qibXNygebkArIdmsbMSp2w0iA9WKI5uUD12AW01qRH+ggrDZxMChnFfPNn/x+KB3ZiuRbV4xewMi7BVB0tNLmdIzTOztCeXsTry7P47z5Bfs8oY089ysCj+0kPl/AGzGlBIC5pWlt/Irqc0qp7e/mVk9iuuyYQwN2hLnpLdSDDiiVFqxVSr/jISKOUwm+b3a/jGk296uyGESZlpGJNrRKglEJgFsDuLGKlb66VxK1GK4gjSRxJmo3AnA46i75gZRiPEOB4FvlCilJ/hspyG42g3YqQUvXko3Ec8NpLUzzwyHjPuTTpMbj1XKv880p1Bb9c7UkwVSyx0x5BpWF23pEkbvoI28JOe6RH+mhNlgkWa8gwonL0PG4+i532aE6ViZtt4nZI4+Icwee/iVfMktvZ2aV3FDy6oWhMzBG3AhACv1zBchxQ5qRaOTwBSqOVMo+PFf5ChcXXzpDdNsjIOx9gz7/8Z8S1Vu9EpKTEn68QtwPCSoPFV04x8uShK/7susFx6h9ewl+qk90xvGJ50eUOVxe95YLB1IUK9XrA3HQd17NJZRxUy0wcA+PfAxBHxrpBYJrMbEegUdTrAa7nEMUhlmURb8ZUmtuQ7pQ0KZVpihMmIERK4TjG6tuxLQ4+OMb3vv8AR16fQclFxneWOHN8AUusMuZjxep6bFsh6THYIq7FjuJqdYVt73t773RhOTYj7zzE9JdewY9i0NCeW8ItZBl82wHzxE6QqRy/gAoirEGH+tkZpB+YBb7b1V9pEFWbxO2Qwt4x0BrLc5j/pyNoNHbaAwG26yCjGOmH6FihpURYFjpekbPKMMZNewgsFl85ib9YwyuZ5tKw2uh1Kms0wVKdb//Sp3j4f/nXbP++xzYcllM5dp7Jv/02wrEIKw2qR8/jzy8z8s5DOLkM/kKFuOUTVpo4hcxN+/e72dy1wWD9SMvu9K5GLaDdDHu2y4NDOcrzDVrNqLdp6k4MQ3TmGkjN4GAW27LQCnJFjyAwBm+2JfD9m2w5ukVIqVZOPHolw6CBbMYjlXYoFFK8/MJFmvWAej0gOLVoPJDs1ekIQSrtEEWSZjNMegy2iGuxoyi/cvKKthJBpbHmdOGV8uz40JOc/6tv0JoqI6wMlutSOTJB5egEg4/dS/+De2nPLGGlXNpzy4TLdWQYrxnvJ1wHLSVRrUn93CyF3aMIIdBKIWwLy7V7v4B2yiVutE0jmCVM8Xk1UhE32rTFMpmRPsLlOiiF2iWZ/9ZRUGBnPJqTC8StgNb0Ei/90u8y/d7HePCXfnTNsJygXKP8ykm00nh9OYa+6z7Ty9AKmfna61iOgwojVBiDgLOf+Qq5bUN3ZCH5rgwG60darpY05ospwkjSXeVcz2ZsW5F2O2Jupo6MFZmMi9IK2za2CkEQU68H5PIeUkosITj08BgXzi5Sq96du1zd6SdYnWLurgHdWQyWLThzqszjT+5kYDBLtdJmebGF78c4rk0cK4ToeC8pTaMekMu7iXPpDXCjPQJXs6O4Wl0h1Ze/5HRhOTaWbeFkUqT686uf0GswyowPUDq4i5mvvoqKZUeSRq/WpjsnCyUjIt0ku2O41zOAhtiPzGvbFiqMVwLAZSYE6k6K1y9XSY2WiBtt5r7+BsFSncxIH/WJWVQkcbJpcjuGiGot2vNVjnz8s7zt136qdzqysyks18HryyPbIeUXjzP8Xfcx/8Ixaicm8QYKeMUsbinLyJOHegN1bodC8pvlrgsGG420hJWxiz/8E49QLKZZLrd694nOdK6BwSyV5TZe2sw27qphRsbztBqRMbGLjQ31zGSVWjW4YzqLr4ZlmY2abcOqBtINU8xaQxDEZrBNKHntxWn23TuIjM1QnyCI0UoTKI3t2FgWlBeaaKWxLCtxLr1ONqtH4Ep2FFe1ix4b6J0u6memCapN2gvLyHbInn/1HnRscvdOJmVqBhcXVoziMi5eKY+KYmRr3SZq1VuqIOLC09/ATqeQkSQ1kKc9t4Q/Y9xx43ZwdWdHpbEcG9kOqZ2apu/+XdTPTBEs1/HnzbQ9J58ht30I0dFIW55NuNxYMyyndmZq5eeW8QgrDWQQMfDwPlozixTuGae4bxvpkb5eX8OdWki+68ZebjTSElbGLs7N1PnhH38E17OpVdpmdkHHxnpoOE/fQJb7Hx6jWEpTKKYZGc3jeY6xU7BMD4GMJIvl5l0TCGDlxH4to5gti566CAFRHPP6S1PYjsX2nX2k0w7ZnMfQSA7bFhSKKfoHsoxtL7G7g9GEAAAgAElEQVT/4NCasZoJ18atGlm5uq6wmtV1hfyuUe77uX+BkgoZmNM3lmDx5VM42RT5XaMoKVl85STLJy4y9/U3kEFI6f7dqFiignhN+nFDFCZgRDHBQpW43iIOQrRUFPaOk9s1etXVK6g2if0A23MYefIQ29/3drJjAzj5DAjI7xnFSrlEjTZhrUVrqkxjaoGF54+iOu6RTja99kU7hW0ZRDgpj+K+bWTHB9c2uN2hheS77mSw2qZ6PV2DtD33DPLz//49fO5PX6NeC3Fdi2zOQylFOusyPJwnl/M4fXyBdjvqNJ0FVJZaxsHTWjWl7A7Hda1ejcSy6TmpdnuC1iMEvfGZSkpkZ0ZCHEuU0ngpi76BHJYwNt5xLOnWDPbfN4xlWYlz6XVwK0dWjn7PIxz9xOeQfkhqoEBqsEh6sMh9/+7DRms/s8SZP/ky3kCe0r07aM0ssvDCMZRUXHzmBaJGu9crAILD//dfcOLTf0tmbAAVx2j5JjcBSqOVREcSZVv4CxVG3vVAJ3d/+TSt9kNwbAYe2odl22TGBkgNFBB1h6jWpDVVJm76qEgi/cDUFoCo1iJutEkPFEkP92Gnvd5ph441RVhrYrnOpZYXcMfaVNx1wWC9TfVqVhukjW4r8m9/4d2XeBR97jOv91xNH3x0nPnZOqdPlAmDmEzGXTOP4G6gFwgsY1IXKInV6aJ2XJt269Idp+PavWAhMIEBjBIrDM3P6dBDI1yYWMb3I0bHC+y5ZwDLWtnKJc6lb45bMbKycWGOV//jH7PwraMd8zlJa7rM4Nv2c/CnP9hr5PKX6lSPTpAaLjHy5KHeglg/PUVUa3YaNgUIgeU6xK2OcsjqXOeNfHeUwi1kqRw9v7Hj2TqEJWjNLNH/sMRybEr37eTiF75F3PKJW75JEWlwS7lOB3RMVG+hIsmFzz/P+D9/lKF33Ef528dNeklA3A5ID5UYeechompr0wcFbRV3XTC4kk31eoO01R5FXT7w4UM8+/RRFuYbaKU5f24Z2xb0d3z6lxZbBMHdpx5Syuzk0fSUQJ5nFn3TZW2+wWYWg+o1nxX70sSRQkrj6OqlXPYfHMJxHQYGc8xM1xkYyq4JBOb9jMXH4Ven16i9EjZmsyyir9RM9sZv/BmLr5zEzWdJD5vvSdwOqJ+b4/lf+CSj3/1gz4PIzqaJ6z6Tz36b/of20ZxcIKqvWJagNAiNloqo2cYrZGlMzN14P44Gf6mGQBg5+FWaPbVSqDA0C7ljcfFvvoWMIoRtm6K1JdBIwkoDy3Nw0h71szNYro1WmgtfeJ7UQIG++3dh3zPOjg8+Qd/9u41CanbppgwK2iruumDQtanuLujr1URX63xd7Wp65oQZi5nOOExdrJJKO2axFAJ5ld9qxzOy066P/51wkpAdnyUv7RCFkiAwjqSptLHdSGdd2s0QpTW2BV7K1AaKpTRKKXbt7Wd8WxGrc1KwbUEm42Lba7dw5QWTItKAbYk1/z5JUXljNmNk5dWayRrnzWJtZ4w6SStlOosXqxBJpB92bte0ZxcRto0MIppTL6KDqHNyWfWL3lUJxRBbvnn+tRSlroKsd3y/XPvqwUUqWrNLVE9cZOmNs6YpzuoEEsc2zXJBDK6FW8gSN320VEitEQgzrjaIULHin/3u/4iXz/Ze+mYNCtoqNsWb6GZxI95EcSTXpIB27Op70xYIh1+d5vnnzmHZghNH5nE9kzapVdvE0WVSRcIsgumMg+s5NGoBcXRrG9McRyClfnMBqPM9FsLUDhzHJp122bajhAYunl8mCmI0pl4iLDP5zLYtXM/mX/+3j/Pyty5eIud94rv38MI3Jnq3K6WZvFBh195+hobXqr2EEIlVxRW4ETWRDEJe+uVPo6XeYI6wYNv73s6RT3wOf76CV8oh/ZDmVBkdS6JmGy0V+d2jbPu+x8zM4fNzRq3jh8h2iPSDqy/Mm23ZIoRJR63vM1j/no6FEBY6VmALhNadlFDneVojUi5ONkVcb4EW2GkXJSW5nSPYrktYb/LY//Zv2PeR773t3Um31JvodmSjFNCbJV9M4Qcxi/MNGo3AjKrs+vM4xsp5/YLrOlYnvWIRtONeR/MtQZg0jus5EMYmIND5Dq67zvUFYgHkCh7pjMvSYgslTZqovNCg3YoYHslTWW7jusa7KQwlSioO3D+CUpp0xr1kTkQ3AO+5Z6B3e63qo2FNIACSovI1cCM70fUF6PWWDAOP7sdOuaCNtUNzqgxKYadcVBSjohhhW8z+42vIMCY1UDA2E0GEiiIuu9KvvvlyqoTrRZuG0Cs/BpAabXesdiPQAlMkW3UtOoqRgWlgE7aNwvRBBOUalucQNX2Of+rzLL16+rLBdzPnRGwFd20w2AxGxgvMz9SJwhjbsgiDCKVUJ53SNeihZ2intcm7S6nIF1PY3sYF2JtFp2aHZXUtNCCTcXE92/gGrRg2YtkCrehdt+NYlPrSVCsB6bTDrr0DDAxmqdd8Tp8sMztdQ6NpNsznTqUcAqWZOL3E+I4ijXpw2QC8+vbDr0735iyvJykqX53rHVm5ugC92pIBIYjqLSasr5MaKNK8OE+41EDH0gSCWCJsC4GDbAW05pcRwsLNpVcav2zLeFbFG6iEVq/VW5WF0NqcHlbLWdfXG5RG+aFJbcUSrYwKTkURdi6F7dqkBoqXbSrbzDkRW0USDK7A/EydoZE8p4/PE0USxxGEwVoJpiU6ensuFXt0+xduBV3JJxosy8JLAZ3xlXGkSKcdwtB0T0eRQiuN6swkEBbc+8AItmURx7rTlV2gUEpz9nS5l9cPAmWG+yA6xWKbWMZMXayQzlxbnvRa1V4JN87qnWp7oYKSEiVXLBm8vs7pTGucdBrLthh4dD8zX3mVqOmjOr0LwrZI9edpTZVRkTQ2ESkXy3VQsUS3JWK9adtqttrR1/iwr7vt0gtySnniaqMTLMwT43aA9EOcXJrstsENpbxX8nN64zf+jHv+q/cS1Vq3/WkhCQZXoFELaLeNA2dXS48whSXLMjvqbq+JWVhXcvXtVrTicXSTWH3StR3L6Pw9m7hzvd1B9gD33j/KyWPztFrR2i+mgNHxIvccGKI812TBbuB1xlTWqn5niI1FrJUpHCOwLOPZpKTCsW1UZ/LbtfBm1F4J18/6naqOJctHzpkO2nbYCwRxO8BOexQPbKN1cYH7fva/ZPx7H+XYbz2NZdu0pss4+TSW5+IVc2YeQSw7xWV6pwO9ak7AJVytwWyr6fzyqqCT7lr9GbQppMt2iNMxu1sv5b1cD4iwBdN//yqN83OkB0u3/WnhrutA3kzSWZfZqSpCCNIZF89zsG0bIQRK0UkL6VUpUYHjWgiLW1IrEJZxBnU9I4PrGgoNDGR56NFx7ntgjFJ/hkzWo1EPUEqbHgJH4DgWqbRNKu0Sx5LzZ5eYmakRhZJ9B4awLIswkLiebSy9pSKVsgET8KRUSKnxUjbbd/XRbl9bOqyr9hJCsDDfoDzfMKovIZJxmJvERt3KhXu2MfjIfqpHzxM12gSVBq2ZRaJ6i8I928wThSCqtRh55wOoKO7dH5RrNC/ME7cChLDIjA0Q1duESzVk0zcpl44XEBulAAVg38ZLjcYETKVM38Hqz9A59WulufjX3ySsNi6R8m7UA6JiyfzzR9GYE9jN6hjfTJKTwRXR2LZFiMkJWZaF1qbpbDW2ZaSmUirs7hg/bVRFN9OyQkAnPWPmDkupGBjM8bZ3bMdxzD/t9l0lJs4uUau0yeVT7NzdT9gpLiulsC1Bq2WM+AqFFO1WxIkjc9z/8BheyvQZFAopolCilEkNqc7J48B9w+zdP8hiufWm0jur5bs3ovZK2JjL7VTzu0dpzW6neXGBYKmKwMIt5akcmaB+dprcjmGcQobjv/00g4/sZ/5bRwlrTZRURI02OpZ4nkNrcmHj+sB62WgnAHRnHWx5uuhKKI32N1igLWGuXyqCaoOLz7zAtvc+1rP7Xj4yQe3UJH65Snb7UG/Ggb9QQfohtueYzuUOt/MQnE0JBkKIDwCfAGzg97XW/9e6+38K+A2g6/r0Sa3172/Ge99M/FbMrr0DnD6+QNCOO+6bGC22s1KA7S74VqcJq1BMMTNZ7Zw4Ndc08uA6vihxrKnXfCwhKPVnSGdc9h0Y7AUCwKiglLnG7hyCbn5fKc3MZI1U2mFsvMjIWIGxbUWOHZ7ljVen2b1vgDhWOI7FE9+9mxNHFohljGPbeCmbvfsHabei60rvbIbaK2FjLtetHFYbNCZmaU4tGgPBlEtYqZPbPoRshdTOTKNj2QskwraY+drrBJUGbjFrbKM7arprQmnT1NV1Kb0T0dpIULVGNn38SNKaKrP48sneKE6lFfWz0zQnFxh/6hG8Ur7nuOoUMmilqJ2d7pn33a7eRTccDIQQNvBbwPuASeBFIcTntdZH1z30z7TWP3+j73cryRdTZLIu73j3Lo69MUcYmrm+zWaEAPIlj6Btdjz5fArLhn0HhkhnXFrNkLHxAhfOLVO9FpvrN/Fd6X7PvYxNynMYGsozNJrjhz7yMP/4D6d7zXZ+EDM/UyeVdmnUgo6jqGRgMIvr2WaCm1IIC2PEB2RzHm97fAfnJ5Y5cN8w73rPXr7zyjTtdsT4jiJTFysorRkcybNYbl1zM1/CrWOjbuVu2kL6EemhEiqM0LFEtgPqZ6fJ7xknv2uE5e+c7f2CZcYGcNIeIeBm0wgg9kP0taY4rkX6eTvTq4EY5ZE3WCRVyhMs13nuv/4/GXj7AUr37sSybdIDRWa++hpTf/8SQ48fJKw0UbGR5JZfPNF7LSvl4hWz1E5N3XYF5c04GbwDOK21PgsghPjPwIeB9cHgjqNb7NRa8/YndlKr+iwttnrdyDt29TF1sYLn2Xgph1YzIp12sWxTY0hnPXbuHUCeXaTZCFf0/h0JqFbXt2EyjWGClGs6omOpcBybMFJ85KOPMXF2ibMnF3j5hUlGxwtksi7tVmgk11KxtNhiZDRv+iSUJpVyKJZW3Bkt26JQSDG2rciDj27jgUfGeymddMZFAO12lKR3blM26lb2FypEtRZ22sP2XNy+XMekzUwN6394X6d5a6WAajk2xQPbac0uIf0QFSvcbNrUCd4KrP5yujay4VNfrBEuF4nqLZZeO0PzwjwjTx7CyWcYevwgy4fPkd0+wD0/+T6+/T99iqjZxkmnsDrdzvWzM6A06dF+yt8+flsVlDcjGGwHLq76+yTwxAaP+5dCiPcAJ4Ff1Fpf3OAxCCE+BnwMYNeuXZtwedfPamuLxcVWpzErRS7v9Tpo+wYynD6+QK3SRghBEEocZTG+vcTF88vEkSLuSDqFZQq3USR7YzYvizC1iG5zWzcVJYRJsWRzXqcmoVgqN+kbyNCoB1QrPi98fYK52TpL5SbNRoDjmEaxQjFFvXNCKC80TRc1cPDQ6CXeQatlnklK585io4lmzckFEDD4tgNUjkwghMDNmxGNYaWBCiMsy2Lgsf3UTk32Akl6uI/s+ABxK8AKI9y+PEGlYdRDV+r8vZtwLBzPM55GQhAFAcKxzShOBbNfex3hOqggImr6lF8+ReXoBaJak7jpE8iaseQOAtx8hvRoP24uQ3Z8sDdO9HYYhrMZwWAjVeH6le4LwJ9qrQMhxH8H/BHwzzd6Ma317wG/B8aOYhOu74bYqNiZybh86ZkTzM/We95F2byxXti2o8jzX5sgnXYY27aLWtWnXg04f26RVjMy8xCkwnYshGVmK5sB8qz5qTmOhW1b5HLmhLG81MJvRXgpm1w+1UsV2bZR/XR37d3BPpmM0zt9hIGZSuY4NplOCstvR2SyLql0mlPH5rn/4TGyWXNcTWSedz7ru5XbCxWmnn2R/O5R6menke2w50EEoEJJerzI0GP3kh4s8eqv/hGVExdQcUxzsozATCsLl+om3eEYIcGG3M6F4mtlTee0hXAsRCDx+gvIdoCSEst1sFIu9bPTpMf6jXwUyO0aZeYfXyOqtSgd2o1sBYS1JkFZIzwH4dlmQA+3V0F5M4LBJLBz1d93ANOrH6C1Xlz1108Dv74J73vL2Ghn/N4fuI/PfeY1wsgoLDI5j4vnlimV0rTbUW/KWl9/lr7+rFH1nFli284SE2cWsSyLWrWNZQlm/BpCmPZ3pTVK0jkNKGxboHRn7oAj8DzHqIBsy6RkY0UcGWsIGZsZw7mcx8WJKvVagOMY7bRSivsfGmN6sko+77F3/yDbd/bh+3GvYLxv/yAIkdQB7hJWdyvLIGTxxRNE1RYjTx4ySqFKAxnGZqEfKfHAL/4I7dkljv/201i2he04VI9dwOvL4ZVyHRO7tBlmKQSFA9toTS+tTRsJNt92YjPpykavdjLvBQIQrt0rCluOTSyEabKzBHHL70lSuz0bApNq1YBsh7iFLCqSxKk2xJ2ehVUKo25BeavtLDYjGLwIHBBC7MWohX4M+InVDxBCjGutZzp//UHg2Ca875YRRZIv/e1xCkWTV+/SbAR89YunyGzQjWtZFsVSmofeto33f+h+/r/ff5HyfAOgYwdtkcm41OsBvd9EbXoXfD8yPQKW6DmIRlFMGBi5J9oohJ79/HEAZiZruJ7ASzl05bF+W3Hy2DyWEGRzHmPbjbvo6oLxvnuH2H9wOKkD3IWsTh2Fyw2K9+4gXG5gOTYHfvqDbH/v20Hrnpld4Z5t2NlFmlML2CkPjab/oX3IdkBZaNpzy8StgNRgkXYQrUhNLev2Th+9WddUIbCEIFiuG8uKyNjXj777IfyFCv5iFdlxdHUyKUbeeQi/XMVOe1i2TdRo4+YzWK4NQhh7D8taOxRHm9d96Zc/vaV2FjccDLTWsRDi54G/x0hL/0BrfUQI8WvAS1rrzwO/IIT4QSAGloCfutH33Uq6ozVXz1gGsyDXawGVpTaua1MopXuDX2AlDz84nOMHfvgBPv8XbxAFMUpp4kgSRZ0B8qs2Vq1maGYIaI3r2Ma7R0AUSmzbBIfBoRz33jfC3GydMyfLZNIu2YLHwKBl5i/4Mb4fgW86ldNZj2NvzLP/4BDZnLemYJzUBu5ermZ0V37l5Jr+hKjeIm6ZlEhYaSJsy+TJhbEnUX5E2DSeVaJjJ62vSUd9B6FBRjGO7aG1ghislIs/v8zQE/dTPTGJimNG/9nDZMcHsGzbNOdZgvRQCWFbhBUzG0VLCVIx+NiB3pjMYKmOU8gw9eyLgLjEzuJW1hM2pc9Aa/0M8My62/7Dqv//FeBXNuO9bgc2Gq3ZaoacPlGmXmsjpabdishkvd6Cuz4Pv3vfAKNjBaqVNvV6SCbrUFlqI6XCyZgAEoaSXN4lCk2hN5v1CPwIxxE0YuOcOjpa4MD9I1i2xch4gZPH5mm3I7IFD9ezGR7JMT1ZxfNs0hkXyxIUSynCMOb0iTIPPjLWO9ImvkB3P1cyultvZrf8xjn8+QoyCFFRTHt2yQgNLIGTSZHfN07l6ATCsbFSHt5AgeaFeTPb4DbFHcjjDRRpTszC1VwCbMv4zIQxyrJID5dQHcfWsN5m5kuvMPbUI+R2DOFm070FPj1segmslMP2738H4WLduMMuj9GcnEPYVmd2hNn9b3vfOzj/51+7JSNNr0TSgXwdrDdbk1Jx+kQZjSKdcdm+s4/52TrtVsDh16bZfc8ghUJqTR6+q1R65nOHiWNJvRobOwvHotSfobrUJpvzKBRTVJbagCCbc0mnHFIZB8tuY9uCnXv7yeZMXtG2LMa2l5ibrvZM8sLAGIuN7ygSBkbFFPoxSmuajTZTFyv09WeSgnFCrz+h25NgpVzk6pnFlpl7bVkOUTugeuw8MogRliBuBYS3YSPVeqKlBlHdKP+uOn5BCPAcdCxxMh5Djx+kdHAnwWKtZ/29/6PvJzs+eMnEs6HHDqABf3rJBAbLorh/nCc+/rOmmLzqZDb7tddv+kjTayEJBtfBerO1etUniiSua+G6NiNjeUbG8tSqPgtzTR5+bBvveNfuS/Lwg8M5fvy/eZzDr8/wtS+eJo4l87OCdjMEAQODWTM3QGkc1yKVcmiriGzOJQiM/Gj94Jy+/gyFQqrnYlpeaCClwm9HpFIO23f18Z2Xp4ljY19x+mSZUl+Gj/7ME0md4C1Otz+hdnoa6YdGPuk5RN3BNZ18uYpiEAKpFMTSzAe4TevFGxJdbU4hxobCtY2EVoPXnyc1UKA5uYDtuViegwwi5l84xsg77mf8fW8nrDRIlfJkxgd6k+euZfbEZo00vVGSYHAdrB+tWVlq4fsRrpti/33DPc1+X3+WOFIUS+nLLrSOa/Po4zt4sNPYNX2xyje/dpaZ6RpRJIlCsysbGMwiLGMfXerP0qhHtFtBr3MYTAG7UEjx3g8e5At/eZjTJxbwWxGtVojrOvQPWFw4t8zAcBatNK1mxK69/fT1p3nhGxPsuWcgCQhvYbpF5hf+h98iarRNQ2PH49/sLlgpwCqF7tkT3WnR4Cp0ynw6VqYGYlk4+QzlF0+gohh/oWKs4tMu9VNTnMx49D2wBzebvmQG8rWkdzZjpOlmcBtbCd7edPsP3v+h+3j7EzvZvqPEg4+O97T6Xa41F9+Vr77rqX383L9/D/vvHWJ4JM/BB0YYHjW/IGEY47o2/f0Ztu0s4noOYSAvcf4cGM6RSjns2tPPwQdGGBrOMzyaQyrFUrmF59rYtkWhmGL7zhKFYoZmM2TyQuVm/KgS7iAyo/3s+hfvJjVcIj3WZ2pjluh10NKxRu+t/dZtLCO9XhTGk6jT/OPmUqgwxi3lCGtNLM9FW4JwqYHbX8DNZ2lenCe7c/i6XEm7QVjYgsbELI3zczQmZhG2WBNYbjbJyeAG6H4FBodzpDMuzXpAoZTp3X+9zVuZrMeP/uTbefbpozSbIYMjOaY6C/X2nX0sLhpPoJ//n7+Hdiu8xPlz4syi8RLabppgin2mSzoMJFEUU6345PLemlNMMmXs7uZaNOzVU5O8+qt/hD+/TDBX6cwsMGmS3qx021qxrAZW5qre0o+zMbbYHC8kC9xSnlR/nsHH7mXx5VP488bOI6q3sDMp007h2Qg0dsYjrDTw5ytkxwevq+h7IyNNN4skGFwniwvN3mItMP5Ap06UGRkrkE47vYHw19u8tb7zOZNxTUHqGjyB1qudslmPBx8dZ+pihdMnyoyOF9hzz8AaC4pETXT3ci0jGaunJvnHj/ya6UxOe9hpj/bsEsKx0aEyu2QhsGwb7Qg0xoTO8hyE5yA7aaUt5UYCgb3SH+EU86QHS3h9eaJai7DWREuFimJjI6M03kCBYL6CiiQqjgmWG8x/+zjFfeM46dR1FX2vd6TpZpEEg+sgimTP9mF1r0GplKLZCHnHu3dT6s/ccPPW9XoCbTRa0rIstu3oY3qyRibrrgkEif3E3cuVRjJ2Nexozav/6x+i/IjMaP/Kk4WgNV3Gcm3slItwbSzHwevL05xZJFqsmel3UiFcBx3GW/AJN4nuSUdAZrSP7OgAUbNN+aUTxics7eEWMohO93G4XEcLjQwi6q/PoMKIsFKneWEOgWD8fY9t7ee5DpKawXXQbTpbPbYRoFDMYDs2pf4Me+4Z3LJi7Gq102r8dsT+g8MUiulkythbhO6gm9WFSQC3lKV2ZpqTv/8M5/7yOfzZCpa3dm+YGiwiLIFWGhlGSD8iWKpTPXYB2fDxhkpYgArjjj3DHb6cODY4FlZ3HojSHYt3CxXFuPmsaa6jM+ZTWzQvzqGlwk57eP0FYzNhwdFP/CVho7WFH+bNk5wMroONms663A659/Vqp25KN5fz+KEfe4RSXzqZMvYWYaNBN2G1wfy3juIvVJGdBrH6uZk1xnUrGAVbaqBIsNxA+qFJF7k2SEXfQ/vM+MswwuvLUT0+2RkIvkn5+83GXmeXYXcGiSsNWmFZNqpT/FWRxELg5rMdDyKJV8zhL1RQsSQ1WKR1sdU7OSk/RLg2xf3bCRaqXPj8N9n/E+/dog/65kmCwXWwURqmy+2Se7/aaMnEduKtwXoNu5KS+W8dBQVuPkOhkzpqTMwRLNWxc2ncrJltEbd8dCzJ7xo2AUUZAzYnn0H6ITqOaU0ukN81Qv3sDCqITbrID1cCwe1SXIZOgFoVCEycw06nOsohgeVY+AsVtNYIjJeQk0ux8weeQAYRcTtA2Ba101MoqYgbPtltg2ilsFwHJ5c2DW1C0Lwwv1Wf9LpIgsF1sL7prMvtlntP5hAkrNawu6Us/3975x4k2XUW9t/X7+mZ2Xns7uxLO9JKlmytLBtbiyTbBBxbdgmFsgzYQQlViAJKIbGrAn8kKKUqoPgnEFdSKQoHI55OisImDsICi8gyxoYklpFs6716rFZa72qlfc3uvKd7uvvkj3O/uad7+jXTM909M9+vqqu777197+k70993zve88sIpli5MkxzIkM7nVgqmZceHoFKhslikGNn+C5dnkVSSsXddy/Tx75MeyVNeKqgMJbt3lKWLM0y/fMa3dy1XGJqcYPHsRcpLy17AJpNQr19yL3AuVk4JWYmKKi/GK/lEZgAEKoVlJJkkNzHCxAfeSXbcJ34Vp+d48xtPUykukxn35ayXLk4zeGiPr9sUIc4xODnR5S/YGVvcyNcb1AwjImZ7N/oajWFfnl/k1F/8Axe+9QKL56+wdO4y5eUSpblFEskkE7cfJZFJkds3yuDkPnITo4wevZr9H3wXiXQqMg0Fc0cRUvkc2bFhkoNZUkMDpHflGZqcIDWcx4lDkklvoUr2iZiJOlgCcfJcIoFoDgUgyQTivJM4MzLInlveTiKdYu71t5g9+SZvPPok4Dj00WPs/6GbSQ3mKC8WmH/j4kqRvuLlOVK78kx+7P1d/4qdYCuDddLKDGMY/cLAvjHSgznGbr6O5dl5Zl55g4H9Y1QKJc4//gIH77iFzMgQu2+5gUN3/iADe0fJ7t7FrusP8b1f/RMW3poC50iOdKYAACAASURBVEjlc0gqSWmpuGISWbo0w/L0Aq7s+ykvnDrP8vyiNxOlE6QGslRKvsQzFReX5G1WSjrhs50llepK0TtXKkfF9tLk9o7iSmXG330do0evZuH0Ba776Q+TSKe48PhxCldmGX/3dStF6Q7/s9s4/dePU5yaZbpQRLJp0gM5bvzU3cy8fKavehy3wpRBB5gZxtgKXH7+dZZnFhh9x2EqpXLUSGWZ1EB2JVkqmc2Q272LI5/4kapEp5t++RM8+5kvcnm5ROHiNJmRQV8mYTjP0vlpClOzuOUS6bEhUtmMb5azVPBJahXvs0ikUyyeu8Ly3AIi4gvgzS42HG8ymyE5kKW0VPB1j0TW3ocgJNHk85EPwS2XSOV9mHhuzwijR6/2Al+E5ZkFDn3kGIVLM+SeGFlRBADp4TyDkxPMnTpPMp/BLft+Bef+73NMfe9EX/U4bkWfrN8Mw9gswoiiRCrJnh98O+WFAnOnz7M0NcPMybMNSx8MTe7j1v/8r/nBz/wrRm48zPC1Bzn4ofcyfN1B8ofGGX3HYdLD+RUBqV2/vCD1xRJFI3SiVpGhbb0eyXyWRDpJbu9YtAIZ8JIqmYD0GkVWwvcel5qwWV19rOBgeXqe0vwiw9cdjM1JQaG4Rs54kQT5/WNkhvLkD+0hMzLYUXmKXmErA8PY5oRCrDg9x8UnXkKSQiKVorJUJJFK8o5Pfbzh7DWZSXP4rts5eMctVeUSFt6a4oX/+hfkJpYozsxHEUZlH4OfEAb2jJK/ai+Lb02RGRlkcalApbBc5bCtx/LcIpnhAXJ7dpEezLBw9pK37UePSqnYOkJJIDU8gCRTJLO+wmhZClSWl33toYrz5SsSCQTBLZcAR2FmgSvPv8bsq2cZvfHqqkJxtQXlls5fobxYJJFN+8xkkZV2lp2Wp+gFtjIwjG2OCrHFC1c4/y2fOZ8ZHSY7OsTg5ATD1x3kxc/+ZcvZq5ZLOPSRY+x57w3k94+TzKZJpFMMX3uQ/FV7yR/cTXpXnlQ+RyKTpnBhmvRw3m8fGYpyGepk6UTmI73O/g+9h/F3X8fIDYdJDQ34FUU24x29jWr/66kyKd9sZ2yYVD6LiPjVSyaFpFLRqsA7iRPJpA8LzWWQdAqWSxSuzFO4Mselp09w46c+vrJaqi0oN3PyLEtTMyzPLpDdOxL3V46+jza972ZPgk7YEGUgIneKyEsickJE7q+zPysiX4z2f1tErtmI6xqG0RoVYsUrsxQuTlMplihemYMETNx+lIE9oxQvz3H5udfWdN6xm65h6Jp9vrfBUpH00ADZ8V3kD0Z+tFKZ0mIhCtNMcPiu26iUyiDBtD4h/pEUEgMZkvksYzcf4cozr3Hum08z9cxJUvkc+QPjjL3rCMlsmuRgrlrwKlHz+tTgAPmDu0lm0uR2j/imPIVlEqkkOIeIkBzy+RSJTMorjFSiyqSVzKQYPnKA4sx81SW0oNw1//xHcOUK4iCRTrNw5qIPqS34plI4Fze972JPgk7o2EwkIkngs8BHgDPAEyLysHPuheCwnwcuO+feJiL3AL8F/FSn1zYMoz2GJvdx7U/fwfLsIpnRIVIDWXITo7EzdB2z12Q2w83/7qcoLxQ4/60XKE7N4pwjmU1z+MduZ/i6g5z7h2cZvmY/uYlRZk+e9f6CcoVKOapj5FxsAnKQn5xg7J3XUJiapTS/BED+4G7KxRIXv32c1PAAyVKFxcWCn8qWnH92vjdxciDL4R+91ZtuyhUuPvEiiUwKJ5BIp9h7+1GKswvMnzpHbu8o5UKR2ZNvesd2xvsz8gf3kEgmmX7x+yyeu7z6izvH2ce+w+iNV1Ne8M7yzEie5ZkFZl9/i6HDEyQHMuQmRrvek6ATNsJncCtwwjl3EkBEvgDcDYTK4G7g16PXXwJ+R0TEue1WCN0w+pf8/nFy47tW9doF1j17HZrcx/t+95e49L1XuPSdVwDH7ltuYPd7rufyc68xd+Is+QO7qZTLXPruK0giQW7/OMUrc5QXC5Gt3Wc2Z8aHede/v4fTX3mc+TMXKC/60hcLb1z0SXK78hRPz7E8v+Tt/NrDWKtpO8iN7yKZTfMDv/YzDOwfj8b1MiDsvuV6dr/neuZOnfMVWhcKLC8s4kplEjkfwZRIJVeyiAuXpilcnl31nbXe09A1+5l431HOf+sFlmcWyI4Ps3D2IksXZ9hz7AYWTl9Y1eymn9kIZXAIOB28PwPc1ugY51xJRKaB3cDF2pOJyH3AfQCTk5MbMDzDMGB9HbXa6YOQzKSZuO0oE7cdbXi9cqGIKzsSqSSuVCY9NLCyAihemWNg/25u/+1PUZye58oLp0gP5cmMxhWBSwtLTL90mkq5EoWnpqiI+ESvaHUxsG+U8fe+jeRgjoH94w3HNXL9VXzoz3+N7/7aH3P5+PcpTPlM60QqyeDhvUhk709k02SDMazcryA6KzMyxMGP3MLS+SuUFgsMTR3g0EdvYfTGq3vSk6ATNkIZ1PPm1M742znGb3TuQeBBgGPHjtnKwTA2CPUd1DZvbzR7bacPQrvXm339LSoV3xfBOcfwkf0kUikGJsYAYfiafey99UZe+9I3qRSWSe6t6RhYcT45DLzdP+ELzCWSCR/OmkyQzKTJjg5TmllYid5ppMx2ve0Q/+SPf4XXvvRNXvrcX1G8PIckhfJicaWnw+iNVzOwf3zV96oNMU0kk+QPeD/JXOIt9v3QzX0fOVSPjVAGZ4DDwfurgLMNjjkjIilgBJjagGsbhrEG2u2o1U4fhHZmvHq91770TU78yaPk9o4ye+KNFaELfqZ4/S/c5YX5yBCJbJrSYiF2wOLDTREhkc9QmS9QXi6tJLY5IBmFd5YWCySSSQqXZloqs2QmzZGf/GEuPfESlaJvUqPXTaRSJDLJuqulfulZvNFshDJ4ArheRI4AbwD3AP+y5piHgXuBbwGfAL5u/gLD6A3tdNQK7eKVUpmlC1coLSyRyucoLxTWFDcfCl1Xdgxfe2DFrFIplhmYGOHQHbcAMHBgnNEbr2b+zAUf8aTnSKdIDedJZlKUE0lcucLy3OJKhJD2aE4NZKkUlkkP59tSZuHqpTTv6zRVCsuk8tmGtv61rrC2Ch0rg8gH8GngUSAJ/JFz7nkR+Q3gSefcw8AfAv9DRE7gVwT3dHpdwzA2D7WLF6fnOP+tF3xtoYjKcokrx0+tyRQSCtCF0xf8LD+RIHdgV5UAHbvpGgYP72Vg33jVTB0RLjz+ApVyhcqFaSSdJJlJ+QQ3ERKpJJnRISSdJJPPgrCizEKy48OrksDW03+4H3oWbzQbkoHsnHsEeKRm268Gr5eAT27EtQzD2Hyyu3fhSuUgSS12pC6eu8yZr/wjRz75wTUJv3YEaKOZemZsiPf/7i/x0h88wrl/eIbF81dw5QpuuUwilyYVVUxNplPc9MufYPql042T0+qE0a6n/3CvexZvNFaOwjCMVYzddA2VSpnlmYWVngcA5cUimdFBJJVYV4mFdgRoM6Wx+5YbuPS9V7jw7eMsvjlFdmKE3O4RcntHGdg3tnLc0sXpKidvFVskCazbmDIwDGMVyWyGq370NqaeetXb7iO7eHIgw8TtRylemd/UEguNlEajcNFatquTdzMxZWAYRl1Gj17N7ltuIDWQXbHda9Zy8fLchs2u28llWCvb1cm7mZgyMAyjLmM3XUNu9y5c2a3E0cPGzq47zWVoxnZ08m4m0s8RnseOHXNPPvlkr4dhGDuWzRTW5UKRJ+//fVzZrTLlSFLazmUwqhGR7zjnjq31c7YyMAyjIZs5uw5zGULqhX8am48pA8MwmrJZIZRhjZ9VbJEeANsJa25jGEZPqK3xU4WFf3YdUwaGYfSEMPwzxMI/e4MpA8MwekJtG8m5U+eYe/0tJCkW/tkDzGdgGEbPsPDP/sGUgWEYPWW71fjZqpiZyDAMwzBlYBiGYZgyMAzDMDBlYBiGYdChMhCRcRF5TEReiZ7HGhxXFpGnosfDnVzTMAzD2Hg6XRncD/ytc+564G+j9/VYdM79QPT4WIfXNAzDMDaYTpXB3cDno9efBz7e4fkMwzCMHtCpMtjnnHsTIHqeaHBcTkSeFJHHRcQUhmEYRp/RMulMRL4G7K+z64E1XGfSOXdWRK4Fvi4izzrnXm1wvfuA+wAmJyfXcAnDMAxjvbRUBs65OxrtE5FzInLAOfemiBwAzjc4x9no+aSIfAN4D1BXGTjnHgQeBN/cpuU3MIwdSQWYB5aBNDCIBQcandDpf8/DwL3R63uBL9ceICJjIpKNXu8BPgC80OF1DWMHs4SfS50BzkXPr0bbDWN9dKoMfhP4iIi8Anwkeo+IHBORP4iOuRF4UkSeBv4O+E3nnCkDw1gXFeA04IBc8HDR9krvhmZsaToqVOecuwR8uM72J4FfiF7/P+DmTq5jGIYyD5TwCiAkjV8ZzAPDtR8yjJZY1VLD6Dua+QOWW3y21X7DqI8pA8PoK5bw5p5SsC0FHMavBlrV+bc+AMb6sPADw+gb2vEHDOKVQ+0KYDnaPtitwRrbDFMGhtE3qD+gdnafjrbP43+yhwHBryL0IdF2+0kb68PMRIaxbjY61r9df0AOuG6Dr23sdEwZGEZb1Ar+JPAGjW3762Et/oAEFjVkbCSmDAyjJUvA96PnCt4kU8AL41DwL+Nt+9exvll66A8IBf9m+wMsm9kwZWAYLagArwGzwbYysckmi1cO0Hmsv/oDTlOdTawrjs0Q0K2il4ydgikDw2jKLDCN/6moMHbEq4MCq4XmemL9w9m5Fv8ts7kz9droJaXTFY6xFTFlYBhNmcULy1Ao6usKMIdXDBmqVwiNqGeSKdKb2bllMxsxpgwMY11UokcRmMErCLX5N7Lt1zPJJPErgCTdn51bNrMRY2tAw2jKMP5nUsavEErE9nwh/gmV8KuIQ9T/WTVKKFNlUjsvC3ML1kIlGsdU9NyscJ1lMxsxtjIwDKBxRM0wsAvvN9BoIo0o0v36uogXwmVW2/n13AlgAb8SyETvdYWRrTOudmfnlejaZ6P3qWhczcxNvYpeMvoRUwaG0TKi5mrgObwABy9kk8S29V342f0SsAi8FX3uWiAffWYeP1NXvwJ4RTAQvdaVRxEv2BPR+3Zm5xr6OhWcV8+9iO918HZW/9x7Eb1k9CumDIwdTjsRNWVgCC8cC3gBmwmOm4qey9GjQCz8b8YrhEusFu7l6Fz6eorYrKMKIUlzdPzF6PhUtE0VUyZ6/1L0XWpXCJbNbHjsL25sI9ZiL1dm8UJTfQEqRCt44ahCUvBmnGG8wC0H1ywE78ELfRXIz+IVgZpswjElo+vmgzEoKbwCeqPF99CIoPCnXIjGoOMjOqZR8xs1d40T+0iMnYatDIxtwgI+OayIF7Ip/Cy3WXjmEnCKuABcMdquYaLqRwht5wKM4H0IJbyiUB8CrBa2eo083px0mdgPoLP+XdFY1VGt/gTBK4nz0ft6s3Y9l24rE/s1HLFSWIrOa+GiRn06mgKIyCdF5HkRqYjIsSbH3SkiL4nICRG5v5NrGsZqFvA2/Vm8cFyKtqmpp95sWM0raorR2bXglUI52vYmfkUQlo1O4WfR+WhfktjGX4sK9ALeryDRQwX1IF7Ia66ChpoWo+vN4n0QjXodq9lJndGl4PvqtXR889E46rGeVZWxneh0ZfAc8BPA7zU6QESSwGfxPZLPAE+IyMPWB9nYGLRcRIXYjg9eoC7gBfYsXhiGNnE1r+Sj12qj15l1CS9IZ4AX8U7k81QL4iywBy+gi9RHFcRCdO3wJ1fCK4rd0f6Z4DPqf0hEn6vghboqsSPESkZNWiPARapXKmreUkURmqJgfVFIxnak0x7IxwFEpNlhtwInnHMno2O/ANwNmDIwNoAwZDNEZ/sFfKRN+K+eIjaVCF4h6EzcBdt1tq+z873EoaO78Db9UzSfRQtxWQkVyBCbcdSWr0rG1TyXo+vPE68gknhfhI7PRft1pVIkVgK6clAFoX6L+ehxkTjKSaOQRogb6lhJip1CN3wGh/D/VcoZ4LYuXNfYEahztx4OP+MeYXWkkDp1IZ6x64xcP6uPEnCh5jOXou1J6puHwjGooM0RrwYcscmmUvNci65atImNHpeiugxGBa+whHgVsEBcTiOJ91m8GYyD6DsPEK+MpvFmMI2KMh/DTqClMhCRrwH76+x6wDn35TauUe+X2vDXIyL3AfcBTE5OtnF6Y2eTjh7qyA1nsWXimX/tZ0r4f8Nlqh3GIfpeo4PUpq9dxtSk1Grm7IgT1sIkslr/QTN0RREep76FLHGdpKui52niexJ+5kz0rCsEXV1owT1dUanZKxxvq1LXVgp7K9NSGTjn7ujwGmfwxkflKmIDZb3rPQg8CHDs2LFWvxBjx6MOWPUDqADUGf4u6s9HBD/7ncULwgzxTLkRqkBqVxDlhp+oRgVsrYAWWiuDeuiqRX0LZXw+QZ5YiDc6b6hcdPWifgmo9ldA68Q8K4W91emG2n4CuF5EjohIBrgHeLgL1zV2BJpFm8LPkEPH6yCr4/d1FqwmmiN4B247wliF5lpJ1owtvJauOjohVEhLeL+GhtdC85956JuovU9akqJRXSX1K5Ra7LfIpK1Ap6GlPy4iZ4D3AV8RkUej7QdF5BEA51wJ+DTwKHAc+HPn3POdDdswQnJ4oZ7AO3X34C2bI3jhfZl4Fj0FXMErgyngZXw0UoH2fg4aabTW8TWKNtosNMqoXeWlJjNNWMsQl6TQyKva0hhqbrvYYv9ai+0ZvaDTaKKHgIfqbD8L3BW8fwR4pJNrGUZzFvGCq9Y/MIZ39mqmMfh/+xH87PkcXgjqzL12FhuacFShrIUkXsCu9XMbQT0/QyOSxGW0D+NNaKocWxXLW2qx30phbwUsA9nYJjQSOCm830Br9AwQO4zDBKwwpLTWjLNe1A7fK9dXq+vqqiEBjEbHH8KvrEJaFcvL0dzfYqWwtwLm6je2Cc0EjhDbscOexWWqwzLVQbxRrNcxvJE0MhMlgucR4nIX9cpoh6WuQ7TU9Z4W+60U9lbAlIGxTWglsGpj5dVxqzP3RmUaOqEfHKeDeD9KjriqKcS5Dxn8d7+CD0mtVyVVnfQaUqsPIXbeN9tvYmYrYGYiY5vQqjZ/hlhZCD4Ov0wcirpdmSOuW6Q5BfqseQ7g700e70zfi78nGqE1SOtS11YKe6tjysDoMpuZmNRKIB0CTuIVgc6M1xMqutXQ8hT6c9ewUc2e3hW9nomO1fsDcXa2moJCBVH7d9NS2MZWxJSB0UW6kZjUSCAt4XsDhAliGj3UD7b9zUZDR/W76t9ARcAMcb6CVk4tE6+yLhMXvBuKXltC2XbC1nBGl2iVuLQRpppGZZjDa6epLhPdTimI7YaW1AB/DzQnAWKFEJbCVrSQ33x0jCWUbSdsZWB0CU1cqp1Jah/hTguiNVt1LAfXLtI6Ln47o8ovQexA14S4cs1zLWFtpWniEiBmGtoO2MrA6BKtEo/Wm5hUwQuml/DOUjX7ZIlXHRoppLPasJLoTkQjqNRMptE/0J4PRR3Qc2xOFJbRC2xlYHSJVolH60lM0tXAHN4sBF6YqdNzhOqGLpoJHPYB2Kno7F+jjFRBqE+hUec23abVTdebWW0VTvsNUwZGlwjzAELBv97EJPUDaEN6bfqiFUuT+BXDALFdfI76JSd2ErXOcsHfI3UWZ2kebhtGYIV5C2vBKpz2I6aKjS7RKnFprf+K6oMIy0jodXRbCR8lM0V1eeudTO1sX/sgQNzyUjuq1TMZhb2iB6mfsaznbeXMtwqn/YStDIwuotVFL+IVQY44fn2tqI9By1XrjDcsL1GM9i1ixdLqofetgu+IpmU7xvAZyYtUm48yxMpEFUcWL+xDc0+R9pz5IRsVSGCsF1MGRhepNQ8s4gXJeswDamrS8glZ4hDJsLlNElMEjQgd6Xq/BvBiYYi4dpMqXIiFucML/qfwmcv6d9CwXa2ASvC503hF0wz7W/UKMxMZXWKjzQPqgwgb1w/glUQGrxy0tWMvK4f2M2F9JvD3MRvs005nqhzUhJTBZy2r4FZfg4buzrB6nqm9DVo5nK3Caa8wZWB0iVYNUtbaACX0QejsXxvNa2G2DHGilBETFqML702O2NymiWdhe081FWkZC61fFOYqJGre16IlLazCab9hZiKjS2xWngFUmyRS+C5nOeAZ1qYIdkJZCqg2oYU+FvAKdYHY1KYO42xwnIaj1jtnouZ9LVrGolFBQZuf9oqOlIGIfBL4deBG4Fbn3JMNjnsdbxwuAyXn3LFOrmtsRTYizyCMTdcuZVDtcFzGR7BAbLZoJuC10Y0qgu2eRBUqPE3QSwITeME/TbwK0FDdBF5BaDvMEVYLe11t6P2svefhzD+BVTjtPzpdGTwH/ATwe20c+0+dcxc7vJ6xZek0z6DW+axmofGa49LECWh54kSz0FatZo483qQEXugtUm0aUSXRj+GOKapNOO1Se+wg3qmbwEcQqfBXhZDF3yetfKoO5mR0nK4cMtH5SsQtRZvN/K3Cab/RaQ/k4wAiO6EMsNEZrfoNNJsV1jqfIc4gnsYrBKk5Xmex6gANyy/k8EIN4pkwrC7O1q9Zyin8d14iFshzTT9RH0cs5B1xeGmYO6BKMh9dc5b475eLPpMjXlGl8OHDGWzmv7Xols/AAV8VEQf8nnPuwS5d1+gr1tsApV6Ru7CMQpFqAabnE/wsVQW+mjJywDXAGeIoGYc3L/Wj8K+lgv8uA9FD6wutpzTEMl6QaxiprjhCU5GurAbx5qTw7zdAnMdR+/e0mf9WoqUyEJGv4T1ytTzgnPtym9f5gHPurIhMAI+JyIvOub9vcL37gPsAJicn2zy9sXVYj3mgnnNZu3eFZh09Nhe8TuNntMXokQTejhdgZbzwcsCFYHza/CXMV+gn57LDj3+JOLEujPBZyzgr+FDQJP6eamN7CZ7DjON6fz8T+tuBlsrAOXdHpxdxzp2Nns+LyEPArUBdZRCtGh4EOHbsWL/8+oyeUs+5rLP+S1Q3YVGzE6w2SQ0Q9+wNFUyR6kzmemUr+ulfUc07Fbydv5Y0a4vOquAVgXY4CxWBKsT9mJlne7PpZiIRGQQSzrnZ6PVHgd/Y7Osa24lGzmeHn/XvwwusWjNFM5OU5iZojX4VfKEwDFccapLqB9R0A7G5TL+/ZgdrJnA7OLwvQH0kuhrQpD6IHcTGdqUjVS8iPy4iZ4D3AV8RkUej7QdF5JHosH3A/xGRp4F/BL7inPvfnVzX2Gk0K3I3iV8hjOPNFYmazw3X2beED0tdwAvBBWK7uTqRawVpWPeoV6iiCjOrNUs4jNJRc1GjYnMhqkjq9Xgo4ZWKlqYIaVSIztiqdBpN9BDwUJ3tZ4G7otcngXd3ch3DWL/zuZYK8H28mSVP+0Xs+sFMFArrSs22NLG/QJ/VlxBWd61HGKKqvgd1HmvobbgisxLU2xEzAhpbiEYz/bUwhW/uvoAXajqD7veaOEniiCk1b6nQ17DOIWIFqdnF4Yy+3ipBVw96L8Pw2hL+XkGcB2IlqLcrVo7C2EFUgDfxgisU/vWav/cTKfx4NZFLZ+5hkT5t8KORUOCVQwGv+JLEfoRw9aArhkYrh1oFstm9rI1eYSsDYweh/Y/Df/uwFEU4o25ma0+02L+R5Ij7Oau/QDN+1akbmnlSxD0itJqomnz0O4bRQvrdQ7OTJutl8XkFEBcS3MwaU0YvsZWBsYPQ0hcaGRQKf0VfN7OxaxZzoY1jQ9aaq5ChuvuYKrE81WGfmn08gJ+VC77ExFQwRg2d1baWJaqjhjQTmWDbSHBNNU0ViZ33GVYrxX43txmNMGVg7CC0naNmJasDVLNuh/DCrlEmbyjMWxXACwmrfGoeQ7vjTUbXWiYW5HPR+UaIE+lq/QeqsLTxjJby1jLf4WpoJHp/hfj+6IqC4HyvEmcsF4gVhob9WgnqrYwpA2MHofkKmp+giVsVYBSfWFUEXieu0hkWqgvDLzWxq1UsvxZ3SxCHr7aDmnWK0eeS+IYymmSmdZl0haJht1qkL00cHTSAVwrTxP0EwvyKWfyKIhU913YoSxJXgtUGQqpMp6LPpLES1Fsb+8sZO4gwX6FAdWXOq/Gz3BxeaeTwQlBXE4oEj0bN4EPUbFO7GtDOYLU/wdBMtUjch1hn7oPB/mV8KQmisWrbT92XJjbbqJNclYyWrR6LxrIbeCexryHM5RgnTmqDuFDeaPTZcXzYr4WVbmVsZWDsMJrlK1TwM94CcWN3ta2HUTeOuHRzKwrEtvUwizlJXAq6HD2Wo2N0tu/qPObxAlzNRer70CijEHUyL+EVnSoDrUMUKp4MfvVQ797UK3mhylA/a/PKrY4pA2MHUq/YmiZSaQE4tYFr2WsVnGrHD804zXwHambJBOeC2A+ghA7tMGms9txq3slGYy0Gn9Vx6jVmiM1F6oRuNHvXWb/eG20kdKWN72hO4+2AKQPDqEqkGiKuBqoCMIVfCaTws/JZYrNMMyGpyqQcnE9XCLqy0LLaGvKpLTwrxCYqNf+Ewl7t/S74vL4O+xhr8bm9VJeaVuo5fmszjB1+FSKs9ieY03i7YMrAMKoSqULhHhavU7PKErETuFX/AF1B6Kxc3y8TO6g10iiFd85qOWk1IUG8CtB8gKXofLuIcyfqlZ0I8wkuA1cBb9C8uVC9RkKKFrOTBp81tjKmDAxjVTlr8PZzLVw3gBe8S9HraWKh3ix3QCOOVLiH6LZBvLDXBDIt+haafbRgnGb5qlLR/AIXfKa2/4LmHGjfh1b1nRplGGtnuHHipDTrXradMGVgGFVmk7AshfoFcsQRRMPEzlst7VCvvLXOnvNUz8TD/Srgw6ikoWDfDPHqI4kX1Oro1RBZra+k4aJ6bheddw+xWUod1M3KRTTLgdAVUm3faWM7YMrAMKr6JYQzNo9bTAAACapJREFUXc3a1f7A6gAuE9vxIZ6Vq4DW91fhhf2rVJuUtOicHqeoY/gIXvCfIvYdaC9iobrv8xhwkepaQ2rbH6s5fzuO3lbHmLN4u2LKwDBW8g+06qaWXVBFENbrfx0/e5/At8rUNpFhVq/OoPdG55giLvegkUVLxKGhtV3atGRGGP0TZg4X8cpiMDpfHq9EktFY1YSjq5O1OHobNRIyZ/F2x5SBYQDV+Qfz+Nm22uRDB682w8nhE7UuEcf4qwLJ4vss689rkjg6R0tGpPFJXppfUGuDrzXX6PnVqT1LrIjUT5Ahtu0Xgs+uxdEbKsZmjmZju2HKwDBWUHv6MF6onqK6SugSsX+gSFzVU8tA5KL3e6n+aa2nMU9te081EdUWl1OfRm2mc6NWoO2wUY2EjK1Ep20vPyMiL4rIMyLykIiMNjjuThF5SUROiMj9nVzTMLqDCtIB4qJttf4EiGsP7QKuBQ5Qf4611sY8oblGE8s09FVzCNRPECqGNLGvo5MmQBvRSMjYSnT6F34MeKdz7l3Ay8B/qD1ARJLAZ4EfBY4C/0JEjnZ4XcPYZOo5SrXsgpZ1UDbDnh7WUVokFv454vBOZREfeTRF7Ki2vgLG2uhIGTjnvuqc0/++x/HhE7XcCpxwzp10zhWBLwB3d3Jdw9h8wpm5IlTXMQqLuW2GPV3NNQeIVx/j0euwH0GaeDUyzepObobRmo30Gfwc8MU62w/hvVHKGeC2RicRkfuA+wAmJyc3cHiGsRZqHalanA586KfWGtpse3oC74cIs3+1ZHXYk0CPLRKXvWaTx2ZsJ1oqAxH5Gr7Qey0POOe+HB3zAP6X8qf1TlFnW8OCLs65B4EHAY4dO7aWtlCGscHozHwK3zsZ/E/mInF0TTfKNtcqJu2LUNv2UmsWZfBhrzrebo3T2Mq0VAbOuTua7ReRe4EfAz7snKsnvM/g/xuVq4CzaxmkYfSWy3gncm3c/Wm8sujGzDuM8JnGK6ihaF+RuP1lijjvoBfjNLYqnUYT3Qn8CvAx59xCg8OeAK4XkSMikgHuAR7u5LqG0T20Vk+tDT4dbZ9f9YnNQyN8DlLtN9CEM/UVhM7lXozT2Ip0OlX4Hfx/52Mi8pSIfA5ARA6KyCMAkYP508CjwHHgz51zz3d4XcPoEq2icnoRtRNGGqkTWxPQtJ9xLRZdZDSnIweyc+5tDbafBe4K3j8CPNLJtQyjN/RrrZ7axLAi3nTU6Cdt0UVGcywD2TCa0otaPdplrFX2b1iBVEtUWE0hY32YMjCMpnS7Vk9tl7HwWs0igqymkNEZpgwMoyXdqtXTqMtYuxFBVlPIWD+mDAyjLVo1hdkIGnUZ0w5n822MoRvjNLYjNmUwjL6hHyOXjJ2CKQPD6Bv6NXLJ2AmYMjCMvqFecTywiCCjG5gyMIy+oV4y2WZWRTWMGHMgG0ZfYRFBRm8wZWAYfYdFBBndx6YbhmEYhikDwzAMw5SBYRiGgSkDwzAMA1MGhmEYBqYMDMMwDEwZGIZhGIDU72HfH4jIBeBUk0P2ABe7NJz10u9j7PfxgY1xo+j3Mfb7+GBrjPHtzrk1J6r0ddKZc25vs/0i8qRz7li3xrMe+n2M/T4+sDFuFP0+xn4fH2ydMa7nc2YmMgzDMEwZGIZhGFtfGTzY6wG0Qb+Psd/HBzbGjaLfx9jv44NtPMa+diAbhmEY3WGrrwwMwzCMDWBLKQMR+YyIvCgiz4jIQyIy2uC4O0XkJRE5ISL3d3mMnxSR50WkIiINow5E5HUReVZEnlqv93+Tx9fLezguIo+JyCvR81iD48rR/XtKRB7u0tia3hcRyYrIF6P93xaRa7oxrjWM72dF5EJw336hm+OLxvBHInJeRJ5rsF9E5Lej7/CMiLy3z8b3QRGZDu7hr3ZzfNEYDovI34nI8ej3/G/rHLO2++ic2zIP4KNAKnr9W8Bv1TkmCbwKXAtkgKeBo10c443A24FvAMeaHPc6sKcH97Dl+PrgHv4n4P7o9f31/s7Rvrku37uW9wX4N8Dnotf3AF/ss/H9LPA73f6/qxnDDwPvBZ5rsP8u4G/wLd5uB77dZ+P7IPDXPb6HB4D3Rq+HgZfr/K3XdB+31MrAOfdV51wpevs4cFWdw24FTjjnTjrnisAXgLu7OMbjzrmXunW9tdLm+Hp6D6NrfT56/Xng4128djPauS/h2L8EfFhEpI/G13Occ38PTDU55G7gvzvP48CoiBzozujaGl/Pcc696Zz7bvR6FjgOHKo5bE33cUspgxp+Dq/1ajkEnA7en2H1TeoHHPBVEfmOiNzX68HU0Ot7uM859yb4f3pgosFxORF5UkQeF5FuKIx27svKMdHEZRrY3YWxVV07otHf7Scjs8GXRORwd4a2Jnr9/9cO7xORp0Xkb0Tkpl4OJDJFvgf4ds2uNd3HvstAFpGvAfvr7HrAOffl6JgHgBLwp/VOUWfbhoZMtTPGNviAc+6siEwAj4nIi9GMpB/G19N7uIbTTEb38Frg6yLyrHPu1Y0ZYV3auS+bfu+a0M61/wr4M+dcQUR+Eb+K+dCmj2xt9PIetsN3gaudc3Michfwl8D1vRiIiAwB/wv4JefcTO3uOh9peB/7Thk45+5otl9E7gV+DPiwiwxjNZwBwtnOVcDZjRth6zG2eY6z0fN5EXkIv8TfEGWwAePr6T0UkXMicsA592a0rD3f4Bx6D0+KyDfws6PNVAbt3Bc95oyIpIARumdyaDk+59yl4O3v431v/cam//91Qih0nXOPiMh/E5E9zrmu1iwSkTReEfypc+4v6hyypvu4pcxEInIn8CvAx5xzCw0OewK4XkSOiEgG78TrSqRJu4jIoIgM62u8Y7xu5EKP6PU9fBi4N3p9L7BqNSMiYyKSjV7vAT4AvLDJ42rnvoRj/wTw9QaTlp6Mr8Zm/DG8rbnfeBj4mSga5nZgWs2G/YCI7Fc/kIjcipejl5p/asPHIMAfAsedc/+lwWFru4+99Iivw4N+Am8Deyp6aNTGQeCRGi/6y/hZ4gNdHuOP4zVyATgHPFo7Rny0x9PR4/lujrGd8fXBPdwN/C3wSvQ8Hm0/BvxB9Pr9wLPRPXwW+PkujW3VfQF+Az9BAcgB/zP6X/1H4Nou37tW4/uP0f/c08DfAe/o5viiMfwZ8CawHP0v/jzwi8AvRvsF+Gz0HZ6lSVRej8b36eAePg68vwf38IfwJp9nAnl4Vyf30TKQDcMwjK1lJjIMwzA2B1MGhmEYhikDwzAMw5SBYRiGgSkDwzAMA1MGhmEYBqYMDMMwDEwZGIZhGMD/B5vXX0TyMhgFAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot the points with colors\n",
    "for k, col in zip(unique_labels, colors):\n",
    "    if k == -1:\n",
    "        # Black used for noise.\n",
    "        col = 'k'\n",
    "\n",
    "    class_member_mask = (labels == k)\n",
    "\n",
    "    # Plot the datapoints that are clustered\n",
    "    xy = X[class_member_mask & core_samples_mask]\n",
    "    plt.scatter(xy[:, 0], xy[:, 1],s=50, c=col, marker=u'o', alpha=0.5)\n",
    "\n",
    "    # Plot the outliers\n",
    "    xy = X[class_member_mask & ~core_samples_mask]\n",
    "    plt.scatter(xy[:, 0], xy[:, 1],s=50, c=col, marker=u'o', alpha=0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Practice\n",
    "To better underestand differences between partitional and density-based clusteitng, try to cluster the above dataset into 3 clusters using k-Means.  \n",
    "Notice: do not generate data again, use the same dataset as above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n",
      "'c' argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'.  Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD8CAYAAACVZ8iyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsvXmQXdd52Pk7d3tbv94b3diaAAFuEAVuEClaFCk7VmIlkeQw8pYZjzJyRmXPaFzj8kzslFQJK2XXJK5y7FLkJfQMJ7IzY0eWZYmaSOMZWQs3UBJIghQJkMS+NHpfXr/9bmf+OO/0u/3wekM3uhvA+VU13nbXh+7vO98upJQYDAaD4dbG2uoLMBgMBsPWY5SBwWAwGIwyMBgMBoNRBgaDwWDAKAODwWAwYJSBwWAwGDDKwGAwGAwYZWAwGAwGjDIwGAwGA+Bs9QUsR39/v9y3b99WX4bBYDDcMLzyyitTUsqBte63rZXBvn37OHbs2FZfhsFgMNwwCCEuXMt+xk1kMBgMBqMMDAaDwWCUgcFgMBgwysBgMBgMGGVgMBgMBowyMBgMBgPbPLX0ZmHmjTOc/8oLlC+OkxseZN+Tj9F7+MBWX5bBYDAssCGWgRDiGSHEhBDizSU+/5AQoiCEON74+Zcbcd4bgZk3zvDm7/4l/myR3J4B/Nkib/7uXzLzxpmtvjSDwWBYYKMsg/8IfAH402W2eV5K+Q836Hw3DOe/8gJedwepnjzAwuP5r7ywLuvAWBsGg2Ej2RDLQEr5HDCzEce62ShfHMfryi16z+vKUb44fs3HTFobwrUZ+X9+wLf/8VO89Cu/ZywOg8FwTWxmAPlRIcTrQohvCiHes9RGQohPCyGOCSGOTU5ObuLlXR9yw4P4hfKi9/xCmdzw4Ir7zrxxhlef+iLPf+p3ePWpLy4Iem1tRPWAqe+/DRK83jzTr542LiiDwXBNbJYyeBW4TUp5H/Dvga8utaGU8mkp5REp5ZGBgTX3Wtp27HvyMfy5EvXZIjKOqc8W8edK7HvysWX3Wy7WoK2NwtsXsVMedsbDSXtEfoDX3cH5r7ywSXdnMBhuFjYlm0hKOZ94/g0hxB8KIfqllFObcf6VuJ7+997DB7j3139m0fHv/KWPLDp+u/MvF2vIDQ/izxYJCmWcziwAUT3A68qt2wVlMBhuTTZFGQghhoBxKaUUQjyMskimN+PcK6FX4F53x6IV+L2//jPXpBCWUixLHWup89cLJfpa9imNTDL5pe9geS5RPcDNZ7A8B4Qgqvn0PXBw1S4og8FgSLIhykAI8efAh4B+IcRl4F8BLoCU8o+BTwC/IoQIgSrw81JKuRHnXi8bme2zWsWSVBjFs6NkdvVddf7ShXH8Qnnh9cxb5xj71qtYaY+O2wapTsxSuTJDUKqSHx5i4JF7sDwXf67Enb/0kbbXZrKPDAbDUmyIMpBS/sIKn38BlXq67ShfHCe3Z3FsIulqOffl73HyC1+lcmWK7K5+7vnMT7P/E0+0PdZqFEurwph4+ST12SJeZ47sUO/C+b3uDvy50sLryaMnwLLo2DuAsCyyQ304mRTCttn9kYcpXxxHWAK/UuO5T/5bAHa8/x7u/pWPAWyo9WMwGG4+bvkKZO1/14Ibmtk+5778PY79xtO4nVkyQ734hRLHfuNpgLYKYSXFAlcrjMyObvxCicLbFxeUgV8o03v/wYXYQfniOFHVp2P/EE4us3AsN5+lOjbDg099kpk3zvDKZ59h/swIXj6LBEa/9zrlkSkyu/qvS62DwWC4ebjllcG+Jx/jzd/9S0AJbr9QXnC1vPyrX8DtzOJ1dTQ+V48nv/DVtspgOcWiaVUY3fcMM/7SW9Qm55BxvOj8yVjDNz70a/iF0qLzVcdniYOQ5z/1OxTPjlKbnsfL57AzHgBCCGrjc5QvTTL80UcX7RvW6ox97Y1FbiPAuJIMhluUW75Rnc728XrylC9P4vXkF9wnlStTuPnsou3dfJbKlfZJUKtJI22tO8gM9tJz737SA91XnR+atQZSQOnCBJWxaWQcUxmdpnJlmq5Dt5HbM0B1Yo7qlSmiMFw4tp1yCes+wKJzVsZmGH/+R9ieu+A2euVzz/DKZ5/Z8LYZS9VKGAyG7cUtbxkAS2b7ZHf14xdKCNumPlUgqgdIKRet7FsDs3s/+iizb55fMo1035OP8crnnmHy5ZOEdR8n5ZEe7Ob9n//MVdegXT+1yTmimk+qL091bIba2CxxFJEZ6iM/PIiwLOVump2nNj6Ll2+mmzopj74HDy6KP8wcPw1A34MHEZZFqifP5MsnQUD/Q3cCTVfSiT98luxQ7zVZCxudqWUwGK4ft7xlsBz3fOanqU/NUzw7ShSESCmJaz5OLs3MG2euKgwrvHuJ47/1Z8wcP91WcM68cYYTf/gs08dPUx6ZJK75IICWvKqZN87wwi//Ht/8iV/nyrdeoXhulMj38QsVoqqP8BxyuwdI9eaZOHqC6vgM3fcM4+SzhKUqUbVOWK3jz5dJD3Zz9698bJH1E9UDhj54mMxg78I5w7pPVPMXXUdYqzP2t69es7WQjI9opWOK4gyG7YmxDJZh/yee4NyXvsPESyeIqj5uZ4b+x99Ldmf/gkDTwq4yNkPhxAWEEMo91LIK1opj6tVTEEtiP6A+F9Fz+HayO/sWjvf2Hz3Lxf/yMmGxSlCsAFCfnsefKyEcG+HahOUavmVheQ5ONs3cyYvs/ND9DLzvbgrvXCLylato5xP3cfevfGxBIenHV5/6Iv5scdG9OilPKaYEs6+fJdXbec2B59UE1A0Gw/bAKIMVcDuy3PmpjyCsphFVGZ3m4tdeJCxXyezso+fQbc3WEGkXf75yleA8/5UXiIKQ6pUprLSHnU0T133GnnuD2z7+Aca+d5wLX32B8oVxglIV4lhZDJYACTKIkLHESrtYQmB3ZCieHcXOeAjbJru7H8ux+eAz/8uC8jn/lRd46/f/6iorpV3QPD3YDRKV5tp4rz49z86feGDR97EWYb6agLrBYNgeGGXQQmsMQHjuouKv6vgMY8+/gdeZw+vK4hfKTBw9QVTzSQ/2LLSFgMWCs3xxnOqVKYRjExYrSAmWbYElmDj6FqWRKSwhlDUQJ/xGLc/jqo90bGpjM2AJ7JRH7IfMvnWO+z/3i4uskKSv/pXPPkNuTz+V0Rn8uRJSgLgwjtfdQe/9B3moUaiWvPddP/kglucu+n7WIsyXy9QyGAzbC6MMEpz78vc4/lt/hgwi0gNdhDWfsFRVK/TblUCbflUFX3vvPwjA5MsnQEBY8wnmyyCh7wH1WVJw5oYHGfve60RBSBzFAAR1H2LJ3NsXsVwHaVtgWRBH7S+wUbQtg5AojEAIAkuQ6uukemWGl//Hz3PuS9+heHaMqFonPdBN193DAMyfGaF4fgxhCYQlkLGk+z37sF1nkdXQGuNYSpivpqJ5NX2ZDAbD9kBsk64QbTly5Ig8duzYppzr3Je/x9HPfB4ZhDj5LG5HFsux6Dp0G6muHOlBlVEz/dpp+h++i9zOfkClaRbevsD82VFSvZ30vGcfnQd24RfKFM+OktvTT1QPsFMuZ770XSxbIGMWFIewLaQfKncQKL99tMz/iWNBGDdfCwGOhduZgzAiDkNlnXTmsNMeTi6Nm89i2RalixN03DaInfaIqj52xqP3vgN4PXkefOqTbU/XTuhDs6I5qSRMlpDBsPUIIV6RUh5Z6363rGWQFHLCcxn/zmvEfoibzyDDiNrUHOn+bqpXpnDS3oKwbA2+Zod6sVMuQx96YFHFsPBcEGB5LumBbvxCGdtziPyQuOZjpz1s1yGs1JGWANuCYAmLIElSEYCyFmJJXK2DEMgwVvqkHmCnPfyZIvXpIvkDOwFVewBgpV2CQnnFGEC7tNtXn/rishXNpg+SwXDjcUsqg1af+sj/e4zq+Cx2NoWMYizHBiAoVUDGDH3ogYX9auMzXGk0jLNTjkr1dG0O/OKHefuPnmXi5ZOAWvF33bV3kcDMDw8RVGrUpwoqTbWu0lSjekBUq1/7DUUxUSWxv20R1X2CeYGMIhW4npxrWilpj7gW4HblKJweoXplmuc/9TurFtytWULaOqpcmaYyNkP1yhT5/TtNbcG2pAhMATUgDfQD+WX32B6s5rpv1HvbHtySdQat+e9xPcDOpZFRjAwj4jBC2BZhsYJwbfY9+diCArE8l577bqc2NkPx9BWstEd2Vz8nf/8rXP6bH2J5DrbnUDw9wuQP32b2xHlGv3ucC197kTgKCQplOm4bJH/bIMJxCIOQOIqWdw2tlSiGSBJHEXEsEZYNQURqoFvVHxRKhLU6VsZj6vsnye7uX1MdQbKKujI2w+TLJ/ALZTI7+5g9fpr50yPEfmBqC7YdReASEACpxuOlxvvbmdVc9416b9uHW1IZtM4ldrtypHrzyDAivaMb4dgExSrCdRYydJIKpD41T/72XXTdPYyXzxJV6iow7Ic4mRR2JoXTkcGfKXLl268xf3qE6ug05ctTYAkyeweozxSJwxBZ8xGtCf4bgZTEQYiT9uh97z4GHz9MVKmT378Tr6uDjv1DxFWfgfcfouuOPW0F91KtJJJtNwpvX1gonOs5dBtxPcDNZ5k7eXHhUkxtwXZhCuUMcFH/aW7j9baYMbUMq7nuG/Xetg+3nDKYeeMMxbOjnP+r5xj97nEqYzN03T2MjGLV3bMnj5NNkxvewaNf+FW67tzDC7/8e5z891/h3Je+y6X/cpTq2AxW2sVOqbTToFBGSokMmz7/9GAPYanaSCOV4NhqxW4J0n1dPPS//jOIYmQYE4eh+nyjCWPiIMTtyVMZmaRyZYre+w+y56PvJ5grM/3aKcZffJPZE+cWdtGCe7mxm8l+TpUr07idOXY8eojMYC9uVw7B4l5IprZgu1Djas+w03h/O7Oa675R7237cEvFDLSAy+7uV1XChRITR9+i+z376Dywe8GnvvPHH1jImtFtoe1MijgMKZ4fU1lAroObSyNsi/p0gXC+iuXaBKUqbkcG23GQAizPg1hipx1SQ71YtsXYd44TV+tqNb6jC3+miIxi8BzwwxXuYg1YAiklky+9heW5uPkMV779CjOvnSG7qw+vu4OwUmX0W68B0HNo/4LgXmk2QzKwnCws67p7mLHnXsfL567qwmrYatIo90mydiRsvL9RXA+//WquezPu7ebmllIGSQHndSpXhur2Od22UdyrT32R2uQcXj6Hk01TuTypisZqPoVTlxcLbiGQSIrnx8ju7EVGMbbnktrRDUFEVAuoTxVw8hnq0wUKb19SPv1iVQWsUy7EkigIr+pVdM1Ygrjqg2tjp1Q66eTRkwjbpjZdVNlNpTpSwNjzb5Ld2Y8/V2Lwg+/l9d/+P0FKvO4Ouu4eJjvU29bdkywsi2o+M8fPEJSqONk04y/8CBlL3O7cguvJBJG3kn6UHx3Un37Y+Nm5QcfXfnuHxX77vaxPIax03UWgDsyilEEOsNnYe7v52RA3kRDiGSHEhBDizSU+F0KIzwshTgsh3hBCPLgR510ryVhBZrCXnR+6n/2feJz87TvbCqnyxXH8uRKVsWkqI1NgWcRhTDhfuXoFL6VqGeGH1KcKDD5+HwPvv4f6xBxR3cfyHKK6T+n8mOp+iiS3dwcyjIh8VRsQVf2NUwTQTEON4kaaqwNIpB8QzMxTnymqtFbAny4Q1gP2fvRRLn39KLbnYqVcpUBePkFlbKatu0e7jGI/YPTbr4GAPR95mM4791C6ME52dz99hw9sWEtsw3rIowSzixKeLusX1Emul99+uevWCsgCuhvbzwExa7+3InAOONl4vLWCzxtlGfxH1FjLP13i848AdzR+HgH+qPG4qay1V46dcqmMzRCHEZYQIARhzVeFXsjFglsAQuJ0pNnzUw/z2B//Gi/88u8xf2qE2A8JSlVkqILMwrGpjE4TzjfiCXEMqygxuGZiSVCuqdYSlgWROpkMApAOCHC782SHepl98zxxGCHjqNH7KIXX08HM8dN03bV3yerj9GAve/7+Iwvf7eg7x3E7s1RGpui6Y4+ZrrZtyHP90i1rKIsgyUb57Ze67qQCclnsLlqrIljKqtHnublTVjdqBvJzQoh9y2zyceBPpSp3flkI0S2E2CmlHN2I86+WlXrltAq5udMjyCiCKCZ2bNUmQgeJWxOApHrTny0tKBfpB/Q+eAdTR08QaWHs2MR1n+rlKXUM29pYa2AJZD1Q/YxarlmGyi3lF0qc+j++iZNJEccx6d5OOvbvpDo2Q2VkisyOHu799Z8BVPWxbro3+f2TXPjq86R68gx98PDCoYNCGTefXRRINllFNztb4bdfSgHNo1b3qxXgSaVC4nEEZWVstOtr+7FZ2US7aTr9AC433ttUlptq1i57ZvaVU2QGe/F68ggkQqjGcgutI5ZAB5+F5zL96inqc2VkLAlLVWTVV79boJRAa0XxZiNptMCIieoB9dki9al5ZBThdmToPLibjtsGyQz1Luq+WjhxgagaYKU9qmOzTL78Nuf/+gUqDYXjduUIipVFKbwmq+hmpx8l/APUL1bQeN1/Hc+ZbpwjSQXlTlpLzcFS2Ujz3Copq5sVQG4nPduuh4UQnwY+DTA8PLzhF7LUVLN22TPCsQmrdbru3LOw3eyJ8yqVNIoXGsctEMekdzYrc/2ZAtXRGWQ92PD72FAaCiHd301lbBphW1RGZ8gf2EVcCxaCwOe+/D3e+Q/P4hfKWK5Nqr9bWTyOjZV2qY7NMPbc6wx+8L1kdvVRvjRB96F9JqvopqY1e6gPJYz1651c3xV0u+ByBchy9Sp/aplrWcqqEdwqKaubpQwu03S+AewBrrTbUEr5NPA0qEZ11/vCzn35e5z8wleZfu0UXk+egYfvoufQfgBywzsonh4hqvnYKZeoHuB2ZJAxhMUysrW7qOvQfffwQrvokW/+cPsrAo0AK+Vguzapvq5GumwFtytHx4Gd1KYKvPjpf0dYrSOjCCmhfHECN5/FcWy8zixOJoWXzzH1g3cY/vgH2PfkB5cdAWrYKjYq/bOdn32azXWh6OBy8n7SKGWQZCUBvlTGUr7xePOnrG6WMngW+IwQ4i9QgePCZscL2nHuy9/j2G88jduZxevuICiUufT1l5n8wTvk9w1hpV2EYzN/9grCskgPdNMxPKimidm2KjLT6soWWEIw+9Z5LNuiNjFHuJ0UQZvxmouIJeVLE3QMDxJWanQMD7H77x1R3VfPjTL96rsIR5DqyVObmluINYTlGrbnkOrNk+rtZOjxw5QvTy409tv/iSc25fYMq2Uj0z+X8rMvtwJvVURZFlsS16KYWoPL51h77KKdUtFpqdczHXf7sCHKQAjx58CHgH4hxGXgX9H4n5BS/jHwDeDvA6dR//P/7Uacd72c/MJXlSLo6sCfr6gJYzKmOjaDROJPF+k+fDuOa1ObLCBsAa5NVKsjbAvh2ImqY4GwLKrjswiJSkUNr2eK0BpZycaSKugrUi6yWMEvlLj49aN03rGLypVpwmJNTVlzXLyuDoK5EjKOkagWHnE9IChVOffl58js6F6oVDZsN65FgLdDC/W4cQy9Ei+jBDFcLdhbFVEZGAU6G/sHwNnGZ5L2yiGpTLT3uXXba62nWCpjaSklsZYA9fZno7KJfmGFzyXwP2zEuTaSypUpMkO9BKUqtfFZLNchjiJkGOFPzeP1deK4Djs/dD8Ac+9e4tLXX0LGqu+P8rULFTuIYqKq6hy6fSdELI/luZTOXGHX33mQzgO7mD9zhcmXTxBLiZ1xif2IYL6M25kj1d9Ffba4MA40lhIZRViuTWZXH6987hkyu/qRfmDaWG8rNiL9Uwt1gSruilEFXzRee7S3OEZQa0GdnRM1Hn1UoVgMVBv79iaOoeMQ86jAcLZxnrnGcbMoxTIG9KByU1oFeFfj9WXWLrxblcT1Kq7bWm6pCmTNzBtnOPGHz1IZnaZ8aQrLtYijCDvlIUKBlbERnoNArZZBdeece+s8MoyxOzLIuZIKIq/kfrmB6LprL2GljteZZey5Nyi8cwkpJZZtkervpjoypQbzlKvYKQcnm+bwZ/8Jl7/+MnJyDq+rg9RAF8VTl5k/PYKTy7DnHzxi2lhvK1aT/rlSTEFbFx0oAW3RLJSxWRy81RZHEaUwHJoKoNbY16LpLnJQSkFn7vgoa8GlqbC0EhGNz2doKqUxYAJVgLY/ce6NFN4bZV1tL27JRnWvfPYZxp97nfRgD1LGhOUaMoiIar5qWLe7T7VsmCxQGZtm9LvHmTr2DsISpHpVmqmdz6gCrptEEWBB8ewo1fEZzv7ld5k7eQG/UCIoV/ELFSzbIrO7H8uziWsBwrI5/Nl/wuFf/znyt+9k3z9+nK67hymeGqEyOoOdTRP5AVPff1tNXjNtrLcJK6V/rqYVtE7DTKFcPFbjWLLxWlseSYtjimZ6pnbbgBLsEUqAlxrvJ9eoFZpuJ/3HFqIsCB+lAGTjPZ2mLYACSoloxZZMD40bx73WSuObsyneLWcZnP/KCwv9huyMh+25FM+MqgyZOCY91Is/W6I+XYAYrGyaqFqneGaE1EA3Oz7wHmZfP0tQrqmK5NnSVt/StSMa/wgQjk1Ub2RNVXwiII5j8EOEJYjDiI69O7A9l8yObt7/+c8Aqn/T9GunmXvrAjJW1pWMYoRt4WbT2CmPwtsXVWDZFJxtA5YKlOoV7WpWvUnrItX40YI46YJKWhw1lCtI9xFK1gbEqJW9bHyWS3xWb3xmJbbRFoG2BpJIlHIRjWvU95lKHG+eZrxhNVZCq6WkFdrNlWF0UyqD5cYuli+OE9V8XN2jaEcPTjZF8ewoTi5DVA8aHUgdUn15onKNuuPg5DI4aY+eQ/tJ9XZRePsihVOXiVyrUcm7xcVja0WAnc8gYoj9AKczh5NJUZuYJZYx+M3gt4wl5XNjhKUqdsqj+4nDFN69zKWvH8Xr7qD/4bsYf/5HVEen6bh9F8K2iOoBmaHehfGapuBsO7FcS4rVxBTaBWhdlBCebjzajfd0sFUrkDRqVd6K3sehqRBClHDX7d11byKZ2Kcd+nPtikoqrwpNt1aEijtYqHjG3S3HKTbeTzbACxrXoI/vN94TwG1LXM+NwU3nJlquDz+o/kR67KPGchyyuwewHItMfxded5buu4fpvH03+duGyO7sZc8/eISoqrqVzp28wNypy9TGZ5FhfGMpAttqzGcWyHqI5Tk4HRm8zizCgjiMlCJoUyZYny0R+QHjL77JS5/+d0y9eoqoHpDb2c/Q4/fh5DKUL02Q3dlLqjeP5djENR8r5eLPlRYqsw3XymY0UmtX0du66m3XOG6IptBuV2Oq3VP+MufWArqOchmVGscMGj/aQlhNlp52W4nG8aYbP9oqqTeO5TS2m6X95LQSTQugiLJEMo19K437cVAxj2lu5OZ2N50yaB1p2Tq9a9+Tj6kB9cUyUbWuxkDOl8nvH6LnvgPK933XMJarjCa9snXSKXofOMDsW+coX5qgPjUPsVw+ZmBbiJS7zAZbgIxVYzwpiWNJWK0TlGuULoyDbWOnvcXbJ/6uLdsimCtRHZsFW+DPFhc6mgqhBvrEdRVP6Dl8O0goj04ThxF+ocz5r7xgupZeM5s11nG1LSXyqADtPY3HCkpI9gEDjccMzbYNWoHo1f5SaIWhhbOPEv5VVMZQRLMlxEpot1Oyo6m+N6/xA81gdbvJaVqh6EC3DnLXGsccRGU+5bjR21TcdMqgdaQlLG6Q1nv4AA/99qcYfPw+Il91Ed35xH089Fufou/+O/ALZbrvGSaq+UQ1f2FlWzw7SvHsWCPQHECw8hAay7Wvblmx1cQgZayC37GyalR6qFQzn+2G9NeXvfC3oNKmVIuOGk42TeyH2CmP6VfeYeLoCYhjOu/YDQJmj58hd9sAffcfZPAD99J7+HbTxnpdbNZYx9W2uW61UnQPH2jOFpgDJmkqrDxKAK9mql9MM6tIr/JDlMjqbBxjpXGxDkoh6W6mfajUU0lT0AcoAV9GpZ2+3bheHSTWcQoa59bBb32M1vPduEHkmy5msJo21b2HD/DYH/9a2/3f/N2/xOvuYOCRexh/8S2qV6ZwunLUZ4tEQUjHngE17QzU70YbD5HwbKRUefthqbqRt7d+BMqisRv1EUiElFieS32uhGVZCM9RxXSxbLTmFo2fRqM+qRRsrTGXoXx5isxQD0joP3I32aFe6rNFZl4/Q+99B5aclmZYC9ezPXQr7WIKrcVedZSg1VZKHfXHECauSa8k3kStnF2U22UtxZitiykbZSWsZpFVp2lJaLI0lVOAsjwkzZX/XGO/VONesihFR2K7EKWQbq4g8k1nGSSHtcs4VuMtV+mvXjTbd3wG4Vjs/sjDdAzvwMmmiCp1ylemm1XHS4QKZBzjpD1VhLYdLAO7IcSTCykpIZLIKFYDbiQQSbJ7Bkh1d+B1dyBcvYKTasSnJYj9ECebQlgWg4+/F4Eg9gO8rhwD7z9EdqgXUMqicmVqWSvNsBZW48tfL3q1/wbwKvB64/U4i11UJZRATtYDeCihWaa5ek5m9RRopoKuB+2+0q6b5Ygb560n3gtpFrgFieuJG9cbNO5DWza6P5HetwNlKe1m8zu0Xl+E3A7CagmOHDkijx07tub9lssmWi2vPvXFBQvjwtdexOvMUrkyTfnyJHEsl3QTWRmPzM5ebM+lfGmSOAiR6xllqQW4EGqlvhyN1trCs5Wg172TbKuRDdfyh5hc9UuJsCwyO3txO3PUJudwMh7+fAUkpAa6qE8WVEZRLk3+9p103z2sOpZmUjgpd5E11s4yACicukxlZIr87TtNZfKaSBZOJVssrKVwarliMn38GCXQ9S+eLgbL0kz51K4pG+V2AVXs1c5K0Smh2u2zXnmjq551CulK26ZR31cv6j4Kjf0kixWTvk6dWGInzpNF3edyrTG2T0sKIcQrUsoja93vpnMTwdJtqtdC+eI4uT2qHbXXlVMxBD/AyaSIgpAobCgD0ZxvIGwbJ52CKKb7PfuoTxUIYlV6v+B2WSuShmumzb7J6md9HbYFUqgJagsZePEyTcQlCIGdTSPDkPpskfz+nQyvfZxFAAAgAElEQVQ9di+Zwd4F4e11d1DMjJLZ2UdcrVObLDD71jnu/9wv0nXnnrZDg+75zE9z6etHF97X7S36H7lnUaaXqUxeDSvVB6zESlW4OiahM3j0ynmOps9eKwMbtcrXbpiIpd1VMavPAFoNOnaglcJyx5WoNhQl1P1XaVoryUI5aBajJffN0FS67QT99Zwat/nclMpgI0jGHrrvGWbi6AmCYhU3n0FUfYQlkLHqx6MLt9xcmq57biOuBRRPjWB3pKlPN/yN7fSAa6tmdivpCEtA1GYjx1EjLB0bN5dWyiaWSAvCig+eo37na2FzVKcFwrKRgfojsjwXO+XSecduKhNz1KfnqU3OMYekXigTlqp03j3M1PdP4qRTKtW24QqqzxaZffM8+z/xBPf++s8sssZ0u+quO/csvF8ZmaL/kXvovlN1MzcxhLWyHuGzUjGZjkloi1engGo/eY2mL91BCVetNFaKi62cbHFtpBvXtZxCmKGZVZRMN23dp/W1Xj3ZqOu/sVtNrAajDJZAj8isTRWYe/sipQvjjTnGEXbaQ0qJaLhY4jBGNlbYYblGfaqAXygTVutNN03c+EXUv4tCYDl2Y/8VrIYwBkuQ3qkCs7La+ENtuKoEqq9Q9co0QbECcUzXHXtI93cS1Xxm3ziLtAWicQphWYRRteHyFeQP7iKsBoTFikrL7e+iOqEG89gdGWrjc5QujGFn0/izRYZ+/AGyQ71XZWm1E+jJ95//1O8sWFsaE0PYLFYKQOvCLAcl3LUwTLp5Sqhf4CKLC7tW06BLcPXqux1LZGVchUS5fVbKKNKVyklhv5pr0K44HROYZLu4ga4XN10AeaPoPXyAvR99lKlj71A6ewUnlyY12E3sB9SnC0SVOpEf6tHHqkC+XKV0YZw4jAjrgXINNUZKknTbSLAcGzefJbN7oLn40n8v7ZCS+lQBGcvGKr+BJRC2xdyb5xCeA45NVPWpTc1Sm54nqvl4PR0LhXFuZw5hNeIEjoWTTUEkCcsV0gPddB3cza4ff4Dbf/ZDWCmX8sUJVa2c9ggKJYoXxrn09ZeojM2suao4Nzy4aCYymFGYm0cyAK1TPydQgr9Is77AQwlO7d/XKZxdjc/naBZe6arj5Vjul7odawkwt/r913NM3SxPb19HxU6S2UYr1XVsRlHg9cMog2WYffM8qZ48XXcOkxnqhTDG6cioDxstm1XhmVq5E8aqt0/NR/qBysZJCu7EUyljwlKN2sSsUhT6b2+pBZZlqYEy9QBsgZ1LYXmOslLCGCkl0g/JDvaoCms/IixV2PHoIQZ+7D0ICXEYE5SrhPUAYVt43R10HtjN3o8+ipNNY6c9uu9pjhqtTsypuHUYEfshMpYI28KfKzP23OsUz46uqap4PZlehmshKZzqNAu3CjRz9j2arSX2ouICHk2LQK9SdOvpPM3YgW7DsByy5WcrWO15tUXQSpj4bKm6js0qCrx+GGWwDOWL4/hzJSpj02q6WaWGpEW4W4JUX1ejP1FEerBbDb4RQmURafePbS1SDJbngi2IKqvME4/iRBZcDBKcvFrJqCZ7kqBSQ1gWu/7uQ+T2DuBk06QHusnt6ie3fxCQyHqAnXJID/Ziew6523ZQvjxJeqCbnnv3kxnsXTiljGNVfzBVwE57eF0dWLaFlDFePkduT/+afP3J1N3y5Um8nrwJHl83WoWT/lPX/n0HlSufrJzVVcXvQWXP5GmuUnTAVVfg6qwbK/H5UiSH0GxXLJR11FrMpoPNughtqb/XpYoCR7hRrAUTM1gGO+VSm55Xwh2IwghZKTbTPCVqxRzFKsvIDymfGyes1bFyaeJyDakleNR4bPyeRZVE7rNtNT9fDcLCyaVxMikQEAcRQqoAts71t1yXysgU5cuT5IYH6X3vfmzXoTZZACSZgW667tpL5x17ePCpTy70dKrPFhcyglI9eaIgJChVsTMpBDF2yqPzzkF2/70jlC9Prvk73YhML8NqWCpgXEe1i0gKvHZCzqIZfHVpCn1dK6BdSRotAPV7SffRdlYCmhDlOtMxAlhc/eyjFGGu7d7tYzK6EV4fN8IQnA2xDIQQPyWEeEcIcVoI8ZttPv+nQohJIcTxxs8/24jzXm9iCU4+A1GMFKJZbGbbCReQJKrVyezoZs9H3gdCBUXdjgw4bcruJY28fvVcOLby4VtXWxytOJ1ZrKynGsBFMWHNR0hI9Xfi9Xaw48fuJbOjm/psEcuxef/nP8MHn/nn7HvyMWZeO4OTSdF55x46hoeQYURY97n4tRd5/lO/w/mvvMDejz66aNX+0G//Er3v2Y/tucSNKW5eb56B991lfP3bnqV67muXzywqKDrbeK2L17RFoV1IesCMpBlP0AVnSXSKqd5nLXGC7UKZq9NLk8/nacYVYLEbrsrV3VjLqO/iercQ2RjWbRkIIWzgD4APo5p7/FAI8ayU8kTLpv9ZSvmZ9Z5vM5F+wO4PP8T0K6conh8lqoDw1IQvYanMISeTxXJthn78AeyUS2qgG7cjTensGMIC3GYaJ46lMoMSmUM6K0nYlrIiGhaH8BxkvdlZVaRd1WBvtkhYrWM1lFOqr4uhH7+f3T/5ILNvnr8qtRNU6maqt1MFuoXAzngE5RoT33uD/O27FnL+L3396FVum64793DiD59l7G9fJdXbSe/9B7A81YX0zl/6yKb8PxiuhXYTzSo020zrFtO6fYTO8kpaFDqzSNcJpFk+2yem2a7iRrAGWlkpBVai3D5aISTrNiKabSuyNDORuluOsX37F22Em+hh4LSU8iyAEOIvgI8DrcrghkPXGgz/w0cBuPh/H6V8eRIhBLk9/fjzFZxsCq8rh91o0zz0oftxUi79R+7m4tdeIAoigtni4hRnHTuIJYTK6lBWuM79VIoIAAucjgxOLkNldBrLtkj3dzH42HuxHHuR8N7/iScWXb+uxD77f30LtzNHUKyQ6sk3rrVIHIT0PXhwobsrXJ3zr/s4Jau6Mzv7Fikbw3akdeaAniFs0xRedZSA76S5qk26O/RsYZ2Xr91Dy7EVgeLNUD4eTQWrLaekG067j3QqaxrVuqKMsiB0m2uL7dq/aCOUwW6av3WgrINH2mz3j4UQjwPvAr8mpbzUZpttha41AOX66bxzD/5skYH3H1oYGD/71nnc7g68nvzCSlk3u+u6ay9+oUwZVA+gKEbGKvOHWBJHsSoak40AnC0QlnIbOR1pZBARVusI1yWuBdiuQxyEBMUKV/72VfoeumPJgi0dA/C6O8js7COYL6ss7zAiqgfEfkj+4K5FAePlcv6Nr/9Go7Vi2acp9LUA09XBWa6uN9BTzNIopaEDqcng8XbheiiCZNBbK8MaStBXUBaTbtGRVJ42qq13ERVz0U31rMbzTGO7k2ynFhawMcpgyUYHCb4O/LmUsi6E+GXgi8BPtD2YEJ8GPg0wPDzcbpNNQ2e/6BVx15172ffkBxfcMZ137OHwb/z8VUJy70cf5eQXvkrx3BhRPaDz7r3UJwqUL02qTJyuHKnuDrK3DTLx0lvENR877WKnU6oYzRbYKY+OfUPUZ+YpnR8nrPsIx8btyGBnUsR1n6lj72K7i+cl6BX8xa+9iO259D14kJ5DtzFx9AR2NoXbkab3voOMfu84XXct/n5NHOBmI1mxfBIltHQAWKeMRixueNdqUWRQwr+TxUrjZqa1VYVOs9WN7PToS+0a6mxsV0J9t+caz7VS1RaVjrd00FSq2yegvBHK4DLqbjR7gCvJDaSU04mXfwL826UOJqV8GngaVKO6Dbi+ddFuRZxssaCH5uhtZt44w6WvH6X3vgMMPX54wXrI7O4jDgKEZZHbu4P0ji7m3x1BWBZORxrb8whKFYTr4HR3EDXiBV5PHmtkijh2cDvU6E0AO+URzVfw55ozmLU1EAUh1bEZ4jCkeGGMoccPs+PRQ0z+8G1mf3Qev1Ahf8cegmJlUfaQiQPczOgVf7IlMzTn+erxlK0WhYMSXno1rIXbUi6jGyFesNI16mI7XV8gEj8OzUC5Zp5mHKWLZk8nt2U7Xaux3HzprWMjsol+CNwhhNgvhPCAnweeTW4ghNiZePkx1DLlhmSlsZqtk9a67tjDzifuY+cT9/OTX/stBh97L733HaA6PkdYrSODEMtxsDMewrKIyjWqE7MEpSql82P4s0XsTGph8locRkgpCWs+lufgdjdT3c5/5QWiIKRw4gLCtrBSHkLA2HNvUJueJ6zU6bxjN8MffZSOXf2N1NjA5PzfEugKY4um4AlotmRuJ4x8lKDTOfhdKOHmoJSLLlBLcqNkES13ncnKZl1pnbzPCGUNdNKc2eyggsVpmkK+NSDdTgElA8pbW8G8bstAShkKIT4D/A1KnT4jpXxLCPGvgWNSymeBXxVCfAz17cwA/3S9590qksIerm62lux2qtG++KTbqTo6TRyE5PYMEBQryCjGznhElToyiLDzWcKajz9XouPALpzGkJy4HhBW6sg4xuvpoHxhglef+iL7nnyM8sVxqlemsFMemaFeypcnsRybqOYz8dJbuJ1Zeu9PBIxvV5bHg099cnO/RMMWkFzxR6jc93b+6mR3U73y14Vm2jcuUAJM1yDooTaw/a0CnTJbX2E7PatBfwc6dgBK4OvvQmda9bJYwaRYPPNBfy+tylO76FbqKnv92ZCiMynlN4BvtLz3LxPP/wXwLzbiXFtNO2Ef1uqMfe0NyhfHKZ4dJaz5C505YbEvPul2Ov+X3yXV34VXyVGbLFCfrSPSLk7KQ8aSqFpFWBbV0Rkyu/uoXJrEzaXJ7lUWieU69D9814J1YmVS1CYLpPq7sIVHbs8A1bEZNZdqpoiTS1N4+yLAVY3mDLcCq+l6mkwtjWj29K/QzELSLbTrjed6Shisf3jN9UQPsm8V0Euhi828xnNtBegUWu1KynP11LMsTUUS0QzU61hBcibFTlbuKnv9MRXIa6R1rGZlbIbx53/UaM8wQFQPmHxZZdV2Hdy9pC9+35OPceGrzxPMl/E6c2SGbIJihcyuPuyUhwwjZBxTGZ8lrvtIP6T3gYNUR6eJyjUyO3rovf/gQjtpUC4f4doLx7QcGyeXXuhh5HV3EFV9Jl8+wcD7D2GnXBMwNrSgU0vrNKd3aaGmu5bqtEm98nVRglZ3NoXtaSFIlCvHQfVnmlthe73SD1FKIYcS6H5j3x6a4dJk0F27425HKdHk8BtoP5PiMps31rQ9pjfRGmlttjZz/DTAQr5+1x17GHj/IaqNqWhL+eJ7Dx/g/s/9IlJCdbKAlXLJ7OpHNtpS2CmPYL6CZQncfBYnnSKu+ux84n7sTIrdf/fIIkXgdeWI6sFVx7RTLnbaY8ePHSKuB6oDacpl5vhp0yTO0IY0zZoEWBwkTtYalFGCSnf41MVrW9mQbiVqwChK8K7kj9cBY30vOZRrLd941C0mtLW1F6UU643HvcAgqtfTPY1HvW3re7A5Y02Xx1gGa6Q13TSqBwx98PCifP3OA7uwUy4ffOafL3us/Z94YlFmkn3IpXx5iuLZKzh9nQTlGsK2SPV3YaVdgkJ5Yaaw7h2k0a6o1mNOv3aawUcPkdupViVTP3xHKRnP4d7/+WdNwNjQQj+qTcVy6JnGWlDqUZk3Aqt1Y7UqtNaeRK2r9vVOPWtN6U26kDYHowyugaTfX89KTrKWfP3W1NWZN87w8q9+gerEHG4ujZPP4nZkiKo+biMFtP+RexZSStulhSbbUEy/dprZ189SnylSPDVCZkcv6R09CASXvn6Urjv3GIVgSJCnmYa61AQzLVBvhDTSayHZoM7m6kpibQHEwNuNfbT77FqLyNY71nT9GDfROtnoHv29hw/w/s9/hoH33cWOD9yLsAR+oURYq5PZ1Yc/V+LQf/+xZVtBJ9NfBx6+G3++zNh3XyeO1LCduB7Q9+BBvO6OhToJg6GJXuW6y2yznd1B60Xfl41y9Vg04wA+KiHSR7mJ5mjGHtY7w2ApF9LmYCyDddLqNmptErfeY4blGv5cCbc716iAfmzh2EudozX9decT93H+r56nNlkg3ddF3wMHyQz2IuPYZBMZWijSnITW2nYiaQlsR6sgg4pn+CttuAK6DYeHcodlUYHnCkop2DRTb7UIraICyrBdisjWilEG10iycVtueHCRkN4I1tMLqDX9NTPYS/ehYaoTBXZ+6P6F9037iVuFIovdD0vVF4yglICL8pHHNAWrDqjaNIVhsu/OVqeUCtS1rtRIbzXH0RaR9t0XaHZxtVFKQgfSdYVymNjnxmzZYdxE18BKVchbTbtZw5ld/VieY0ZO3nKsZhyj3qZEUxDWUBW2uuK4A7VC1oriNppTwbaDhSBZnyIQiUet3GZpptY6NBWA7kyqazCqjZ/LqGyl7dbIb3UYZXANtLacSPXkt5X/vV0cw3Yd7v/cf21GTt5yLDWOcarNNrrKVheRBTQtiCrKfZQDhmgWoenjbbUo0ZbLSliotFA9eErfb3L/pEKAxYVqNup78VBuqSpKOeh9I5RyvfHcr8ZNdA0s13JiO7BcHKN15oHhZmepcYwFmm6jeZTQt2kKPQsl/H2U0MuhMmrmURZEHmUtFBrH7Gs8r9O0FDbTYljtuSRqxQ9KkWmLYmHYCM2urnr+cydNCyGFUgQ5moN8bJoKVB9zBBV8vnEwyuAaaK1Chu3nfzfzBwyK1olndZTQTvbAqaMEWrKjqW7Wpn3kRZqr/whlGXSiAqslVEaNQ9OdYtMcl7mdSGZBaReXS3P2sb5+rRA8FsdX9He5v/H6Za4e86nbWNxYbLVtd0Oy0emkBsP1Q3cr1cJOt4vooOk2yqKEe7Kjqc6a6aKZVtnqOtGuIl2Q1YVq2KatCj0HYDsgWp5rKyhmsY9ff1dZ1HeklamkWXtRp9lZVDeqS9La4vrGwCiDa0C7YYz/3bD90cVMMWr2cZWrfetZFrde7gPuRbk5kn39oelf17n3sHjwe7qx/3IZRlvR5jppoWhXUDJOoLOi9BCbCio2cjuLi8xobKetKt28LmwcQ1cO775+t3KdMG6ia8S4YQw3FjFq5a6Fu57QlaLZeXN/m/30vN+I5opaCz0dTNVWxWWaM5aTvveky0n74ZNod03y+Lp99PVwM2mh3qog9LVbNAf6DNK0lrQlkOws2kHTzaQ7nN7GjRYvAGMZGAy3AMmMomzi/TLNzqT9bfbTVkUHTSugi6ZwdFncnE4rgGLi8zTNAKtWKK2WQWvAWXB10Hujab2GpNLSraZb6wX05LckOgurBxVoB5VJtPnDadaLUQYGw01PUoilaE7o0sHQ5Qao5IG7UW6jvsZ7OVTLhMMoJeHRFPJJodrJ1SLGYXk3kbYmLJpBbr29xepdTJnEc73qbyV5LO02SlpLrR1D23UW1S2qyzTdRdXG6/W0pth8jJvIYLjpac0ogqbbZrUs1ZVTu0Z0IDbpi8+jBKRueQ3NNtjtSApnC2WFVBvHXUtxm943YLECaXVRtT6PaGZFuVzdMbRdZ9EKzRkHySC7j7KobpzWFMYyMBhuepIZRTWUwAtRwmq9zdU8mhPR0ijBmGq81nUJyZGRyykC0bimnsajHiaj21/oCujl0C4mHRDXCiAp+HVbCS28k602lqsebje3IEVzqlnSggm50VpTbIgyEEL8lBDiHSHEaSHEb7b5PCWE+M+Nz78vhNi3Eec1GAyrISnEdBvm5PD21orktbCbZgaNTDzuQa2O48bxe1l6UIt2L+nt+lApsBONY+xBKYgsi0dstqKDu90oBajHXGq3lRbW2lVlt/wkhblN++8k2Vm0H6UEJmjGXkjc8+YOp1kv63YTCSFs4A+AD6PSCX4ohHhWSnkisdkvAbNSyoNCiJ8H/i3wc+s9t8FgWC3azaMrkpMumfWsYHXWzAhXZ9NUUYI5ea7WQfTJz3SAexplGSRdMX2Nc4jGNmHjua54FijBm0YpqDwqLXQEpQB1JXEaJbz1/joIrq9FZzyVWb75nu7n5NFUgDr91KIZX9i84TTrZSNiBg8Dp6WUZwGEEH8BfBxIKoOPA081nn8Z+IIQQkgpt1t5osFwk9MufrDeFewg7VMpk+eqo5RDMk6R7HgqUEqkQrMKWLuztDWgq6H9xnGTxV4eKk6gFQE0g9+taMUyS9NdpHET19QaLE6SzNByGtet3VG6fce1DrrZGjbCTbSbZlQFlHXQWnGxsI2UUveE7aMNQohPCyGOCSGOTU6uNH7PYDCsjdaK5OVSSzfyXDqALFCCW6ecglrpD6AUSg0l5OdpulskzYK5Os2gdBK9Ml8NrRlSOm6QVASwfIC9NUOrB1Wk1oPKstr84TTrZSOUQbtcr3ZVJStto96U8mkp5REp5ZGBgYF2mxgMhmtmqeHtSwmuIipnXrdfWEugOXkuXa3bSTNgm6YZGNbrxzRKcSRjA8luou3aVGuhHLC22If2/+9CubN0wNlqXFPnMvtu/QD7jWYj3ESXUf/jmj3AlSW2uSyEcFD23MwGnNtgMKyZ1Q5v137xZFO7SyyvPJY7l3YZeSi3im4FnTxePzBGs/unDv5qxaAnjGm0nz9q/LTGPlYz2KcfZXUk4xQrWUtbP8B+o9kIy+CHwB1CiP1CCA/4eeDZlm2eBT7ZeP4J4NsmXmAwbHdWMwthtSRdRh5K8OaBO1gsnPM0x0dqJdDZOLeOGej5CckiNx1TSK7MVzPYR59zLdbSte6zvVm3ZSClDIUQnwH+BvW/8YyU8i0hxL8GjkkpnwX+d+DPhBCnURbBz6/3vAaD4XqjM4/qqJW8rie4lk6kWngmV+k7aS88dYhRt8QOaQa8bZotMDT6ulwWr+aTyozEY7tCsNVaS+vdZ/uyIRXIUspvAN9oee9fJp7XgJ/ZiHMZDIbNQvvvdSM6XUwVolbX10t4LqU4aLwXN97XMwN03UQykwjaD/a5sQrBNhPTjsJgMCxBP6Az+pJZNlmuf5uFpRTHWs55PdJob15MOwqDwbAEeZopoDrfvxOlDG6E1fVmptHe+BjLwGAwLEOeq1fXATfG6notcQqDUQYGg2EZNiOFcjXpn9fKzRXkvZ4YN5HBYFiG651Cudr0T8P1xlgGBoNhBa7n6not6Z+G64mxDAwGwxay1CjJGyFAfXNhlIHBYNhCbr4ePzcqRhkYDIYtxKR/bheMMjAYDFvIzdfj50bFBJANBsMWY9I/twPGMjAYDAaDUQYGg8FgMMrAYDAYDBhlYDAYDAaMMjAYDAYDRhkYDAaDgXUqAyFErxDi/xNCnGo89iyxXSSEON74aZ2PbDAYDIYtZr2WwW8CfyulvAP428brdlSllPc3fj62znMaDAaDYYNZrzL4OPDFxvMvAj+9zuMZDAaDYQtYrzIYlFKOAjQedyyxXVoIcUwI8bIQwigMg8Fg2Gas2I5CCPEtYKjNR59dw3mGpZRXhBC3A98WQvxISnlmifN9Gvg0wPDw8BpOYTAYDIZrZUVlIKX8yaU+E0KMCyF2SilHhRA7gYkljnGl8XhWCPFd4AGgrTKQUj4NPA1w5MgRueIdGAwGg2HdrLdR3bPAJ4F/03j8WusGjQyjipSyLoToBz4A/M46z2sw3LJcPD/DK0cvMTVZpn8gx0OP7mV4X+9WX5bhBme9MYN/A3xYCHEK+HDjNUKII0KI/62xzT3AMSHE68B3gH8jpTyxzvMaDLckF8/P8M2vnqBc8unry1Iu+Xzzqye4eH5mqy/NcIOzLstASjkN/J027x8D/lnj+UvAe9dzHoPBoHjl6CVyuRS5Dg9g4fGVo5eMdWBYF6YC2WC4gZiaLJPNuovey2ZdpibLW3RFhpsFM9zGYNhmLBcT6B/IUS75CxYBQKUS0D+Q26rLNdwkGMvAYNhGrBQTeOjRvZTLdcolHxlLyiWfcrnOQ4/u3eIrN9zoGGVgMGwjkjEBYQlyHR65XIpXjl4CYHhfLx/56UPkOjympyvkOjw+8tOHTLzAsG6Mm8hg2EZMTZbp68sueq81JjC8r9cIf8OGY5SBwXANXK9cfxMTMGwVxk1kMKyR65nrb2IChq3CWAYGwwq0WgGFuep1y/XXMYHk+R7/8IHr5hYy1cwGjVEGBsMyXDw/w1/9p+MUCjUCP2Lk0hyFmSqHH9oFNF05G5nrv1kxAW3h5HKpRRaOCUjfmhhlYDAsw7e/+S5jV+bJ5jxyOY8giKjXQ069Pcn7Hr1tYbv1+PW3anVuqpkNSYwyMBiW4cw7U2SyLq5nA+B6Npmsw8Roke+/cJ58V4revhyWLXj8wweWPVY7oQ9s2ep8NZlLhlsHowwMhmWQSJJ5FrVqQLUa4jgWnV1pivN1KuWAj//svcsK76VcMinP2bLVuclcMiQxysBgWIaDdw5w8q1x6vWAWjWgOF9HSsnAjjz3vFfNfCqXfEYuzqvm7EuwlEvm7bfGePB9izOFrnV1vlZ300OP7uWbXz2xcM5KJaBcrq9o4RhuTowyMNzyLCdEf+Lv38nlS7NcuTyPAEDgOII4lhTmqnR1Zwj8kFe+P7asEF7KJSMQVCrBulfn1xIM3uzMJcP2xigDwy3NSkJ0eF8ve4Z7CENJGMRUynXSGZdUyuHK5QJIePP1UerVgLmZClLC0efP8Yn/6j4e/sD+hfP0D+QYG5lndrZCtRyQybn09GQ5cFc/5XKd4nyNmekyxUIdx7H42M/eu6Z7+E9/8kOKhTr5rhS7dncBMDoyz5/8/ks8+MjeJa0EU81s0JiiM8MtzUq9gACCIObw/bt46JG9HH5gN5YQSKncQ+++PcH8XJUwUtGFIIiYHC/yH//oB/zgxXMLx9g93MnpU5NUyj7pjEOl7HP61CSHDg/ywPt2c/niHKV5n86uNLv3dvPaD0dWVcSmldn8XI2OvIfvR5z40Rgn3hgFCbGUZgCOYVUYy8BwU3At6ZkXz8/w6vcvEUtJrsNj154uurozVxw9ACYAACAASURBVPnsk4HWrp4MB+8a4Py5GSwhmCvU8FIOQgjq9RAhBI5jU6+HfO1LbzK0u4vhfb2MXJzn4J0DzM5UqJQDsjmP3Xu7VawBuOfeoUWuonLJX1UQWSuzzu40gR/heTaT5RqBLylXAlIpmzCIFhScsQIMS7Euy0AI8TNCiLeEELEQ4sgy2/2UEOIdIcRpIcRvruecBkMrP3jxHP/h917kBy9eYHKsyNjI/IorYb2itl0L17UpFeu88vJFjj5/jjeOX8F1m38arS0iHNdm5+5O/rv/6cfo6k4DyiIQQgCSMIyIwpi52Srf/ua7AFw4N8PMdLmhCFx27elicCjP1GR5yYE158/O8Nd//jp/8vmX+Os/f73t/eh9d+3pwvcjSvM1atWIKIoASTrjcvqdSQI/NCmjhmVZr5voTeBJ4LmlNhBC2MAfAB8BDgG/IIQ4tM7zGgyAEurPfulNhICOTg8/UFXCcSgXuXpa0Svqfft7qVR85marSCmpVwOKhRonXh/l9377O/z1n78OsGTb6IN3DiBBKYAool6LiKPGSaTkzddH+cGL55gcL1GtBGSzLoEfcertScbHivQP5OgfyDE+WuTkm2O8+v1LnHxzjLOnppieKK3Y/6h/IEelEtDVneGOuweo1UJAYts2/QMddORTOJ7NhfOzSwalL55fWekYbn7WOwP5JNBYES3Jw8BpKeXZxrZ/AXwcOLGecxsMoIR6GMZ0dqVBgNcoDpudqeA0nrdDZ/cIyyOddqjXQgI/olbzsSyBl3KYnSozsKNjIaD8j37hvkXHuHh+BgQIAWEoQTY/syxBGMa4rsW3vvEue/Z2M3JpDj+IcF2LIIg4d2aaTMZl9EqBC+dmyaTVn+PEeBG/HpLNeVw8P7PgvtL3m3T1JNNDOzvTpFIuvX0WCHUNSEBKSiX/qmZ3F8/P8O1vvsubr4/SkU9x274e05LiFmYzYga7geQS7TLwyCac13ALMDVZJt+VIgiiZpWwa1Gcr3PPe4eWjCUk4wBxDF3daaYmy8iGQK9VQ0YuFxgY7KCnL8e3v/EuXT2ZheOk0hYvffc8YRiTSjmU5v1F1xVLiRAC17GZna5y96FBMlmXKyMFquWAMIiYnanwzokJajUf17EozFWRUiClREooFX3OnZrmwtkZevuyuK5NtRoCkEpbHD92hdnpKtmcS/9AjmrVI9+Voqc3q851uUClHOC4Fvfet5Phfb0L38eFczNMjpeI45iORqzizLtTHLxrwMQXblFWVAZCiG8BQ20++qyU8murOEc7s0G2eU+f79PApwGGh4dXcXjDrUz/QI4wiLh8cQ5ACcxKgGULdg93Lpk2mlxRZ7IOoyPzxGGMjCVYAiEkAsG7JyfYPdzFhbOzdPdkyXemKMxWOP3uFJ2dadIZh4nx4tUXJkFYUCkHxDLmh0cvku9Kkc+nqFdDJsfLWDY4rsAvRNSqARKlCJIEQYQVCybGSrgpmx07OjjzziTvnpyguy9Ld0+aajXg/LkZfu6/eYCh3XfzV//pOBfPz1Ccr+PXI+IoZnqixLlT03hpm9v291Eu1REC5mZqDOzIkc4qUXBlpMDdhwZNfOEWZEVlIKX8yXWe4zKQtE/3AFeWOd/TwNMAR44cWVJpGAzQdJPsGe5mZrrMfKGG41h8/GfvZeTi/JKtHv7RL9y3UHCV60gRhTFYqv1EvR6BAMcWVMoB77w1gePYFOdrVCp1Al8pjSAICcOIOFr8a2pZICVUSgG1akj/QI4gCJmbibh0fhYhBLGMyXgeczNVLFsQxRIhoEUXICXIWBJLiSttJieLXPr/2zvzILuutLD/vru9rV9v6la31Fo9tmzJHmu8zGIMBAYCM05qBgNDYSpkKEhck4Qs/MVQU5Wi+IcipJIqloRhCFUkQ8w2KEyCJzNDBkIKbCNvsuWRLVlSa+mWWr2+futdT/44971+3Xq9SN3qRTq/KpXee/e8e78+ap3vnG+9NIdCMTNVpV4LQAkKxZ//6bf5mX/5FNWqT2lOV1kN/BixIAgirlyaIUng8sVZLEsYHOrCy1iUSg2yeRfXtahXw45Jb2uJ1jLlsHc2m5FncBJ4QEQOi4gH/Djw1U14ruEeoJlFO7y3m93D3Xz0Ow/x/M89zUeePrxslE5z13vgUD/PPnecH/6J4/T05mjUogXnLxDHiihKSBKwHUnNUEKjEZIkCYEfU6sGNy3gSbKwqPcP5CgUM1i2RaMeEceKRj3EsoQgTFBKEUcJKmHRs5feD6DRCKlXI5JEoRIIg4RqOaBWbVCZ97l0YYb//O//H9VKwO7hIrYteBkb17HxGzFxev/AjwnDmIlrZXI5V78PYsIgwXatm5rprKWZz51s+GPYHNblMxCRZ4FfBwaBPxeRN5VSPygie4HfUUo9o5SKRORnga8DNvC7Sql31i25wZDSzKJt7ky//tV3GRgs4DrWolIPpbk6oxdmiMKYEy+c4omn9nN9rMSf/dFpSnN17QxO79m+SxfRC7II2LZg26IXcJc0eqcz2ZzDrl1dIOA3hFotIG4qF4QwjvAbYNkrBmAAWhbh5pNDvORUcm1sHtezGd5bJAwTHEc7q5XS8uufS5vA/EbE9FSV/oECCqiWfR4+voePP3Nk0Y5+LaWuTTnsnc96o4lOACc6fD4OPNP2/kXgxfU8y2BYiU5lJWZnaiRKMTTcTRhGvHt6AkR46NhuqpWAP/5vbzAxXsbL2Ni2RTbnUq+GgF6g83mXajkgk7OJwoQ4VoRBRBjqnXyjEWjfQLt5R6ca4HoWw3u7CcOYOE4oz2vzUpM4VljpuTyJFbYj2hyUsCxLFcFSJA0eCoOY2ek6rmsRx0r7QdBRf0nS1HBascWxwm/EHH24/yYl0GQtpa5NOeydjylHYbgr6FRWYnCoyK6BAoUujwtnp8kVPI4+PERvf55Cl0d53qda8cnnPbyMTSbj4KTJZkmsCIME17OIQr1CV+brVCshSoHrCXGU2vQV2Lb2FVgiiEA26zKwW0csXR8vU634N8ncvvDvHuqi2J1tKYgmS9+vhMjC6aFeC+kqZrRPo3kqsLRfwnFSmR2b/Qf6ePJjB+jpyy27g2/mMrSz1K+wljGG7Y1RBoa7guX8A2GY8OxzxzlwXz+PfmgvPX251vUoTBARwjDGyzhUyn5r4bcdvcXX5hUhk7VbpiLHFlSy2LSjWBiXzTlkMg6j52fw/ZA4Tpb1BzSZm21on0TbbZvmqTWRjrVswfUsMlkby7bIZFxsu1lv1cJ1LQrFLJmMQ09vlsMP7Fp1B780A7taCW7yK6xljGF7Y2oTGe4KVmvU0um641p092apVHzKJV9HFKXEkWrZ2C3LIpt1qVd1dm8UKWBJ9I+CWCniWHH8ib10FbO8+reXaPiRTv5ahXotpF4Lac/fXM0stAi1ILflCUc/OIyX0Y1z8nmXiWtlLpyfolbRfotdgwXue2CAnt4c1Uqw4g5+LaWuTTnsnY9RBoa7gtUatXS6XuzOYEmW8aslfD9sLb62vWB7dxwd5VOrBgv29hSlaDmdk0SbdAS4cG4apRS1eojnOQwOFbh8cW5NP8ctKYBlCIOEq5fnGNzdxfSkrodkWToHo6uYxXEt9uztQSWKU6+PUSn7PHJ8D5dHZ1bsfbDawm7KYe9sZGmSy3biySefVK+++upWi2HYIawW5968PnphhnotIF/wsC3h1GtjRHECCiwbbNsi8Ndqn4FsTodutkcfNZ25jm0xcqCX0fObG2JpWYLjWgzuLuA4NtNTOru6q+hx+P4BLp6bIggTevtzHDzch+s6VKs+j314hLHL8yZXYAcjIq8ppZYtHLrs94wyMGwW2yEpqT3qKAwiTp+6xtRktWWqV0ujg1ZBBGxnwcksAgioNl3ieRZhmkuw2TgOFLqyJEmik+mAw/fvolL28TIOxx8fAXTY7dkzN5iZrtHfn8dxLKrVAN+PGNnfw0MPDxvFsEO4XWVgzESGTeF22jJu9PNfe+kKr79yhThOcFyL6ckqvq/zBNoX/9UUgZf2CGgmlzUVQVOjtC/6IhAEW6AFUqIIyvONNJTUwvOcNPS0RrE7A2hFcO7dSSplHxTMTFUJQq04LEu4dGEG17V59eXLDA51cfBwv1EMdyEmmsiwKaylo9idoj07tuGHzM3VmZqsUq8HqCTpWDxrOUS07b1TlE+npDBkbUllKxf+XR9Jmt3sujrSyfVsHMdqKcLxqyU8zyaOEuI4Jop1uQ2VKCwR/EbM+bOTiEC14pvs4rsUczIwbAp3OilpJRNUuyJq1EPCIAKl8wSaDmDb1vV9VgsB1Qti0PHacicKS6ClO9KktLV+dyPx/Yju3ixBEJPN6bDbaiWgWglwXRvLEuI4wbbtBVuZCIgi8PXY6ckaIkJff95kF99lGGVg2BRWC/28HZoK4Mzp64xfLZHNuPQP5ImCeJEJqqmISrN1Aj9GJQsJWCptQ2BZglJCwsomHcuSm6KKVmTp0C100dmWMF9qkM24OK6FiHB5dIYoTHA9myPHdvPmq2NpuQoWVVGNk4QojMnltJnpyqXZ1snCcHdgzESGTWGjk5Kapp/rY/OMXy2hEkWjEVAp+zd1Omtmx46PlXDTRVAnnC3cL4oTojhBVvkfoXMM1o5St/6d22U1U5PjalPQfLlBox6SzTmUSz5KKRxH6OsvMDhUaPVTsGzBdW1UAnaaCl3szuJ6+hRRq3Y+IRl2JuZkYNgUNjopqWn6uTw6g0oUmYxDnCgajZDu3tyiTmdPPLWfr3z5Ta5cnKVeDxeFgC5CbY655k7QrKUUhrr15lKyOadVJ2l+to4CyvM+uZxLLu8xX/KJwpihPT0kiQ6JDaOYwI9wfAvXsyh2Z/GyNkGgw2hzee9mQdbAdogqM9yMUQaGTWMjk5Kapp96NcTL2MSJwrJE1xNKO53t3dfDiRdO8e4717l0cUZ3CVsmaqjZU2AnousNCZZl4bi67LZKaJmzHEfS+kkWUZQQK4XrWIgl1OshUZSQL3j09OX4+DNH+NbXznL+vSksWzhwqI/JCd2LudEIiWNFT1+WoaEiwyPdHeVZabHf6qgyw/IYM5FhR9I0/eQKemcbx9qm7boW9VpIFOu2khfOTXL+7BSBHy+21y85FTi2hedZeJ6NZTeTtu5giM8Gok1RCdVqg2o5JI7Uoo5pUaQI/CiNENI/umULliWtXs1JkjB6QS/UjmPz2If3MbK/l9ELM3T3ZMlkHfIFjyRJqMz7nH9/itJc/aaIotX6GmxlVJlhZczJwLBpbKR5oFleoq9PLzj5gket4iNpE/gDh/pwHJt33rqmO5F1MAm1E0dJWtjNAV9XLI1CdgxLo6CWmrvE0n4R/UaPt+30u4nCsizqtYDB3V0tJ//sTE03vwlijhzdzYVzU5RLPpmczYee3Ifj2Hzly2/S158njBIGBguUZusr9jUwpa63L+ZkYNgUNroTVqvD2Ug3A7u76O3Lcf9Du/nO7/0Az//c02RzHjPTVaIwQaE6Zv82y0OLQJzG1cdxgoiF7cgdjf3fbATdlKdZURXRykEp7W/J5VzyBW9R5ddaNSSX9njo6c2RzbkMj3RT7M7R25cnimKuj88zenGm9W/6zqlrOnS3jfbF3pS63r6Yk4FhU7gTnbBW8kEMDBY4/94kKIVtWR17DNu2he2kPYYThVIJSaKTsey00qjjSmp20cpjzSWlNwHdo2D5dplNHNciX3ApdGWYna7hN2K8jEUc6+9mijafSntGt4f/5gsutao+dQHUqyG2I+QLWmGMXy2Ry7vEYbJg8unOcGl0lt7+hd1/+2K/WkFBw9ZhlIFhU7gT5oGVzE4jB7qpVnx8P2qFSeqtvl7YHUfo7s3iug7ZnEOSJJRLPuX5RisvTNKjgViCilX6fmudzJa1IFuz/WaLDuKJBfm8Qy7vkc25FLu1QggjhedZ9PXn2T1UZHikh+GRHr7y5Td5/2yDMIhRSocAj+zrRSUK27WoVwMOfWAXoE8OjmORKywsIwcP9/HOW9e16a7DYm9KXW9f1tsD+TPALwJHgY8opTpWlRORUaAMxEB0O0WUDDubjU46WykqBeCNk2McPNzPe2cm9M4+0a0ovYxDNufSVcy0irSBziqOwphX/maUIIh1RnKiFhWYW9pveLMRgb6BPCP7e5m+UWV2po5ItFBhtVNmcwLl+QDfj/EbEaEfMbSnGy/r4Hk2Rx8ZploJeO2lKzzx1H4StdAi0/VsCl0eubzL9HSNQ4f7dciuY6MShZM66w/dt7CQu67DI8f3UOjyVux9YBb/7cd6TwangR8GvriGsd+rlJpa5/MMO5T1mgeWngJWclQCFAoZdg8VyWQdzp65QRQpsjmbI0eHaNRDknTX2y5LJuPQ15/nxkQ5PU1Y2EoRpSYm2No8BMe1qMwHfOYnH+PLXzpJkigyWYfxq6W0rpC6aXzGc6jVA11iO/GJogS73KDPzVOvaiXSPKG99tIV8nkPvxERhQmeZ9PXn6enL8ezzx0HFv87HLqvn+mpKo6rlUNzHk2Y6M5kXcpAKXUGFo7TBsNyrMc80OkU8M6paxx7dBhYOGm0m52aJql9B/sodusGNrMzdYb3dreynr/1tbO8cXICheL+I4NcGyvhOBYCRGGMZVlpCCatfsbVarBllqIk0W06DxzqZ3CoyEPHhnT0FDBxvYzfCEniZnc27SARW+kEsiRpnXLmSz71Wsjg7i5g4YR26eIMN66XyWQcncDWoezE0l39UiVtTD47l83yGSjgGyKigC8qpX57uYEi8jzwPMCBAwc2STzDZnC75oFOzufVHJXtJqmevhyOa3P0g96iHa7vRzz48BD5vMuFc1O8/96UtoHnPYIgIkx7JLuurkcUhfGWugziCIpF7bxtN7sNDBb0Iu45+H5EkiyYtCrlsCWzSmj5FeJIMTtT4/r4PJYtfPff/wBnTl/HsrR5CMD1bMIwXrHshDH53D2sGloqIn8hIqc7/Pn0LTznaaXU48AngX8hIt+93ECl1G8rpZ5USj05ODh4C48w3K10anZ/8HAflXLnWkdrqYPUrmDmSw0uXZzBdS3dvF4pHNemuzuDbVscObqbbM5dUy/jO4llQRgmnHjhlHaQpz9jed6n2J0hTpK0blD6hQ4OZZ1oph3ocaw49+4N6rWQ1166QqMeMjdbZ+zyHJMTZSplf11lJww7i1VPBkqp71/vQ5RS4+nfN0TkBPAR4K/Xe1/DvUEn5/NqjsrVTFLt0U3jYyWSRJFL7eWeZ+M3IqqVAMsSbkxUqNeDLa/SKZaQydhcH5/ntZcvYztCed5nvtSg2J2huztL/2CB8SslqtXgptwKrQQsslldp8hxbTzP4cDBPq6PzzM5USGbdRAR/CAmDBscuq9/kYPYcPdyx81EIlIALKVUOX39A8Av3ennGu4elnM+r+SoXM180a5gmvWNojAmm3MY2N1FoxZy9fIc2bxDrabNJBvhPBaL22p/6XkWjmvRaERcfH8aBGqzAb39OcIgxnVtPS+VAMuydOawHy2qmKqAXM4linRkkVIQBBHjY3NcODdNEisqFZ+B3V3s2l2gXguZm6nfdmVZw85ivaGlzwK/DgwCfy4ibyqlflBE9gK/o5R6BhgCTqROZgf470qp/71OuQ33EBsdm355dIbSXJ3Tp67RVcwQhhF+I6ZeC7FsYW62Tpwuon4jJpcTvIxzS72Rl+OmMhGr3FPSTmm5gqdrCyndpKaZb5AveCSxojzfIJ93mS81SGJFoxGRzTl4ipa/I44V5bIuWW1bguNY2LbF6VPXQUGu4OLXIyrzPnGk6O3LkSt4xidwj7DeaKITwIkOn48Dz6SvLwDH1/Mcg2GjHJXtkUkPPzrMuTOTTE5U0mxkXfW0nShMqCV6t207Fskt9jN2HFncz2DJwr9sdzRLXyt0eQRBrHMMdhWYma6lRfl0GQmAQtHT1UR7swR+mUQUjqP7Hff0ZKnVfGam6ngZu6V8kliX3pibrQPavRBFCbmCS3dvDs+zOXCof5Fpbuk8mjLUdxcmA9lwT/GtF89ybWyeONQLn+Na7NnXQ2mujt+Isey09Rl6QY2iBKWgXg/JZl2iyL+lvge3WxZb5zmAl3V46OEhJq6VEVsQC+JQ6d4EUcLUjQrZnEtPX5YDh/o5+sFhpiarOI5wbWyeWjUk8BPyeZc40cqp0OUxN1MnSrOXBQjjBNXQ1U6nJiooIJN1+MxPPnaTbKYM9d2JKVRnuGe4PDrDO6eugYJs3iEIYiYnqumuWOFl7HRHbekqn6ki0CWiY6IoxvNsMlkH21lbZFG7f2DN6TjpuHxem4Acx+YDRwbS++l+BV1dHrZr4TciZqZqJHHMu6evM3phhsmJMo16xNFHhrn/yABhqOtWZ3MOlm0xNVnFcgTPc3BsiyTtdWBZFo5j62gqx8JaRmBThvruxJwMDPcMr710hUJ3BkEnSnqejZexmJutk8s5xGn9oSiKUYkgolqF6VRMGrZpYVngeXlqVZ9q1SfuEGTkZayFMhGA6wpiWa12m2JBFC5zalBgOcKDDw/hODazMzUGvSLHHx+hNFunXguZna1Rmm0QhhFxpJi4VuHYB/cwtKfIxLUy75+dBGB2uqbDSBNFd08WgPJ8g9CPKXZnyeYcZmfq2JaQKEXfLt1D+v4HB3Fcu2MhQVOG+u7EnAwM9wxTk1UOHu4jCGLCQCeQ5Qsevh+zd18PXppsFccKRaoIBGxb735RcP+RAbyMS6MRppVOJc34XfwsQSdv2Y5gO4JlW+waKNDbn0NEWg5qkQ4nBoF8zmVivEw26zA4VOSf/qvv4NnnjhNGCUN7ihx9ZJiPfdchvuvj9zO0p0ihK8PwSDdiCcMj3dz/wCCzMzVmZ+r09ufoKmYQS/AyNvmcS6x0T+Su7izHHh1GIdi2bu5z/4OD9PTlll3gTRnquxNzMjDcMzTDSR94aJDxqyVq1ZBM1uXRx/cwsr8Px7XJF2qMXZnTjmQBSyCXd/EyNvVaxMXz0ziOTeCHxIk229i2tDJ+RQDRyWGuZ9E/0EUmYzM9VaM838BxLd1kppkV3NaPuVkm23Ftdu3uQizh0ujsooJ6nXIuyvM+xZ7Mop91aE8Rp60QXRTGjI+VKM02QATHFvJdHntGunFdh97+HPsO9DK8d6GVZfsC3+4wdl2L6akqQ8Pdpgz1XYRRBoZ7hma+QqGQ4aFjQx3zFS6PzvDF//g3zEzVyOYcqpUgLfGgiOMYvyF0D2XJd3k4jnDx3PTiWH61oBB6enO4ro3j2AwMFJidqVGZD8hkLJSCMEhayqCpCBDo7cuRyTgEQUSlEjByoJsTL5zquBBPXCtTLjcIw5gzp69TLGYol33KJa0gjj+xhzdOjlEoZNgz0q0L/BUzHNk/yMx0nXfeus4jx/fw6R97hDdOjnUsPb3UYVyrhVgiRGHM9HRoahLdJYjayjKMq/Dkk0+qV1/tWBXbYLgtVguJPPHCKa6PzfPtt6+lfQOgVglIEshkbTJZXf76gYcGqddCTr12lXrtZqdBvsuhf1fXolNILu9w5eIcwyNFvIyjyz/M1KnVQuIoZnikm77+HNfHy7qUhlJksi49PVlGDvQyNFykVguZnCjT159ndrbO9I0KfbtyTE/ViMI4LU2RxXEt9h/sw7KFxz48wtjleV5/5Qq2a3HocD89fTlgoYbTs88dX3ZuTrxw6qbTSPv3DNsLEXntdtoEmJOB4Z5itXyFqckqQ3uKRFHM2TM3SBLId3lU5n0AentzDAwVGL9a4trVElGU4DgQRbSigCxp7vpVyxwkaNORWKo1ULeRdAn8mFrVZ/dQkevj82SyDo16hFKKMIwIgpirl+fI5V16enNAkUKXR09fjmzG4fr4PKXZOrVqCKL7DDzx6AF6enNUKwFjl+d59rnjLcevtNVYavcLLDc3xmF8b2CUgcHAwonh8oUZxq7McehwP0989ADjYyXKJR/XtXngod3kci7vvzeJ49nEaYN5EcF2dAiooHMLBKFeDTjzzgS5nIvt6PfZnEepVCeb1SeDIFUaRx8eotidpVRqMD1ZxfVsenqzTE/VCMOIQpfH+NUSPb0Ljt1K2efalRLVqo/tWNhpNze/sVBddWn/4dtpMLTRjYkM2xOjDAz3PO028fuO7OLtN8b5u5cu0dWVId/lMTBU4Ls+fh9vnBxj9OIMjqsXfRS6ymnagziOYiwLXM9haE8Xvh+Ty7vUq0GrMFwu71IuNSg1IgSwHR3bH0UJs7N1Hv3QXl4/eVVXaRVdkyjw09pDVR3B01yIJyfK1OoBjmtj26Kd01FCHMe89Yb2EzTNQrD2BkNLzUUjB7p54+TYqt8z7GxMaKnhnmdREpUIjmNjWdIqUGeJMDzSwyd/6BhxmBCGMUmswJKbmtEnCWQ8i737esl4DgcP6yilnt4cvWmZhyCI6e3N0berwPDebh59fITBoSL1WkCtFpIvuIShvnGzdHatFpDLO4vKcecLXtpjQZujmkolScBvRK3TyOxMjcujM60aT4Uuj+npGoUu76as4aZirFaCVnbxGyfHeOzDIyt+z7DzMScDwz3PonLWV0upPT5LoxZx/PGRVo/gZ587zuMf3U+1EnB5dIYoTHDdhZIVkIaIIli28IEHBxi9MIPn2a2GMaTKplD0OPbBPS0ZVKLIFzyqVZ++/jxXLs0Spt3WDh7uZ3a2TqFLK6xm5M7Bw/3cuK6dzWGY4GVscrFLFCZYli6ud+gDu3CcheSx1XwmnRoJAS2/g+HuxSgDwz1Pu028Vg3J512CMCZX0A112u3uTVPL/FyDKIqxbQvXs3E9izjSNYOSRPHJHzoGwK//8v+lUMykzuCEKIjZNZCnXPIXyVCrhRw83M8TT+3ntZeu4PsRtWqgy0yHMYNDXfT25SjN1vn6V99tmW9Gz3cRR/Pk8i4g1KohPf1ZHjm+LXzNEgAACgdJREFUN3U2a0WzVmevcRbfuxhlYLjnabel55r9CxQtW3u7s7Rpavnyl04yPVkFFNmci+vaxLHC9Wz27utp7b4fPr6H0YszNGoRuYLLocP91GshY1fmOsb0t+/c230ZYRhx5u3rIMJDx3a3zDdPf+9hvv3WBOffm0KRMDRcZHiku6UIlsq/GsZZfO9ilIHhnqe9X0KhK0OtGrJvfy9KKU69Pkal7PPI8T1cHp0BtCkll/foKnqt1ppRmBCEMZmMw/c/c6R1748/c6S1oDcXfsuJ+NSPPcLY5fkV+zO0m2zOnJ4hX/BQwLXxeY4+Mgxo881P/bOPtb7TbvO/HWfvWp3MhrsPk3RmMCzh8ugM33rxLO+cukahO8PBw324rsPkRJlEqUXZv++dmSDwYxSwe3cX/+BHjvGRpw/fdL/bqf3/pV/721ZewGuvXCGfd1EoGrWIxz+6n7mZGufPTXHwvl2L7rveXgOmV8HOxiSdGQwbxIFDOkP3+JP7FplL3j/bAOC++/VnwyPdFHuyq2bi3m5jnnaTTb7gUplvUK3qAnlvnrzKfLlBT2+uY0+B9SzeG9VIyLCzWFdoqYj8qoi8KyJvicgJEeldZtwnROQ9EXlfRD6/nmcaDJvB1GRVx/q3EQYxUbi409mddK4+8dR+qlWfaiWgq+gxPVXD9yN6+7JMT1epVQJ2DeRNTwHDhrDePINvAo8opR4FzgK/sHSAiNjAbwKfBI4Bz4nIsXU+12C4o3Qq0+x6No67+L/MnXSutucF3LheYddggaGhIkoJSaxwHYuz357kzOnrlObqJurHsC7W2wP5G21vXwZ+tMOwjwDvp72QEZE/AD4NfHs9zzYY7iSdHKk9PVkSpW7bOXs7NE027XWFSnN1Xp2ttepfh0HMuXcnbypBbTDcChuZgfzTwNc6fD4CtJ9dr6afGQzblk7Zuj/yjz7EZ37ysS3JxG0/qYxfLdHdnU37JOsyFCIwdnmOJ57af8dlMdydrHoyEJG/AIY7XPqCUurP0jFfACLg9zvdosNny4YwicjzwPMABw4cWE08g+GO0dyVN6NrmsleWxFd035SqVYCvIxDoStDNufQqEXk8x65gnZsN3sfmEggw62w7tBSEfks8Dng+5RStQ7XnwJ+USn1g+n7XwBQSv3yavc2oaWGraY98avdNLQVtXmaSmm5vgRRGOMH0baQ1bB1bEloqYh8Avh54O91UgQpJ4EHROQwMAb8OPAT63muwbBZLFerp1Oj+DtN86TSPCU4ro1KVGvRz2ScbSOrYeexXp/BbwBF4Jsi8qaI/BaAiOwVkRcBlFIR8LPA14EzwB8ppd5Z53MNhk2hU4jpVkftLFd9NAyTbSerYeew3mii+5f5fBx4pu39i8CL63mWwbAVbNdaPZ0Sw7arrIadgelnYDCsQHvil0rUon4C242dJKth+2FqExkMq7CZtXpMXSHDerldB7JRBgbDNmE7RS4Zdi63qwyMmchg2CYsar9p6g0ZNhmjDAyGbcJ2jFwy3DsYZWAwbBM6Fccz0UCGzcIoA4Nhm2CigQxbiVEGBsM2YblkMuM8NmwGptOZwbCNMF3GDFuFORkYDAaDwSgDg8FgMBhlYDAYDAaMMjAYDAYDRhkYDAaDAaMMDAaDwYBRBgaDwWBgm1ctFZFJ4FKHSwPA1CaLc6sYGTcGI+PGYGTcGHaCjA8qpYq3+qVtnXSmlBrs9LmIvHo7JVo3EyPjxmBk3BiMjBvDTpHxdr5nzEQGg8FgMMrAYDAYDDtXGfz2VguwBoyMG4ORcWMwMm4Md62M29qBbDAYDIbNYaeeDAwGg8GwgewIZSAivyoi74rIWyJyQkR6lxn3CRF5T0TeF5HPb7KMnxGRd0QkEZFlow1EZFRE3haRN2/X678JMm7lPPaLyDdF5Fz6d98y4+J0Dt8Uka9ukmwrzouIZETkD9Prr4jIoc2Q6xZl/CkRmWybu3+yyfL9rojcEJHTy1wXEfm1VP63ROTxzZRvjTJ+j4iU2ubw326BjPtF5C9F5Ez6f/pfdxhza3OplNr2f4AfAJz09a8Av9JhjA2cB+4DPOAUcGwTZTwKPAj8FfDkCuNGgYEtmsdVZdwG8/jvgM+nrz/f6d86vVbZ5LlbdV6Afw78Vvr6x4E/3IYy/hTwG1vx+5c+/7uBx4HTy1x/BvgaIMDHgFe2oYzfA/yvrZrDVIY9wOPp6yJwtsO/9S3N5Y44GSilvqGUitK3LwP7Ogz7CPC+UuqCUioA/gD49CbKeEYp9d5mPe92WKOMWzqP6bN+L339e8APbeKzV2It89Iu+58A3yciss1k3FKUUn8NzKww5NPAf1Wal4FeEdmzOdJp1iDjlqOUuqaUej19XQbOACNLht3SXO4IZbCEn0Zru6WMAFfa3l/l5snZDijgGyLymog8v9XCdGCr53FIKXUN9C88sHuZcVkReVVEXhaRzVAYa5mX1ph081ICdm2CbDc9P2W5f7sfSc0GfyIi263B8lb//q2Vp0TklIh8TUQe3kpBUnPkY8ArSy7d0lxumwxkEfkLYLjDpS8opf4sHfMFIAJ+v9MtOny2oaFSa5FxDTytlBoXkd3AN0Xk3XQnsl1k3NJ5vIXbHEjn8T7gWyLytlLq/MZI2JG1zMsdn7tVWMvz/yfwglLKF5HPoU8yH7/jkq2drZ7DtfA6cFApVRGRZ4D/ATywFYKISBfwFeDfKKXml17u8JVl53LbKAOl1PevdF1EPgv8Q+D7VGoQW8JVoH2Xsw8Y3zgJV5dxjfcYT/++ISIn0Ef7DVMGGyDjls6jiEyIyB6l1LX0SHtjmXs05/GCiPwVemd0J5XBWualOeaqiDhAD5trblhVRqXUdNvbL6F9cNuJO/77t17aF12l1Isi8p9EZEAptak1i0TERSuC31dK/WmHIbc0lzvCTCQinwB+HviUUqq2zLCTwAMiclhEPLQDb1OiTNaKiBREpNh8jXaMd4xY2EK2eh6/Cnw2ff1Z4KbTjIj0iUgmfT0APA18+w7LtZZ5aZf9R4FvLbNx2TIZl9iMP4W2NW8nvgr84zQS5mNAqWk23C6IyHDTFyQiH0Gvo9Mrf2vDZRDgvwBnlFL/YZlhtzaXW+kRvwXP+fto29eb6Z9mxMZe4MUl3vOz6B3iFzZZxmfRmtgHJoCvL5URHeVxKv3zznaUcRvM4y7g/wDn0r/708+fBH4nff0dwNvpPL4N/MwmyXbTvAC/hN6kAGSBP05/X/8OuG8z526NMv5y+rt3CvhL4KFNlu8F4BoQpr+LPwN8Dvhcel2A30zlf5sVIvO2UMafbZvDl4Hv2AIZvxNt8nmrbV18Zj1zaTKQDQaDwbAzzEQGg8FguLMYZWAwGAwGowwMBoPBYJSBwWAwGDDKwGAwGAwYZWAwGAwGjDIwGAwGA0YZGAwGgwH4/yfYa1wkH7MgAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# write your code here\n",
    "\n",
    "\n",
    "from sklearn.cluster import KMeans \n",
    "k = 3\n",
    "k_means3 = KMeans(init = \"k-means++\", n_clusters = k, n_init = 12)\n",
    "k_means3.fit(X)\n",
    "fig = plt.figure(figsize=(6, 4))\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "for k, col in zip(range(k), colors):\n",
    "    my_members = (k_means3.labels_ == k)\n",
    "    plt.scatter(X[my_members, 0], X[my_members, 1],  c=col, marker=u'o', alpha=0.5)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Double-click __here__ for the solution.\n",
    "\n",
    "<!-- Your answer is below:\n",
    "\n",
    "\n",
    "from sklearn.cluster import KMeans \n",
    "k = 3\n",
    "k_means3 = KMeans(init = \"k-means++\", n_clusters = k, n_init = 12)\n",
    "k_means3.fit(X)\n",
    "fig = plt.figure(figsize=(6, 4))\n",
    "ax = fig.add_subplot(1, 1, 1)\n",
    "for k, col in zip(range(k), colors):\n",
    "    my_members = (k_means3.labels_ == k)\n",
    "    plt.scatter(X[my_members, 0], X[my_members, 1],  c=col, marker=u'o', alpha=0.5)\n",
    "plt.show()\n",
    "\n",
    "\n",
    "-->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "\n",
    "<h1 id=\"weather_station_clustering\" align=\"center\"> Weather Station Clustering using DBSCAN & scikit-learn </h1>\n",
    "<hr>\n",
    "\n",
    "DBSCAN is specially very good for tasks like class identification on a spatial context. The wonderful attribute of DBSCAN algorithm is that it can find out any arbitrary shape cluster without getting affected by noise. For example, this following example cluster the location of weather stations in Canada.\n",
    "<br>\n",
    "DBSCAN can be used here, for instance, to find the group of stations which show the same weather condition. As you can see, it not only finds different arbitrary shaped clusters, can find the denser part of data-centered samples by ignoring less-dense areas or noises.\n",
    "\n",
    "let's start playing with the data. We will be working according to the following workflow: </font>\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### About the dataset\n",
    "\n",
    "\t\t\n",
    "<h4 align = \"center\">\n",
    "Environment Canada    \n",
    "Monthly Values for July - 2015\t\n",
    "</h4>\n",
    "<html>\n",
    "<head>\n",
    "<style>\n",
    "table {\n",
    "    font-family: arial, sans-serif;\n",
    "    border-collapse: collapse;\n",
    "    width: 100%;\n",
    "}\n",
    "\n",
    "td, th {\n",
    "    border: 1px solid #dddddd;\n",
    "    text-align: left;\n",
    "    padding: 8px;\n",
    "}\n",
    "\n",
    "tr:nth-child(even) {\n",
    "    background-color: #dddddd;\n",
    "}\n",
    "</style>\n",
    "</head>\n",
    "<body>\n",
    "\n",
    "<table>\n",
    "  <tr>\n",
    "    <th>Name in the table</th>\n",
    "    <th>Meaning</th>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"green\"><strong>Stn_Name</font></td>\n",
    "    <td><font color = \"green\"><strong>Station Name</font</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"green\"><strong>Lat</font></td>\n",
    "    <td><font color = \"green\"><strong>Latitude (North+, degrees)</font></td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"green\"><strong>Long</font></td>\n",
    "    <td><font color = \"green\"><strong>Longitude (West - , degrees)</font></td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>Prov</td>\n",
    "    <td>Province</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>Tm</td>\n",
    "    <td>Mean Temperature (°C)</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwTm</td>\n",
    "    <td>Days without Valid Mean Temperature</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>D</td>\n",
    "    <td>Mean Temperature difference from Normal (1981-2010) (°C)</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"black\">Tx</font></td>\n",
    "    <td><font color = \"black\">Highest Monthly Maximum Temperature (°C)</font></td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwTx</td>\n",
    "    <td>Days without Valid Maximum Temperature</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"black\">Tn</font></td>\n",
    "    <td><font color = \"black\">Lowest Monthly Minimum Temperature (°C)</font></td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwTn</td>\n",
    "    <td>Days without Valid Minimum Temperature</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>S</td>\n",
    "    <td>Snowfall (cm)</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwS</td>\n",
    "    <td>Days without Valid Snowfall</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>S%N</td>\n",
    "    <td>Percent of Normal (1981-2010) Snowfall</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td><font color = \"green\"><strong>P</font></td>\n",
    "    <td><font color = \"green\"><strong>Total Precipitation (mm)</font></td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwP</td>\n",
    "    <td>Days without Valid Precipitation</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>P%N</td>\n",
    "    <td>Percent of Normal (1981-2010) Precipitation</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>S_G</td>\n",
    "    <td>Snow on the ground at the end of the month (cm)</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>Pd</td>\n",
    "    <td>Number of days with Precipitation 1.0 mm or more</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>BS</td>\n",
    "    <td>Bright Sunshine (hours)</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>DwBS</td>\n",
    "    <td>Days without Valid Bright Sunshine</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>BS%</td>\n",
    "    <td>Percent of Normal (1981-2010) Bright Sunshine</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>HDD</td>\n",
    "    <td>Degree Days below 18 °C</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>CDD</td>\n",
    "    <td>Degree Days above 18 °C</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>Stn_No</td>\n",
    "    <td>Climate station identifier (first 3 digits indicate   drainage basin, last 4 characters are for sorting alphabetically).</td>\n",
    "  </tr>\n",
    "  <tr>\n",
    "    <td>NA</td>\n",
    "    <td>Not Available</td>\n",
    "  </tr>\n",
    "\n",
    "\n",
    "</table>\n",
    "\n",
    "</body>\n",
    "</html>\n",
    "\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1-Download data\n",
    "<div id=\"download_data\">\n",
    "    To download the data, we will use <b>!wget</b> to download it from IBM Object Storage.<br> \n",
    "    <b>Did you know?</b> When it comes to Machine Learning, you will likely be working with large datasets. As a business, where can you host your data? IBM is offering a unique opportunity for businesses, with 10 Tb of IBM Cloud Object Storage: <a href=\"http://cocl.us/ML0101EN-IBM-Offer-CC\">Sign up now for free</a>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--2019-06-18 16:06:07--  https://s3-api.us-geo.objectstorage.softlayer.net/cf-courses-data/CognitiveClass/ML0101ENv3/labs/weather-stations20140101-20141231.csv\n",
      "Resolving s3-api.us-geo.objectstorage.softlayer.net (s3-api.us-geo.objectstorage.softlayer.net)... 67.228.254.193\n",
      "Connecting to s3-api.us-geo.objectstorage.softlayer.net (s3-api.us-geo.objectstorage.softlayer.net)|67.228.254.193|:443... connected.\n",
      "HTTP request sent, awaiting response... 200 OK\n",
      "Length: 129821 (127K) [text/csv]\n",
      "Saving to: ‘weather-stations20140101-20141231.csv’\n",
      "\n",
      "weather-stations201 100%[===================>] 126.78K  --.-KB/s    in 0.06s   \n",
      "\n",
      "2019-06-18 16:06:07 (1.93 MB/s) - ‘weather-stations20140101-20141231.csv’ saved [129821/129821]\n",
      "\n"
     ]
    }
   ],
   "source": [
    "!wget -O weather-stations20140101-20141231.csv https://s3-api.us-geo.objectstorage.softlayer.net/cf-courses-data/CognitiveClass/ML0101ENv3/labs/weather-stations20140101-20141231.csv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2- Load the dataset\n",
    "<div id=\"load_dataset\">\n",
    "We will import the .csv then we creates the columns for year, month and day.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Stn_Name</th>\n",
       "      <th>Lat</th>\n",
       "      <th>Long</th>\n",
       "      <th>Prov</th>\n",
       "      <th>Tm</th>\n",
       "      <th>DwTm</th>\n",
       "      <th>D</th>\n",
       "      <th>Tx</th>\n",
       "      <th>DwTx</th>\n",
       "      <th>Tn</th>\n",
       "      <th>...</th>\n",
       "      <th>DwP</th>\n",
       "      <th>P%N</th>\n",
       "      <th>S_G</th>\n",
       "      <th>Pd</th>\n",
       "      <th>BS</th>\n",
       "      <th>DwBS</th>\n",
       "      <th>BS%</th>\n",
       "      <th>HDD</th>\n",
       "      <th>CDD</th>\n",
       "      <th>Stn_No</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>CHEMAINUS</td>\n",
       "      <td>48.935</td>\n",
       "      <td>-123.742</td>\n",
       "      <td>BC</td>\n",
       "      <td>8.2</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>13.5</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.0</td>\n",
       "      <td>12.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>273.3</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1011500</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>COWICHAN LAKE FORESTRY</td>\n",
       "      <td>48.824</td>\n",
       "      <td>-124.133</td>\n",
       "      <td>BC</td>\n",
       "      <td>7.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>15.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>104.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>12.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>307.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012040</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>LAKE COWICHAN</td>\n",
       "      <td>48.829</td>\n",
       "      <td>-124.052</td>\n",
       "      <td>BC</td>\n",
       "      <td>6.8</td>\n",
       "      <td>13.0</td>\n",
       "      <td>2.8</td>\n",
       "      <td>16.0</td>\n",
       "      <td>9.0</td>\n",
       "      <td>-2.5</td>\n",
       "      <td>...</td>\n",
       "      <td>9.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>11.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>168.1</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012055</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>DISCOVERY ISLAND</td>\n",
       "      <td>48.425</td>\n",
       "      <td>-123.226</td>\n",
       "      <td>BC</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>12.5</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>1012475</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>DUNCAN KELVIN CREEK</td>\n",
       "      <td>48.735</td>\n",
       "      <td>-123.728</td>\n",
       "      <td>BC</td>\n",
       "      <td>7.7</td>\n",
       "      <td>2.0</td>\n",
       "      <td>3.4</td>\n",
       "      <td>14.5</td>\n",
       "      <td>2.0</td>\n",
       "      <td>-1.0</td>\n",
       "      <td>...</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>11.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>267.7</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012573</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 25 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                 Stn_Name     Lat     Long Prov   Tm  DwTm    D    Tx  DwTx  \\\n",
       "0               CHEMAINUS  48.935 -123.742   BC  8.2   0.0  NaN  13.5   0.0   \n",
       "1  COWICHAN LAKE FORESTRY  48.824 -124.133   BC  7.0   0.0  3.0  15.0   0.0   \n",
       "2           LAKE COWICHAN  48.829 -124.052   BC  6.8  13.0  2.8  16.0   9.0   \n",
       "3        DISCOVERY ISLAND  48.425 -123.226   BC  NaN   NaN  NaN  12.5   0.0   \n",
       "4     DUNCAN KELVIN CREEK  48.735 -123.728   BC  7.7   2.0  3.4  14.5   2.0   \n",
       "\n",
       "    Tn  ...  DwP    P%N  S_G    Pd  BS  DwBS  BS%    HDD  CDD   Stn_No  \n",
       "0  1.0  ...  0.0    NaN  0.0  12.0 NaN   NaN  NaN  273.3  0.0  1011500  \n",
       "1 -3.0  ...  0.0  104.0  0.0  12.0 NaN   NaN  NaN  307.0  0.0  1012040  \n",
       "2 -2.5  ...  9.0    NaN  NaN  11.0 NaN   NaN  NaN  168.1  0.0  1012055  \n",
       "3  NaN  ...  NaN    NaN  NaN   NaN NaN   NaN  NaN    NaN  NaN  1012475  \n",
       "4 -1.0  ...  2.0    NaN  NaN  11.0 NaN   NaN  NaN  267.7  0.0  1012573  \n",
       "\n",
       "[5 rows x 25 columns]"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import csv\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "filename='weather-stations20140101-20141231.csv'\n",
    "\n",
    "#Read csv\n",
    "pdf = pd.read_csv(filename)\n",
    "pdf.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3-Cleaning\n",
    "<div id=\"cleaning\">\n",
    "Lets remove rows that don't have any value in the <b>Tm</b> field.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Stn_Name</th>\n",
       "      <th>Lat</th>\n",
       "      <th>Long</th>\n",
       "      <th>Prov</th>\n",
       "      <th>Tm</th>\n",
       "      <th>DwTm</th>\n",
       "      <th>D</th>\n",
       "      <th>Tx</th>\n",
       "      <th>DwTx</th>\n",
       "      <th>Tn</th>\n",
       "      <th>...</th>\n",
       "      <th>DwP</th>\n",
       "      <th>P%N</th>\n",
       "      <th>S_G</th>\n",
       "      <th>Pd</th>\n",
       "      <th>BS</th>\n",
       "      <th>DwBS</th>\n",
       "      <th>BS%</th>\n",
       "      <th>HDD</th>\n",
       "      <th>CDD</th>\n",
       "      <th>Stn_No</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>CHEMAINUS</td>\n",
       "      <td>48.935</td>\n",
       "      <td>-123.742</td>\n",
       "      <td>BC</td>\n",
       "      <td>8.2</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>13.5</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.0</td>\n",
       "      <td>12.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>273.3</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1011500</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>COWICHAN LAKE FORESTRY</td>\n",
       "      <td>48.824</td>\n",
       "      <td>-124.133</td>\n",
       "      <td>BC</td>\n",
       "      <td>7.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>3.0</td>\n",
       "      <td>15.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>-3.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>104.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>12.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>307.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012040</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>LAKE COWICHAN</td>\n",
       "      <td>48.829</td>\n",
       "      <td>-124.052</td>\n",
       "      <td>BC</td>\n",
       "      <td>6.8</td>\n",
       "      <td>13.0</td>\n",
       "      <td>2.8</td>\n",
       "      <td>16.0</td>\n",
       "      <td>9.0</td>\n",
       "      <td>-2.5</td>\n",
       "      <td>...</td>\n",
       "      <td>9.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>11.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>168.1</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012055</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>DUNCAN KELVIN CREEK</td>\n",
       "      <td>48.735</td>\n",
       "      <td>-123.728</td>\n",
       "      <td>BC</td>\n",
       "      <td>7.7</td>\n",
       "      <td>2.0</td>\n",
       "      <td>3.4</td>\n",
       "      <td>14.5</td>\n",
       "      <td>2.0</td>\n",
       "      <td>-1.0</td>\n",
       "      <td>...</td>\n",
       "      <td>2.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>11.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>267.7</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012573</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>ESQUIMALT HARBOUR</td>\n",
       "      <td>48.432</td>\n",
       "      <td>-123.439</td>\n",
       "      <td>BC</td>\n",
       "      <td>8.8</td>\n",
       "      <td>0.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>13.1</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.9</td>\n",
       "      <td>...</td>\n",
       "      <td>8.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>12.0</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>258.6</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1012710</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 25 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                 Stn_Name     Lat     Long Prov   Tm  DwTm    D    Tx  DwTx  \\\n",
       "0               CHEMAINUS  48.935 -123.742   BC  8.2   0.0  NaN  13.5   0.0   \n",
       "1  COWICHAN LAKE FORESTRY  48.824 -124.133   BC  7.0   0.0  3.0  15.0   0.0   \n",
       "2           LAKE COWICHAN  48.829 -124.052   BC  6.8  13.0  2.8  16.0   9.0   \n",
       "3     DUNCAN KELVIN CREEK  48.735 -123.728   BC  7.7   2.0  3.4  14.5   2.0   \n",
       "4       ESQUIMALT HARBOUR  48.432 -123.439   BC  8.8   0.0  NaN  13.1   0.0   \n",
       "\n",
       "    Tn  ...  DwP    P%N  S_G    Pd  BS  DwBS  BS%    HDD  CDD   Stn_No  \n",
       "0  1.0  ...  0.0    NaN  0.0  12.0 NaN   NaN  NaN  273.3  0.0  1011500  \n",
       "1 -3.0  ...  0.0  104.0  0.0  12.0 NaN   NaN  NaN  307.0  0.0  1012040  \n",
       "2 -2.5  ...  9.0    NaN  NaN  11.0 NaN   NaN  NaN  168.1  0.0  1012055  \n",
       "3 -1.0  ...  2.0    NaN  NaN  11.0 NaN   NaN  NaN  267.7  0.0  1012573  \n",
       "4  1.9  ...  8.0    NaN  NaN  12.0 NaN   NaN  NaN  258.6  0.0  1012710  \n",
       "\n",
       "[5 rows x 25 columns]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pdf = pdf[pd.notnull(pdf[\"Tm\"])]\n",
    "pdf = pdf.reset_index(drop=True)\n",
    "pdf.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4-Visualization\n",
    "<div id=\"visualization\">\n",
    "Visualization of stations on map using basemap package. The matplotlib basemap toolkit is a library for plotting 2D data on maps in Python. Basemap does not do any plotting on it’s own, but provides the facilities to transform coordinates to a map projections. <br>\n",
    "\n",
    "Please notice that the size of each data points represents the average of maximum temperature for each station in a year.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from mpl_toolkits.basemap import Basemap\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import rcParams\n",
    "%matplotlib inline\n",
    "rcParams['figure.figsize'] = (14,10)\n",
    "\n",
    "llon=-140\n",
    "ulon=-50\n",
    "llat=40\n",
    "ulat=65\n",
    "\n",
    "pdf = pdf[(pdf['Long'] > llon) & (pdf['Long'] < ulon) & (pdf['Lat'] > llat) &(pdf['Lat'] < ulat)]\n",
    "\n",
    "my_map = Basemap(projection='merc',\n",
    "            resolution = 'l', area_thresh = 1000.0,\n",
    "            llcrnrlon=llon, llcrnrlat=llat, #min longitude (llcrnrlon) and latitude (llcrnrlat)\n",
    "            urcrnrlon=ulon, urcrnrlat=ulat) #max longitude (urcrnrlon) and latitude (urcrnrlat)\n",
    "\n",
    "my_map.drawcoastlines()\n",
    "my_map.drawcountries()\n",
    "# my_map.drawmapboundary()\n",
    "my_map.fillcontinents(color = 'white', alpha = 0.3)\n",
    "my_map.shadedrelief()\n",
    "\n",
    "# To collect data based on stations        \n",
    "\n",
    "xs,ys = my_map(np.asarray(pdf.Long), np.asarray(pdf.Lat))\n",
    "pdf['xm']= xs.tolist()\n",
    "pdf['ym'] =ys.tolist()\n",
    "\n",
    "#Visualization1\n",
    "for index,row in pdf.iterrows():\n",
    "#   x,y = my_map(row.Long, row.Lat)\n",
    "   my_map.plot(row.xm, row.ym,markerfacecolor =([1,0,0]),  marker='o', markersize= 5, alpha = 0.75)\n",
    "#plt.text(x,y,stn)\n",
    "plt.show()\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5- Clustering of stations based on their location i.e. Lat & Lon\n",
    "<div id=\"clustering\">\n",
    "    <b>DBSCAN</b> form sklearn library can runs DBSCAN clustering from vector array or distance matrix.<br>\n",
    "    In our case, we pass it the Numpy array Clus_dataSet to find core samples of high density and expands clusters from them. \n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cluster import DBSCAN\n",
    "import sklearn.utils\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "sklearn.utils.check_random_state(1000)\n",
    "Clus_dataSet = pdf[['xm','ym']]\n",
    "Clus_dataSet = np.nan_to_num(Clus_dataSet)\n",
    "Clus_dataSet = StandardScaler().fit_transform(Clus_dataSet)\n",
    "\n",
    "# Compute DBSCAN\n",
    "db = DBSCAN(eps=0.15, min_samples=10).fit(Clus_dataSet)\n",
    "core_samples_mask = np.zeros_like(db.labels_, dtype=bool)\n",
    "core_samples_mask[db.core_sample_indices_] = True\n",
    "labels = db.labels_\n",
    "pdf[\"Clus_Db\"]=labels\n",
    "\n",
    "realClusterNum=len(set(labels)) - (1 if -1 in labels else 0)\n",
    "clusterNum = len(set(labels)) \n",
    "\n",
    "\n",
    "# A sample of clusters\n",
    "pdf[[\"Stn_Name\",\"Tx\",\"Tm\",\"Clus_Db\"]].head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see for outliers, the cluster label is -1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "set(labels)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 6- Visualization of clusters based on location\n",
    "<div id=\"visualize_cluster\">\n",
    "Now, we can visualize the clusters using basemap:\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from mpl_toolkits.basemap import Basemap\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import rcParams\n",
    "%matplotlib inline\n",
    "rcParams['figure.figsize'] = (14,10)\n",
    "\n",
    "my_map = Basemap(projection='merc',\n",
    "            resolution = 'l', area_thresh = 1000.0,\n",
    "            llcrnrlon=llon, llcrnrlat=llat, #min longitude (llcrnrlon) and latitude (llcrnrlat)\n",
    "            urcrnrlon=ulon, urcrnrlat=ulat) #max longitude (urcrnrlon) and latitude (urcrnrlat)\n",
    "\n",
    "my_map.drawcoastlines()\n",
    "my_map.drawcountries()\n",
    "#my_map.drawmapboundary()\n",
    "my_map.fillcontinents(color = 'white', alpha = 0.3)\n",
    "my_map.shadedrelief()\n",
    "\n",
    "# To create a color map\n",
    "colors = plt.get_cmap('jet')(np.linspace(0.0, 1.0, clusterNum))\n",
    "\n",
    "\n",
    "\n",
    "#Visualization1\n",
    "for clust_number in set(labels):\n",
    "    c=(([0.4,0.4,0.4]) if clust_number == -1 else colors[np.int(clust_number)])\n",
    "    clust_set = pdf[pdf.Clus_Db == clust_number]                    \n",
    "    my_map.scatter(clust_set.xm, clust_set.ym, color =c,  marker='o', s= 20, alpha = 0.85)\n",
    "    if clust_number != -1:\n",
    "        cenx=np.mean(clust_set.xm) \n",
    "        ceny=np.mean(clust_set.ym) \n",
    "        plt.text(cenx,ceny,str(clust_number), fontsize=25, color='red',)\n",
    "        print (\"Cluster \"+str(clust_number)+', Avg Temp: '+ str(np.mean(clust_set.Tm)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7- Clustering of stations based on their location, mean, max, and min Temperature\n",
    "<div id=\"clustering_location_mean_max_min_temperature\">\n",
    "In this section we re-run DBSCAN, but this time on a 5-dimensional dataset:\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "from sklearn.cluster import DBSCAN\n",
    "import sklearn.utils\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "sklearn.utils.check_random_state(1000)\n",
    "Clus_dataSet = pdf[['xm','ym','Tx','Tm','Tn']]\n",
    "Clus_dataSet = np.nan_to_num(Clus_dataSet)\n",
    "Clus_dataSet = StandardScaler().fit_transform(Clus_dataSet)\n",
    "\n",
    "# Compute DBSCAN\n",
    "db = DBSCAN(eps=0.3, min_samples=10).fit(Clus_dataSet)\n",
    "core_samples_mask = np.zeros_like(db.labels_, dtype=bool)\n",
    "core_samples_mask[db.core_sample_indices_] = True\n",
    "labels = db.labels_\n",
    "pdf[\"Clus_Db\"]=labels\n",
    "\n",
    "realClusterNum=len(set(labels)) - (1 if -1 in labels else 0)\n",
    "clusterNum = len(set(labels)) \n",
    "\n",
    "\n",
    "# A sample of clusters\n",
    "pdf[[\"Stn_Name\",\"Tx\",\"Tm\",\"Clus_Db\"]].head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 8- Visualization of clusters based on location and Temperature\n",
    "<div id=\"visualization_location_temperature\">\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from mpl_toolkits.basemap import Basemap\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import rcParams\n",
    "%matplotlib inline\n",
    "rcParams['figure.figsize'] = (14,10)\n",
    "\n",
    "my_map = Basemap(projection='merc',\n",
    "            resolution = 'l', area_thresh = 1000.0,\n",
    "            llcrnrlon=llon, llcrnrlat=llat, #min longitude (llcrnrlon) and latitude (llcrnrlat)\n",
    "            urcrnrlon=ulon, urcrnrlat=ulat) #max longitude (urcrnrlon) and latitude (urcrnrlat)\n",
    "\n",
    "my_map.drawcoastlines()\n",
    "my_map.drawcountries()\n",
    "#my_map.drawmapboundary()\n",
    "my_map.fillcontinents(color = 'white', alpha = 0.3)\n",
    "my_map.shadedrelief()\n",
    "\n",
    "# To create a color map\n",
    "colors = plt.get_cmap('jet')(np.linspace(0.0, 1.0, clusterNum))\n",
    "\n",
    "\n",
    "\n",
    "#Visualization1\n",
    "for clust_number in set(labels):\n",
    "    c=(([0.4,0.4,0.4]) if clust_number == -1 else colors[np.int(clust_number)])\n",
    "    clust_set = pdf[pdf.Clus_Db == clust_number]                    \n",
    "    my_map.scatter(clust_set.xm, clust_set.ym, color =c,  marker='o', s= 20, alpha = 0.85)\n",
    "    if clust_number != -1:\n",
    "        cenx=np.mean(clust_set.xm) \n",
    "        ceny=np.mean(clust_set.ym) \n",
    "        plt.text(cenx,ceny,str(clust_number), fontsize=25, color='red',)\n",
    "        print (\"Cluster \"+str(clust_number)+', Avg Temp: '+ str(np.mean(clust_set.Tm)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h2>Want to learn more?</h2>\n",
    "\n",
    "IBM SPSS Modeler is a comprehensive analytics platform that has many machine learning algorithms. It has been designed to bring predictive intelligence to decisions made by individuals, by groups, by systems – by your enterprise as a whole. A free trial is available through this course, available here: <a href=\"http://cocl.us/ML0101EN-SPSSModeler\">SPSS Modeler</a>\n",
    "\n",
    "Also, you can use Watson Studio to run these notebooks faster with bigger datasets. Watson Studio is IBM's leading cloud solution for data scientists, built by data scientists. With Jupyter notebooks, RStudio, Apache Spark and popular libraries pre-packaged in the cloud, Watson Studio enables data scientists to collaborate on their projects without having to install anything. Join the fast-growing community of Watson Studio users today with a free account at <a href=\"https://cocl.us/ML0101EN_DSX\">Watson Studio</a>\n",
    "\n",
    "<h3>Thanks for completing this lesson!</h3>\n",
    "\n",
    "<h4>Author:  <a href=\"https://ca.linkedin.com/in/saeedaghabozorgi\">Saeed Aghabozorgi</a></h4>\n",
    "<p><a href=\"https://ca.linkedin.com/in/saeedaghabozorgi\">Saeed Aghabozorgi</a>, PhD is a Data Scientist in IBM with a track record of developing enterprise level applications that substantially increases clients’ ability to turn data into actionable knowledge. He is a researcher in data mining field and expert in developing advanced analytic methods like machine learning and statistical modelling on large datasets.</p>\n",
    "\n",
    "<hr>\n",
    "\n",
    "<p>Copyright &copy; 2018 <a href=\"https://cocl.us/DX0108EN_CC\">Cognitive Class</a>. This notebook and its source code are released under the terms of the <a href=\"https://bigdatauniversity.com/mit-license/\">MIT License</a>.</p>"
   ]
  }
 ],
 "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.6.8"
  },
  "widgets": {
   "state": {},
   "version": "1.1.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}