Commit 091f5843 authored by Kelvin Rodriguez's avatar Kelvin Rodriguez Committed by GitHub
Browse files

Merge pull request #163 from jlaura/coupled

Coupled Decomposition (closes #142)
parents d01a584f bcd6e8cb
Loading
Loading
Loading
Loading
+278 −8
Original line number Diff line number Diff line
@@ -3,6 +3,7 @@ from collections import MutableMapping

import numpy as np
import pandas as pd
from scipy.spatial.distance import cdist

from autocnet.utils import utils
from autocnet.matcher import health
@@ -10,9 +11,9 @@ from autocnet.matcher import outlier_detector as od
from autocnet.matcher import suppression_funcs as spf
from autocnet.matcher import subpixel as sp
from autocnet.matcher.feature import FlannMatcher
from autocnet.transformation.decompose import coupled_decomposition
from autocnet.transformation.transformations import FundamentalMatrix, Homography
from autocnet.vis.graph_view import plot_edge
from autocnet.vis.graph_view import plot_node
from autocnet.vis.graph_view import plot_edge, plot_node, plot_edge_decomposition
from autocnet.cg import cg


@@ -91,32 +92,298 @@ class Edge(dict, MutableMapping):
    def health(self):
        return self._health.health

    def match(self, k=2):
    def decompose_and_match(self, k=2, maxiteration=3, size=18, buf_dist=3,**kwargs):
        """
        Similar to match, this method first decomposed the image into
        $4^{maxiteration}$ subimages and applys matching between each sub-image.

        This method is potential slower than the standard match due to the
        overhead in matching, but can be significantly more accurate.  The
        increase in accuracy is a function of the total image size.  Suggested
        values for maxiteration are provided below.

        Parameters
        ----------
        k : int
            The number of neighbors to find

        method : {'coupled', 'whole'}
                 whether to utilize coupled decomposition
                 or match the whole image

        maxiteration : int
                       When using coupled decomposition, the number of recursive
                       divisions to apply.  The total number of resultant
                       sub-images will be 4 ** maxiteration.  Approximate values:

                        | Number of megapixels | maxiteration |
                        |----------------------|--------------|
                        | m < 10               |1-2|
                        | 10 < m < 30          | 3 |
                        | 30 < m < 100         | 4 |
                        | 100 < m < 1000       | 5 |
                        | m > 1000             | 6 |

        size : int
               When using coupled decomposition, the total number of points
               to check in each sub-image to try and find a match.
               Selection of this number is a balance between seeking a
               representative mid-point and computational cost.

        buf_dist : int
                   When using coupled decomposition, the distance from the edge of
                   the (sub)image a point must be in order to be used as a
                   partioning point.  The smaller the distance, the more likely
                   percision errors can results in erroneous partitions.
        """
        def mono_matches(a, b, aidx=None, bidx=None):
            """
            Apply the FLANN match_features

            Parameters
            ----------
            a : object
                A node object

            b : object
                A node object

            aidx : iterable
                   An index for the descriptors to subset

            bidx : iterable
                   An index for the descriptors to subset
            """
            # Subset if requested
            if aidx is not None:
                ad = a.descriptors[aidx]
            else:
                ad = a.descriptors

            if bidx is not None:
                bd = b.descriptors[bidx]
            else:
                bd = b.descriptors

            # Load, train, and match
            fl.add(ad, a.node_id, index=aidx)
            fl.train()
            matches = fl.query(bd, b.node_id, k, index=bidx)
            self._add_matches(matches)
            fl.clear()

        def func(group):
            ratio = 0.8
            res = [False] * len(group)
            if len(res) == 1:
                return [single]
            if group.iloc[0] < group.iloc[1] * ratio:
                res[0] = True
            return res

        # Grab the original image arrays
        sdata = self.source.get_array()
        ddata = self.destination.get_array()

        ssize = sdata.shape
        dsize = ddata.shape

        # Grab all the available candidate keypoints
        skp = self.source.get_keypoints()
        dkp = self.destination.get_keypoints()

        # Set up the membership arrays
        self.smembership = np.zeros(sdata.shape, dtype=np.int16)
        self.dmembership = np.zeros(ddata.shape, dtype=np.int16)
        self.smembership[:] = -1
        self.dmembership[:] = -1
        pcounter = 0

        # FLANN Matcher
        fl= FlannMatcher()

        for k in range(maxiteration):
            partitions = np.unique(self.smembership)
            for p in partitions:
                sy_part, sx_part = np.where(self.smembership == p)
                dy_part, dx_part = np.where(self.dmembership == p)

                # Get the source extent
                minsy = np.min(sy_part)
                maxsy = np.max(sy_part) + 1
                minsx = np.min(sx_part)
                maxsx = np.max(sx_part) + 1

                # Get the destination extent
                mindy = np.min(dy_part)
                maxdy = np.max(dy_part) + 1
                mindx = np.min(dx_part)
                maxdx = np.max(dx_part) + 1

                # Clip the sub image from the full images
                asub = sdata[minsy:maxsy, minsx:maxsx]
                bsub = ddata[mindy:maxdy, mindx:maxdx]

                # Utilize the FLANN matcher to find a match to approximate a center
                fl.add(self.destination.descriptors, self.destination.node_id)
                fl.train()

                scounter = 0
                decompose = False
                while True:
                    sub_skp = skp.query('x >= {} and x <= {} and y >= {} and y <= {}'.format(minsx, maxsx, minsy, maxsy))
                    # Check the size to ensure a valid return
                    if len(sub_skp) == 0:
                        break # No valid keypoints in this (sub)image
                    if size > len(sub_skp):
                        size = len(sub_skp)
                    candidate_idx = np.random.choice(sub_skp.index, size=size, replace=False)
                    candidates = self.source.descriptors[candidate_idx]
                    matches = fl.query(candidates, self.source.node_id, k=3, index=candidate_idx)

                    # Apply Lowe's ratio test to try to find a 'good' starting point
                    mask = matches.groupby('source_idx')['distance'].transform(func).astype('bool')
                    candidate_matches = matches[mask]
                    match_idx = candidate_matches['source_idx']

                    # Extract those matches that pass the ratio check
                    sub_skp = skp.iloc[match_idx]

                    # Check that valid points remain
                    if len(sub_skp) == 0:
                        break

                    # Locate the candidate closest to the middle of all of the matches
                    smx, smy = sub_skp[['x', 'y']].mean()
                    mid = np.array([[smx, smy]])
                    dists = cdist(mid, sub_skp[['x', 'y']])
                    closest = sub_skp.iloc[np.argmin(dists)]
                    closest_idx = closest.name
                    soriginx, soriginy = closest[['x', 'y']]

                    # Grab the corresponding point in the destination
                    q = candidate_matches.query('source_idx == {}'.format(closest.name))
                    dest_idx = q['destination_idx'].iat[0]
                    doriginx = dkp.at[dest_idx, 'x']
                    doriginy = dkp.at[dest_idx, 'y']

                    if mindy + buf_dist <= doriginy <= maxdy - buf_dist\
                     and mindx + 3 <= doriginx <= maxdx - 3:
                        # Point is good to split on
                        decompose = True
                        break
                    else:
                        scounter += 1
                        if scounter >= maxiteration:
                            break

                # Clear the Flann matcher for reuse
                fl.clear()

                # Check that the identified match falls within the (sub)image
                # This catches most bad matches that have passed the ratio check
                if not (buf_dist <= doriginx - mindx <= bsub.shape[1] - buf_dist) or not\
                       (buf_dist <= doriginy - mindy <= bsub.shape[0] - buf_dist):
                       decompose = False

                if decompose:
                    # Apply coupled decomposition, shifting the origin to the sub-image
                    s_submembership, d_submembership = coupled_decomposition(asub, bsub,
                                                                         sorigin=(soriginx - minsx, soriginy - minsy),
                                                                         dorigin=(doriginx - mindx, doriginy - mindy),
                                                                         **kwargs)

                    # Shift the returned membership counters to a set of unique numbers
                    s_submembership += pcounter
                    d_submembership += pcounter

                    # And assign membership
                    self.smembership[minsy:maxsy,
                                minsx:maxsx] = s_submembership
                    self.dmembership[mindy:maxdy,
                                mindx:maxdx] = d_submembership
                    pcounter += 4
        
        # Now match the decomposed segments to one another
        for p in np.unique(self.smembership):
            sy_part, sx_part = np.where(self.smembership == p)
            dy_part, dx_part = np.where(self.dmembership == p)

            # Get the source extent
            minsy = np.min(sy_part)
            maxsy = np.max(sy_part) + 1
            minsx = np.min(sx_part)
            maxsx = np.max(sx_part) + 1

            # Get the destination extent
            mindy = np.min(dy_part)
            maxdy = np.max(dy_part) + 1
            mindx = np.min(dx_part)
            maxdx = np.max(dx_part) + 1

            # Get the indices of the candidate keypoints within those regions / variables are pulled before decomp.
            sidx = skp.query('x >= {} and x <= {} and y >= {} and y <= {}'.format(minsx, maxsx, minsy, maxsy)).index
            didx = dkp.query('x >= {} and x <= {} and y >= {} and y <= {}'.format(mindx, maxdx, mindy, maxdy)).index
            # If the candidates < k, OpenCV throws an error
            if len(sidx) >= k and len(didx) >=k:
                mono_matches(self.source, self.destination, sidx, didx)
                mono_matches(self.destination, self.source, didx, sidx)

    def match(self, k=2, **kwargs):
        """
        Given two sets of descriptors, utilize a FLANN (Approximate Nearest
        Neighbor KDTree) matcher to find the k nearest matches.  Nearness is
        the euclidean distance between descriptors.

        The matches are then added as an attribute to the edge object.

        Parameters
        ----------
        k : int
            The number of neighbors to find
        """
        def mono_matches(a, b, aidx=None, bidx=None):
            """
            Apply the FLANN match_features

        Returns
        -------
            Parameters
            ----------
            a : object
                A node object

            b : object
                A node object

            aidx : iterable
                   An index for the descriptors to subset

            bidx : iterable
                   An index for the descriptors to subset
            """
        def mono_matches(a, b):
            fl.add(a.descriptors, a.node_id)
            # Subset if requested
            if aidx is not None:
                ad = a.descriptors[aidx]
            else:
                ad = a.descriptors

            if bidx is not None:
                bd = b.descriptors[bidx]
            else:
                bd = b.descriptors

            # Load, train, and match
            fl.add(ad, a.node_id, index=aidx)
            fl.train()
            self._add_matches(fl.query(b.descriptors, b.node_id, k))
            matches = fl.query(bd, b.node_id, k, index=bidx)
            self._add_matches(matches)
            fl.clear()

        fl = FlannMatcher()
        mono_matches(self.source, self.destination)
        mono_matches(self.destination, self.source)



    def _add_matches(self, matches):
        """
        Given a dataframe of matches, either append to an existing
@@ -404,6 +671,9 @@ class Edge(dict, MutableMapping):
        # Else, plot the whole edge
        return plot_edge(self, ax=ax, clean_keys=clean_keys, **kwargs)

    def plot_decomposition(self, *args, **kwargs): #pragma: no cover
        return plot_edge_decomposition(self, *args, **kwargs)

    def clean(self, clean_keys, pid=None):
        """
        Given a list of clean keys and a provenance id compute the
+11 −0
Original line number Diff line number Diff line
@@ -296,6 +296,17 @@ class CandidateGraph(nx.Graph):
        """
        self.apply_func_to_edges('match', *args, **kwargs)

    def decompose_and_match_features(self, *args, **kwargs):
        """
        For all edges in the graph, apply coupled decomposition followed by
        feature matching.

        See Also
        --------
        autocnet.graph.edge.Edge.decompose_and_match
        """
        self.apply_func_to_edges('decompose_and_match', *args, **kwargs)

    def compute_clusters(self, func=markov_cluster.mcl, *args, **kwargs):
        """
        Apply some graph clustering algorithm to compute a subset of the global
+14 −1
Original line number Diff line number Diff line
@@ -154,7 +154,7 @@ class Node(dict, MutableMapping):

        return self.coverage_area

    def get_array(self, band=1):
    def get_byte_array(self, band=1):
        """
        Get a band as a 32-bit numpy array

@@ -167,6 +167,19 @@ class Node(dict, MutableMapping):
        array = self.geodata.read_array(band=band)
        return bytescale(array)

    def get_array(self, band=1):
        """
        Get a band as a 32-bit numpy array

        Parameters
        ----------
        band : int
               The band to read, default 1
        """

        array = self.geodata.read_array(band=band)
        return array

    def get_keypoints(self, index=None):
        """
        Return the keypoints for the node.  If index is passed, return
+6 −1
Original line number Diff line number Diff line
@@ -25,10 +25,15 @@ class TestNode(unittest.TestCase):
    def test_get_handle(self):
        self.assertIsInstance(self.node.geodata, GeoDataset)

    def test_get_byte_array(self):
        image = self.node.get_byte_array()
        self.assertEqual((1012, 1012), image.shape)
        self.assertEqual(np.uint8, image.dtype)

    def test_get_array(self):
        image = self.node.get_array()
        self.assertEqual((1012, 1012), image.shape)
        self.assertEqual(np.uint8, image.dtype)
        self.assertEqual(np.float32, image.dtype)

    def test_extract_features(self):
        image = self.node.get_array()
+25 −11
Original line number Diff line number Diff line
@@ -28,9 +28,10 @@ class FlannMatcher(object):
    def __init__(self, flann_parameters=DEFAULT_FLANN_PARAMETERS):
        self._flann_matcher = cv2.FlannBasedMatcher(flann_parameters, {})
        self.nid_lookup = {}
        self.search_idx = {}
        self.node_counter = 0

    def add(self, descriptor, nid):
    def add(self, descriptor, nid, index=None):
        """
        Add a set of descriptors to the matcher and add the image
        index key to the image_indices attribute
@@ -46,6 +47,10 @@ class FlannMatcher(object):
        self._flann_matcher.add([descriptor])
        self.nid_lookup[self.node_counter] = nid
        self.node_counter += 1
        if index is not None:
            self.search_idx = dict((i, j) for i, j in enumerate(index))
        else:
            self.search_idx = dict((i,i) for i in range(len(descriptor)))

    def clear(self):
        """
@@ -55,6 +60,7 @@ class FlannMatcher(object):
        self._flann_matcher.clear()
        self.nid_lookup = {}
        self.node_counter = 0
        self.search_idx = {}

    def train(self):
        """
@@ -62,7 +68,7 @@ class FlannMatcher(object):
        """
        self._flann_matcher.train()

    def query(self, descriptor, query_image, k=3):
    def query(self, descriptor, query_image, k=3, index=None):
        """

        Parameters
@@ -76,6 +82,10 @@ class FlannMatcher(object):
        k : int
            The number of nearest neighbors to search for

        index : iterable
                An iterable of observation indices to utilize for
                the input descriptors

        Returns
        -------
        matched : dataframe
@@ -86,22 +96,26 @@ class FlannMatcher(object):

        matches = self._flann_matcher.knnMatch(descriptor, k=k)
        matched = []
        for m in matches:
            for i in m:
        for i, m in enumerate(matches):
            for j in m:
                if index is not None:
                    qid = index[i]
                else:
                    qid = j.queryIdx
                source = query_image
                destination = self.nid_lookup[i.imgIdx]
                destination = self.nid_lookup[j.imgIdx]
                if source < destination:
                    matched.append((query_image,
                                    i.queryIdx,
                                    qid,
                                    destination,
                                    i.trainIdx,
                                    i.distance))
                                    self.search_idx[j.trainIdx],
                                    j.distance))
                elif source > destination:
                    matched.append((destination,
                                    i.trainIdx,
                                    self.search_idx[j.trainIdx],
                                    query_image,
                                    i.queryIdx,
                                    i.distance))
                                    qid,
                                    j.distance))
                else:
                    warnings.warn('Likely self neighbor in query!')
        return pd.DataFrame(matched, columns=['source_image', 'source_idx',
Loading