diff --git a/README.md b/README.md index cd8724b..4e48216 100644 --- a/README.md +++ b/README.md @@ -1,41 +1 @@ -# Week 12 Deliverables (E8) - Due 4/12/16 -For this week make sure that you have completed the following: - -* Fork Assignment 9 to your own github repository. - * You can access assignment 9 - [HERE](https://github.com/Geospatial-Python/assignment_09) -* Clone the repository locally - -## Deliverables -1. Create a GUI application that can be launched from a script named `view.py`. - The GUI should include: - * A single window (QMainWindow or QWidget) - * A `File` menu with at least one entry - `Quit`. The quit item should - exit the program. - * The code should be presented using the following conventions: - a. Include an `if '__name__' == __main__:` block - b. Write your GUI in a class, inhereting from QMainWindow, QWidget, or -another appropriate parent. You should have at least an `__init__` method and -a `init_ui` method. - c. Call a function called `main` when the application is launched and -an instance of your GUI class is created - * A central widget. This can be an empty widget, a text box, a graphic - - we will replace this next week with something else, so the important part -for now is having something in the window. - - Note: For the above requirements, the linked readings offer step-by-step -examples of how to achieve this. -1. If you wrote a script for last week demonstrating your point pattern - analysis code, please move that code into an iPython notebook. Extend the -notebook to: - * Plot the point pattern. To do this, please write a function that takes - two vectors (x, y) as arguments and then generates a plot. - * Plot the results of the G-function. You can see an example of this in a - previous week's readings. In addition to plotting the observed function, -plot the upper and lower confidence thresholds in a different line color. For -example, if the observed value is red, then the upper and lower thresholds -(with p=0.05 maybe) could be blue. - - Hint: To get plots to automatically display in a notebook add `%pylab -inline` to the first cell (where your other imports are) -1. Update any other support code as necessary. +# assignment_07 \ No newline at end of file diff --git a/__init__.py b/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/analytics.py b/analytics.py new file mode 100755 index 0000000..ed40fb9 --- /dev/null +++ b/analytics.py @@ -0,0 +1,176 @@ +import point +import math +import statistics +import random + +def mean_center(points): + """ + Given a set of points, compute the mean center + + Parameters + ---------- + points : list + A list of points in the form (x,y) + + Returns + ------- + x : float + Mean x coordinate + + y : float + Mean y coordinate + """ + total = len(points) + y = 0 + x = 0 + for point in points: + x += point[0] + y += point[1] + + x = x/total + y = y/total + return x, y + +def minimum_bounding_rectangle(points): + """ + Given a set of points, compute the minimum bounding rectangle. + + Parameters + ---------- + points : list + A list of points in the form (x,y) + + Returns + ------- + : list + Corners of the MBR in the form [xmin, ymin, xmax, ymax] + """ + x_list = [] + y_list = [] + + for p in points: + x_list.append(p[0]) + y_list.append(p[1]) + + mbr = [0,0,0,0] + mbr[0] = min(x_list) + mbr[1] = min(y_list) + mbr[2] = max(x_list) + mbr[3] = max(y_list) + + return mbr + + +def mbr_area(mbr): + """ + Compute the area of a minimum bounding rectangle + """ + width = mbr[2] - mbr[0] + length = mbr[3] - mbr[1] + area = width * length + + return area + + +def expected_distance(area, n): + """ + Compute the expected mean distance given + some study area. + + This makes lots of assumptions and is not + necessarily how you would want to compute + this. This is just an example of the full + analysis pipe, e.g. compute the mean distance + and the expected mean distance. + + Parameters + ---------- + area : float + The area of the study area + + n : int + The number of points + """ + + expected = 0.5 * (math.sqrt( area / n )) + return expected + + +""" +Below are the functions that you created last week. +Your syntax might have been different (which is awesome), +but the functionality is identical. No need to touch +these unless you are interested in another way of solving +the assignment +""" + +def manhattan_distance(a, b): + """ + Compute the Manhattan distance between two points + + Parameters + ---------- + a : tuple + A point in the form (x,y) + + b : tuple + A point in the form (x,y) + + Returns + ------- + distance : float + The Manhattan distance between the two points + """ + distance = abs(a[0] - b[0]) + abs(a[1] - b[1]) + return distance +def create_random_marked_points(n, marks = None): + point_list = [] + rand = random.Random() + for i in range(n): + rand_x = round(rand.uniform(0,1),2) + rand_y = round(rand.uniform(0,1),2) + if marks is None: + point_list.append(point.Point(rand_x, rand_y)) + else: + rand_mark = random.choice(marks) + point_list.append(point.Point(rand_x, rand_y, rand_mark)) + return point_list + +def euclidean_distance(a, b): + distance = math.sqrt((a.x - b.x)**2 + (a.y - b.y)**2) + return distance + +def average_nearest_neighbor_distance(points, mark = None): + new_points = [] + if mark is None: + new_points = points + else: + for point in points: + if point.mark is mark: + new_points.append(point) + + dists = [] + for num1, point in enumerate(new_points): + dists.append(None) + for num2, point2 in enumerate(new_points): + if num1 is not num2: + new_dist = euclidean_distance(point, point2) + if dists[num1] == None: + dists[num1] = new_dist + elif dists[num1] > new_dist: + dists[num1] = new_dist + + return sum(dists)/len(points) + +def permutations(p=99, n=100, marks=None): + neighbor_perms = [] + for i in range(p): + neighbor_perms.append(average_nearest_neighbor_distance(create_random_marked_points(n), + marks)) + return neighbor_perms + +def compute_critical(perms): + return max(perms), min(perms) + +def check_significant(lower, upper, observed): + return(lower <= observed or observed <= upper) diff --git a/io_geojson.py b/io_geojson.py new file mode 100644 index 0000000..c62228e --- /dev/null +++ b/io_geojson.py @@ -0,0 +1,78 @@ +import json + +def read_geojson(input_file): + """ + Read a geojson file + + Parameters + ---------- + input_file : str + The PATH to the data to be read + + Returns + ------- + gj : dict + An in memory version of the geojson + """ + # Please use the python json module (imported above) + # to solve this one. + + with open(input_file, 'r') as f: + gj = json.load(f) + return gj + + +def find_largest_city(gj): + """ + Iterate through a geojson feature collection and + find the largest city. Assume that the key + to access the maximum population is 'pop_max'. + + Parameters + ---------- + gj : dict + A GeoJSON file read in as a Python dictionary + + Returns + ------- + city : str + The largest city + + population : int + The population of the largest city + """ + max_population = 0 + city = None + features_list = gj.get('features') + x = 0 + + for f in features_list: + if f['properties']['pop_max'] > max_population: + max_population = f['properties']['pop_max'] + city = f['properties']['name'] + + return city, max_population + + +def write_your_own(gj): + """ + Here you will write your own code to find + some attribute in the supplied geojson file. + + Take a look at the attributes available and pick + something interesting that you might like to find + or summarize. This is totally up to you. + + Do not forget to write the accompanying test in + tests.py! + """ + #find the largest city west of the Mississippi River + + largest_western_city = None + features_list = gj.get('features') + for f in features_list: + if f['properties']['longitude'] < -95.202: + largest_western_city = f['properties']['longitude'] + + + return largest_western_city diff --git a/plotbook.ipynb b/plotbook.ipynb new file mode 100644 index 0000000..76a2c55 --- /dev/null +++ b/plotbook.ipynb @@ -0,0 +1,259 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 1, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "%matplotlib inline\n", + "from point import Point, PointPattern\n", + "import pysal as ps\n", + "from pylab import *\n", + "import matplotlib\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "\n", + "shapefile = ps.open(ps.examples.get_path('new_haven_merged.shp'))\n", + "shapefile\n" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "DataTable: /home/qstin/anaconda3/lib/python3.5/site-packages/pysal/examples/newHaven/new_haven_merged.dbf" + ] + }, + "execution_count": 2, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "dbf = ps.open(ps.examples.get_path('new_haven_merged.dbf'))\n", + "dbf" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "zipped_pp = PointPattern()\n", + "\n", + "for geometry, attributes in zip(shapefile, dbf):\n", + " zipped_pp.add_point(Point(geometry[0], geometry[1], attributes[1]))\n" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "This average nearest neighbor for animal bites is: 0.01261354238016035\n" + ] + } + ], + "source": [ + "nn = zipped_pp.average_nearest_neighor_KDTree('animal-bites')\n", + "print(\"This average nearest neighbor for animal bites is: {}\".format(nn))" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "crit_max, crit_min = PointPattern.get_critical(zipped_pp.generate_realizations())\n", + "def is_sig(n):\n", + " if n < crit_min or n > crit_max:\n", + " return \"is significant!\"\n", + " else:\n", + " return \"is not significant.\"" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "The nearest neighbor of animal bites is significant!\n" + ] + } + ], + "source": [ + "print(\"The nearest neighbor of animal bites {}\".format(is_sig(nn)))" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "The neareset neighbor of all data is significant!\n" + ] + } + ], + "source": [ + "nn_all = zipped_pp.average_nearest_neighor_KDTree()\n", + "print(\"The neareset neighbor of all data {}\".format(is_sig(nn_all)))" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "'Begin assignment 9'" + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "\"\"\"Begin assignment 9\"\"\"" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "new_array = []\n", + "for item in zipped_pp.points:\n", + " new_array.append([item.x, item.y])\n", + "nparray = np.array(new_array)" + ] + }, + { + "cell_type": "code", + "execution_count": 26, + "metadata": { + "collapsed": false + }, + "outputs": [], + "source": [ + "def plot_points(x, y, **kwargs):\n", + " plot(x, y, **kwargs)\n" + ] + }, + { + "cell_type": "code", + "execution_count": 32, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAEACAYAAACd2SCPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXecFeW9/99zOOyeXX8hNxrbjawo7iIirCtiiQIrlkiR\nuhQpFroBxRa9SRQREEGvLdfEEgFhC9swSkyiuQliieK6sDSjLGIBc6M3JrnGyNZzvr8/vs/szJkz\nZxttifN5vc6LZeozzzzz7cUSEQIECBAgQIC2InS4BxAgQIAAAY5MBAwkQIAAAQK0CwEDCRAgQIAA\n7ULAQAIECBAgQLsQMJAAAQIECNAuBAwkQIAAAQK0C61mIJZlhSzLqrIsa535f55lWTssy4palnV2\nknNSLct6y5y33bKsuz37b7As612zb+n+PUqAAAECBDiUCLfh2HnAO0AX8//twCjgyWQniEidZVkX\ni8g+y7I6AX+wLOs3IlJhWdbFwJVAbxFptCzr2+18hgABAgQIcBjQKg3EsqyTgCHA0/Y2EdkpIrsA\nq7lzRWSf+TMVZVh25uJsYKmINJrjPm/b0AMECBAgwOFEa01YDwM/wCH+rYZt+gI+Bf5bRN42u7KA\nAZZlbbQs62XLss5p67UDBAgQIMDhQ4sMxLKsocBnIrIF1Taa1Ti8EJGYiOQAJwHnWZZ1htkVBr4l\nIucDtwOlbRp5gAABAgQ4rGiND+RCYLhlWUOANOAblmWtFpGr23IjEfmHZVkvA1cAfwQ+AZ41+962\nLCtmWdYxIvJX93mWZQXFugIECBCgHRCRNgn8bUWLGoiI/EhEMkTkVGACsN6HefgO0rKsb1uW9U3z\ndxpwGfCe2f0cMMjsywI6e5mHawxH7O/uu+8+7GP4uo7/SB57MP7D/zvSx38o0O48EMuyRlqWtRc4\nH3jBsqzfmO0nWpb1gjnsROBly7K2AG8BL4nIr82+FcCplmVtB4qANmk0AQIECBDg8KItYbyIyCvA\nK+bv51AtwnvMn4Fh5u/tgG+OiIg0AFPaON4AAQIECNBBEGSiH2Tk5uYe7iHsF47k8R/JY4dg/Icb\nR/r4DwWsQ2Uray8sy5KOPsYAAQIE6GiwLAs53E70AAECBAgQwA8BAwkQIECAAO1CwEACBAgQIEC7\nEDCQAAECBAjQLgQMJECAAAECtAsBAwkQIECAAO1CwEACBAgQIEC7EDCQAAECBAjQLgQMJECAAAEC\ntAsBAwkQIECAAO1CwEACBAgQIEC7EDCQAAECBAjQLrSpnHuAry/q6+tZtmwZAHfccQcpKSmHeUQB\nAgQ43Ag0kAAtYsmSnxKJjGH+/B7Mn9+DSGQMS5b89HAPK0CAAIcZQTn3AM2ivr6eSGQMIs/jyBsx\nLGsEtbVrA00kQIAOiqCce4DDjmXLliEyhfilEkJkcpNJK0CAAF9PBAwkQIAAAQK0C4EJ62uOWCzG\npk2b2LlzJ5mZmQDs2rWLHj160LdvXxobGwMTVoAARyAOhQkriML6GqOq6h3Gj1/G++/XIzICeB94\nETgby/otmZlQXHwHixdfwZ13jkBkMgCWlc/ixYMD5hEgwNccgQbyNUVjYyMnnzyM//mfz4AngH6o\nhhEDbgIeAm4mOxs2b36UxsbGIIw3QIAjCIdCAwkYyNcQVVXvcPHFd/LFF3lAZ+AFs+cOoBewFugG\nfEQkspvXX7+Evn37Hp7BBggQoF0IorACHHDEYjGuueZxvvhiLTAJGAc8A3wT1URi5shG4BUaGzfT\n2Nh4WMYaIECAjo2AgXzNsGZNOdu3X4A3LBdyga7AJqAAWAZcRGPjaHJzl7FmzbpDPtYAAQJ0bAQm\nrK8RYrEYPXtOpLp6JDDBs3ctUA1sAFKB53BHXaWl5fGPf5QSDgdxFwECHAkITFgBDiiqqqrYuzcP\n+AOOqQrz9wbgFdT3MRmvhlJTM57i4uJDNNIAAQIcCQjEya8ZLCsEzEb9H8PQJfBL4BMgD/i3wze4\nAAECHFFotQZiWVbIsqwqy7LWmf/nWZa1w7KsqGVZZyc5J9WyrLfMedsty7rb55hbLcuKWZZ1dPsf\nI0BrkJOTQ1bWBqAnUAh0B9ajpqsXgbnAOCwrH6+GkpZWwoQJXrNXgAABvs5oiwYyD3gH6GL+vx0Y\nBTyZ7AQRqbMs62IR2WdZVifgD5Zl/UZEKgAsyzoJuAz4uF2j/xdDLBajqqqK+vp6fvvb3xIKhQ5o\nzkUoFGLFillMnXoT1dUDATjppBhwOp988isAMjM3MG7cFSxenEdNzXgAIpFili+/LvB/BAgQIA6t\ncqIbQr8SuBe4RUSGu/a9DNwqIptbuEY68CpwvYi8bbaVAQuBdUBfEfmbz3lfCyd6VdU7TJ36JNu3\nn0w0+nvgQuAULGsNixdfwY9+NKeJwYBqE6FQ+1xY3uvo/eOv29jY2OTzmDBhQsA8AgQ4wtBhEgkN\nob8XTRa4tS0MxFKj+ybUXvJTEfmh2T4cyBWRWyzL+pCvGQNxE+hx48Zx3nm3sWXLNOAe4CrgI3Ta\nfohl3cnGjYuYNWsF1dW5AGRlbWDFilnk5PQ6rGMPmEuAAB0THYKBWJY1FBgsInMty8pFmcWVrv2t\n1UC6oLGhc4EPgZeBy0TkS8NAzhGRv/qcd0QykNraWm6++WYAHn74YSKRSNO+NWvWMW3aM9TUqE8h\nNfVxYrFZNDT8DrUUPg7sBgSIAP054YS3+fTTYtyhtWeddRObNj3Sbk2kLbC1lhdffIXFi1+jtvYq\nANLSilm+/Fquump4C1cIECDAoURHYSBL0LjORiAN+AbwrIhcbfa3ioGYY+8CvgJ+C/wO2AdYwEnA\nn4BzReR/PefI3Xc7vvfc3Fxyc3Nb+XiHB3Pn3sVPf1qJ4y76B3PmnMNjjy2isbGRLl3GUVPzDOpC\nArgTeAM4DfipOe9qs2818D9Y1s2mL4eDzp2LefPNzINeZsQ2r1VXD2TfvjrgTTSSqxdBjkiAAB0D\nGzZsYMOGDU3/v+eeew4/A4k72LIG4m/Cuk1ENvkc/22gQUS+sCwrDXgJWCoiv/Yc9yFwtoj83eca\nR5QGUltbS1paf+AEwCb4+cCn1NS8Rnl5OVOmLAeOIp5J7EJDaytQl5CjacBw4GzUXeRGIW++2Z3z\nzz//oD1PLBajb9+b2LLlEc+YbgLsbSXk5zcwefLkgzaOAAECtA0dOpHQsqyRlmXtBc4HXrAs6zdm\n+4mWZdnV+U4EXrYsawvwFvCSl3kYCKqJHPG48cYbgeOB59E6U+PM38dz44038vnnn6PMY51r/zog\nE/gJylS8ZUamoDEM3uS/F3nppZcO6vNUVVUZv4t3TAOBqoNyT7tHyaZNm4jFYi2fECBAgMOCNtkc\nROQVNF0ZEXkO9Wl4j/kzmqGGiGxHReeWrntqW8bRkaHRTD/AjwlUVf0nZWVlaORz4v5OnV4kGk12\n5ZNQqX+g+f8GIIdQ6KsDOPr2wM4RKT0gV3PMZbkAZGWtOmzBAgECBGgeQSmTA4zBgwe3ax/AUUcd\nhZqzvJpGPqq1PISWGukGPAysp7HxqP0Zbotwkg+9Y3oOqCYSGcPy5dceEP9HLBZj6tQn2bLlEfbt\nG82+faPZsuURpk598oBoIo2NjRQUFFBQUBBXYXjfvn3k5eWRl5fHvn379vs+AQJ8XRAUUzzAqK+v\nJzV1NH5+jLq6Z/nb3/7GiSdO993/y1/O5sorV6OxBbb/ZDXKPOrM/6eg1r4NwEzgh2zcuIjzzss5\naM/kdqIDnHbay4wdm0G3bv9+QMN4N23axIABH7Nv3+i47enpa3n11W77FSzgjXyLRIq5887+bNhQ\nwe9+9yVuf9S4cRmUlPys3fcKEKAjIGhpewQiJSWFe+8dnLQF7N69n6FRzMOJZxIf8fbbb6P1qIYD\nN5t9a9GGT4+iEdCnmO2PYpu+Ro68mT/9af1BC+fNyenFpk2PuJINf3JIQocPFBobG5k69Rlqa8ux\nmXZtbR533jkF2Itqc33NvjxKS4ezcuU+0tPTgfjEy969e1Naqua6IAcmQDIcqKTfDg8R6dA/HeKR\nh7q6Olm4cKEsXLhQampqpLKyUioqKqRr1xECxQJfCYwxv68ECmXmzJkCIwWiAmJ+dQKXCcwU+L5n\nnwgUSzg8VyorKyUajUplZWXT30caotGonHXWDZ5n1G378zyLFi0VKPTMmwgUCNwnsFbgBoEdTXM6\naNAgM5fb5KyzbpD09LXSufPDYlnDzLWKJCXlSlm8+KEjcq4DHDxs3ryjac2kp6+Vs866QTZv3iEi\nckC/UTeNqaurS9hvaOfBpc8H+wb7PcAjlIHYcC+mSOR+gQWGcCUSsyeeeELgJIEpAuUCPxIYJlBi\nflcKjDB/Vwo0CAyXSGSNFBSUJl20RxKc+SqX9PRyyc6eu1/PEY1GJSNjhJlP75yXCVSY3yqBcWZO\niyUUOl/S0sokLW2UYWgNAlPNsZVmW1RgjPTu/X0pKCg5Yhl3gAOHZEJQRsYkueee+yQ7e+4B+Ubv\nvfcxgSECc81viNx772NxxxwKBhL4QA4S7FDUyZOfpLr6KdQ8sgm4H23Y9AzxPpBrGDHiSyCV55/v\nDMxBuwI+6zluBDAWzVAvBuro3bsboVCIrVsfdR3bSFbWZAoKbqVv374dUoVOpuYfSPX/7bffpn//\nD6ir+wNO3groXF6LvotBQCfUVFgL/AX4FfAuWlLmdLRf/NHASHP+BmAWsAVYh2WNJC0thaysV3jw\nwXFMnjwOgG3btnH00Ud/PcwZAZL68WANGniyhmTVJFq77tXPOgg4DqcxXDHwv9TVrW8qvnoofCCH\nXcNo6ccRqIFs3rxDsrPnSkrKjQJr4iQROE/gIWMyKTe/uQIPSrdu3eXpp58WGChwqdE0vFJzscBC\nlwYyVmbO/L6kp691HbPDXPM+SUlZut9S/IGAV3VvTs0/UNi8eYdkZY038/i8wCgzfwUCkwSmJ0iK\nOm+XC7wukC9QarZN8zn2BoHVRiuxt28TGGrumS9wqkBPiUSKj2jNMEDrUFlZab7FqPlG7e/0eg8t\n0DUUiSyT/Pz8OFOp3zpxfz933XWXwHCf9Thc5s+f33QOgQZy5Gkg2jb2WqqrvwlkACejyYI2qoC7\nUGlkq9mWjUq2Y9C+HGnmuB97zgUoQcueDAR+D7xPOLyZxsbJwIPmmGvRupcXm/+/TFbWF7z77jOH\nRfr15nZkZr7Mvn1fsGvXDag0phFkB7K2l5NB/xCaP2OhzvJitHlWFRqwMNZzpt3ad68Z2y+BK9G5\nHOM5thTVamajiaDvAnaB6bfQwgt2oEQx+l6GHdIaZgEOLWKxGKeffi27dtnfXwwoQunABTjf8zto\nPthFRCIWlrWGmpp7gN72lZrWydat78Z9P507384XXyzBjzaMGVNGeXk5EGggR6QGUlFRIZY1xWUj\nT7SHwtlGSi02v6ECj7n2TxGYJTDY59xRRqJ5TOL9I0ONtLvanB9/nmVNkYqKikM+H8lswipBlYnb\ngZ2eXi6VlZUH5L6OJLhD1FZsO9ErzN8l4u9YLxe4ycy9fdwogfs94y8RuMLMd5HAxQKjzfZS8z78\n392BfM4AHQfRaFQqKirkhBOGmG/Uq9nObYYuNAiMN+tTt6enl0tFRYXP97NcklknVq5c2TQeDoEG\nEsQgHmDs3LkTkWE4ds5ZqAR8IdBIZubLDBo0kCef7I3a0f8P9XPYTaNCqMRbDpyJaiaTzL5CYDp2\nGRMtkWLfxw7//RlwI95Md5Gh7Ny5k379+h3YB24GsViMoqIi3nuvf8J4VDI/BQ2fHQnchEj/A3p/\nkRgq5S0A9tijQrW4tWio9ATi/SIvo9rHC8TPbR5aT3Qn8ATwXfRdvALMQAtL2j6oTcB1JD7zeFQT\nSUsY68Hq9RJoOYcGtpb93nv9qa29BrgF/fZ7oe8+F4iitOBE4CKc9WFrI6NRn1u+OVfpSWIpoavR\nbyYP99oNhwuZPPnZg/OASRAwkHYi2Yfao0cPLGsXjtWtF2rmWAJY7N37V5555jSUUN2AEitvx8EG\n4DagH0q0ClEieC/KJBahBNiPKK8nPmtcYVlRevToEbftQPT1SDYPmzZtZ+LEe/noo+Opr7/Q50y3\nZh0CBnDSSWXk5Kxp8xj8kJOTQ9euD1BdPRZlUvlo9eNdqPnvFrT6sbs3/AvA34BrSJzbq8yx3gCI\nUWhC5xWec5JZDmJkZb1CTs6opi1VVe9w3XVPsHNnBgA9ejzD8uUzgfqmZ2kNI2hPGZiA4ew/3BUU\n4oWOeTi17aJoS6Q8OnW6kU6dulNfD/qtPkl8gMcYYB6ZmdCjx9UojXAjTDicTSg0mvp6basQiRSz\nYsX0Q56XFKyWdqCq6h369r2JAQM+ZsCAj+nb9yaqqt4B9CPs3LmMRCJeCfyZ2tq11NVdBpyDJhA+\nhxKKAvOrR7ULO+s6jBK0JSgR9JYU2WR+9vZ/N+fHH5eZuSEuk3vNmnV06TKOKVNSmDIlhS5dxrFm\nzbqkz1xfX8+iRYtYtGgR9fX1xGIxCgvL6NlzJgMGfBQ3D5s2bad//7uprh5LfX1/VOr2jvsVbN+H\noo66un+ydeu7vvdva4HFUCjE/Pmjsax6dJnbmuDHaLHLR1Ap8BaUSb+Oam4L8Cf+IWAoqh16mcvl\naHKojRxUQHgb592oLfzUU9fz9NMz4iLOJkxYxtatUFvbndra7mzdCt/97m307/9hwvpKhvaUgWlu\nHQdoPZIXHL0AtTJ8hH67ygjOPNOiZ8/X0DVRhWon8eda1gXk5X0HUJ+h9/s588wv+fLLUvLzG8jP\nb+DLL8sOT0+eg20j298fHcwH0lKyW2VlpaSmPiqJUVYPuezotg20QeCH4kTt2L6M2b72TbjF2NFX\nG5v7XONDWGv+HiSaePignHDCWIlESiQSKZHs7DlxER0NDQ2u/AbnGdLSRklDQ0PCM99772Mmgc4Z\n47HHfk8gN85mC1Hp02eOSZZ0b98mTgRUkWg+yzYfG3GDb9JgeyO2otGoZGfPdY2jQWC+aASW+9kr\nzLjcY/GLzqowc+19N2WiOST2OTtEI7wKzHxNNnN1taSmFsdFxcX7zNz3myJOdFfLyZSOzyd+bMn8\nLQcrafNfHdFoVDZu3CgLFiyQBQsWyBtvvCGrVq2SSKQ0ybqodL3T4dKr1yzZvHlH05qORJaJnz/D\nsgolErlf0tPXSlbWDMnKmtLm3CiCRELpcAyksrJS0tLsheEklNnheBUVFZKeXi7xYXy203WZa3uJ\naGa5XzjeUMMI3NuubCKySljG+Jw3wewfKj/+8SJZsGBBU5ZqQ0OD5OfnS35+vqxcudKzaO0xLZFV\nq1bFPW9dXZ1hHn5jtAnkXEM0dxhCWiKJ2d12OGy+wC9EM+6XCCwWDVt+UCCaQPD2l9DZH2pq6iOG\nUC+W+DDqUtHkzPE+DCBflOlNN0xgqahD3usgnWTez2Dz7H4MwR0GHJWsrBkSjUYlPz/fl4DoffOb\n/t+S472iosIkqlbG3ds+zxtG7R9umjj/ARTRaFQKCkrk+OMvNu/aFqaGC1wqljXB5517122xnHzy\nyU3voaGhQSoqKiQryy9E3C3EqCBUUVHRpmTVQ8FAAh9IG/Heex9QW/tLnISyR4EItbWXMmOGRUbG\nUzQ27kFt47bJKIZlvUAo9Gei0b1oeF8YDRedQaIl8UI0tHc+EKJz50JisWOJRu3zQsBEElXmkag6\n3Jt7792IXSBw0aKxiPwfjY1zAAiHi3BCfG0HXi7QnbvuWkfv3n3p3bsHxcXFPPfcc2gnRO+9rsVx\ngufhhA67E6XUOa7mIgtN0JsATEXL03dHbcMZqIo+j2j0NLTasCKZeaC6eiBVVVUtFljMyenF228/\nRK9es6mufsbcZzLqQI+hjTFHAMvNeMejfo4U1Nz1OrAUeNo870mojXog2jTsN6hTNBPogSYiev0h\nITRZcRPq1wqxa1cumzZtMj6z91FZyY2ouV5y2D6sjz76E2Vle6mruxA1l6xCTXY9ycp6Bciib9+b\n4nwjt93Wn2j0z+j7yTVXXEU02h33/AdQf96ECY/w/vsZ6HspB7abvWuBcYiEgXlY1gWINAL/DfwH\n8evgT/z5z90ZMOBjwPFRFRffzNSpN1FdPZBYLEpd3S8RcZ+r6yUUCh307qNtxsHmUPv7owNpIIkm\nkWTmjulme7lAmaSkjJCVKwslJcWrbazySJ87jNRSYqT7K2X06AlSU1MjGRmXGWlZjLToZ0opErjH\nSNTeEMFcgY3ihBGOEtVyEqX7rl3HSyQyyozjXkke7lrp+v8i1/i8x1UYSW2b+XuyJJYEUZOe14yW\nzDQTiZRIfn6+NDQ0tFhbKD6k9wZR00KhqBb1kGiC14UC/cVJ0rTnZIv4h+ROM/M51zz3QnGSB1vW\nKKBAFi1aKg0NDZKR4aexuJMcEzWuoqLnjRnSDke2tUA93rKmSJ8+c5oS1LzX79NnjnnHrTNjfl1R\nWblNUlIuN/NbYtbKKIFHxdGyHzT/Vkgkskwee8wuM+L9BhPXkf1eba0wPz/fWDji1097NEMCE5Z0\nKAaihMhdUykZIbeJZoVAvqSm3icLFizwISwN4piw/HNGMjImS58+cyQ11bal28fO8CE6kwXOkORM\nqUQcs1KpaJHGEtf5lWbME+OIlz+T9I51iYdAugnnpaIM7L9EzXhTxPHd2ONRhhSJlCaYsJRp2wyn\nTpTxXijh8AOSljZKIpGyZn0jjtnRL/Z+qCgzKBclxJPMeOz5sBmD3zuuFGVG4824RokylZZ9GjBX\nunadLNnZcyU19VGxrCliWQWSklIsmZnTmrV5J/NhuZ8vEilpMnn4M+BlkpKS+Fze+T9YOBIKf0aj\nUcnMHJfkfdpzbX93cwUqmwj9qadeJPG5XheKXw28lJTiA2qydeNQMJDAhNUGvPfebmpqkrYM9OBD\n1PyRS13dqTzxxEq0ppIbYdRMMQrtDOyODQcIsWfPMPbseRYNB+6Fhgaeh5qD5uGYH15GM6wHuM73\nCxHMQ80W/VFTjkW8GavRHPsuTgz7bNQcNtT8Px+4z3XNmLnGDuJNazHgF+YZZ6Eq/VFoVJbXzKU5\nIN4w0q1b36Wmps6YebYBd5ux3EhjYz6NjdOwzYlbtoxk6tTELO/4kF537P0TaH2rp1zbJ5ixRlAz\n319o3pRUj8bvp5hxLQQ64/9uXkXzUTYAs/nkk63s3dsNrSjQHTiFbt1+wjvvFBAKhVzhtY/GPU9x\ncbHpa5KszXBfQqFOLYbkhkKJ0WaHIozXDjfeuXMAIh+SkfEARUU/pm/f3i2ffAhgmwY//PBDPvro\nOHRteud6APZcq8nyeSCbrKxV5OSMYvfu1/jb3/7WZHJ66qmnuPzyvyTcq76+Ma65WSgUYsWKWU0m\nLYDMzA2sWDG7Y4ZYH2wOtb8/OogG4kjCrTFheY9TadfPGR0OXyHwPXGcz8m0GbfEY6vTXkd9mZGA\nR7n2NRc11CAwJ8kzeKWgpwT6imoSdgSYHWU2x2yzo8EKzD3miGo6FWb7YtEIMr/xjBdokMzMqbJx\n48YmJ6MjjTW4nss9TjszX7clU/ULCkrEstyRVjdI8qiqAnG0hWT3taPo7GCCYtFqvQ8I5Jh3ulRU\n45ojar6bJvHRaYUCF4mjHQ6TUGhWixpAcse7rRXFm0b8JNrs7Lkec2z7Jd22wBnPNjOHawXKJS1t\nlFRWbjto920t4k2Da0SrDfiZcAtd76BAOne+tdnoqFWrVom/JjM5IXBF5MBoaAQaSMdBVVUV772X\nieZzXItK41HgQzIyruXzz0cAcNJJL1JT82f27rUTiGyECYcvwbJGUl+vmeWRSDHXXXcGjz/eF5VC\nX8ObXar5Eo+gkqst8Rxv9oVwHPW4tl1rrtMPlWy9qCMcfofGxhLg7ziahQ2vhHUnWnX2drN/NSpN\n/7v5/yfAjVjWs3TqlEVj4wc4Urug2tj/AX8jFOpGYlpCnXmmPN5/P40LLviA1NTOZGQ8zscf9zHj\nKSI+a9wep53dPRGoorFxF42NJyU88VVX5fHAA/PYunU8/rH3boRd+8LofI4y948Bv0ZzPa5DNZD/\nRB3qz6MJoHegmpytcebg5IvsROfUrpH0Mk4+bx6x2PCmPBtbA8nOzmbr1q1Nf2dmZtK58xIaGrxr\n5TlgCJmZP2PFilsAXbe33dafBx6Yx65duYBKtCtXzgY45JLupk2bePfdC4Gf49aMa2pGMXHiTN59\n96nDJmnX19dz7bUPU19/P+4GY/rzVix40/w9mszMlyksnNWKqtffRtfpGHOtDYB/8maHdJj74WBz\nqP390UE0kMQaV7a/YHKTxGxLCxpSmehQTk8vl40bNzaF0zY0NLjs2duMhDLeSLN2/ojtGHWkS5Vo\nr/CRZqa4zmkQuMpIwt7jxsry5cvlD3/4gwwaNCiJNFtopOcC8XciDxW4U9QfoRpAJFIis2fPNpLb\n1T7nXCqpqd6GWVFJSblS4A3x04ScOW/OOX2fOGG5hRKJjJKioucT3qETe79UHA3Or07XZM82EdWc\nekqXLueasboDAGaY9zdZNNfHDudNNtbJoprHLeINu4VimT17TlPeS2pqiaSm5krnzrdJSsrDTT6f\nTp3s/CHbxj5S4CGJRJZJRUVFXO5MWlqZZGaOk0WL7pOKioo4ifZQ+iKc6sh2E6/4+YlESg5bCPHm\nzTvkhBP8QtCj5juaLPF5XTsEiqVr1+Et5mS89VaVHH98nvmWykU10VLxCxo5kOAQaCCHnUG0OMAO\nxUCKEha9ZRUmFClMZjbIyprm+wEvXvygpKbmmgW2UOAC8SboKXOwcym2SSh0nliWbT4pcy1qmyja\nTl3bVFBujhsqsEzS08vlrLNukIqKLUkcsnaUka3O+xHDsaKMTO9nM8hOnfpL8g6AiyQtbZSkp5dJ\nenq5ZGZONYmX/uY2yyo0c5HMlDTUfOCtiyayE8GOOeYyc117fsrM+AZLYlfIqMAk6dRpbJKEMdvM\nmCcqAHw3yTXcZsjhoow2sRviCSecZ46xAyDKRSPs7PdpX7NONEDBZuK67uJNf841LKuoycxiM46X\nX35ZTj9hpYbSAAAgAElEQVT9dOnZs6f8/e9/P2jfj/NN2EUDExt8WVahFBSUHNQx+DHLxOhK9zen\nbRH0XbnNxSJQJAsWLGiW8S5c+Kj4C2CTJTV1hK+gc6CeK2AgHYiBVFZWSiSSGF6XLGrFlgDT0kol\nJWWppKTkSkrKw3HRQvEZ1uWSkTFJFi1aJqtXrzVE3ZYuh4n6GJaJzURSUx+RFStWSDh8n/hJscqE\nCl0L1l78djdDaSI4BQW/MCGdRYaIThQNb50scKL4M5AHRRmIHb00VI4//jKJRqPSr993k5xTLNBX\nUlOLTA+ESpN4uVaSMZBIpFSyssZLenq5hMM/Mn4knZeUlOEyatRkUd9K4hzk5+cnfS/6Lu2oqwfF\n8UE1mLmeZO5TKMpML5XTThvricJz7nXUUWdLJNJN1Fe0ROBW0YSzInMdtzYpkpil7PhUUlPXSDzD\n8WNCDtGNRJbFRWvFJwn6CTJTJDt7rlhWgXnfU8z7vkJOPfW8g6KJxEeD2ZUJEn2H2dlzD8r9m6tm\noN+2VzCwGUjyCEkYJpFISdLov40bN4tqmn5rpkDefPPNg/pcAQPpQAykPeF1lZXbJDNzmpGiy8Wt\nFidzYjrbGyS+IY07lFUzhvPz831DNMPhQpk/f34SYhefv2E7nWtqaiQcvkgSNZ9ZkihB+ce0w1Cp\nq6uT+fPnSzKHoTbWuVIWL37YM68Nvh9pWtooqajY0iRh1dXVSX5+vixadJ+cdtokScx83ybJsuqT\nl5afJirN2/NdJGqWWiWqYX1PunYdmTSnQmP+15hxDBH4D9c82cEDXpOYN4+mWOC7Mnv2Hc0y1MT3\nV9ZUceCNN96QiooKVy6B3zWiSUqn2Gtukpx55vUHvOlVYjhxsTitmx2z0MHIhG9N+aFEBmILW/b/\nbU3OLsfzPdfcJtKBaDQqxx13riQz10GBr4BzIJ8rYCAdiIGItK1fd3JipdsikWW+5pBk270Sa6Kp\nQlq1z0uk7Q82eWRPoahZx21vv8j86z22uKl0in+k1mWiZrEK6dRpiNTV1cXNa2rqw6JMplDc5jq/\njzMxIs5+vhFiaw7u9/Pll1/KOeecI6GQXz/6UlFGdL/5TTLPuEOgSEaOHNl0//g1UCbKMLx1vexS\nNMkj8vR+bmadL089tdK1bpqru+Wsg5SUyw1D+E+ByWJZRRKJlLk0WO81vITR/Z5LxDbHHehorMTv\nYaM4uUPxPTD2l4EkL9sS/8zuMi+JwpzdN8Z9TlQ0sq63qMCQb36JPV7uuus+gZNENTy/7/DK/e7P\n09JzHQoG0gEDizsucnJ6sWnTI7z6ajdefbUbmzc/mrRUdvIKnQPRKKC2wbLqiUR+T3r6WrKz57Fi\nxSzC4TArVszirLNuIi2tjEjkfjIzr+Lpp2fE7UtPX0t6ejkpKaOJL50S46STXiInJ6eZO6egkSNn\nc/rpCzj66GVor/bOSc8Ih8N8//vnAx+g5Vrs39FoNNkeotGjycgYSlXVO+Tk9OIPf1jKhRf+Esva\nhZY2ORUtE9O7qWyJjaqqKlP6/GIS53cyWg6/B1u3PszUqU9y2WUT+cY3JlBZ2Z9YzK/Sbr35dTe/\nb6HlPJ4w2zs1VQF2r4HBg4vREvru/AUtq3/MMblEIv+J5hDMRnNd1ppnGmN+e8327RxzzPNMm3Z1\nUx5AdvZqLOs5vFVYLWsVqanvk55eTmrqaOrrj0dkBVp2ZRUiV1Fbm0dNTTmRSLHvNTRCzIsUMz6N\n5PPO+f7Cfq6srJloNN/t5l6paF7Ruzil7ptbj83Dr8Lwe+/tbnFsK1fOJitrNpZVhPaLedSU/Il5\njn0NzdkqRecsBRhHff3mpmOuv/5OFi16HXgA7Uo5A+f9lwNX0r17lyMjyqolHGwOtb8/OpAG0hYk\nkw5sCa9lE1b89szMabJq1SrZuHFjQlG1ysptkpU1TSKR0ibnuC152xLYxo0b5aST8sSrFWRlTZFo\nNNqK7OZiiUTul6ysKXLaadMksaKtSt7z5y+S7Oy5kp6+Vjp3LhDoKlrk0e/aU+X4478n/foNMlL7\nGlGV/0KBH5nj6iQcHi+zZ89uKgq5cOFCCYcniL9tuVBUi1CTVqdOS821bXPSHEnUBtxlXuyxTRG4\nWdSEpTkkWVkz4jROjTjz9/XMnDlT8vPzXdpk1Ei1iYXzLGuoVFRsiVs/WryvVLKypkl6eplEIiWS\nmTm1yZyXn58vKSn3iZpT/M1d6ellMmbMtRKvPY6Q5JnV+WKX4GirJuCW+hsaGuSNN96Q2bNny4IF\nC5o0TRGRxx9fLXC+keRLxTE77r/pbH9zXuwIyvz8fNm4caNkZtpO9DLznqdK9+6TTPBK1KydfFFT\np2rUNTU1El/GZIe5xmpRjeu7ctxxlx4QE2FgwvoXZiDJXq5doyjeiR5vEnNvj0RKJRIZKampj0pq\n6qMJpTuS2eXdH4cTPpmYfOgmFEVFzxtnujuyy65fNV7sSJ/s7Lly/fW3i2Vd4SJMQ0TrN02ReP/N\nl6IJiIkRbHqfxaJmq/jEMiV6M8XbtjcU6i9OkldiYmZiomdfc3yFqNP/YjPGYlFGYxPXxGgoDUSw\nnagatn3aaVOboruUWPj7gmpqanzWQKUkY3rJ7OHJImwqKyslJeVGcYIiEhlIWlqJdO5sVw+238cW\ngctdc1Bo3tsO885HCjS0yZntBIyUSefOPxDL+q5olJmTjLd48X+5IpLKJd5ndYNAgcycOXN/Prlm\nTToFBaWtNj/HXyv+e1GmvUbgeXEiFO02DJ0EkMSSPlFR5nGJ/Nu/DTqgCZPNmdUDBnIEMxCRRHt5\nZuZUKSgobVUcvi0NKeFvkOQRNeN9neVu+25zNnWvpNnQ0CCLFi2VjIwREg6rXV2Jv11LynF01tXV\nyYIFCyQcPl2gxnxkD4vDCB4V9SUMSUI4i0SlN786VVFR6d9bOt2deb7NfLhF5mfXsXLfY7ZohJEt\nidsMoULgGlEmNduM3z3HhaIh1bbz1K7dNVm6dh3a9JHOmXOnxEv4Q2TOnDub3mFBQYlkZY030Xj3\niV89JCiQp59+uk1rK75Ok38AwgknXGGInXv+3KHEdi7TXNc1igQuldWri1s9DiezfK55B37a5hBR\nxuWn+ZQILJHZs2e3aQ68aI2vo7U5L8mupe8w3/OMj4lX0NFt9jNqUAcMEG/Pm/bm4ezevVvC4bCE\nw2HZtWtXxw/jRY27VcA68/88tPhRFDg7yTmpwFvmvO3A3a5996OGzy2ocbBLkmu0elI7IvYnUSt+\nEScLc/V3utsfTUshne4+A+7KtnV1dZKR4Wem0tBkm+nEO983SnzBxxtEGcQa33trmOs4UXOGn7nP\nNs+4t3mr2hbJzJkz5fjj3cXqHKkxHF5lPmhvMUs7usyWki8SNS/Z2soV4iQo+pvebPNMTU2NzJw5\nUwYNGiTLly838xhvVszKmiYDB14m8X1H7OuNlxkzZrR5/VRWbpNw+BJRwv2IKDMpFCgQy5oss2bN\n9TxzssiuYvMedpg5zJNjjrlAVqwoaHHdagRTiThlWrz3tH9LxZ95lotttqypqWm6bl1dnSxcuLAp\nKKM1OJCFCJNdq0+fORIO93U9Y534a8LDBDaLk8dj5xjdIpHIGqmoqJCCghI57bSx0rnzXElJubHJ\nMtESQqFurrVra+bdEo47FAykLU70eWgFOhvb0foOrzTjX6kDLhaRHOAsYLBlWeea3b8FeonIWWij\n6h+2YSxHDOySBC2XOWgfLOsUMjJewusotZ2RsViMWGw3ysNnoK/xfuA+unYdR01NHbm5e7nwwj/Q\npcs4LrpI26j27v19PvtMe1e4ngYYQNeua5M4OotxemHY5UJ6AoLTUnat+U1B2/rWo702Gnyu17r5\nOvHEE4EuaHmJ7Titaz+isfHnWNZVnmvVouVX1qEFI98w51yMliO52czRr9CgB+8cXMJnn+WyYEEP\nIpExXHPNLeTn/4X162cxbdr/Iz09jwsu+CHV1U9RWzuWffvGUF39FDt2HGXuPc81DzcCtb5ro6rq\nHc4+ex4XXbSeiy76PaefPoPCwrKmFrV9+/bmjTcepHPnPea6Z6DlTMKkpv6T8847Hygkfm2Izwxa\naAmWGPoOx/LXv97E1KlrOPfcnzfb6lb746xCy7SEgM9RmdLvHp18tgvwPOeck8IDDzxAfX09S5b8\nlEhkDPPn92D+fJ3jJUt+2nRGsvbGtqPeCRxxAk7a+u0lu9Yzz1zPjBn9XUcuQ9eyXzDHbLRcyxhg\nEvAC8D61tasZPPhWJk/exPvv19HQ0J/6+gvZtu1/GDHizmZbNn/wwQfEYr3QtTvO/NYRi/Xigw8+\naNMzHhC0hsugXXT+G6UI6zz7XiaJBuI5Lh0tJNXPZ99IID/JeS1y5CMVLWkn8VKQvwbh9oP4+VI0\nYaxI1GE5RVRDKBXId/WD8L92Yr5AYraw43yvE+0saCdb2tJuVJzeFrZmUCzquJ4kTj+LwaKdCt1S\nnO38dm+z72VfZ4jJsykx9/GaSSqM9OeWFG3nt/9z631XCtydRGp2h9LWSaIfxC8EVExJ/kXi1pD0\n7wJ59dVXE959VpbtxLXNZ3MFhkpm5ri43JiNGzebSgZLxdEmiiUzc6qEwz8Up51wsrI0k5rZN1Sg\nxleKjy+1v1bUvzJE4rVQ20x2ncS3CXBf/2xxpOnBvuOwrGFSV1fXtKbVN1gqZ5wxSwYNulTGjBkj\nX331Vau+q/39RhsaGlzfTrJy/8WSvDX1KlFfk1+1guGycePGpOMJh8NJ7xcOh+OOpaOYsIAyVIMY\n2FYGgiOO/gO4L8kx64CJSfa19l0fUWhtn2+3HyU19RGJREY1mUXcDjPvQk+eh2KXw3A7dFsqI+Kc\n7+dcLSp6XuBMUdOSl+HZNb6miOPzGC7+ORxDRHuWLxKNFrrTRfyKDXE51xDINYbouTO87egp73XH\nmn9tW/Wl0pzz2ckcHy7qSPcmdbrn1S9/xj/XIhwuFA0a8N6vUBYuXBg3pxUVyRpvTRa4RSxrWFMw\nRVbWNOnc+UHRSgWTxHZSW9YU6drVzkmxcxaq4krJZGfPkXvuuU9CoWTVoIsFFvpGZcWbR+eKE2G2\nQ9RUN1kcAWGowDnmHdrJeMNEI7JsH1+lqE/KfxwLFiyQrCx3H5zrxWvKGTfu+gP8pfqjqOh5SUkZ\nKRpd5cd4h4h/OR/bBFsgarpL3O9dC250NAbSYjVey7KGAp+JyBbLsnJRXbQtGk4MyLEsqwvwnGVZ\nZ4jIH13X/zHQICJFya6xYMGCpr9zc3PJzc1tyxA6HGKxGFOnPsmWLU410uS9LDTvQGPyu5GdPaep\nMqu7T0QoFCInJ4eqqiqqqqqIxWKmyqpXtZ4I3AvcSkuvUiviPsgnn4wFnCqu7vHFYjGi0Ro6depM\nNJqKmqrmAV3RKrP/AfzSHF2FWkH7oeYWv7F9iLaO3YRWsT0BNS81MnToQH71qxDxbXOvwmmbG0Ir\nC9tVhO3rnosquVebeyxDKwrf3czTn4fmbFzuOhe038c0mjev5ZjxxFfLPeaYX/DZZ3U4VZMxf5fQ\nteuouCu89NLL6PvZa7bYbWqHAL9C5Hlqa/Ua1dWZwI+BfzPPpdtFRmJZ15KdfSvV1QOIxT7g2GNX\nMW3aQE49dR89e/akb9+fAPD444P49NNmHqlZhNB8F7uXTA/UNFnkes4JZj5SUZNlN1JTq6mrOwU1\nsz2OmjTtKrV2u2Rnnl944Xfs2vV9s20fmlO0znVMHqWlw1m5ch/p6entfZhWYfz4YSxd+lu2bRNg\nKzAcNWWBvoP30JwPbxXfEjSHpCTptU8++eSk+3bu3En37jeSWLU7n1WrVsXRykOCljgMsAR9Ux8A\nfwb+Cax27X+ZVpiwzLF3Abe4/n8t8AcgtZlzWisUdEj4qcAtRYu0B16NJitrvG/HOZVOl0lij5Hk\nGe3JzAGOicw2IY0SNUPNFTWT2V0U3ZFRlaJtcpNJ/uPFKQp5maipa4o4bX6HippKvM9km5TcvTxs\nafAGc844Ue2mSKCPxJtbvFraRnHa7yYLFW4QNXWdJSrlu4+x26CWC5SJZQ2V0aOniWo13hDa++PM\nFl999ZV06eIn1d4gajIp9dmeLMquTFavLjZtcwvNb4rAA3F5LW+8USkaSOAtZTNMmjNhOVquW6Nt\nrnLybIHzxNEcHjJ/F5p3N9xsKzJzYwc1jBTVHO06YcmqSBfLmDFj2vUNtQWJ3/DHAqkCiFOpYbRo\nMIZdD22kaPivPa9+lalHtlid13GiO5F/h8uJ3raDk5uw+iY5/tvAN83faWgTiSHm/1eg4ugxLdyz\nFa+zYyKZmepAM5CGhgbJyvImqDWYMul+RKjUfPDbmswZzZnH/JDcRDZMEv0Wdq0r2x5+lSRvYuU2\nEa0Wv8Q7/Xh+YY5dIWrWyhXNOZks2tipRBxT1DbXve0wzAbRekYnijIq+2O8UpSJrRXNE0lWVXiq\nOLkOdujmHeIk611q7nOfKOPaLPFhxE4IbefORU3vXcOCkxWjLDPP6Q6Jts1wyXJBSiUz028Oda77\n9Jkjq1cXGzOY7SeZLFpgcoiEQrOaXQtOA6YicZhtcwxkhijjbpD4/CLv2OywYJvpXS3arCvXbE9m\nciuSM87ofdDL0zefKLxB4NuizKPUrIH+ouVmCgUmytFH95dw+EdmDenai0RGtro6rzuMd/fu3b7H\ndGgGgur1e4Eao5n8xmw/EXjB/N0b2IyG6m4Dfuy61i40VGaz+f0syT1bNaEdDc2FFDZXw6qtiz4+\nSdD9EVZKp043S0rKCEN47IJ128SypkhaWqlkZ8+VysptcRnE+xsnrx+DN/T2YVFGYCcJXiEqQY6V\n5nufLJF4adt9jxxxugGWiErVuaK+mG2ivpZFEp+gaNeuWmrmqFi0BPtgUcn4LInPVE+W+PeMxGcb\n23M+ROBGUQZiZ9XnijK7ZDWoyiQjY4REo1FXFvMqUeYV9RxbaMbsTU4sk2RaZFbW+GZqq5WI+kwG\ni1fzOOGEsfL66683uxbiC2FWilO/zPZz+Wl2/UQ1DFuQKZH4BE57bH7hzpNEE0OjEh8U4Q7OGCJQ\nkFCR4UDAm22fvNZcjahQ4t03TWCJfOMbA5u+tY0bN8qqVaua+gMdSHQ4BnI4fkcqA2lJy2hLYUY3\n/Bex23zhTnwrkU6dBskJJwyVSGSZpKWVSNeuE2X27Btl48aN+yWhtZ6B2I5f78c0RWC5IbKVnv3l\n5pmSJSAWiyZl+ZmWhgiMl+7dr5Fjjx0jTr+MLRLfW2OuqJS9ynzwswW+I/FRV/5EORw+U5JL2MeZ\n63kJ58Yk5xTIokVLRUTk3HMHitO4yB7jDtd9h0hq6hpJSXnIlLUvE6fcv+28vkHs3iZ9+syRgoKS\nJAyk0BA529Edn4nfGm3Yfw1EJRJZKueee5l5F27N7mKBbEleSsXeViLKQBPnSiP97OMeExVE7GdI\nDKo4UEUh/awJRUXPu0y47p48C8S/8kK5wCKZOXPOfo+nNQgYyL8wAxFpe7ihX/+Qzp1/II7Zxz8j\nORIZKffcc78J9dSy42lp2rWvvSGPyU1YV0q8SaIiycdUaIiq3zUmi2optrkiWdip3we6RKCnHH98\nriGshQLfF38/xVBRSbnInNdDEon8DnEqBOeb55vgc5zNQCaa+z7m2f6M+CWcpaSMkIaGBlPBeJDP\ns9rRa3lyzz33NdXBSk21GXW+aBKhzRy1snDnzrc2NS/zb5aUJ/5htfo+0tPL2sBAbA3gddEoqkHS\nqdMCSUkpMPb6U0Wj674vWuMs2buzBYlp4oSDu39loua9VeL4oKZ6nqFBVAPUJlsHqrpvc9YEd7WB\n9PRy6dQpQzRR1TaVusd/YasTI/cXAQM5ghnIgcyKbe56Tte0582HkyzEb6AP8RopffrMaTGUOBmT\n8SvVsnjxw6aYYrmkpNiht37EoEBUQn1QlEAXSChUKOHwFRIKTRH4gcAZolLrZHHyGMaLFjn0agqV\nosEBiwVON3PxC4mvVzRC4CbXh10oTgkU22fiZ5oaKWoOGyjKhDZKchPNHHPtYeIwLLuu1g/FKTev\nUvn1198mIiL9+l0i/oS1RODSJkYj4m2A5M63cXJL3ISzoKDU5PS4fRzeOXQT8grJyprm6zS3iw1W\nVFQY39sU88yjRZ3H7r4oeaICxJWiZWOGmveTjPkuMsc+IMmDG+zrTxEnHNs+xtbA7OoCoyQc/tF+\nMxBtUX2/eDVlP2EQTpb40OJR4jjOh8qCBQ/v11jagoCBHMEMRKRt/UNaQkvVffXDWZyEKBRKsphz\nldSSM7iW8lX8mIu7AnBq6gjx1yJGiUYUzRZ1Zt9kCJD9jNsMsa0zY1wo8Kao1DzYdU1v29ehogyk\nUOLrFXmL340SZVKDPWMrl8Re478w29x1vqaI2wHqRNjYknSxGXNUOnUaJqpRlZt7DxLoLnCsrF+/\n3vg+kneug2xZvPihuPl1tAp/M5v7PVZUVJjKyLni+Dpsn4XfmhgsBQWlce+4oKBEMjJGmvGsERgi\nXbsOlhNOGCnKQOx+3zeYeZjr2jZJHMYcNcd781tGiub93CPK4O0AiHJxAhO8fVfGuebMfx7sBMT2\norJymynps8b1fP5mvt27d4t/TshwgWEyevSUdo+jPQgYyBHOQEQOXFZs8wzEJlhuwuos4M6dB0t8\nUT03sYivHOqVqvZXiyoqel46dRrgIbYjzLjdBM3NWOokWd9slVILzUc8Xfwl1csMkbKl0+b6qXud\nnSLajvYs0UQ127l+rsTb7u05d5oKJb6PPIlERoqaXeyijWtEpew1Yleq/fa3T5Hk9cIGSUrKvZKW\nVhbHvDV4YrookX7QzG+BgFZvXrVqjaxatUruuWeJnHbaWEmM7PIz/6j5KDV1RFPFWKerZpF5thni\ndBLMl+SamNt0Y9dDsxn5cHEaWE0y62CeQC/R8Oc3xEnatCO1/NZvqThResmTQtvb+a+yclszLQ4a\nEr4Dy0qe5Acn7ncHwrbiUDCQoKHUQcaBqoWVk5NDVtYGEpsDvYImrkUJhwcQ37xoLXA1P/7xxaSl\nlfqcW4QmOvkjWVOstjQbyso6hc6djzFjOc38yoFfoIlx/cyRm9AEwLFool8O3mY+CjG/XsBMNIHL\nm5B4LVp5x0Yx8Qld9nFT0JpY7sZB69AAwUvR9KeL0TJtx6JJjfY1ctAUpolo3aMwzvvIBvJ56qnL\nOfPMGrRB1SXmWgVocuIE8/sVn39uX9dbL+xq4DvU159OTc0rbNlyIRMmLKKxsZGcnF4UFMwkNfUT\n4DjgR0AWINTWCtdcU8I11/w3d999Mu+/Px6tRvQnM/Z3gO+jdcomool8JeZ+9dTVLWL69J/z9ttb\n6d//bnbtegqRq9CaTk8A30TL4PU054TMs7vf4VbPXIMmGj6DvvuJaB2nO9AqR+cC84Hjzb/l5joL\n0Hed7NvJQtfBf+PfKKt9iMViTJx4LzU1k0hcNxeSmTmJp5+eQVVVFZs2baK+vh6RY5u54v8xYULy\nb+1IRYuZ6AE6BuzibuPHz2TXrovRjN5XUIYBlrWOTp1G0djYC82CVgKfljacYcO6c9JJXZk+3Z0t\nW4x+DO6Pwy7CGJ8V3V7EYjFGj/4PamttAtsX/ciL0eKHr6LZ1i+ijAXgj2jRPUEJ3RjiM243Av+L\nZqCH8C/QlwIchWYE57UwyjPM2OYBFwJPAv+JEq74LGclmo3meJvgX4Nmh8eA36AVf4YxduwpTJo0\niZkzn/e5zghgOk5t0huBlWjxRvvdxYCv0CpCO9F0qxjV1aPp1Ws2xcU3EwqF6NQpExjtmp/VKLOw\n0Ex6+76jzH3noEzgCRzC/zZwJ5rd/XsgQnX1APLyfkBNzSwSCWguTldN+z09abYD/BrnfYMyytFo\ndro3M/vnaKUC9/zMQ4tiTjHHvIEydXf29Xa0IsEktHtlCfo9xGdop6WVMGFCKW1FVVUVe/b0xa9a\nQ0oK3HXXaKZP/7kRsODb316E5lzbay4+S3zp0mWEw/965PZf74n+hZGT04s//vEJuncfwZ49/w9d\nqO8Bj3PaaSmkp7/O1q1j0cXbF4jRo8cqcnLymD17Nlpl1q56axOmawiFvkckkkZm5gZWrHBKlajW\ns4otW0bi/iBay2Teeust9uz5Ak0XiqHVSJ9BJcuLUIL5DMo8eqBlTBpRSfl3Znw3AQNQ4rMKTS06\nDyUyA1DiEz8+1SguQxnUcJTZFJL4YRegkr62JVVmMA0l7j8kkXBONNe2y5r0RMuHnGau9zGQRSQy\nlV27Xueccy5Mcp3JKLG93Wz7NRABrnRdu9CMI4QSZ4cZVFfnMXXqTbz99kOe91OFlpCx0Na8fvf9\nLlp5F1RjANUCZ6Nr6QHgLmKxGJ9+2hP/cjdRtBPDOGAxyrDczGokyhRrzLYzzd//BfzAdR27YrMf\ng/p3dK1OBsab8c1D1009jhZrnzsafdd55niIRIpZvvy6dhNuyzoFLasSv75OPvm3PPhgOlu3Os+8\nZ09XlMllkFjW5B0uv3xZu8bQ0RGYsI4whMNhnnvufrKzjyMS+ZBIZDd9+kBJyR2sXDm7hVLWYfSD\nvAJd6JcB59KjxxLfHu+hUIjbb7+UtLQ8VMIrIRIZw+23X9qiOW7NmnUMHLgYJSRd0Y//v9AP/yqU\nGM8FjkFNEK+ifcj7oJL3CtR08QhaF+k5lNC76yMJDtEtNdceg5pAdqEmpufMcUejH3aJ+V1p5iHF\nflqz/wnggmae7EmU8RSa5+hvxmQztaHU1o5jy5aHePfdfc1cZzBK9EajTPQbwN9RxlmCagzfQRnW\nALxEtrp6IFu3bo0rOR6J/B7/svg2wuic78Epef+x+fsT13ExMjJ+Syh0Ibbm496ngkAENYPFUKbk\nZQIXAP9j7nkn2hv8NZSgtlRe3osQ8D30PX+MVlUa73PPa1CT2GvMnr2BL78s46qrhrfi+onIycmh\nR56uGn0AACAASURBVI9X8fYz1zLzlxsrgPv+fVGBqBPwLPoeGlAmdyXXXfdEs2Xaj1RY6mvpuLAs\nSzr6GA8HYrFYkx8iJyeniaAn2/6Pf/yDb35zImr+OQ5H0l1NOPwuDQ27fe/Rt+9NbNnyEI5NO5uz\nzroloeijG42NjXTpMo6amnLiJf5rUC3C3nYDavt+G/1Qf45Kn4KadP4KPAa8jtq413mudxPwEFpA\n8R7gfZSwvY+aNcZ5RrYKJWBZaFEFr026EJWoF6BE3Dv+0ai0/qm5x3RUQt+AmrN2oNpIX1S6f9dc\n53nPdYYDJ6ME2GbYRcCHWFYVXbq8zxdfnGiuH0O1rDtQjUffbWrqLl57rTv9+vVreuexWIzp01ex\nbRskmrBiWNa19OyZzh//uJd4s1EMldxrSU2dwumnv8Hy5TOZOPEBqqvtuRmKaocvo5psL3PeJei7\nm+iZy3KU8QN8hDJ20Pf5Io55qpDE9zoPZWhlOEaStajQk4+arXb73LMIyAQ+ID+/gcmTJ7M/qKp6\nh6lTn2Tnzv6IfEjXrptYs+ZOoJ4BAz5m377RnjPmoYKEd1xriUR28/rrl9C3b18OFSzLQkTaVPy2\nrQhMWEcobOd8a7d36dKFnj338e67x+G1yTc2DmfPnj1kZGTEneM40d32bNi58yIuu+wyvvWtb7F6\n9eqEyqfFxcXU1PhJiMNwKuVuN7/voFL8z3Eq6oJK4ENRovwiTqMq9/UuBAaZY+9HiXI/HOK1DyW6\nYdTpHUbNL/8F3IIyGK9/5TQch7JjDlHG809zj52otGn7ELqiVXr/iUqcNlLMuEegmh8oARwMXE98\nFeHOwOeIfMEXX/QmntGOQRlUBHXEC3V1a5g06ZuUlNxBTk6vpnf+zDPpTJiwjOrqL1GGPQzLitK1\n629YvHgYzz1Xyh//6NcAaQy9ej3MqlVZ5OTYzxxBNTJw/DIbzZzac3Yc2sjL69t4BWW4VThmsHeA\najMXz6JM+Edo0MP3zPv5JWrCnIdDnmIoE+5i5m0fuia893wR6N5uv4cX8ZWwTyUn5zZCoRCxWMzX\ntKsCw3n7fd8jCQED+Rrhgw/ewF3uW6HRSFlZWdTW1rbiKu9QU7OO9etnACHWrh3HuHEZlJT8zHWf\nvai06EUj+qHFUOn4dygxieJvC78alTx/hRKZmOeYBtQpPB7bfKRM6WbUlLAG+AtK9J9BiVcjSoxP\nRwnsYJR4v4ra2XegvpjVKOF/HfUNnAw8hRLEdDMOt/N4DMp4HkJNg3awQDlwNvAT1GH/LI7Z7CKc\nSLhXUdPMEtRn5OcXyMKJWhvFrl3zGDFiCWVlN9KvXz9CoRCnnPIdQqG36datlltvvZUuXero2bMn\noVBvpk//Odu3/9NzbRth+vfv18SINm3axCefuJm2LUDsMXMQARaZ5+6FMsOB5phfmLkLoYx7Fap1\nPYnDMCeYORph5lNQIWCumZ+fo10Ww2j59+tQ31kI1Vj7eO65ATiLlJT/YPnymxP8Hsk085bgJ5DZ\nAS1Tp95kyuTHOPbY59m7tzP+zHQDPXqQpIvnkY2AgQRIikQnegyVSN1MKL4HQywWIz//fdT+7W4l\nqx+SSrDr0TDZMNrHYjpK0L0Iuf69AjUL2QS0ETV/LEAJCmi4bIHrvHEokbkV1SZmoX3N/oyaOnqi\nfSgaUYLUHWUgJ5nzP8VhXF6iHiOeIGLuMRxlnh+gGsksVKsZRaJJTVBTzUi01uh9+BN3cCK/3HOT\ny9691Zx//kKyso6hsfFDPvigC+oDghtuKOD00/fRs+e3WL9+G198scVcJw+/SKGHHy5Pcm/vmO1n\nv9XMS3zkn76rbq5xzkA1pxtc97SZ7zWohvIaapr6C8qkOqFMdRfqhH/BnPuBOfbvqO/sBHQ9ncNJ\nJ73E7t0vkZJiM2iFbYqqrh5ILBYjI+NxiormkZPTq1mmUltby8033wzAww8/TCQSadqXk9OLp56a\nRl7eD/j005588sk30Gi8s1BBJNc8wy85+ui/sXLlAwelpfVhx8FONNnfH0d4ImFHwscffyzJ2pZ+\n/PHHvue4s+k7dVosybqs2T0YnDIb7gxxd1+QBklsA7ox6bjgK1fSmJ1c+LA4meF2rail4k2KjE/s\ns/9eJE6/kDvF29FOq/x+JVomY7ErSc0ux2Inkvn333DGaVceHizastXuiuhNSLMzse3Cg/4Z1Zpw\n5014LDH3WiCa4HeFz3nfEydp8SLRRL07RBMo7eq+Q2TOnDvj3nvzHS3vFa1soOXo/cvrzBG7qKOO\nq4/r/TT3jHNFa2YNNz87E/1iM5dTRBMx7XIlQ0XLpOv7u/fex1r1HCkpIyUzc3LSCgtaVj9+bbjn\nyEkyLDfv+gpxkh2jZm7yBRbLjBkzDtQn3CZwCBIJAyf61wypqd2pr++JO8wwJeU96uoSneg2bPX/\n9ttvZ/36WSRK0iWMGVNGeXk5mzZtYsCAj9i3bxSOI/lDVOLPRjUYdwivLZVdBfwfTmLaatRs9DNU\n4h2BmgY+Qs1Kv0T9GLYG8DbqWPU6xtei0nBf83c1msjXC5XCvQ7c4ahJ7TpUM3kT1ZDuwXGqv4P6\nPMb4zIXtyLU1pRiqxYDjwwiZZ5iN40QvQc16r5s5exKNvqrj+ON/wTe+kc7777vnazvqP7Df46Oo\nU36SZzyF5pw3UXNPrRnjIKABy/oVK1feyDXXjMULW3J/550LaGiwx3yJmYfhqDbzornmCDNf5WhU\n1mZ0ro9HEw+/Z45djWoqH6M+khiO5vJHNJruWFQ7zDXzuB3tRXeVucdG19zZwRSPAGBZI/jnP0so\nLS3lzTffBGD16kuprfU+XznqpLe7W8bo0+dGbr31XBoaGpg+/VkSAw2GU1NTTkpKCj17zqS6+inX\n/kZUy/QGTIxjxYqhXHfddQnze7BxKJzoAQP5GmLPnj1kZWUBUF1dneA8T4Z9+/Zx1FHj8CO6X31V\n2mTC6tnzWqqrv4kSrAY0we4L1OH6JA4RfgI431xjHWoeeQM1/UxEiZyghGs6Skg+RHNGuuMQIfta\nXxDPlOKJi/pJ/o4SsTkoMU9khurzuNt1jQmoc/RjHJv7etRM9yMcWz8os3CPAZwIotVoDsMpaA6I\n+5gS1HdwnLlHDA1l/oi7755COBxm9epdvP/+xYhEUWf8i65r3IMSXb/neRidV/95iUTy+PLLUt98\niVgsxpw5c3jiiZPR/BE307avNQ/1ZT2EMiw7K/8G4HNzXl+UMT6BMpUslAG4ExCfQwMUHsKynkTk\nGbPdzr/xG3+IeCGhBMu6A5EzzXl70LlPjIxSgeMSc946YDnKgKPomr0Dh8HrXM6evYHp06dz0UUf\n+DClB8z4B5v/vwz8kTfeuJcLLmguNPzg4FAwkH9Bo1yAlpCRkUFtbS21tbWtZh4A6enpjBtnJ0o5\n+RTjxmV4IrEiqESch0qNz6ARNP1RaXMTSuweRRnCYpTw3WPOWYDa1l8z+x9FtZfJqJ3ZDdse/yj6\nwV+HErE1qG+hK05Ule2baQleX8MINDflEZRQdUMl4FRUyt6NMo7xOH293RCzrT/KMF/3eYZ8lGG4\n73EasJ177nmPu+7KYO/ez7n44l+gRNhOMrRxCuon8MvZONfnmQai7yJEbe14iouL/WciFOLBBx80\nY06W+HcRylhG4bzfd1B/0yj0Xd5kjn8UJdCrzXM8gr6nk1GGcxpwLSJ7UIEhH83ZSTb+RIikoQxh\nIsqo3/SZl9+gzPwpdL2tRBn4OJw1+6TnvPh5SUQG6vPqZn4/AWaxa9cu32v8KyBwogdoE0pKfsbK\nlfu4+mrNI1m9uiyOeVRVVXmidzB/D0Ed3NPMtlXox5aJOtRfQ01Cj6AS7SSz70nPCLJR89Fac42u\nOEStF0oI3gYeRAn8qyiBOxOVNv+EMoQ8kpWd0EgpL94y//bFkbrdUrFdgmMXas6w82Z6o0TzfDQa\n6y9oGOow1IGMGcc3zfNfi36WG1AmFTbj7UZtbQnr1482c7TXM76JZuzzcCT6l4HPUALcPKLRaNJ9\nzz77LBpufQtK5L1oNM/3a5yotiJUIOhtjhmJMpFBqHkrC2W2trnOHrOFzs2FaAKinT+UDHbI8Cic\nUi7u7P8QOo/XmrEJGh13Jqpp2gLJB2YsvVznDcAJO3cCDVJSUkxwySji147XLGuu9K/oPDcIGEiA\nNiM9PZ3y8tZE7LghqKRn+wZG4piUzkCJ5E2olL4JzWTuh4Z4zkMJlF3/axoaqpuBRhy57csh1Nw0\nHpUCc8z5oISgH1CHag5pJJadsMtv2GU+slFmZWckD0DNFH7Z17momWscavb6k3nmG1BTXjGqbd2G\nZuivRwnXWDT6ax7K8AYRnwRooxTVpPqizMKdhxBCie8n5tkwY6lBbfleRukmukWItFSa5hLUFHal\nz7XeNPd+xrV9NPFmJpsgP4rW26oyz+6NZBuJRlfZWfpTkoz/OTTaa76ZrzIcE2F8FJYyhaHo/P0N\nzYq/jfhowgme8YKuk9+imutq5sw5uykSyxvG27Xri3z00Wc0xBUCaH8triMF/7qsMcBhQU5ODpmZ\nL5NoMngRdzKiU+7iU7TkSE/04/0QlcZfNOf1Qon+XtRM8yhK9B9FP/Cx6EfuvZ9dpdi+T3+c5R5G\nJc8yVFp9GpX416KEfRRqltqNMoOrUWL/kBnDp/jLXgJUoA7aPPMsz+OYRe5EGcjbqOY0FS1hko1T\nrmST+b+X2LtzCELEV+4tQYltg/m7K1pxoBOaeGeh0neB+V2JMtcyVOvp2Wy9qAkTJpCWVoyaJj9C\npfky1FR4jRm7t7SHn5mpHnWm2z6jteYY73lXYJvXnHcxD6dK8Y3YvhKd293AvTz99BRSU49CGYV7\nPWw383K2mY/t+Jvi/n97Zx5eRXU28N+5XG5uqKLFVq0CApoAooSY0qooRLsJKAQIEHYEylJURKtY\nLKKIC62tWrVFgUAIJCEQK2jrWomohYaEsKVAUEHAr7RfPxdQst/z/fHOcOduZCErnN/z3CdzZ+bM\nPTOTOe+cd7VnHFjts0lKymP69FxKSnJ4/vlHT+5tBxlu2tSZDz64nL17l5KWNjMk7c+yZRPPyCSK\nNsaIbqh3Vq9ey7hxOWidhDyImUiMRiqBA+9q2rX7AyNH3sSf/rQHeQt8BXk73Iy8+fVDBvcbCW8g\nvgIZ2F5EVFW20d3p4WRHpjvVT/Yb/lbrd2wPn4kEvpnaBtsbkfQoCcjAuR5/RmN7v9HIG/QIRIUT\nhQg4Z8BhJaIi+yEiJKOQWZYbv1dShnXeFyNqukWIKqjc2sd2YvAhAud+Wrduy+TJV7N48XZE/fdX\nAqPZK5FB/iAy+NpquhFER6dw7Fh4I7pNZuYGJk9eQUnJIUTY34WofKKse2A7NDgJvu4/RWw3tqdc\nNmIPGRXUzmkUX4cYwt9HbDkg97cHImj/yEUXKdavfw6328fSpVksXvwWkoxxrHX+WYh9w/b+uixC\nf1cjQv9yIBu3+0q2bBlaq/QjlZWVJ+1JKSkpTSo8jBHd0CLp1q0LUVFDgReQiOIxiF7b9uACGVDW\n8M03s/n737/hiiu+i6hehiEDx0fIG2YnZGBeTngDsQ9xI7VnB2mI/aMUGVwrkQH5E2QwGoCoizKt\nz9PIwLHL6me4+iL9kIHF/t31+NVf6cgAPwTxoDqKCJx/ERh09wwyYNnC5QCiMz8fETonkIGyMzKw\n7kKCHIdb55lJ69ZDGT68A34nhrXAo7hcVaSlTeZPf1pISck6Lr74Hus8nechiRS7dTuf6OixiJ2i\nNV7vyBq9JY8aNYhjx7J58cWpiJpspXUN/oIM5LmE3p/liGvuGkQVdQmB99FDYC0Wu93r1l/bLfgu\nYA2tW+ewYsUllJX9mby84aSna/LynuEvf/kT06cv54YbDrB4cbzVn58gs62X8Uewx1v9jNTfLVa7\nlUAmV111vNbR4263m7FjxzJ27NgzeuZhc+afoaHRETVWKrt2nUf4mhol1vqbKCtLYceOEcTGTiQm\n5p/s398df/TycGRQV0iMyEDkLbqz1f5zZNBViOD4OXDc0Q5ER34BIpy2I95U38GfYv44Ymiegxh9\ng43TIIIgF/geoobqiQiCg/ijqccguvoNSIoWEEHSGb+B2LarxCNCaQciDDojA+UQZJBNQ1yVrwTO\nQWYqPXC5bqS4+AhbtjzKX//6Ch988AF9+vShf//+eDwefD4fXq+XSy+9lKNHw9VJgXPO8XLsWLbj\nLXltjQc6t9vN1KlTWbt2I++884nV1wvwe5eNR1RlbmTG0A7JMOBCBC2IvWoY8r+QjaRusVOSaGSW\n2g1RSz2B3PeJQBku1ziuuuoqPB4PvXv3PplIUhJ+Om0Xtj3jVUQgODMuT0OcATo5+uvBn8qmBzCR\nmJgJpKY+FJCktKCggH379tG1a9fTLhB3pmBUWIYGYcaMO1m8OJLaqRjxsIrBtou0aZPDu++2Z/z4\nJVaAVoG1Xzdk9rAJMVy7ESNoOaExH7chA1rw+gmITeA6RI3x56DtfRBDdzhDqh1c+AUirJ5AVCOV\n1rkFZ+yd6Pj9IiRf1A+RGUeitV8uog67EVEnua1jrUHsGCMQw/CV+FN53GS13UiHDkdo1649e/fG\nUF6eD/QnKqo13bpt4v77f8yDD67kwIFSQuN1+jNhwvd46aWXQtJ91IaCggJ+8IMd+HyPIQPxJGTG\n9R6iRnQhHk5xyOzBSQ6SKmYF/jgcZzBhMTLDOYy8EPwGCZBcRnT0z3j//S4BKiUJXA2XGddWg32M\n2MiCg/5uIimpI6+80g2ZrcWf3O71ZrNpU2d69xaHj8LCIlJSnmX//pvQuhVKvUZMDGRlzQkof9Dc\naAwVVpOnKqnug0ll0iJZsCA4XYk/7YnUyA5MgWHXYrdTp3g8d1jtI6W9GKdDU3s8oSWVSLh0Jlka\nOlj75Ae1neHoqzMFS5aGH1lt0jUs05L+I13765uHO798x/cyLak3gvs/2Np2h+NaZGm4R0u6EDuF\nxgAN47W/7nqFlhQfFWGvi6TXKNOSquRW65gZ1nF+ryOl/KgNCxc+Zx3Prkdfm3tkp5ZJi/D/YW+3\nj7lKS+qbtTo2dmRADXKt7dQ5a8McJ9s6zjpr+U4tqVVW64svTtZ5edt1VVWVjosLTcXirHUeaR+4\nQ8fF3RHSn+YEjZDKpMkFRLUdNAKkRVJWVqaVujXMgzdIu90DNew85UOblpZmDUCRck4FD9Raw+M6\nfK6utdaAN9oatHKsAWW3tX2hDsw1ZecyukVLHit7MO9rHWeNFqGyKsxvrdb+XFvaWg7XpwwtgiDb\n+uRpSNIiPIKv2XDrGEO05P9aowPzczk/mda23VpySg3XkgOrIuCYSt2qy8rKTuO+Ou9LfoR7tCro\nWlRpETrZVt8GhzlXpyCyc5et0Eqt1qtWrQnpT0VFhSU0w+XsKtOB+cbytde7SOfl5Z1s78z11qbN\nOh0Xd0dATix/brdQQefxPKnz8/Pr8HQ0Do0hQIwSz9AgeDweFi68BaUG44xanzGjO1u2PEGvXksi\nVk50uVyMHTsWybbyFKJyCMZODW/jQ6LMw0Udv4Ho49MRlYkdjPci/nxXBwl0Ff0tEkndBbFRJCHq\nsVetY9yPGF2Df2sjotdfZ31mExqXYO97PeJqu8H67S4EpgLHWh5m7fMIov6yzz1c9cEq/O7PzyHe\nW3cSaO50ofVYFi2qfZnVRYsWoXVwPZFIeBCHhnWIo8I45Dp3QuwNCrH7ZCL/H7MQG4V9bI24eMfS\ns+dmRo0KrW+/Y8cOfL5E/C7NtpvvZYjXVhXwGHLN4unW7UiACszvjtspbFXOU1FeXnlGVhmsDcaI\nbmgw5s6dyS9/+fOTA9WcOX8+qXv3F+qB+PhnQwySLpeLrKw5TJjwR3btsuMqnPr8tcgAOt5al4nY\nH2IQm4ddrjYDsUF0IXRg7oPo4Uci8RufIgOWHa2uEOHzIf608+A3hk9FPMRuQx6lDYgQuwWxbbyN\nGOhXIELL2f8t+F2Jk6zP8QhX0o2kTl9i9fV3iIfb6wSWdrWF5RWOc7wswjFPF7vOR1LQsrMv7yOe\nUDtwu99i7twraN06jrVr09i16xha5+B3W85EhLOziFQGcBMxMc+wfHmgQdv+3/H5fLRq1R6pI2IH\nJyrknn8fsZddAQzloovKSU19Ouz/WiRX3fj4eDp0+BP79w8LOrdcxMD/01pfuTMJY0Q3NGt8Ph8j\nRkwkJ+dz/BHj6Yjb61s4S+2K91VvZBBPQLylihHj7o2EZqrNRGYXTyKDwjrEYOs0uNpxINcghvQD\nBBrDOyIC5hhiQP4V4hr8mfXb3yAJBTsjxvePkViXqUhQm81qZMD9H2t78O8/gxj/ixGHgGIkBcwB\nnEWVWre+nM6dt3DkyHB8vipKS/+MCKZgg/ptPPTQtcyb96tauZuWl5fj9Q5D6/UEZg0+gMRnDEdc\nhANjcdq0yWHTpk4kJCSwdetW+vY9GJSMcBfi4WaXun0TuIvY2OcoKlp8so+FhUXcfvti9u3rC0DX\nrpsoKfmK4uIV1nGcMT7O6/d7lBpCaWlOrR0IJK7pZbS+DZkxvgYooqIS+PDDGxq1TG1taFZxIEop\nl1KqUCm1wfqerJTarZSqUkpdE6FNlFLqH1a7XUqp+Y5t31ZKvaWU2qeUelMpdd7pn46hJVBaWsqM\nGTOYMWNGtVUQZSaSSlSUG4k8L0NmI1Pxl9pNQN74z0FmGn3xl0ZVSG31NYSqm9Yib6nXITOJbMR9\nNqAHyADdFUnBYcdz2GqwXMS9dyYixO5H1CWSSVdUSJcgHlmrreP8ComfcCayjEK8jqYgb/JrEHWM\nU62jkQDB8fTs2ZmYmCL8sTKdgKfp0eMTiopWsWlTJ5YsKad168GIq2qgKhF6smBBV849dwSZmRuo\nKbZqUlyqtyNCMxvxmhoMpNOu3RJEiNqqIB+xse+djKlwuVy4XMFuxlfj8fTjvPNewO3eT3T0z4iL\nW0pWllQX9Pl8bN26laSkRezY8SylpcMpLR3Ojh3PUlamiIubRatW9yL3Mlwcz446q+1GjRpGz57f\nQQIMKxAVWSrdu398RlYZrBU1NZYgytxVwAbre1dEX/AucM0p2rWx/rZC5u0/sL4vAu63lucAT0Zo\nX592JUMTU12hnkg89tjzWgoKDddwsw4sHhXJC2iw9XlCiweS7V21Tovn0zOWAfheLQb7bKtPTgO7\nbcxN12L4rs5Q/GcthvYnrM/PNKAjF8z6S5j+l2kpaJWnAw37A63zGaDz8rZXawAOLKZUpqXw1DXW\nsr8f0dFDdEVFRa3u49KlS7U4HwR7tK3SjzzyRC36FehIUVFRofPz83V+fv5Jpwr7PKOintTiTGAX\n+ZLfVmq1Tk1dqeE6616FGrxl3yy9YMGCWp2nTXXXujlCc/HCQhy330bm7huCtm08lQBx7NcGyAd6\nW9/3AhdZyxcDeyO0q9+ramgySkpKIg6kJSUlEdtVVVXpdu1u0uKWmqmlIt5PHMfJjzBwZGlYag3+\n2UEDzzfWMZKsPlU4tjldZO3BfYuOVI1R2pRomGIN0NO1VP5L0CL0WunILs1of7VG57bZWjzD1lr7\njdaw3TrPx3V6evrJaxM84DrZtm23jo2dbB3jzoj9sI9XU9LT061jBh8rU6enp9eoXzUZkCsqKqz+\nV1n30fkikK1hsobf6gsuuEb73XXDeXZV1NnzzKa6c2puNIYAqany82ngPiTLXa1QStkZ4i4HXtBa\n2wWsL9Ra/9uSEEeVUhfW9tiGloXUlx5PoH66EOjD3XffzeLFi8O227x5M59/3gbxblqCvMd8Cwnw\nm4DYFSLVNXkesT3swJ8mZR5SMGoKYnRthUQnJ1pt0hCV0OOInr+ntb+bQC+pSmvfC5F3oxutYx1G\nIq4PIrr8U6uho6K+pqysu2ONfV0mIbYTELuPC0nx4udUBmAQL6OiosX06DGW4uLuEferLSkpKUyd\nOoKSEjtNiPQ7OjqblJTsGvWrOkcKCeB7lOJi2wEhBlE7OgNBhwGD+Pzzrsi1t5NM3oiom/4K9EKp\nJBYu7H/S/uE0xIerhx6O6s7pbKRaAaKUGgj8W2u9XSmVSHVPQxBaax8Qr5RqC7yilLpSa/3PcLtG\nOsbDDz98cjkxMZHExMTadMHQLHEmGLyctWv/wbRpRVx9ddeQZHRPP/00IgSWEBglPtNq/1NEt+8c\n3H3Wus2Id9YQRHcdhwzOtlF5K+LuusLRNgkRdHuAZ+jY8Xfs2/cK9967gD/+0U7//hmivR2NPBI+\npMbFYUQwvWv1LQYZ2FYSrvZITEwMjzxyOxMmDKGiwq6atxqXq5L27d/m0KExBHv/eL1HSEm575RX\nN3iAzMqax8SJf2Lnzm0h/ahLynG3282yZROZPDmZkpKRAHi9WSxbdnuNjfLhBmS73z6fj8mTV1hC\n7wAirPchLtTBNo6JaP0CIiySkP+RQqASl+sj5s27grlzXz4pPOxSvcXFiQDExqaRmjqtWUeV14Tc\n3Fxyc3Mb90erm6Igr2GHkFe1fyH1Rlc6ttdIhWXtOw+4x1reQ6AKa0+ENvU2pTM0LX4VVrCKKF9D\nnv7Wt27VXm+Stu0j0dFDdEbGej116lQt9oRwwWqrLZXR01rsA1nWJ0nD8w61jb3uIh0YQZ6nI9s2\nJNDP631SL1++XHfsOM7q+xYdGHgYqCqRILaKoG2drXO3+zdAQ4eT16aiokIvXbpU33zzzXrKlCm6\npKREb9u2W1922XirL1kaxujWrX+iMzLWh1zbL774Qnfv3l1369ZNP//8Szo2drKlHsrRvXrdqbdt\n262rqqr0woVPa4/Hf5283qSwx6spFRUVOj09Xaenp9fajhKMX62Vo93u31vXK1P7gyhn6/DR/3bA\n4e3Wftka1mqlBupVq/4c8Bunsr+0BLVUbaC52EBO7izuDOFsIAkR9v8OcJ61HI0kNBpgfV8Ef7a1\nYwAAIABJREFUzLGWjRH9LGHEiBlaIqPXaH/akBzrM0YH6tXFwHv8+HGrTTg7R6Y1GFdZA/tdWnTl\nzihkZ3qSdC0Gc9vmECmie7WG5VrsE+kartKw0ur3T/Sp03CEi5LPss4RDeiYmJiA65KRsd6KqA4U\nnlVVVXrLli16wYIFOi0tLewgfdVVP9bBjgmSPiX8AFmfg3594R/YK7QI9ZFhhPAdWmwe4QR3poYf\nakjVcItu3/42nZe3PeR38vPzdZs2oS8idiqdM4lmLUCQueJhJLXqv4DXrfXfA16zlq9G6lduB3YC\nDzqO1Q5JW7oPceg/P8JvNsS1NTQB/kHiAy1v1eEMnsEDhBh4f/GLB3V4A/wYDU9Zx8rW/vxLdnqN\nSB5adn4j57Jz+wgt+Zr6aojXktZktJbZyiIdftZSnQC5Vi9c+FzIdYmUjqMm3lFffPFFhOsyUMNX\nLWaAzM/P11FRzzruYzhvuGzrOg7RIvRtb7o/a3GwsAXoIL1gwe8i/o4RIE0kQJriYwTImYP/4bVd\nbMPNKNYGDb5+D6HNmwv0Oef0twb2VZbw+I32q8XyNDxpDf6PaH/OqEh5mh63BqRELW/sq6zPcEsQ\nJerAN/shGtafQujYb9B2f5zb+uthw0aEvS7i0RTeOyotLS3A86esrEwvWLBAL1iwQJeVlenu3btH\nbAvfbzEDZF5enlYqOPlisPC384zt1OLdtlDLzCNUgHo8g8MKX6PCqt+PSWViaALsMrNl1ewXaOC9\n9tpreOed+Vx77dtIOvf/IBHpfRFD+qVIkFwHxFvnAyQPUjiiEAP6aiSNiddq80ckSFEhqTSC65kk\nI+m/pyPG+YGIgfd1JFp9OJLi4jbE2N4KqYj3XzIyXqnuwgShmTfvZf77329RWVmF1/sgx461QjzP\nYP78YVxwwRc1OI4dyFdd3fOGJ1zFPp/Px549e5CAx2AD+Q2IE2cCco0rkYqI30a8rPYT6Nkn7crL\nR5GVlcXYsWMDft/lcjnqmfcDICYml9TU6aa+Rx0wqUwMjUZg8R8f4etpDEK8lzoQFZXJuHHt+eST\nPVx++eX84Q9/oKioiBtuOEBp6VDEXfMpZIDOQgZ7Zy6lHyNR4HbBqeAaIT5gLhIxHVzqdhWSSytc\nPZMKxCtsLaLBfRHR1t6AeILtQ6oJYu37KY89dhtz584Me10qKytp23YEJSWB10KpwUi+qF1W/1YS\nmqbDzvkVnKpkEDAUpbz07Pl3li+f0aheRuEExerVrzBpUirl5aMBFx5PBlOmXMmHHx5n376OlJa2\nJ7S87RrkvL9AXgyGIGaklShVitbdkf+X0PuUnl4RIkBs6uLG29Iw9UCMCuuMwxlA5nbPsdQPTn32\nTt2x4xg9bdps3arVLZZKaY2llvqBnjZtrkMF8aCWSO8MHV5nvtZSedgBZ7Y3U38NXbTUmXAGIzrT\ndkeqV5Glxahu17uYaanNbglSi1RouE4PHpxco+A1vxFdvKM8nls13K/9Tga/0eHtLln6u9+9Sod6\neF2vY2JG6FWr1jS6aiacQ8CCBc9GSO8/UIvDg+29FkktaNdA8W/zeofoJUuW6NatB4W0q0t0/ZkG\njaDCMjMQQ6Njv/3t2bOHqVOjKCnpYm2RqnDR0dlUVaVTXr6e0Gp//8fkyTG88cYnfPaZj9A3b2dF\nQbsqXQKi+rgH+C7i9DcFKXeah6jAbF+QHyNv9C8iPiLB1QuHILOMd5C3/2G43Wu45pp08vJc+BM+\npnHzzefzt79l1Pi62G/tPp+PkpISpk9/H3n7tmNxDyKBc07WsGBBMXfeeSfXX389Wmuef/55zj//\n/CZ5sy4vL+fcc0dQXv4ygdetH3JvQvsvs7TuSFLKnch9cSN5xuyEjGuQrLoJAW3T0yto1aotkyev\nCIhHSU29nVGjBjXAGbYcGmMGYmwghkbHGUCm1KcEDgpQVfUJ5eXh6mIMBjqybNl9iO1jTph9+iFB\nZPFIFt4h+IMWb0DUWb+0fvNKJEPvMCSyvT8S7pSJqNZeQ2weI63jpyGC5gr8adR9XHLJX9i8+TVO\nnDjB4MGDAVi/fj3nnHNO2PMvLS3lrrvu4uOPP2b06NFMmDABt9uN2+2mR494Jk16kZ07L0RsLfb5\n2SnThxCo5lrFnDmSYfaf/wwXn9t4FBYWkZQ0m/LyKYTelwEEJrO08SGBgl2RjEnXIYLjZgLVdeHw\nceDAAR588EGGDx9QpzrvhtOkoac4p/vBqLDOWCJ5xHTsOFiHr/bndJO9XoeP38jSEnRoVz08VblV\nO1lhSpAqy+kdVmGprO7QkozQVpXZSRnH6EcfXVTjc5ZkknZVwzVavL+u1jNmzNZlZWWO6xEuPmW3\n1W/JQ6XU6ZWmrU/89zKS6i9Dh4/hGGDdszs0/FxHKtUr+wV7to3V0dHZJwMlDYFgVFhGhXWm408r\n0Q8Qj5glS6bQt+/DlJbmELkuRnskFjXYCH8rokZagxhbr0XSsAUbWTOQNCMgnjx2GpECpLDU0KD9\n1yHJE25BZi+FgA+v93GOH6/ZG29paSnR0cn41W5FwGKk+NEnwAaUGobWvyJUHWef320MGnQu3/9+\nD+bMmVPr2hYNRUFBAX37fsqJE4MI5xwhDgELgWXILFEj6WPs9PU+ZCboQ2aKW5EZpwtRbf0YWIrc\n20rg74jqsQKIIS5uFdu2hebTOpsxKizDGU+kpHqpqbczbtxAqqrGIi637yH6cKzlZ5ABZCqierIH\nJA+SE+llZDB6gvCqk9bI4PQWImBO9ojw1fXeQ9yGtyE2EICVvPTSpBqrSwKTSfoQtdpU/Akie6D1\ncqTQ1J/wJwbsg524MTm5C2vXhk862TxwI7Yqv+rP48lg/vxbSEt7geLivkgi7lzEgy7OaudCPLDe\nR+qu2PnGtiHVIK8mKqqU22/PZenSb1FZWY5cQw/wAjt2HCczM4cxY5xFqgwNjRHXhibHtokkJCSc\nfIMcNWoQpaWvMmnS35GKgdcjA4+zwNKFSLbaYiRb7aOIjcAuCuVGBp+/E1pM6lUks+2biK3D3u4C\nJiMzmXRkJjPSOoYd57EU0duP5uDButodCpE3bTtB5FD8tc8PIjaeHsDvOe+8xQwd+jLffPNKsxUe\n8fHxxMbm4nchzgbK6NhxGcePr+WBB2YQHR2F3JvLEMEZF3SUI8hM7K+Ibep+ZKa4BPDRqpWb6667\njqqq/yDCdxQy21kBtOeRR94462uUNzZGgBiaLZLx9QXee+95pAZ4R8SwasdtbEbeWDcjs4hfI//S\nscib+w5kALoUiftYgwQOjkcC/SqQjL4gginHOv4jyOCUj2TqHQQ8gAz6RxCh8wDgrtWAJVmFV1p9\n9yGqs0RCDc7jgUTatMkhLm42Gzc+T05ODm3atKnxbzU2doBer15306ZNDm3arCcubiuvvPI0Ho+H\nwsJC9u+/CSk5PBq5b85rV4qo8yYQej36AgXExLxn6d5vCbNPIocPx56cyRoaB6PCMjR7+va9luTk\nq1m37iEkgM+NvKkfRQb+t5GBqQ+iapqGqDbsGtudkQHrMyS/523ITMK2d8xBBq/9SNq2ddYvb0Xy\nfzprT8yyliUde8eOwYFvkfF6vcycGc8LLwxEvI22AN3C7tup0/+ybl2nsHUymis1qfEhuAhUz20D\ndiP1T8JRAcxnz54qfv5zLzITDEfDxswZQjFGdEOL4cSJE4wbN44vv/ySlJQU7rrrPEpLPyDU0HwX\n8D+EGtiTkYj1SchMAPyGapDZShSBwiXYmL4GUbNsAW6nY8ccDhxIq/Eg7/P56NXrLnbt+gHwXyTd\nSnA/B/LFF5mcf/75NTpmSyAwC4GzINdA5IXgVWtdOMeB8fjrtdgCPDgifxY9e2oKC//QYgRuQ2OM\n6AaDgzZt2pCTkwPIgPTUU6MpLnZWxMNadhFYXMpePxIZgHYis5C+iHrLznF1GAk8PBV2FcKdgIdD\nh75hy5Yt/O1vfwOo1jOqsLCQjz++GTHS341UOpyK5PJyAcvp0aPijBIeEC4HlebSS1/nvPO+Ij//\nHvz3yjkzARH4C/APVYmIwJiOeHPZwudr9u4tYc2a1876AMLGxIhqQ4vE5XLx0ENDUaq8li0PInmr\nRgNdEFfabyFG8b/gN7jHIYbgrfh19T5k5nED/gHtM/r0eZSHHurKQw91xesdzMyZD9TANmKrcZYh\nwqOY1q1/Q27ur9m9+51anlPLwFZxbdrUiU2bOrN37xIuu6x90F49kBnIJ0gtuwXI/XLSBUl6uQpJ\nyDkLeJny8teZPHkFlZWVDXsihpMYFZahxeLz+bjmmlns2BGszhgHfEVompNkoA3+9CD2+gmIYPlf\nJPliW+BcRL1ShRjNExCvremIN1gnRMgMIzTdyQR69jyPFStCExiGqnJ8QAGxsb+jqGjVWRdB/eGH\nW7nhhkcIvVdDEdfs4Gtrq7gKCIzfsTl1EsWzCaPCMhhOgcvlYvny6ZZapC8+n48LL1zPwoWDOHDg\nKPPn2/XLQVQhNyLVk8Ol2SgHornggpeoqKjk2DFn9t5RiJppMfLI7AbeABZa24KPl8TOnR2ZNOlF\nCgqeCdDJR0onvnz5vLNOePh8PmbOXAlchXi62fdqNeKVdQcyu7gWj8eNy5VFaen1SIxPOqHCw9DY\nmBmIocUTKTV3Wto6Jk9+gqqq64Dr8HiWAjMoLw+OSs9GYhNWExV1KUp1obQ0OCDNTswYjxh177TW\nf0zoQLYOKMXr/R82bbqJ3r1717jPZxOrV69l3LhX0ToJ8bT6IyLguyH1WWKRitl96NhxOWvWLGD/\n/oMsXPgKn346mLKyLIJnKNHRyRw7ln3WCeNwNMYMxAgQwxmNsy7FiBEjuOyy8Rw9mkGoCusYcDte\nbzHQOYwAWYNET29H0qP8AJnVVCCBiMGqstFIIae3yMqa3ai1OFoCPp+P7t2nUlz8EqEqqusQ1WF3\nxB17CXADXq+iW7f3Wbr050A5b7zxPgsXbqK0dCTgw+PJIDV1MmPGJDXBGTU/jADBCBBD/VFZWUl0\n9AAqK7sig5QHf4qU7kAyl1/eBpcriv37lxA4sCUh7r1RiNfPf5CYknLgayQ+JQpRvzgNvz569bo7\nRJV1tiO5sw5y4kRwevd1wHJO5dbrvJ5bt+4gOXk+R49eR6tWXeja9X1SU6cZgU3jCBDzH204a8jK\nyqKycgoy2P8D0bM/jXj+uIDRfPPNF4wdG4Pkc1qLxIYkIaVxxyGpM9KAixBD/UgkrmQFMM/67vQa\nclFc3M9ESNcQpcrp2NH+Voi4WgfamPbuvZHCwkJ8Ph9Tpy7j0KGXKS+fQ0nJMLZvH0dKyqPGE6uR\nMALEcJZxGElPch0QjRSZKrK2KT7//EY6dboEmXWsQ+I+xhPqSqoRz6FkZGbyJnA1rVs3+AmcIXjQ\nejXBOcp69tzMwoXjiI0djcfzFjLDC6SsrOKkDam4OBF/ZuO7gcMUFw+lR4/pFBYWhbQ11C9GgBjO\nGqqqvIhRdiUSaDgMUY+8iAiK92jVqgvdu3fnggsOIzMPe3bipBAJRAz2vhpG27bLCB4UY2PfIz4+\nviFOqUXi8/mYMmUJJSWPIIN+DrAOj2coJ058xfTpHo4cGc5FF21BvN2CE2G+EXxE5B7aSSlHUFz8\nEpMmvWiSKzYwxlXBcFZQWVnJ5MnPALMJHfj7INlff03XrktISEgmMfEScnI0kdO7h1ORVNKvXzt2\n7x7N4cPDUMpFTEwuqanTjf3DQUFBAXv3dkBmF79Hkl5CeXky+/d3RfKawYkTHREb1SzEBRtgEx5P\ngrV/OW3bZnHiRAfCJaW0VYd29UtD/WP+qw1nBenp6VRV3UL4hHtVREXFExf3Eqmp03C5XNx7711I\nWnHwp9dYh9hEViBG3uA340xee+1ajhwZTocOb/LSSz62bXvWGHQdFBYWMXbsi5SWXobkGrsHcdlN\nwF+jxSYBKSalrX0/BeDCC3NJTn6O668/wNGj1yElio3No0lo6JKHp/vBlLQ11APDhg3TkBm2XOo5\n5/TXeXl5uqqq6uT+VVVVum3bG60SsqutMq0jrL8V1vohVtnZdVZJ1p2O40uJV+cxz3YilTCWa1ah\nlRoXsk2plJB1Hs/goPK2ZVZp4MD9zvbrTyOUtDUzEMNZQffu3ZGAwZ/j1LtDEqNGdaF3794hEeOj\nR/dC8lS5rf1tu8dspODUKCRVfCckO+zVSIK/Qoz3VSiBRm8bUSHGxIwmJia4RQESoR64f3n5KGy1\nl8z8diHliVPweNZYdVRmnZxNGhoOc3UNZzyZmRt46qkiRI/+W0THXowIlOM8//zvw7abOHEMYrBN\nBjKRuJG1SBBhMaI2SbA+5lGqK16vYvXq+8jKmuMoSJVDTMxTREWdyq3N9rz6FKkr/zWDBr3Gpk2d\njOqwsWjoKc7pfjAqLMNpUFFRoaOjhzjUG1Ua8jT01XCLfuyx5yO2raqq0h06jLTUU+uszy/0BRf8\nUM+fv1C3bh2qNjEqrMhUVFTo2NjJp1Q1VVVV6fz8fJ2fn68rKirCqry83iGW2ir8toqKiiY+0+YB\nzUmFpZRyKaUKlVIbrO/JSqndSqkqpVTYEmFKqfZKqXeVUkVKqV1Kqbsc2+KUUputY+Yppb5/mrLQ\nYAghKyuLkhJnbRAXMgP5BcuXj2Tu3JkR27pcLtavn0fPnhqPZz8ez3569oS3317Gww8/yIMPXouo\nsbKRmckwpLTuWjp0GGdUKA4KCnbRo8dYDhxoh1ITUSoDr3dtiKrJ5XKRkJBAQkICbrc7qEyurZqa\nyAUX3Iqk1Xd6xhVSWtqbjIyMJjrLs4/auPHOQuaMba3vu4AhiAN2JCqBe7TW25VS5wAFSqm3tNZ7\ngd8A87XWbyml+iO6hZtqfQYGQx2pScK9+PgeFBb+IWziwwEDfsLDD+9B7CAuRIgUAg+zdu1DRoVi\nUVCwixtvnE9JyRjEC24j0IqOHXPYsmUF69ato6iokJSUlJB7EqlMbkXFl0yYYHvBFSHDUCJwOfPm\nbeDqqxPM9W8MajJNAdojhacTgQ1B2zYC19TwOK8AP7KWXweGW8ujgFUR2tTblM5w9hGqwhJVR3T0\n6as6qqqqdGzsuCAV10wdGzvOqK4s5BqFqq3gTu12P2B5VK3RsEZHRw/RGRnra3TMLVu2OLyxQlVZ\nRn3YvFRYTwP3IQ7ZdUIp1QnohSQhAnFleUopdQiZjfyqrsc2GCLhdrtZtmwi0dHJSEbdNXi9w1i2\nbOJpp/x2uVxkZc0hLg683o/xej+mZ0/IyppjVFcWhYWFHDr0M0I9r26gsnIX5eUvI/nFRlBSsq7a\nioKFhUUkJNzNzTd/hs93FaK0cKqy5PjGA65xqPYJUkoNBP6tRQ2VSPhIrOqOcQ7iMzlLa/21tXqG\n9f0VpVQykAr8JFz7hx9++ORyYmIiiYmJte2C4Sxm1KhBDB8+4GRa95SUtfVWLyI+vgfbtj171tf2\niERpaSkVFUuAd5H3UA+i5nsfSWoZOPCXlIwkKysrbEVBn8/HpEkvOqo5DkVKDH/cwGfRMsjNzSU3\nN7dRf7PadO5KqceRO12JZJ87F3hZaz3e2r4RuFdrvS1CezdSMOF1rfWzjvVfaq3Pd3z/Smt9Xpj2\nuro+GgyG5scdd8zjhRcKkWSUIDnIjiHm1D1IbfOal6SVFPCfcuLEUMdaH0pNROsVREr5frbSLEra\naq3nAnOtDvVDhMX4oN1O1clU4J9O4WHxmVKqn9b6PaXUjxDHeoPBcAZQWlpqCQ9nrfNkJDBwIFJt\ncD6S1NJZUXANKSnZtfglFx5PApddNpUjR/oDmPxjjUid5/FKqSQkL/Z3gNeUUtu11v2VUt8Dlmit\nb1VK2VnqdimlChEbylyt9RtIkelnlVKtkMIMU0/3ZAwGQ/Ng9uzZyMwj2PYxDpmBRAF2Nt5+gEap\nFSxZMiWiejE+Pp7Y2DS2bw9MbNm9+8ds3bqYHTt2WPs9a4RHI1ErAaK1fg9Jj4nW+hXEqyp4n38h\neR7QWn8ItIpwrA8BE/thMJx1/B+iuroaScEu9iOvdzzdunWI2MrlcpGaOo1Jk+6muLgf4J9tuN1u\nk3W3CTAlbQ0GQ71TWlpqeb45VVg+RIX1K+AoEnjpp02bHDZt6lStILCLSYFxWjgVzcIGYjAYDLXF\n6/Uyc2Y8L7wwCFFbAaykXbvjfP11HuXlG4GOiHCRYltSeGtItce2o9UNTY+ZgRgMhgbjxIkTjB8v\nPjcrV67E4/HQpctoDh8+D7gZ0XC/RuvWR9m8+XckJASXDjbUFTMDMRgMLZbCwiImTXqR4mJx1e3T\n5wHuued6jhzxIqlH/N5ZlZUT8flKAaOiakmYGYjBYKh3fD4fCQl3O4L+AHxcfPGPOHp0BhJ97kTi\nP3r0iLeETiIAsbG5pKZOM3mt6oCZgRgMhhZJpOJRR4/eDFSF7K9UFTExMUGR5rB9exKTJpmgwOaK\nuSMGg6FB8Pl8YdZegeRRDawn36HDG7hcrrBCx+S1ar6YGYjBYKh34uLiUOoxJF+V0433r8D9SDDh\njQAo9SqPPz7IzDBaIEaAGAyGemfHjh34fIn4I81BKj98DfRAasgXAj569jyfUaMkJiRcpHlN3XsN\njY8RIAaDoUFo1ao9cAd2pDn8gaio50LyVi1f7s9bFSnS3MxOmifGC8tgMNQ7kbywevW6m61bf+/I\nWxXqphvsxgsYt9460BheWEaAGAyGeqGystJRcyWFXbv2WS65/tnE8uXTw7rklpeXs2jRIgDmzJmD\nx+PB5/ORmbmOBQte5tChocCndOxYQEbGgybgsAYYAYIRIAZDSyAzcwOTJ6+gpCQFgOjoLJYtm8jI\nkbeGzB7sGYbtpZWamsPixbux64Yolc4vftGLDz74kp07r0fr/yA5XEcDiujo1bz//iNGiFSDESAY\nAWIwNHcqKytp21ZK0gbW9kjm2LHsgPTsdnT63r19KSurQOu/Ap8B7xCadHEdUsHwbiRrr9OwPpU9\ne14y6qxT0BgCxFx9g8FwWmRlZVkzj/DlaW2cJWlLS5PRehSQBrQPOqJdN+S3iAE+MeTYhw791MSG\nNAOMADEYDI1CpOh0qVBYO2FgZh7NA3MXDAbDaZGSkkJ0dBbB0eVSnjYlYN/w0emVIW0hHbgPSfWe\nG7I9NnbTSQ8tQ9Nh4kAMBsNp4Xa7WbZsIpMnJ1NSMhIArzeLZctuD7B/RI5Oz7SW/XVDkpMvZf/+\n+9i/P5Gqqi5oPRQYhdvdipiY90xsSDPBGNENBkO9EOzGG1zbvKCggD59PqSs7CP80em5QGfgTZQq\nweX6X1yuebRu7SUmZiP33deXbt26EBcXd8rYEUMoxgsLI0AMhjOFgoIC+vb9lBMnkvDbPOLxetex\neHEJjz/+PsXFL+GcncTFzWLbtmeNwKgDxgvLYDCcMcTHxxMbm2t9S7A+0K3bB1x55ZUcOTKAYAP7\nzp3XkZm5rlH7aag5RoAYDIZGweVykZo6jV697qZNmxzatMkhLm4WqanTrBlGqKZBaw+PPJITwfhu\naGqMCstgMDQq4UrW+nw+unefGqLCgqlAX1at8jJmTHAVQ8OpMDYQjAAxGM4WVq9ey7hxr6L1bUAF\nUniqN5BPhw6fc/DgBmMLqQVGgGAEiMFwtuDz+ejadQofffQN8EvERuICfCg1kX/840569+7dtJ1s\nQRgjusFgOGtwuVxMmNAVGIzMPOzhyYXWA9m3b1/Tdc4QFiNADAZDs+FnP7sZpUJfmpWqomvXrk3Q\nI8OpMALEYDA0GxISEoiJ2Uhw6pKYmFwSEhKaqluGCNRYgCilXEqpQqXUBut7slJqt1KqSil1TYQ2\n7ZVS7yqlipRSu5RSdwVtv1Mptcfa9uTpnYrBYGjpuFwusrJmERc3C683m6ioLGJiUnjooR/j8/ko\nKCigoKDAuPU2E2psRFdKzUasWm211oOUUl2R14QXgV9qrbeFaXMxcLHWertS6hygABistd6rlEoE\n5gIDtNaVSqnvaK3/G+YYxohuMJxlSDXCHBYseJMjR/rj82mUysDnS6RVq/bExuaSmjotbHVDg9Bs\njOhKqfbAAGCpvU5rvU9rvR+I2EGt9VGt9XZr+WtgD3CptXkG8KTWutLaHiI8DAbD2ctTT0lqkxMn\nhlFamkxJyTrKyj7ixIkktm9/hkmTXjQzkSampiqsp5HcynWeCiilOgG9gH9Yq2KBvkqpLUqpjUqp\n79f12AaD4cwicu2QfkgeLRfFxf1MUakmploBopQaCPzbmkkoTjHjOMUxzkHqU86yZiIgqeS/rbW+\nFrgfyK7tcQ0Gg8HQdNSkHkgfYJBSagAQDZyrlFqptR5fkx9QSrkR4ZGutV7v2HQYeBlAa71VKeVT\nSl2gtf6/4GM8/PDDJ5cTExNJTEysyU8bDIYWiiReTGP79iTkPbcSyABWAa8hRaXeIz5+SFN2s1mR\nm5tLbm5uo/5mrSLRlVL9gHu11oMc6zYiRvSCCG1WAv/VWt8TtH4qcKnWer5SKhZ4W2t9WZj2xohu\nMJyFFBYWMWnSi+zefS6VlbuBMdaW1XTo4GL9+oXGiH4Kml0qE6cAUUolAc8B3wG+BLZrrfsrpb4H\nLNFa36qU6gNsAnYh9hMNzNVav6GUag2kInaRMuu474X5TSNADIazlNLSUtq2HUlFxZ9xJln0epM5\nfjw7pGiVwU+zEyBNgREgBsPZQXCW3h079pCUNJtDh6YAwZl415CeXsHYsWMbvZ8thcYQIEZ8GwyG\nJsEpMMDDlClLLM8riIlZQUnJVxw6ZAREc8bMQAwGQ6Nj2zeKixPR2gdkUFKyDqeaSqmJaL0USEH8\ncPzboqOTOXbMqLBOhZmBGAyGMwZ7xuHz+ZgyJY2dOydgp2sXA3lgzIfWAxHz6UQgGRgJgMeTwbJl\nk43waAaYO2AwGOoVO2fVvn376Nq1KwkJCRQU7CI5+T6OHu2Oz9eKysqvEE9+kBCw5JAyUTcSAAAI\nTklEQVTjKFVlzU4GIYkwMujYcRX797+Gx+NptPMxRMYIEIPBUG8UFhaRkrKI/ftB61tR6iMuuuj3\nHD36H2AaMttYBWzAP+MYhBjJh+JUU8XE5BIdvZn9+w8BEBOzleXLnzbCoxlhbCAGg6Fe8Pl8XHPN\nLHbsAHiWwNrms6x1hcBBYFhQ62fweN7F7R4PKGJiclm+fDpxcd1D6qcbaoaxgRgMhhZDYWEhe/a0\nB2IIzWGViAgPCJcNSakLWbp0NFde2RmA+PhnTwoLUwek+WIEiMFgqBd8Ph/l5ZXV7BUPpAF2ihIA\nHz17bmbMmGfNDKOFYQSIwWCoR/YCRwkWEJCLX631c2AQXu8EXC7XSXWVER4tDyNADAZDveByuYiK\n6k1ZWT7iejsQqAJycLu/orJyLSBuuEuXTuLKKyX1nVNdZWhZGAFiMBjqhfj4eLp3T2P79lTE3rEP\niCEu7hLy8jLJzpaKDSkpOSaG4wzBeGEZDIZ6wx9h3g/gpHrKZM1tfEwyRYwAMRhaGsFJEY16qmkw\nAgQjQAwGg6EuNIYAMa8GBoPBYKgTRoAYDAaDoU4YAWIwGAyGOmEEiMFgMBjqhBEgBoPBYKgTRoAY\nDAaDoU4YAWIwGAyGOmEEiMFgMBjqhBEgBoPBYKgTRoAYDAaDoU4YAWIwGAyGOmEEiMFgMBjqhBEg\nBoPBYKgTRoAYDAaDoU7UWIAopVxKqUKl1Abre7JSardSqkopdU2ENu2VUu8qpYqUUruUUneF2ede\npZRPKdWu7qdhMBgMhsamNjOQWUCR4/suYAjw3inaVAL3aK17ANcBM5VS3eyNSqn2wE+AT2vRjxZF\nbm5uU3fhtGjJ/W/JfQfT/6ampfe/MaiRALEG+gHAUnud1nqf1no/ELFgidb6qNZ6u7X8NbAHuNSx\ny9PAfXXod4uhpf8TtuT+t+S+g+l/U9PS+98Y1HQGYg/0dS4NqJTqBPQC/mF9HwQc1lrvqusxDQaD\nwdB0VCtAlFIDgX9bMwnFKWYcpzjGOcA6YJbW+mulVDQwF5jv3K22xzUYDAZD01FtTXSl1OPAWMSe\nEQ2cC7ystR5vbd8I3Ku13hahvRt4DXhda/2ste4q4B3gBCI42gOfAT/QWv8nqL0piG4wGAx1oKFr\nolcrQAJ2VqofIiwGOdZtBH6ptS6I0GYl8F+t9T2nOO4B4Bqt9Rc17ozBYDAYmpQ6x4EopZKUUoeB\na4HXlFKvW+u/p5R6zVruA4wBbrZcgLcppW4JcziNUWEZDAZDi6JWMxCDwWAwGGyaJBJdKZVlzUa2\nKaUOKKW2Wet7WzMV+5MUof23lVJvKaX2KaXeVEqd59j2K6XUfqXUHqXUTxu5/z9WSuUrpXYopbYq\npW6K0L6nUurv1n7rLScDlFJupdQKpdROK/jygZbU/6Btu63tnpbUf2t7R6XUcaVURLVrc+x/Tds3\n1/5b2xr0+a2HvscppTZb41OeUqq3tb6lPLvB/f++Y1vtn12tdZN+gKeAX1vLXsBlLV8M/Nv+HtRm\nEXC/tTwHeNJavhIoBNxAJ+AjrFlWI/U/DrjYWu4BHInQJg+4wVqeCCywlkcBGdZyNHAA6NiC+t8K\n2AFcZX3/dku6/o7ta4E1SBBsY/7/n+71r1H7Ztz/Rn1+69j3N4GfWsv9gY3Wckt5diP1v07PboP+\nc9XwIhwCLg+zvjPwL8ILkL3ARdbyxcBea/kBYI5jv9eBHzZF/61t/wVah1n/hWO5PVBkLacA662b\neYF1nue3oP73B1Y2h/+fuvTf+j4YeUF5iMYRIPXa/5q0b679b+znt459fx0Ybi2PAlZZyy3l2Y3U\n/zo9u02aTFEpdSNwVGv9sWPdD5RSuxFpOF1r7QvT9EKt9b9Bot2BC631lwKHHft9RmDke70Srv+O\nbcnANq11RZimRUoCKQFGIA8RSKzMCURwHgSe0lp/We8d9/exvvsfa7V9w5pON2iWgfruv6VKuR94\nhEZw6miA61/T9vVCA/S/0Z7f0+j7bOAppdQh4DfAr6z1LeXZjdT/Oj277jr1vgYopd4GLnKuQryt\nHtRav2qtGwVkOttprfOAq5RSXYGVSqnXtdbl1fxcvXsC1LX/VtsewBNInq9wTAKeU0rNAzYA9vn9\nEIm3uRh5i3lfKfWO1vpgC+m/G+gDfB8oBf6mlMrXWm9sIf2fDzyttT6hlLJ/s040Uf9r2r5Z9/90\naeC+z0ACol+xBupUa9+W8uxG6n/dnt2GnGJVM/1qBRwFLjnFPn9D4kOC1+8hUIW1R4efAr9BA02B\nI/UfeZvaB1xbw+PEAFus5eeBMY5ty4DkFtT/kcByx7ZfI3FDLaX/m4BPrM8XiBrgFy2l/3Vp35z6\n31jP7+n0Hfgy3PeW8uyeov91enYb7B+sBhfhFiwDjmNdJ6CVtXwZcARoF6btIvsfjfBGdA9iQ2kw\nI1yE/p8HbAeSqmn7XeuvC0gDJljf7weWWcvfQrIfX9UC+j/R+n4+kI84Q7iBt4H+LaX/QfvMpwFt\nIA14/att34z73yjP72n2vQjoZy3/CNhqLbeUZzdS/+v07DboP1k1J7IcmBq0biywG9hmncxtjm1L\nsGYjQDskFco+4C0cxipEp/cRMkv5aSP3/0HguNX/Quvvd8L0/y6r73uBxx3tvwVkW9dgdwMPYPXe\nf2vbaKvvO4EnWlr/HcdqaAHSEP8/Edu3hP5b2xr8+T3NvvdBxqZCYDMQb61vKc9u2P5b22r97JpA\nQoPBYDDUCVPS1mAwGAx1wggQg8FgMNQJI0AMBoPBUCeMADEYDAZDnTACxGAwGAx1wggQg8FgMNQJ\nI0AMBoPBUCeMADEYDAZDnfh/1uk09dhqOL8AAAAASUVORK5CYII=\n", + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "zipped_pp.plot_points = plot_points\n", + "zipperita = Point(nparray[:,0], nparray[:,1])\n", + "\n", + "zipped_pp.plot_points(zipperita.x, zipperita.y, marker='o', linewidth=0)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": true + }, + "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "collapsed": true + }, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.5.1" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} diff --git a/point.py b/point.py new file mode 100644 index 0000000..17adb82 --- /dev/null +++ b/point.py @@ -0,0 +1,166 @@ +import analytics + +import math +import random + +import numpy as np +import scipy.spatial as ss +import pysal as ps + +class Point: + + def __init__(self, x = 0, y = 0, mark = ''): + self.x = x + self.y = y + self.mark = mark + + def __add__(self, other): + return Point(self.x + other.x, self.y + other.y) + + def __eq__(self, other): + return self.x == other.x and self.y == other.y + + def __radd__(self, other): + return Point(self.x + other, self.y + other) + + def check_coincident(self, b): + + if b.x == self.x and b.y == self.y: + return True + + def shift_point(self, x_shift, y_shift): + + self.x += x_shift + self.y += y_shift + + def to_array(self): + return [self.x, self.y] + +class PointPattern(object): + + def __init__(self): + self.points = [] + + def __len__(self): + count = 0 + for item in self.points: + count += 1 + return count + + def add_point(self, point): + self.points.append(point) + + def remove_point(self, index): + del(self.points[index]) + + def average_nearest_neighbor(self, mark=None): + return analytics.average_nearest_neighbor_distance(self.points, mark) + + def average_nearest_neighor_KDTree(self, mark=None): + nn_distances = [] + points_to_measure = [] + ps2m = [] + + if mark is not None: + points_to_measure = self.points_by_mark(mark) + else: + points_to_measure = self.points + + for p in points_to_measure: + ps2m.append((p.x, p.y)) + + treeKD = ss.KDTree(ps2m) + + for p in ps2m: + nearest_neighbor_d, nearest_neighbor = treeKD.query(p, k=2) + nn_distances.append(nearest_neighbor_d[1]) + + nn_distances = np.array(nn_distances) + return np.mean(nn_distances) + + def average_nearest_neighbor_ndarray(self, mark=None): + points_arrays = [] + + for item in self.points: + if item.mark is mark: + points_arrays.append(item.to_array()) + stacked_array = np.array(points_arrays) + + dists = [] + for num1, point in enumerate(stacked_array): + dists.append(None) + for num2, point2 in enumerate(stacked_array): + if num1 is not num2: + new_dist = ss.distance.euclidean(point, point2) + if dists[num1] == None: + dists[num1] = new_dist + elif dists[num1] > new_dist: + dists[num1] = new_dist + return np.mean(dists) + + def count_coincident(self): + counted = [] + count = 0 + for index, point in enumerate(self.points): + for index2, point2 in enumerate(self.points): + if index != index2: + if point.check_coincident(point2) is True: + count += 1 + return count + + def list_marks(self): + marks = [] + for point in self.points: + if point.mark not in marks and point.mark is not None: + marks.append(point.mark) + return marks + + def points_by_mark(self, mark=None): + points_to_return = [] + for point in self.points: + if point.mark == mark: + points_to_return.append(point) + return points_to_return + + def generate_random_points(self, n = None, marks = None): + if n is None: + n = len(self.points) + point_list = analytics.create_random_marked_points(n, marks) + return point_list + + def numpy_point_generator(self, low=0, high=1, n = 100): + nparray = np.random.uniform(low, high, (n, 2)) + point_list = [] + marks = ['black','white','ying','yang'] + for index in range(len(nparray)): + point_list.append(Point(nparray[index][0], + nparray[index][1], random.choice(marks))) + return point_list + + def generate_realizations(self, p = 99, marks = None): + neighbor_perms = [] + for i in range(p): + neighbor_perms.append( + analytics.average_nearest_neighbor_distance( + self.generate_random_points(100))) + return neighbor_perms + + def get_critical(neighbor_perms): + return max(neighbor_perms), min(neighbor_perms) + + def comupte_g(self, nsteps): + ds = np.linspace(0, 1, nsteps) + dist_counts = [] + for i, d in enumerate(ds): + min_dist = None + for n in range(nsteps): + if n != i: + if min_dist is None or min_dist > d: + min_dist = d + dist_counts.append(min_dist) + dist_counts = np.array(dist_counts) + return np.sum(dist_counts)/nsteps + + + + diff --git a/tests/__init__.py b/tests/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/tests/functional_test.py b/tests/functional_test.py new file mode 100644 index 0000000..16ee979 --- /dev/null +++ b/tests/functional_test.py @@ -0,0 +1,75 @@ +from .. import analytics + +import random +import unittest + +from .. import point + + +class TestFunctionalPointPattern(unittest.TestCase): + + def setUp(self): + random.seed(12345) + i = 0 + self.points = [] + marks = ['ying', 'yang', 'black', 'white'] + rand_marks = [] + for mark in range(100): + rand_marks.append(random.choice(marks)) + while i < 100: + seed = point.Point(round(random.random(),2), + round(random.random(),2), rand_marks[i]) + self.points.append(seed) + n_additional = random.randint(5,10) + i += 1 + c = random.choice([0,1]) + if c: + for j in range(n_additional): + x_offset = random.randint(0,10) / 100 + y_offset = random.randint(0,10) / 100 + pt = point.Point(round(seed.x + x_offset, 2), round(seed.y + + y_offset,2), random.choice(marks)) + self.points.append(pt) + i += 1 + if i == 100: + break + if i == 100: + break + + def test_point_pattern(self): + """ + This test checks that the code can compute an observed mean + nearest neighbor distance and then use Monte Carlo simulation to + generate some number of permutations. A permutation is the mean + nearest neighbor distance computed using a random realization of + the point process. + """ + random.seed() # Reset the random number generator using system time + # I do not know where you have moved + # avarege_nearest_neighbor_distance, so update the point_pattern module + observed_avg = analytics.average_nearest_neighbor_distance(self.points) + self.assertAlmostEqual(0.02779598180193161, observed_avg, 3) + + # Again, update the point_pattern module name for where you have placed the point_pattern module + # Also update the create_random function name for whatever you named the function to generate + # random points + marks = ['ying', 'yang', 'black', 'white'] + rand_points = analytics.create_random_marked_points(100, marks) + self.assertEqual(100, len(rand_points)) + + # As above, update the module and function name. + permutations = analytics.permutations(99) + self.assertEqual(len(permutations), 99) + self.assertNotEqual(permutations[0], permutations[1]) + + # As above, update the module and function name. + lower, upper = analytics.compute_critical(permutations) + self.assertTrue(lower > 0.03) + self.assertTrue(upper < 0.07) + self.assertTrue(observed_avg < lower or observed_avg > upper) + + # As above, update the module and function name. + significant = analytics.check_significant(lower, upper, observed_avg) + self.assertTrue(significant) + + self.assertTrue(True) diff --git a/tests/point_test.py b/tests/point_test.py new file mode 100644 index 0000000..a58ce5a --- /dev/null +++ b/tests/point_test.py @@ -0,0 +1,136 @@ +from .. import point +from .. import analytics +from .. import utils + + +import unittest +import random + +class Test_Point(unittest.TestCase): + + def test_points(self): + + """set up random x and y values for point instantiation""" + random.seed(12345) + rand_tup = (random.randint(0,10),random.randint(0,10)) + x_val = rand_tup[0] + y_val = rand_tup[1] + new_tup = (random.randint(0,10),random.randint(0,10)) + x_new = new_tup[0] + y_new = new_tup[1] + + """check that point is created properly""" + rand_point = point.Point(x_val, y_val, "random_mark") #(6,0,"random_mark") + self.assertEqual(rand_point.x, 6) + self.assertEqual(rand_point.y, 0) + + """check for whether the points are checking coincidence properly""" + new_point = point.Point(x_new, y_new, "different_mark") #(4,5,"different_mark") + self.assertTrue(rand_point.check_coincident(rand_point)) + self.assertFalse(rand_point.check_coincident(new_point)) + + """check whether points can be shifted""" + rand_point.shift_point(2,2) + self.assertEqual(rand_point.x, 8) + self.assertEqual(rand_point.y, 2) + + def test_marks(self): + + random.seed(12345) + marks = ['ying', 'yang', 'black', 'white'] + marked_points =[] + for i in range(20): + marked_points.append(point.Point(0,0, random.choice(marks))) + + ying_count = 0 + yang_count = 0 + white_count = 0 + black_count = 0 + for i in marked_points: + if i.mark == 'ying': + ying_count += 1 + if i.mark == 'yang': + yang_count += 1 + if i.mark == 'black': + black_count += 1 + else: + white_count += 1 + self.assertEqual(ying_count, 3) + self.assertEqual(yang_count, 7) + self.assertEqual(black_count, 6) + self.assertEqual(white_count, 14) + + def test_nearest_neighbor(self): + random.seed(12345) + marks = ['ying', 'yang', 'black', 'white'] + point_list = analytics.create_random_marked_points(20, marks) + + points_with_mark = analytics.average_nearest_neighbor_distance(point_list, + marks) + points_without_mark = analytics.average_nearest_neighbor_distance(point_list) + + self.assertNotEqual(points_with_mark, 0.2, 5) + self.assertAlmostEqual(points_without_mark, 0.11982627009007044, 1) + + + def test_magic_methods(self): + """check whether add and eq magic methods work""" + new_point = point.Point(2,4) + point_to_add = point.Point(3,6) + added = new_point + point_to_add + #This asserts that the add method worked as well as eq method + self.assertEqual(added, point.Point(5,10)) + + """check whether reverse adding works on random points""" + point_list = analytics.create_random_marked_points(20) + self.assertTrue(type(random.choice(point_list) + + random.choice(point_list)), point.Point) + +class Point_Pattern_Test(unittest.TestCase): + + def setUp(self): + self.point_pattern = point.PointPattern() + self.point_pattern.add_point(point.Point(3, 6, 'ying')) + self.point_pattern.add_point(point.Point(9, 6, 'yang')) + self.point_pattern.add_point(point.Point(3, 9, 'black')) + self.point_pattern.add_point(point.Point(9, 6)) + + def test_remove_point(self): + length_after_removal = len(self.point_pattern) - 1 + self.point_pattern.remove_point(0) + self.assertEqual(length_after_removal, len(self.point_pattern)) + + def test_coincident(self): + self.assertEqual(self.point_pattern.count_coincident(), 2) + + def test_list_marks(self): + self.assertEqual(len(self.point_pattern.list_marks()), 4) + + def test_points_by_mark(self): + self.assertEqual(len(self.point_pattern.points_by_mark('yang')), 1) + + def test_point_generation(self): + self.assertEqual(len(self.point_pattern.generate_random_points(10)), 10) + self.assertEqual(len(self.point_pattern.generate_random_points()), 4) + + def test_g(self): + self.assertEqual(self.point_pattern.comupte_g(100), 0.5) + + def test_nnd_kdtree(self): + self.assertTrue(True) + self.point_pattern.points = self.point_pattern.generate_random_points(100, + marks=['black','ying','yang']) + kd_result = self.point_pattern.average_nearest_neighor_KDTree(mark='black') + self.assertAlmostEqual(kd_result, 0, 0) + + def test_nnd_ndarray(self): + self.point_pattern.points = self.point_pattern.generate_random_points(100,['black', 'ying']) + nd_result = self.point_pattern.average_nearest_neighbor_ndarray(mark = + 'black') + self.assertAlmostEqual(nd_result, 0, 0) + + def test_pont_generator_with_numpy(self): + self.point_pattern = point.PointPattern() + self.point_pattern.points = self.point_pattern.numpy_point_generator(10,20) + self.assertEqual(len(self.point_pattern.points), 100) + diff --git a/utils.py b/utils.py new file mode 100644 index 0000000..541c803 --- /dev/null +++ b/utils.py @@ -0,0 +1,114 @@ + +import random + +def shift_point(point, x_shift, y_shift): + """ + Shift a point by some amount in the x and y directions + + Parameters + ---------- + point : tuple + in the form (x,y) + + x_shift : int or float + distance to shift in the x direction + + y_shift : int or float + distance to shift in the y direction + + Returns + ------- + new_x : int or float + shited x coordinate + + new_y : int or float + shifted y coordinate + + Note that the new_x new_y elements are returned as a tuple + + Example + ------- + >>> point = (0,0) + >>> shift_point(point, 1, 2) + (1,2) + """ + x = getx(point) + y = gety(point) + + x += x_shift + y += y_shift + + return x, y + + +def check_coincident(a, b): + """ + Check whether two points are coincident + Parameters + ---------- + a : tuple + A point in the form (x,y) + + b : tuple + A point in the form (x,y) + + Returns + ------- + equal : bool + Whether the points are equal + """ + return a == b + +def check_in(point, point_list): + """ + Check whether point is in the point list + + Parameters + ---------- + point : tuple + In the form (x,y) + + point_list : list + in the form [point, point_1, point_2, ..., point_n] + """ + return point in point_list + + +def getx(point): + """ + A simple method to return the x coordinate of + an tuple in the form(x,y). We will look at + sequences in a coming lesson. + + Parameters + ---------- + point : tuple + in the form (x,y) + + Returns + ------- + : int or float + x coordinate + """ + return point.x + + +def gety(point): + """ + A simple method to return the x coordinate of + an tuple in the form(x,y). We will look at + sequences in a coming lesson. + + Parameters + ---------- + point : tuple + in the form (x,y) + + Returns + ------- + : int or float + y coordinate + """ + return point.y + + diff --git a/view.py b/view.py new file mode 100644 index 0000000..9ec4b27 --- /dev/null +++ b/view.py @@ -0,0 +1,36 @@ +import sys +from PyQt4 import QtGui, QtCore + +class GUI(QtGui.QMainWindow): + + def __init__(self): + super(GUI, self).__init__() + + self.initUI() + + def initUI(self): + + exitAction = QtGui.QAction(QtGui.QIcon('exit.png'), '&Exit', + self) + exitAction.setShortcut('Ctrl+Q') + exitAction.setStatusTip('Exit application') + exitAction.triggered.connect(QtGui.qApp.quit) + + self.statusBar() + + menubar = self.menuBar() + fileMenu = menubar.addMenu('&File') + fileMenu.addAction(exitAction) + + self.setGeometry(300,300,300,200) + self.setWindowTitle('Menubar') + self.show() + +def main(): + + app = QtGui.QApplication(sys.argv) + ex = GUI() + sys.exit(app.exec_()) + +if __name__ == '__main__': + main()