def color_labels(labels, distance_transform = False):
    '''Color a labels matrix so that no adjacent labels have the same color
    
    distance_transform - if true, distance transform the labels to find out
         which objects are closest to each other.
         
    Create a label coloring matrix which assigns a color (1-n) to each pixel
    in the labels matrix such that all pixels similarly labeled are similarly
    colored and so that no similiarly colored, 8-connected pixels have
    different labels.
    
    You can use this function to partition the labels matrix into groups
    of objects that are not touching; you can then operate on masks
    and be assured that the pixels from one object won't interfere with
    pixels in another.
    
    returns the color matrix
    '''
    if distance_transform:
        i,j = scind.distance_transform_edt(labels == 0, return_distances=False,
                                           return_indices = True)
        dt_labels = labels[i,j]
    else:
        dt_labels = labels
    # Get the neighbors for each object
    v_count, v_index, v_neighbor = find_neighbors(dt_labels)
    # Quickly get rid of labels with no neighbors. Greedily assign
    # all of these a color of 1
    v_color = np.zeros(len(v_count)+1,int) # the color per object - zero is uncolored
    zero_count = (v_count==0)
    if np.all(zero_count):
        # can assign all objects the same color
        return (labels!=0).astype(int)
    v_color[1:][zero_count] = 1
    v_count = v_count[~zero_count]
    v_index = v_index[~zero_count]
    v_label = np.argwhere(~zero_count).transpose()[0]+1
    # If you process the most connected labels first and use a greedy
    # algorithm to preferentially assign a label to an existing color,
    # you'll get a coloring that uses 1+max(connections) at most.
    #
    # Welsh, "An upper bound for the chromatic number of a graph and
    # its application to timetabling problems", The Computer Journal, 10(1)
    # p 85 (1967)
    #
    sort_order = np.lexsort([-v_count])
    v_count = v_count[sort_order]
    v_index = v_index[sort_order]
    v_label = v_label[sort_order]
    for i in range(len(v_count)):
        neighbors = v_neighbor[v_index[i]:v_index[i]+v_count[i]]
        colors = np.unique(v_color[neighbors])
        if colors[0] == 0:
            if len(colors) == 1:
                # only one color and it's zero. All neighbors are unlabeled
                v_color[v_label[i]] = 1
                continue
            else:
                colors = colors[1:]
        # The colors of neighbors will be ordered, so there are two cases:
        # * all colors up to X appear - colors == np.arange(1,len(colors)+1)
        # * some color is missing - the color after the first missing will
        #   be mislabeled: colors[i] != np.arange(1, len(colors)+1)
        crange = np.arange(1,len(colors)+1)
        misses = crange[colors != crange]
        if len(misses):
            color = misses[0]
        else:
            color = len(colors)+1
        v_color[v_label[i]] = color
    return v_color[labels]