started tests (#37)

* started tests

* updated travis stuff

* updated travis stuff

* python 3 error

* skip hexarray tests

* updated tests

* fixed tests

* fixed py3 error

.travis.yml

@@ -20,7 +20,8 @@ addons:
       - g++-4.8
 
 env:
-  - CONDA_PACKAGES="numpy scipy cython pandas bcolz pytest future xarray dask"
+  - CONDA_PACKAGES="numpy scipy cython pandas bcolz pytest future xarray dask scikit-image"
+    PIP_PACKAGES="tqdm codecov pytest-cov"
 # need to use include
 #  - CC="gcc-4.8"
 
@@ -43,7 +44,7 @@ install:
   - conda create -q -n test-environment python=$TRAVIS_PYTHON_VERSION $CONDA_PACKAGES
   - source activate test-environment
   - gcc -v
-  - pip install codecov pytest-cov
+  - pip install $PIP_PACKAGES
   - CC=gcc-4.8 CXX=g++-4.8 python setup.py build_ext --inplace
   - pip install -e .
 

floater/rclv.py

@@ -1,3 +1,5 @@
+from __future__ import print_function
+
 import numpy as np
 import xarray as xr
 from skimage.measure import find_contours, points_in_poly, grid_points_in_poly
@@ -53,6 +55,38 @@ def point_in_contour(con, ji):
     return points_in_poly(np.array([i, j])[None], con[:,::-1])[0]
 
 
+def contour_area(con):
+    """Calculate the area, convex hull area, and convexity deficiency
+    of a polygon contour.
+
+    Parameters
+    ----------
+    con : arraylike
+        A 2D array of vertices with shape (N,2) that follows the scikit
+        image conventions (con[:,0] are j indices)
+
+    Returns
+    -------
+    region_area : float
+    hull_area : float
+    convexity_deficiency : float
+    """
+    # reshape the data to x, y order
+    con_points = con[:,::-1]
+
+    # calculate area of polygon
+    region_area = abs(polygon_area(con_points))
+
+    # find convex hull
+    hull = qhull.ConvexHull(con_points)
+    #hull_points = np.array([con_points[pt] for pt in hull.vertices])
+    hull_area = hull.volume
+
+    cd = (hull_area - region_area ) / region_area
+
+    return region_area, hull_area, cd
+
+
 def find_contour_around_maximum(data, ji, level, border_j=(5,5),
         border_i=(5,5), max_footprint=None):
     j,i = ji
@@ -65,7 +99,7 @@ def find_contour_around_maximum(data, ji, level, border_j=(5,5),
     grow_down, grow_up, grow_left, grow_right = 4*(False,)
 
     while target_con is None:
-        
+
         footprint_area = sum(border_j) * sum(border_i)
         if max_footprint and footprint_area > max_footprint:
             raise ValueError('Footprint exceeded max_footprint')
@@ -112,38 +146,6 @@ def find_contour_around_maximum(data, ji, level, border_j=(5,5),
     return target_con, region_data, border_j, border_i
 
 
-def contour_area(con):
-    """Calculate the area, convex hull area, and convexity deficiency
-    of a polygon contour.
-
-    Parameters
-    ----------
-    con : arraylike
-        A 2D array of vertices with shape (N,2) that follows the scikit
-        image conventions (con[:,0] are j indices)
-
-    Returns
-    -------
-    region_area : float
-    hull_area : float
-    convexity_deficiency : fload
-    """
-    # reshape the data to x, y order
-    con_points = con[:,::-1]
-
-    # calculate area of polygon
-    region_area = abs(polygon_area(con_points))
-
-    # find convex hull
-    hull = qhull.ConvexHull(con_points)
-    #hull_points = np.array([con_points[pt] for pt in hull.vertices])
-    hull_area = hull.volume
-
-    cd = (hull_area - region_area ) / region_area
-
-    return region_area, hull_area, cd
-
-
 def convex_contour_around_maximum(data, ji, step, border=5,
                                   convex_def=0.01, verbose=False,
                                   max_footprint=None):
@@ -195,7 +197,7 @@ def convex_contour_around_maximum(data, ji, step, border=5,
 
     for level in contour_levels:
         if verbose:
-            print ('  level: %g border: ' % level) + repr(border_j) + repr(border_i)
+            print('  level: %g border: ' % level) + repr(border_j) + repr(border_i)
 
         try:
             # try to get a contour
@@ -204,7 +206,7 @@ def convex_contour_around_maximum(data, ji, step, border=5,
                 max_footprint=max_footprint)
         except ValueError as ve:
             if verbose:
-                print ve
+                print(ve)
             break
 
         # get the convexity deficiency
@@ -222,7 +224,7 @@ def convex_contour_around_maximum(data, ji, step, border=5,
             # re-center the previous contour to be referenced to the
             # absolute position
             if verbose:
-                print("  moving on to next contour level, region_data.shape: " + 
+                print("  moving on to next contour level, region_data.shape: " +
                         repr(region_data.shape))
             contour[:, 0] += (j-border_j[0])
             contour[:, 1] += (i-border_i[0])
@@ -274,8 +276,12 @@ def find_convex_contours(data, min_distance=5, min_area=100.,
         pool = ThreadPool()
         map_function = pool.imap_unordered
     else:
-        from itertools import imap
-        map_function = imap
+        try:
+            from itertools import imap
+            map_function = imap
+        except ImportError:
+            # must be python 3
+            map_function = map
 
     plm = peak_local_max(data, min_distance=min_distance)
 
@@ -334,6 +340,6 @@ def fill_in_contour(contour, value=1):
                                                        contour_rel)
 
     for n, con in enumerate(contours):
-        fill_in_contour(con, n)
+        fill_in_contour(con, n+1)
 
     return data

floater/test/test_hexarray.py

@@ -1,9 +1,13 @@
 from __future__ import print_function
 
+import pytest
+pytestmark = pytest.mark.skip(reason="hexarray module is being deprecated")
+
 import unittest
 import numpy as np
 from floater import hexgrid
 
+
 def _get_test_array(shape=(3,3), nonzero=[(1,1)], dtype=np.float64):
     a = np.zeros(shape, dtype)
     for (j,i) in nonzero:

floater/test/test_rclv.py

@@ -0,0 +1,122 @@
+import numpy as np
+import pytest
+from floater import rclv
+
+@pytest.fixture()
+def sample_data_and_maximum():
+    ny, nx = 100, 100
+    x, y = np.meshgrid(np.arange(nx), np.arange(ny))
+    x = (x-5)*2*np.pi/80
+    y = (y - nx/2)*2*np.pi/100
+    # Kelvin's cats eyes flow
+    a = 0.8
+    psi = np.log(np.cosh(y) + a*np.cos(x)) - np.log(1 + a)
+    # want the extremum to be positive, need to reverse sign
+    psi = -psi
+    # max located at psi[50, 45] = 2.1972245773362196
+    ji = (50,45)
+    return psi, ji, psi[ji]
+
+@pytest.fixture()
+def square_verts():
+    return np.array([[0,0], [1,0], [1,1], [0,1], [0,0]])
+
+
+def test_polygon_area(square_verts):
+    assert rclv.polygon_area(square_verts) == 1.0
+    # try without the last vertex
+    assert rclv.polygon_area(square_verts[:-1]) == 1.0
+
+
+def test_get_local_region():
+    # create some data
+    n = 10
+    x, y = np.meshgrid(np.arange(n), np.arange(n))
+    (j,i), x_reg = rclv.get_local_region(x, (2,2), border_j=(2,2), border_i=(2,2))
+    assert x_reg.shape == (5,5)
+    assert x_reg[j,i] == 0
+    assert x_reg[j,0] == 2
+    assert x_reg[j,-1] == -2
+
+    with pytest.raises(ValueError) as ve:
+        (j,i), x_reg = rclv.get_local_region(x, (2,2), border_j=(3,2), border_i=(2,2))
+    with pytest.raises(ValueError) as ve:
+        (j,i), x_reg = rclv.get_local_region(x, (2,2), border_j=(2,7), border_i=(2,2))
+    with pytest.raises(ValueError) as ve:
+        (j,i), x_reg = rclv.get_local_region(x, (2,2), border_j=(2,2), border_i=(3,2))
+    with pytest.raises(ValueError) as ve:
+        (j,i), x_reg = rclv.get_local_region(x, (2,2), border_j=(2,2), border_i=(2,7))
+
+
+def test_is_contour_closed(square_verts):
+    assert rclv.is_contour_closed(square_verts)
+    assert not rclv.is_contour_closed(square_verts[:-1])
+
+
+def test_point_in_contour(square_verts):
+    assert rclv.point_in_contour(square_verts, (0.5, 0.5))
+    assert not rclv.point_in_contour(square_verts, (1.5, 0.5))
+
+
+def test_contour_area(square_verts):
+    region_area, hull_area, convex_def = rclv.contour_area(square_verts)
+    assert region_area == 1.0
+    assert hull_area == 1.0
+    assert convex_def == 0.0
+
+
+def test_contour_around_maximum(sample_data_and_maximum):
+    psi, ji, psi_max = sample_data_and_maximum
+
+    # we should get an error if the contour intersects the domain boundary
+    with pytest.raises(ValueError):
+        _ = rclv.find_contour_around_maximum(psi, ji, psi_max + 0.1)
+
+    con, region_data, border_i, border_j = rclv.find_contour_around_maximum(
+                                                            psi, ji, psi_max/2)
+
+    # region data should be normalized to have the center point 0
+    assert region_data[border_i[1], border_j[0]] == 0.0
+    assert region_data.shape == (sum(border_j)+1, sum(border_j)+1)
+
+    # the contour should be closed
+    assert rclv.is_contour_closed(con)
+
+    # check size against reference solution
+    region_area, hull_area, convex_def = rclv.contour_area(con)
+    np.testing.assert_allclose(region_area, 575.02954788959767)
+    np.testing.assert_allclose(hull_area, 575.0296629815823)
+    assert convex_def == (hull_area - region_area) / region_area
+
+
+def test_convex_contour_around_maximum(sample_data_and_maximum):
+    psi, ji, psi_max = sample_data_and_maximum
+
+    # step determines how precise the contour identification is
+    step = 0.001
+    con, area = rclv.convex_contour_around_maximum(psi, ji, step)
+
+    # check against reference solution
+    np.testing.assert_allclose(area, 2693.8731123245125)
+    assert len(con) == 261
+
+    # for this specific psi, contour should be symmetric around maximum
+    assert tuple(con[:-1].mean(axis=0).astype('int')) == ji
+
+
+def test_find_convex_contours(sample_data_and_maximum):
+    psi, ji, psi_max = sample_data_and_maximum
+    res =list(rclv.find_convex_contours(psi, step=0.001))
+
+    assert len(res) == 1
+
+    ji_found, con, area = res[0]
+    assert tuple(ji_found) == ji
+    assert len(con) == 261
+    np.testing.assert_allclose(area, 2693.8731123245125)
+
+    # also test the "filling in" function
+    labels = rclv.label_points_in_contours(psi.shape, [con])
+    assert labels.max() == 1
+    assert labels.min() == 0
+    assert labels.sum() == 2693