Commit 09ba2960 authored by jlaura's avatar jlaura Committed by Kelvin Rodriguez
Browse files

Exposes point funcs (#344) 

* Exposes point funcs

* Docstrings

* Linspace supports 2D

* Fixes to properly use linspace

* cleaner func

* Adds docs
parent 142718b0

autocnet/cg/cg.py

+85 −41
Original line number Diff line number Diff line 
from math import isclose

import warnings



import pandas as pd
@@ -210,12 +211,11 @@ def single_centroid(geom): 
    Returns

    -------



    valid : list

     : list

            in the form [(x,y)]

    """

    x, y = geom.centroid.xy

    valid = [(x[0],y[0])]

    return valid

    return [(x[0],y[0])]



def nearest(pt, search):

    """
@@ -238,6 +238,51 @@ def nearest(pt, search): 
    """

    return np.argmin(np.sum((search - pt)**2, axis=1))



def create_points_along_line(p1, p2, npts):

    """

    Compute a set of nodes equally spaced between

    two points, not including the end points.



    Parameters

    ----------

    p1 : iterable

         in the form (x,y)

    

    p2 : iterable

         in the form(x,y)



    npts : int

           The number of nodes to be returned



    Returns

    -------

     : ndarray

       (n,2) array of nodes

    """

    # npts +2 since the endpoints are included in linspace

    # but this func clips them

    return np.linspace(p1, p2, npts+2)[1:-1]



def xy_in_polygon(x,y, geom):

    """

    Returns true is an x,y pair is contained within

    the geom. 



    Parameters

    ----------

    x : Number

        The x coordinate

    

    y : Number

        The y coordinate



    Returns

    -------

     : bool

       True if the point is contained within the geom.

    """

    return geom.contains(Point(x, y))



def distribute_points(geom, nspts, ewpts):

    """

    This is a decision tree that attempts to perform a
@@ -275,44 +320,29 @@ def distribute_points(geom, nspts, ewpts): 
    # Find the points nearest the ul and ur

    ul_actual = geom_coords[nearest(ul, geom_coords)]

    ur_actual = geom_coords[nearest(ur, geom_coords)]

    dist = np.sqrt((ul_actual[1] - ur_actual[1])**2 + (ul_actual[0] - ur_actual[0])**2)

    m = (ul_actual[1]-ur_actual[1])/(ul_actual[0]-ur_actual[0])

    b = (ul_actual[1] - ul_actual[0] * m)

    newtop = []

    xnodes = np.linspace(ul_actual[0], ur_actual[0], num=ewpts+2)

    for x in xnodes[1:-1]:

        newtop.append((x, m*x+b))

    newtop = create_points_along_line(ul_actual, ur_actual, ewpts)



    # Find the points nearest the ll and lr



    ll_actual = geom_coords[nearest(ll, geom_coords)]

    lr_actual = geom_coords[nearest(lr, geom_coords)]

    dist = np.sqrt((ll_actual[1] - lr_actual[1])**2 + (ll_actual[0] - lr_actual[0])**2)

    m = (ll_actual[1]- lr_actual[1])/(ll_actual[0]-lr_actual[0])

    b = (ll_actual[1] - ll_actual[0] * m)

    newbot = []

    xnodes = np.linspace(ll_actual[0], lr_actual[0], num=ewpts+2)

    for x in xnodes[1:-1]:

        newbot.append((x, m*x+b))

    newpts = []

    newbot = create_points_along_line(ll_actual, lr_actual, ewpts)



    points = []

    for i in range(len(newtop)):

        top = newtop[i]

        bot = newbot[i]

        # Compute the line between top and bottom

        m = (top[1] - bot[1]) / (top[0] - bot[0])

        b = (top[1] - top[0] * m)

        xnodes = np.linspace(bot[0], top[0], nspts+2)

        for x in xnodes[1:-1]:

            newpts.append((x, m*x+b))

    valid = []

        

        line_of_points = create_points_along_line(top, bot, nspts)

        points.append(line_of_points)



    points = np.vstack(points)

    # Perform a spatial intersection check to eject points that are not valid

    for p in newpts:

        pt = Point(p[0], p[1])

        if geom.contains(pt):

            valid.append(p)

    valid = [p for p in points if xy_in_polygon(p[0], p[1], geom)]

    return valid



def distribute_points_in_geom(geom):

def distribute_points_in_geom(geom, 

                              nspts_func=lambda x: int(round(x,1)*10), 

                              ewpts_func=lambda x: int(round(x,1)*5)):

    """

    Given a geometry, attempt a basic classification of the shape.

    RIght now, this simply attempts to determine if the bounding box
@@ -320,11 +350,26 @@ def distribute_points_in_geom(geom): 
    is made, the algorithm places points in the geometry and returns

    a list of valid (intersecting) points.



    The kwargs for this algorithm take a function that expects a number

    as an input and returns an integer number of points to place. The 

    input number is the distance between the top/bottom or left/right

    sides of the geometry. 



    This algorithm does not know anything about the units being used

    so the caller is responsible for acocunting for units (if appropriate)

    in the passed funcs.



    Parameters

    ----------

    geom : shapely.geom object

           The geometry object



    nspts_func : obj

                 Function taking a Number and returning an int

    

    ewpts_func : obj

                 Function taking a Number and returning an int



    Returns

    -------

    valid : list
@@ -333,7 +378,7 @@ def distribute_points_in_geom(geom): 
    """

    coords = list(zip(*geom.envelope.exterior.xy))

    short = np.inf

    long = -np.inf

    lng = -np.inf

    shortid = 0

    longid = 0

    for i, p in enumerate(coords[:-1]):
@@ -341,21 +386,20 @@ def distribute_points_in_geom(geom): 
        if d < short:

            short = d

            shortid = i

        if d > long:

            long = d

        if d > lng:

            lng = d

            longid = i

    ratio = short/long

    ratio = short/lng

    ns = False

    ew = False

    valid = []



    # The polygons should be encoded with a lower left origin in counter-clockwise direction.

    # Therefore, if the 'bottom' is the short edge it should be id 0 and modulo 2 == 0.

    if shortid % 2 == 0:

        # Also if the geom is a perfect square

        ns = True

    elif longid % 2 == 0:

        ew = True



    # Decision Tree

    if ratio < 0.16 and geom.area < 0.01:

        # Class: Slivers - ignore.
@@ -365,16 +409,16 @@ def distribute_points_in_geom(geom): 
        valid = single_centroid(geom)

    elif ns==True:

        # Class, north/south poly, multi-point

        nspts = int(round(long, 1) * 10)

        ewpts = max(int(round(short, 1) * 5), 1)

        nspts = nspts_func(lng)

        ewpts = ewpts_func(short)

        if nspts == 1 and ewpts == 1:

            valid = single_centroid(geom)

        else:

            valid = distribute_points(geom, nspts, ewpts)

    elif ew == True:

        # Since this is an LS, we should place these diagonally from the 'lower left' to the 'upper right'

        nspts = max(int(round(short, 1) * 5), 1)

        ewpts = int(round(long, 1) * 10)

        nspts = ewpts_func(short)

        ewpts = nspts_func(lng)

        if nspts == 1 and ewpts == 1:

            valid = single_centroid(geom)

        else:

autocnet/cg/tests/test_cg.py

+41 −34
Original line number Diff line number Diff line 
import os

import sys

import unittest

sys.path.insert(0, os.path.abspath('..'))



import numpy as np

import pandas as pd



from .. import cg

from autocnet.cg import cg

from osgeo import ogr

from shapely.geometry import Polygon

from unittest.mock import Mock, MagicMock
@@ -17,20 +15,24 @@ from autocnet.graph.network import CandidateGraph 
from autocnet.graph.edge import Edge

from autocnet.utils.utils import array_to_poly



import pytest



class TestArea(unittest.TestCase):



    def setUp(self):

@pytest.fixture

def pts():

    seed = np.random.RandomState(12345)

        self.pts = seed.rand(25, 2)

    return seed.rand(25, 2)



@pytest.fixture

def nspoly():

    return Polygon([(0,0),(.2,0),(.2,1), (0,1), (0,0)])



    def test_area_single(self):

def test_area_single(pts):

    total_area = 1.0

        ratio = cg.convex_hull_ratio(self.pts, total_area)

    ratio = cg.convex_hull_ratio(pts, total_area)



        self.assertAlmostEqual(0.7566490, ratio, 5)

    assert pytest.approx(0.7566490, 5) == ratio



    def test_overlap(self):

def test_overlap():

    wkt1 = "POLYGON ((0 40, 40 40, 40 0, 0 0, 0 40))"

    wkt2 = "POLYGON ((20 60, 60 60, 60 20, 20 20, 20 60))"
@@ -39,37 +41,42 @@ class TestArea(unittest.TestCase): 


    info = cg.two_poly_overlap(poly1, poly2)



        self.assertEqual(info[1], 400)

        self.assertAlmostEqual(info[0], 14.285714285)

    assert info[1] == 400

    assert pytest.approx(info[0]) == 14.285714285



    def test_geom_mask(self):

def test_geom_mask():

    my_gdf = pd.DataFrame(columns=['x', 'y'], data=[(0, 0), (2, 2)])

    my_poly = Polygon([(1, 1), (3, 1), (3, 3), (1, 3)])

    mask = cg.geom_mask(my_gdf, my_poly)

        self.assertFalse(mask[0])

        self.assertTrue(mask[1])

    assert mask[0] == False

    assert mask[1] == True



    def test_compute_voronoi(self):

        keypoints = pd.DataFrame({'x': (15, 18, 18, 12, 12), 'y': (6, 10, 15, 15, 10)})

        intersection = Polygon([(10, 5), (20, 5), (20, 20), (10, 20)])

@pytest.fixture

def keypoints():

    return pd.DataFrame({'x': (15, 18, 18, 12, 12), 'y': (6, 10, 15, 15, 10)})



def test_voronoi_keypoints(keypoints):

    voronoi_gdf = cg.compute_voronoi(keypoints)

        self.assertAlmostEqual(voronoi_gdf.weight[0], 12.0)

        self.assertAlmostEqual(voronoi_gdf.weight[1], 13.5)

        self.assertAlmostEqual(voronoi_gdf.weight[2], 7.5)

        self.assertAlmostEqual(voronoi_gdf.weight[3], 7.5)

        self.assertAlmostEqual(voronoi_gdf.weight[4], 13.5)

    for i, v in enumerate([12.0, 13.5, 7.5, 7.5, 13.5]):

        assert pytest.approx(voronoi_gdf.weight[i]) == v



def test_voronoi_keypoints_with_geom(keypoints):

    voronoi_gdf = cg.compute_voronoi(keypoints, geometry=True)

        self.assertAlmostEqual(voronoi_gdf.geometry[0].area, 12.0)

        self.assertAlmostEqual(voronoi_gdf.geometry[1].area, 13.5)

        self.assertAlmostEqual(voronoi_gdf.geometry[2].area, 7.5)

        self.assertAlmostEqual(voronoi_gdf.geometry[3].area, 7.5)

        self.assertAlmostEqual(voronoi_gdf.geometry[4].area, 13.5)

    for i, v in enumerate([12.0, 13.5, 7.5, 7.5, 13.5]):

        assert pytest.approx(voronoi_gdf.geometry.area[i]) == v



def test_voronoi_keypoint_intersection(keypoints):

    intersection = Polygon([(10, 5), (20, 5), (20, 20), (10, 20)])

    voronoi_inter_gdf = cg.compute_voronoi(keypoints, intersection)

        self.assertAlmostEqual(voronoi_inter_gdf.weight[0], 22.5)

        self.assertAlmostEqual(voronoi_inter_gdf.weight[1], 26.25)

        self.assertAlmostEqual(voronoi_inter_gdf.weight[2], 37.5)

        self.assertAlmostEqual(voronoi_inter_gdf.weight[3], 37.5)

        self.assertAlmostEqual(voronoi_inter_gdf.weight[4], 26.25)
 
    for i, v in enumerate([22.5, 26.25, 37.5, 37.5, 26.25]):

        assert pytest.approx(voronoi_inter_gdf.weight[i]) == v



@pytest.mark.parametrize("polygon, nexpected",[

    (Polygon([(0,0), (.2,0), (.2,1), (0,1), (0,0)]), 10),

    (Polygon([(0,0), (1,0), (1,.2), (0,.2), (0,0)]), 10),

    (Polygon([(0,0), (.2, .1), (.2,1.1), (-0.1, 1), (0,0)]), 11)

],

    ids=['vertical', 'horizontal', 'verticalskewed'])

def test_points_in_geom(polygon, nexpected):

    pts = cg.distribute_points_in_geom(polygon)

    assert len(pts) == nexpected
 
 No newline at end of file