Search code examples
pythonnumpycomputational-geometry

sorting a complicated collection of 2d euclidian points in in clockwise/counterclockwise fashion to form a closed ring


This looks like a repeated question but I tried the solution that already exists and none seems to work so far for me. .

this solution gives a hint but it works only for a regular geometry. I have a rather complicated geometry from which I extract boundary points which are unsorted.

Below is a picture of the geometry and of the boundary vertices that I extract from the geometry. geometry_triangulation enter image description here

The (x,y) coordinates of points as in the image are :

import numpy as np

pts = np.array([[  30.        ,   -6.25      ],
                [  30.        ,   -8.10127917],
                [   0.        ,   -6.25      ],
                [  34.14082772,   -6.75584268],
                [  36.49784598,  -10.        ],
                [  44.43561524,  -10.        ],
                [ 100.        ,  -10.        ],
                [ 100.        ,   10.        ],
                [  84.1244615 ,  -10.        ],
                [  84.1244615 ,   10.        ],
                [  36.49784598,   10.        ],
                [  34.14082772,    6.75584268],
                [  44.43561524,   10.        ],
                [  30.        ,    8.10127917],
                [  30.        ,    6.25      ],
                [   0.        ,    6.25      ],
                [ -30.        ,    6.25      ],
                [ -30.        ,    8.10127917],
                [ -32.92183092,    9.05063958],
                [ -35.84366185,   10.        ],
                [ -51.88274638,   10.        ],
                [-100.        ,   10.        ],
                [-100.        ,  -10.        ],
                [ -83.96091546,   10.        ],
                [ -83.96091546,  -10.        ],
                [ -35.84366185,  -10.        ],
                [ -51.88274638,  -10.        ],
                [ -32.92183092,   -9.05063958],
                [ -30.        ,   -8.10127917],
                [ -30.        ,   -6.25      ],
                [ -67.92183092,   10.        ],
                [ -67.92183092,  -10.        ],
                [  68.24892299,   10.        ],
                [  52.37338449,   10.        ],
                [  68.24892299,  -10.        ],
                [  52.37338449,  -10.        ]])

In the boundary vertex data, we can see that the points are unordered. Is there a way the points could be ordered clockwise/counter clockwise so that these points form a closed ring when connected successively?

My goal is to create a polygon or a linear ring as described here and later find if an arbitrary euclidean point lies inside the polyogn/ring

Update: the approach of computing angle between centroid of the pts and the individual euclidean point in pts also doesn't work. Here is a sample code of what I tried :

def sort_circular(pts):
    cent = coords.mean(axis=0)
    idx = list(np.arange(0, len(pts)+1, dtype=int))
    angle = []
    for i, cc in enumerate(coords):
        dx,dy = cc[0] - center[0], cc[1]-center[1]
        angle.append(math.degrees(math.atan2(float(dy), float(dx))))
    #simultaneously sort angle and indices
    _, idx_sorted = (list(t) for t in zip(*sorted(zip(angle, idx))))
    pts_sorted = pts[idx_sorted]
    return pts_sorted

The result from this is still not the way i expect it to be (image below): enter image description here


Solution

  • Method 1:

    Define a center point, compute the angle between every coordinate and the center point, then order them by angle:

    import pandas as pd
    
    # Define function to compute angle between vectors
    import math
    
    def clockwiseangle_and_distance(point, origin = [0,0], refvec = [1,0]):
        # Vector between point and the origin: v = p - o
        vector = [point[0]-origin[0], point[1]-origin[1]]
        # Length of vector: ||v||
        lenvector = math.hypot(vector[0], vector[1])
        # If length is zero there is no angle
        if lenvector == 0:
            return -math.pi, 0
        # Normalize vector: v/||v||
        normalized = [vector[0]/lenvector, vector[1]/lenvector]
        dotprod  = normalized[0]*refvec[0] + normalized[1]*refvec[1]     # x1*x2 + y1*y2
        diffprod = refvec[1]*normalized[0] - refvec[0]*normalized[1]     # x1*y2 - y1*x2
        angle = math.atan2(diffprod, dotprod)
        # Negative angles represent counter-clockwise angles so we need to subtract them 
        # from 2*pi (360 degrees)
        if angle < 0:
            return 2*math.pi+angle, lenvector
        # I return first the angle because that's the primary sorting criterium
        # but if two vectors have the same angle then the shorter distance should come first.
        return angle, lenvector
    
    import pandas as pd
    
    # Compute the center point
    center = pts.mean(axis=0)
    
    angle = []
    for i in range(len(pts)):
        ang, dist = clockwiseangle_and_distance(pts[i,:] - center, origin=[0,0], refvec=[1,0])
        angle.append(ang)
    
    df = pd.DataFrame(pts)
    df['angle'] = np.degrees(angle)
    
    df = df.sort_values(by='angle')
    df['clockwise_order'] = np.arange(len(df))
    
    import matplotlib.pyplot as plt
    
    # Create plot to show the ordering of the points
    plt.figure()
    df.plot(kind='scatter', x=0, y=1, s=100, alpha=0.5)
    plt.title('Points by clockwise order')
    
    for idx, row in df.iterrows():
        plt.gca().annotate('{:.0f}'.format(row['clockwise_order']), (row[0], row[1]),
                ha='center', va='center_baseline', fontsize=6, color='k', fontweight='bold')
    
    plt.gca().annotate('Center', center,
            ha='center', va='center')
    

    Plot of points by clockwise orientation

    If this clockwise ordering doesn't give you what you want, try Method 2.

    Method 2:

    To sort the points for the given geometry in clockwise order such that they form a closed ring, you can do the following:

    1. Divide the data set in quadrants
    2. Choose a center point such that the remaining points of the quadrant lie on about the arc of a circle centered at the center point
    3. Order each quadrant by clockwise angle
    4. Place each quadrant in clockwise order
    # Compute the center point
    center = pts.mean(axis=0)
    
    df = pd.DataFrame(pts)
    
    # Group points into quadrants
    df['quadrant'] = 0
    df.loc[(df[0] > center[0]) & (df[1] > center[1]), 'quadrant'] = 0
    df.loc[(df[0] > center[0]) & (df[1] < center[1]), 'quadrant'] = 1
    df.loc[(df[0] < center[0]) & (df[1] < center[1]), 'quadrant'] = 2
    df.loc[(df[0] < center[0]) & (df[1] > center[1]), 'quadrant'] = 3
    
    quadrant = {}
    for i in range(4):
        quadrant[i] = df[df.quadrant == i]
    
    # Intelligently choose the quadrant centers
    x = 35
    y = 5
    subcenter = [[ x,  y],
                 [ x, -y],
                 [-x, -y],
                 [-x,  y]]
    
    # Compute the angle between each quadrant and respective center point
    angle = {}
    points = {}
    df_sub = {}
    for j in range(len(quadrant)):
        angle[j] = []
        points[j] = quadrant[j][[0,1]]
        for i in range(len(points[j])):
            ang, dist = clockwiseangle_and_distance(points[j].values[i,:] - subcenter[j], origin=[0,0], refvec=[1,0])
            angle[j].append(ang)
    
        df_sub[j] = quadrant[j]
        df_sub[j]['angle'] = np.degrees(angle[j])
        df_sub[j] = df_sub[j].sort_values(by='angle')
    
    # Combine the data frames
    df = pd.concat(df_sub)
    df['clockwise_order'] = np.arange(len(df))
    
    # Plot the points by clockwise order
    import matplotlib.pyplot as plt
    
    # Create plot to show the ordering of the points
    fig, axis = plt.subplots()
    df[[0,1]].plot(x=0, y=1, ax=axis, c='lightblue', legend=False, clip_on=False)
    df.plot(kind='scatter', x=0, y=1, s=100, ax=axis, c='lightblue', clip_on=False)
    plt.title('Points by quadrant in clockwise order')
    plt.axis('off')
    
    for idx, row in df.iterrows():
        plt.gca().annotate('{:.0f}'.format(row['clockwise_order']), (row[0], row[1]),
                ha='center', va='center_baseline', fontsize=6, color='k', fontweight='bold')
    
    plt.gca().annotate('Center', center,
            ha='center', va='center')
    
    for i in range(len(subcenter)):
        plt.scatter(subcenter[i][0], subcenter[i][1], alpha=0.5, s=80, marker='s')
        plt.gca().annotate('Quadrant \n'+str(i)+'\n', subcenter[i],
            ha='center', va='center_baseline', color='k', fontsize=8)
    

    Plots by quadrant in clockwise order

    # Plot with axis equally-spaced
    df2 = df[[0,1]].reset_index(drop=True)
    df2.loc[len(df2),:] = df2.loc[0,:]
    df2.plot(x=0, y=1, c='k', legend=False, clip_on=False)
    plt.axis('equal')
    plt.axis('off')
    

    Plot with axis equally-spaced

    If this doesn't give you what you want, you may have to order the coordinates by hand.