X-Git-Url: https://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=src%2FGeometry.py;fp=src%2FGeometry.py;h=89a6cdb929e87bc787574d62396c0fdea21df372;hb=27ce45a6f4e2844f3d7bc5d570f3cb06e087bb5a;hp=0000000000000000000000000000000000000000;hpb=39329a27446673a022c7b0525af916ad94fa70e6;p=dead%2Fcensus-tools.git diff --git a/src/Geometry.py b/src/Geometry.py new file mode 100644 index 0000000..89a6cdb --- /dev/null +++ b/src/Geometry.py @@ -0,0 +1,262 @@ +""" +Geometry classes and functions. We don't want to rely completely on +some huge third-party library, nor do we want implement everything +from scratch. We use Shapely for the hard stuff (Geos integration), +and wrap its classes with our own. +""" + +from math import sqrt +import os +import site +import sys + +# Basically, add '../src' and '../lib/Shapely' to our path. +# Needed for the imports that follow. +site.addsitedir(os.path.dirname(os.path.abspath(sys.argv[0])) + '/../src') +site.addsitedir(os.path.dirname(os.path.abspath(sys.argv[0])) + '/../lib/Shapely') + +import shapely.geometry +import shapely.wkt + + + +class TwoVector: + """ + Represents a two-dimensional vector. + """ + + def __init__(self, x=0, y=0): + self.x = x + self.y = y + + + def __hash__(self): + """ + We have to define __hash__ in order to override __eq__. + Just return the hash of our ordered pair. + """ + return (self.x, self.y).__hash__() + + + def __eq__(self, v): + """ + We would like to override equality to compare against our + underlying ordered pair. + """ + if (v == None): + return False + + if isinstance (v, TwoVector): + return ( (self.x == v.x) and (self.y == v.y) ) + + # If it's not a TwoVector, assume that it's a tuple. + if (len(v) > 1): + return (self.x == v[0]) and (self.y == v[1]) + + # If v is neither a TwoVector nor a tuple, we ain't it. + return False + + + def __str__(self): + """ + Print the contents of our ordered pair. + """ + return "TwoVector<%d, %d>" % (self.x, self.y) + + + def to_tuple(self): + """ + Returns the two-tuple underlying this TwoVector. + """ + return (self.x, self.y) + + + def perpendicular(self): + """ + Return a vector perpendicular to ourself. + """ + return TwoVector(-self.y, self.x) + + + def dot_product(self, other_vector): + """ + Return the dot product of this vector with other_vector. + """ + return (self.x * other_vector.x) + (self.y * other_vector.y) + + + def add(self, other_vector): + """ + Add other_vector to this vector, and return a new TwoVector + containing the result. + """ + return TwoVector(self.x + other_vector.x, self.y + other_vector.y) + + + +class Polygon: + """ + Wraps shapely.geometry.Polygon. + """ + + def __init__(self, points=None): + """ + We can initialize from either a set of two-tuples, or a set of + TwoVectors. + """ + if (points == None): + return + + # Is the first element of points a TwoVector? + if ( (len(points) > 0) and isinstance(points[0], TwoVector) ): + # If so, convert points to a list of two-tuples. + points = map ( (lambda p: p.to_tuple()), points ) + + self._shapely_polygon = shapely.geometry.Polygon(points) + + + def __eq__(self, p): + """ + Override equality to test whether or not our boundary is equal + to the passed polygon's boundary. + """ + if (p == None): + return False + + # Shapely polygon exteriors know how to compare themselves. + ext1 = self._shapely_polygon.exterior + + if isinstance(p, Polygon): + ext2 = p._shapely_polygon.exterior + return ext1.equals(ext2) + elif isinstance(p, shapely.geometry.Polygon): + return ext1.equals(p) + else: + return False + + + def __hash__(self): + """ + Needed when we implement __eq__. Delegate it to our + shapely.geometry.Polygon. + """ + return self._shapely_polygon.__hash__() + + + def __str__(self): + """ + Write out the well-known text representation of this polygon. + """ + return "Geometry.Polygon<%s>" % self._shapely_polygon.to_wkt() + + + @classmethod + def from_shapely(self, shapely_polygon): + """ + Create a Geometry.Polygon from the shapely.geometry.Polygon. + """ + p = Polygon() + p._shapely_polygon = shapely_polygon + return p + + + def coords(self): + """ + An alias to retrieve our coordinates. + """ + tuples = self._shapely_polygon.exterior.coords + vectors = map ( (lambda v: TwoVector(v[0], v[1])), tuples) + return vectors + + + def union(self, other_polygon): + """ + Return a new Polygon representing the union of this polygon + and other_polygon. + """ + # First, get the union of self and other_polygon as a + # shapely.geometry.Polygon. + u = self._shapely_polygon.union(other_polygon._shapely_polygon) + + # Then, create a new Geometry.Polygon from that union. + return Polygon.from_shapely(u) + + + def translate(self, v): + """ + Translate the self polygon by the TwoVector v. + """ + regular_coords = self.coords() + f = (lambda p: p.add(v)) + translated_coords = map(f, regular_coords) + return Polygon(translated_coords) + + + + def drag_vertices(self, v): + """ + The vector v is that vector along which we are 'dragging' the + polygon. If we are dragging from p1 = (x1,y1) to p2 = (x2,y2), + then v = (x2-x1, y2-y1). The 'drag' vertices of a polygon are + those two vertices which are most/least in the direction of + the v_perp. + """ + v_perp = v.perpendicular() + + current_min_vertex = self.coords()[0] + current_max_vertex = self.coords()[0] + + for w in self.coords(): + if (w.dot_product(v_perp) < current_min_vertex.dot_product(v_perp)): + current_min_vertex = w + + if (w.dot_product(v_perp) > current_max_vertex.dot_product(v_perp)): + current_max_vertex = w + + return [current_min_vertex, current_max_vertex] + + + @staticmethod + def from_wkt(text): + """ + Create a Polygon from a Well-Known Text string. + """ + # The WKT string must contain the word "polygon" at least. + # Ex: POLYGON((1 1,5 1,5 5,1 5,1 1),(2 2, 3 2, 3 3, 2 3,2 2)) + if (text.lower().find('polygon') == -1): + return None + + return shapely.wkt.loads(text) + + + def drag_rectangle(self, v): + """ + When we dtag a polygon in the plane, we drag it from its + starting position to some terminal position. Since there are + two 'drag vertices' on the polygon, we have two points + corresponding to the drag vertices at both the initial + location and the terminal location. From these four points, we + can form a rectangle. + """ + initial_dv = self.drag_vertices(v) + + terminal_dv = [0, 0] # Just initializing it. + terminal_dv[0] = initial_dv[0].add(v) + terminal_dv[1] = initial_dv[1].add(v) + + return Polygon( (initial_dv[1], + initial_dv[0], + terminal_dv[0], + terminal_dv[1]) ) + + + def drag(self, v): + """ + Drag this polygon along the line v = (x2-x1, y2-y1), and + return the new polygon consisting of the union of all points + that I've contained along the way. + """ + terminal_p = self.translate(v) + dv_rect = self.drag_rectangle(v) + + return self.union(dv_rect).union(terminal_p)