How to find a minimum angle in polygon?

650
3
08-16-2019 12:23 AM
KK2
by
New Contributor III

Is it possible to calculate the minimum angle in a polygon? I have several polygons that have sharp edges, as in the attached picture, and I would like to indicate them automatically. For example, can I create a column in the attribute table with minimum angle values in a polygon?

0 Kudos
3 Replies
JohannesBierer
Regular Contributor

To find sharp edges you can convert the polygons to raster, boundary clean, reconvert to polygon, buffer with a appropriate distance and clip?

0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

PS, you picture is for polygons, so moving to Python since this isn't a Spatial Analyst thread.

Feed polygon shapes to 'angles_poly'.  Numpy required

Polygons can be converted to numpy arrays using.  This is the implementation for a whole featureclass converted to a numpy object array

def poly2array(polys):
    """Convert polyline or polygon shapes to arrays for use in the Geo class.

    Parameters
    ----------
    polys : tuple, list
        Polyline or polygons in a list/tuple
    """
    def _p2p_(poly):
        """Convert a single ``poly`` shape to numpy arrays or object"""
        sub = []
        for arr in poly:
            pnts = [[pt.X, pt.Y] if pt else null_pnt for pt in arr]
            sub.append(np.asarray(pnts))
        return sub
    # ----
    if not isinstance(polys, (list, tuple)):
        polys = [polys]
    out = []
    for poly in polys:
        out.append(_p2p_(poly))
    return out‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

In a script, a search cursor can be used to get a polygon, and it can be converted using _p2p_ in the above..

    def _p2p_(poly):
        """Convert a single ``poly`` shape to numpy arrays or object"""
        sub = []
        for arr in poly:
            pnts = [[pt.X, pt.Y] if pt else null_pnt for pt in arr]
            sub.append(np.asarray(pnts))
        return sub‍‍‍‍‍‍‍
def angles_poly(a=None, inside=True, in_deg=True):
    """Sequential 3 point angles from a poly* shape

    Parameters
    ----------
    a : array
        an array of points, derived from a polygon/polyline geometry
    inside : boolean
        determine inside angles, outside if False
    in_deg : boolean
        convert to degrees from radians

    Sample data: the letter C

    >>> a = np.array([[ 0, 0], [ 0, 100], [100, 100], [100,  80],
                      [ 20,  80], [ 20, 20], [100, 20], [100, 0], [ 0, 0]])
    >>> angles_poly(a)  # array([ 90.,  90.,  90., 270., 270.,  90.,  90.])
    """
    if len(a) < 2:
        return None
    if len(a) == 2:
        ba = a[1] - a[0]
        return np.arctan2(*ba[::-1])
    dx, dy = a[0] - a[-1]
    if np.allclose(dx, dy):  # closed loop
        a = a[:-1]
    a0 = np.roll(a, -1, 0)
    a1 = a
    a2 = np.roll(a, 1, 0)
    ba = a1 - a0
    bc = a1 - a2
    cr = np.cross(ba, bc)
    dt = np.einsum('ij,ij->i', ba, bc)
    ang = np.arctan2(cr, dt)
    two_pi = np.pi*2.
    if inside:
        ang = np.where(ang < 0, ang + two_pi, ang)
    else:
        ang = np.where(ang > 0, two_pi - ang, ang)
    if in_deg:
        angles = np.degrees(ang)
    return angles

# ---- demo --- a polygon "C"
import numpy as np

a = np.array([[ 0, 0], [ 0, 100], [100, 100], [100,  80],
              [ 20,  80], [ 20, 20], [100, 20], [100, 0], [ 0, 0]])
angles = angles_poly(a, inside=True, in_deg=True)

min(angles)
90.0

angles
array([270., 270., 270., 270.,  90.,  90., 270., 270.])
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

So the code would entail

open a search cursor

convert each shape to an array using _p2p_

call angles_poly

take its minimum

KK2
by
New Contributor III

Thank you for the answer. I have found an easy solution in QGIS: Vector > Geometry Tools >Check Geometries > Minimum angle between segments (deg). After its finished you can export a point shapefile indicating places where polygons have lower angle than a specified threshold, and the ID of the original polygons.