ESDA Tabular Functions

.travis.yml

@@ -18,7 +18,7 @@ before_install:
   - conda create -y -q -n test-env python=$TRAVIS_PYTHON_VERSION
   - source activate test-env
   - chmod +x ./.travis_testing.sh
-  - if [[ $TRAVIS_PYTHON_VERSION == 3* ]]; then rm -rf pysal/contrib; 2to3 -nw pysal/ > /dev/null; fi
+  - if [[ $TRAVIS_PYTHON_VERSION == 3* ]]; then 2to3 -nw pysal/ > /dev/null; fi
 
 install:
   - conda install --yes numpy scipy nose pip

.travis_testing.sh

@@ -1,5 +1,5 @@
 if [[ $TRAVIS_PYTHON_VERSION == 3* ]]; then
-    cd pysal; python -m 'nose';
+    nosetests --with-coverage --cover-package=pysal --exclude-dir=pysal/contrib
 else
     nosetests --with-coverage --cover-package=pysal; 
 fi

pysal/__init__.py

@@ -40,6 +40,19 @@
 import pysal.cg
 import pysal.core
 
+try:
+    import pandas
+    from pysal.contrib import pdio
+    pysal.common.pandas = pandas
+except ImportError:
+    print('Pandas adapters not loaded')
+    pysal.common.pandas = None
+
+# Load the IOHandlers
+from pysal.core import IOHandlers
+# Assign pysal.open to dispatcher
+open = pysal.core.FileIO.FileIO
+
 from pysal.version import version
 #from pysal.version import stable_release_date
 #import urllib2, json
@@ -84,19 +97,6 @@
 import pysal.examples
 from pysal.network.network import Network, NetworkG, NetworkK, NetworkF
 
-try:
-    import pandas
-    from pysal.contrib import pdio
-    pysal.common.pandas = pandas
-except ImportError:
-    print('Pandas adapters not loaded')
-    pysal.common.pandas = None
-
-# Load the IOHandlers
-from pysal.core import IOHandlers
-# Assign pysal.open to dispatcher
-open = pysal.core.FileIO.FileIO
-
 #__all__=[]
 #import esda,weights
 #__all__+=esda.__all__

pysal/cg/comparators.py

@@ -0,0 +1,153 @@
+import numpy as np
+from pysal.common import RTOL
+from .shapes import Point, Chain, Polygon, Geometry
+
+#### This approach is invalid :(
+#### OGC has a different definition of equality, so we can't be ogc-compliant
+#### by focusing on vertex set relationships. For example:
+
+# from shapely import geometry as geom
+# p1 = geom.Polygon([(0,0),(0,1),(1,0)])
+# p2 = geom.Polygon([(0,0),(0,1),(.5,.5),(1,0)])
+# p1.equals(p2) #is true
+# no relationship on point set (without simplifying boundaries) will make these
+# approaches work. 
+
+ATOL = 1e-12
+
+def _point_equal(a,b):
+    """
+    Test that a point is exactly equal to another point, using numpy array
+    comparison
+    """
+    return np.array_equal(a._Point__loc, b._Point__loc)
+
+def _point_almost_equal(a,b, rtol=RTOL, atol=ATOL):
+    """
+    test that one point is equal to another point, up to a specified relative or
+    absolute tolerance.
+    """
+    return np.allclose(a._Point__loc, b._Point__loc, 
+                       rtol=rtol, atol=atol)
+
+def _chain_equal(a,b):
+    """
+    Test that two chains are equal. This considers reversed-incident chains,
+    chains with the same pointset but in reverse order, as different chains. 
+    """
+    for a_part, b_part in zip(a.parts, b.parts):
+        for a_seg, b_seg in zip(a_part, b_part):
+            if not np.array_equal(a_seg, b_seg):
+                return False
+    return True 
+
+def _chain_almost_equal(a,b, rtol=RTOL, atol=ATOL):
+    """
+    Test that two chains are equal, up to a specified relative or absolute
+    tolerance. This considers reversed-incident chains,
+    chains with the same pointset but in reverse order, as different chains. 
+    """
+    for a_part, b_part in zip(a.parts, b.parts):
+        for a_seg, b_seg in zip(a_part, b_part):
+            if not np.allclose(a_seg, b_seg, 
+                               rtol=RTOL, atol=ATOL):
+                return False
+    return True
+
+def _poly_exactly_equal(a,b):
+    """
+    Test that two polygons are equal. This is a straightfoward linear comparison
+    of parts and holes. That is, it is assumed that the parts and holes are
+    conformal.
+
+    Thus, shapes with the same parts or holes, but listed in a different order,
+    will be considered different. Solving this will require some way to relate
+    ring/hole sets
+    """
+    for a_hole, b_hole in zip(a.holes, b.holes):
+        if not np.array_equal(a_hole, b_hole):
+            return False
+    for a_part, b_part in zip(a.parts, b.parts):
+        if not np.array_equal(a_part, b_part):
+            return False
+    return True
+
+def _poly_exactly_coincident(a,b):
+    """
+    Check that two polygons have coincident boundaries 
+
+    Thus, this returns True when two polygons have the same holes & parts in any
+    order with potentially-reversed path directions
+    """
+    n_holes = len(a.holes)
+    n_parts = len(a.parts)
+    if n_holes != len(b.holes):
+        return False
+    if n_parts != len(b.parts):
+        return False
+    b_in_a = [None]*n_holes
+    a_in_b = [None]*n_holes
+    for i, a_hole in a.holes:
+        for j, b_hole in b.holes:
+            i_j = coincident(a_hole, b_hole)
+            if i_j:
+                b_in_a[j] = True
+                a_in_b[i] = True
+                break
+        if not a_in_b[i]:
+            return False
+    if any(in_b):
+        return False
+    for b_hole in b.holes:
+        in_a = [not coincident(b_hole, a_hole) for a_hole in a.holes]
+    if any(in_a):
+        return False
+    return True
+
+def _coincident(a,b):
+    """
+    Check if two vertex lists are equal, either forwards or backwards. Thus,
+    this checks for equality of the path or equality when one path is reversed. 
+    """
+    return np.array_equal(a, b) or np.array_equal(np.flipud(a),b)
+
+def _almost_coincident(a,b, rtol=RTOL, atol=ATOL):
+    """
+    Check if two vertex lists are equal to within a given relative and absolute
+    tolernance, either forwards or backwards. Thus,
+    this checks for equality of the path or equality when one path is reversed. 
+    """
+    return (np.allclose(a, b, rtol=RTOL, atol=ATOL) 
+            or np.allclose(np.flipud(a),b, rtol=RTOL, atol=ATOL))
+
+def _poly_almost_equal(a,b, rtol=RTOL, atol=ATOL):
+    is_equal = True
+    for a_hole, b_hole in zip(a.holes, b.holes):
+        is_equal &= np.allclose(a_hole, b_hole,
+                                rtol=rtol, atol=atol)
+    for a_part, b_part in zip(a.parts, b.parts):
+        is_equal &= np.allclose(a_part, b_part, 
+                                rtol=rtol, atol=atol)
+    return is_equal
+
+_EQ_MAP = {Point:_point_equal, 
+           Polygon:_poly_equal,
+           Chain: _chain_equal}
+_AEQ_MAP = {Point:_point_almost_equal, 
+            Polygon:_poly_almost_equal,
+            Chain: _chain_almost_equal}
+
+def equal(a,b):
+    ta, tb = type(a), type(b)
+    if not all([isinstance(a, tb), isinstance(b, ta)]):
+        return False
+    return _EQ_MAP[ta](a,b)
+
+def almost_equal(a,b, rtol=RTOL, atol=ATOL):
+    ta, tb = type(a), type(b)
+    if not all([isinstance(a, tb), isinstance(b, ta)]):
+        return False
+    return _AEQ_MAP[ta](a,b)
+
+def is_shape(a):
+    return isinstance(a, Geometry)

pysal/cg/shapes.py

@@ -10,6 +10,7 @@
 import math
 from warnings import warn
 from sphere import arcdist
+import numpy as np
 
 __all__ = ['Point', 'LineSegment', 'Line', 'Ray', 'Chain', 'Polygon',
            'Rectangle', 'asShape']
@@ -37,8 +38,15 @@ def asShape(obj):
         raise NotImplementedError(
             "%s is not supported at this time." % geo_type)
 
+class Geometry(object):
+    """
+    A base class to help implement is_geometry and make geometric types
+    extendable.
+    """
+    def __init__(self):
+        pass
 
-class Point(object):
+class Point(Geometry):
     """
     Geometric class for point objects.
 
@@ -345,7 +353,7 @@ def __str__(self):
         return "POINT ({} {})".format(*self.__loc)
 
 
-class LineSegment(object):
+class LineSegment(Geometry):
     """
     Geometric representation of line segment objects.
 
@@ -776,7 +784,7 @@ def line(self):
         return self._line
 
 
-class VerticalLine:
+class VerticalLine(Geometry):
     """
     Geometric representation of verticle line objects.
 
@@ -853,7 +861,7 @@ def y(self, x):
         return float('nan')
 
 
-class Line:
+class Line(Geometry):
     """
     Geometric representation of line objects.
 
@@ -982,7 +990,7 @@ def __init__(self, origin, second_p):
         self.p = second_p
 
 
-class Chain(object):
+class Chain(Geometry):
     """
     Geometric representation of a chain, also known as a polyline.
 
@@ -1018,7 +1026,7 @@ def __init__(self, vertices):
         else:
             self._vertices = [vertices]
         self._reset_props()
-
+    
     @classmethod
     def __from_geo_interface__(cls, geo):
         if geo['type'].lower() == 'linestring':
@@ -1180,7 +1188,7 @@ def segments(self):
         return [[LineSegment(a, b) for (a, b) in zip(part[:-1], part[1:])] for part in self._vertices]
 
 
-class Ring(object):
+class Ring(Geometry):
     """
     Geometric representation of a Linear Ring
 
@@ -1385,7 +1393,7 @@ def contains_point(self, point):
 
 
 
-class Polygon(object):
+class Polygon(Geometry):
     """
     Geometric representation of polygon objects.
 
@@ -1756,7 +1764,7 @@ def contains_point(self, point):
 
 
 
-class Rectangle:
+class Rectangle(Geometry):
     """
     Geometric representation of rectangle objects.
 

pysal/contrib/geotable/ops/_accessors.py

@@ -1,7 +1,7 @@
 import functools as _f
 
 __all__ = [ 'area', 'bbox', 'bounding_box', 'centroid', 'holes', 'len', 
-               'parts', 'perimeter', 'segments', 'vertices']
+            'parts', 'perimeter', 'segments', 'vertices']
 
 def get_attr(df, geom_col='geometry', inplace=False, attr=None):
     outval = df[geom_col].apply(lambda x: x.__getattribute__(attr))