reorder x/y edge coordinates

order_edges.py

# this should let ya surf the line.
# take edge coordinates (1 px thick) and reorder 
# so first-last goes from one end of the line to the other

import numpy as np
from skimage import morphology
from scipy.spatial import distance
from collections import deque as dq

def order_edges(edges):
    # convert edge masks to ordered x-y coords
    
    # check if we got a list or what
    if isinstance(edges, list):
        edges = np.array(edges)
    
    # check if we're starting w/ coords or a labeled matrix of ints
    
    if edges.shape[1]==2:
        # got coords. only built to handle one edge at a time this way for now
        uq_edges = [1]
        got_image = False
        
        
    else:
        # got an imagelike thing
        got_image = True #see?
        
        uq_edges = np.unique(edges)
        uq_edges = uq_edges[uq_edges>0]

        # if we only get one edge, try to label in case they just forgot
        if len(uq_edges)==1:
            edges = morphology.label(edges)
            uq_edges = np.unique(edges)
            uq_edges = uq_edges[uq_edges>0]
        elif len(uq_edges)==0:
            Exception('nothing here you dadgum rascal')
    
        
    edges_xy = []

    for e in uq_edges:
        if got_image:
            e_pts = np.column_stack(np.where(edges==e))
        else:
            e_pts = edges
           
        dists = distance.squareform(distance.pdist(e_pts))
        # make tri...nary connectedness graph, max dist is ~1.4 for diagonal pix
        # convert to 3 and 2 so singly-connected points are always <4
        dists[dists > 1.5] = 0
        dists[dists >1]    = 3
        dists[dists == 1]  = 2

        # check if we have easy edges - 
        # any point w/ sum of connectivity weights <4 can only have single connection
        dists_sum = np.sum(dists, axis=1)
        ends = np.where(dists_sum<4)[0]
        if len(ends)>0:
            pt_i = ends[0]
            first_i = ends[0]
            got_end = True
        else:
            # but if the edge has forex. a horiz & diag neighbor that'll fail
            # so just pick zero and get going.
            pt_i = 0
            first_i = 0
            got_end = False

        # walk through our dist graph, gathering points as we go
        
        # this tight lil bundle will get reused a bit...
        # we are making a new list of points, and don't want to double-count points
        # but we can't pop from e_pts directly, because then the indices from the
        # dist mat will be wrong. Instead we make a list of all the indices and pop
        # from that. But since popping will change the index/value parity, we
        # have to double index inds.pop(inds.index(etc.))
        inds = range(len(e_pts))
        new_pts = dq() # use a dq so we can quickly leftappend
        forwards = True
        
        # inital conditions n we're off to the races.
        new_pts.append(e_pts[inds.pop(inds.index(pt_i))])
        while True:
            # get dict of connected points and distances
            # filtered by whether the index hasn't been added yet
            connected_pts = {k: dists[pt_i,k] for k in np.where(dists[pt_i,:])[0] if k in inds}

            # if we get nothing, we're either done or we have to go back to the first pt
            if len(connected_pts) == 0:
                if got_end:
                    # still have points left, go back to first and go backwards
                    # not conditioning on len(inds) because we might drop 
                    # some degenerate diag points along the way
                    # if we accidentally started at an end it won't hurt much
                    pt_i = first_i
                    forwards = False
                    got_end = False
                    continue
                else:
                    # got to the end lets get outta here
                    break

            # find point with min distance (take horiz/vert points before diags)
            pt_i = min(connected_pts, key=connected_pts.get)
            if forwards:
                new_pts.append(e_pts[inds.pop(inds.index(pt_i))])
            else:
                new_pts.appendleft(e_pts[inds.pop(inds.index(pt_i))])
        
        # finish'er up
        new_pts = np.array(new_pts)
        # if we're only doing one, return here
        if len(uq_edges)==1:
            return new_pts
    
        # otherwise stash it and head back for more
        pts_xy.append(new_pts)
        
    return pts_xy