Commit 71452dca authored by jay's avatar jay
Browse files

Merge branch 'subpixel' into refactor

parents fb5eb22d f74f628a
Loading
Loading
Loading
Loading
+1 −1
Original line number Diff line number Diff line
@@ -39,7 +39,7 @@ before_install:
  # Install dependencies
  - conda config --add channels conda-forge
  - conda config --set ssl_verify false
  - conda install -c conda-forge vlfeat scipy networkx numexpr cython matplotlib pillow runipy geopandas opencv numpy libgdal gdal
  - conda install -c conda-forge vlfeat scipy networkx numexpr cython matplotlib pillow runipy geopandas opencv numpy libgdal gdal scikit-image
  - conda install -c menpo cyvlfeat
  - conda install -c usgs-astrogeology plio
  
+9 −5
Original line number Diff line number Diff line
@@ -38,9 +38,13 @@ We suggest using Anaconda Python to install Autocnet within a virtual environmen
   * ``conda create -n <your_environment_name> python=3 && source activate <your_environment_name>``
   
   Note, that you might want to specify either ``python=3.5`` or ``python=3.6``, depending on your requirements. Both are currently supported by autocnet.
#. Bring up a command line and add three channels to your conda config (``~/condarc``):
#. Make the newly created environment the active one:

   * ``conda config --add channels conda-forge``
   * ``conda config --add channels menpo``
   * ``conda config --add channels usgs-astrogeology``
#. Finally, install autocnet: ``conda install -c usgs-astrogeology autocnet``
   * ``conda activate <your_environment_name>`` (or ``source activate`` on an older conda system)
   
#. Bring up a command line and add three channels to your conda environment-specific config file:
   
   * ``conda config --env --add channels conda-forge``
   * ``conda config --env --add channels menpo``
   * ``conda config --env --add channels usgs-astrogeology``
#. Finally, install autocnet: ``conda install autocnet``
+2 −2
Original line number Diff line number Diff line
@@ -491,8 +491,8 @@ class Edge(dict, MutableMapping):
            raise AttributeError('This edge does not yet have any matches computed.')

        matches, mask = self.clean(clean_keys)
        domain = self.source.geodata.raster_size

        rs = self.source.geodata.raster_size
        domain = [0, 0, rs[0], rs[1]]
        # Massage the dataframe into the correct structure
        coords = self.source.get_keypoint_coordinates()
        merged = matches.merge(coords, left_on=['source_idx'], right_index=True)
+73 −86
Original line number Diff line number Diff line
@@ -45,133 +45,116 @@ def distance_ratio(edge, matches, ratio=0.8, single=False):
    return mask


def spatial_suppression(df, domain, min_radius=1.5, k=250, error_k=0.1):
def spatial_suppression(df, bounds, xkey='lon', ykey='lat', k=60, error_k=0.05, nsteps=250):
    """
    Spatial suppression using disc based method.
    Apply the spatial suppression algorithm over an arbitrary domain for all of the spatial
    data in the provided data frame.

    Attributes
    Parameters
    ----------
    df : dataframe
         Input dataframe used for suppressing

    mask : series
           pandas boolean series
    df : object
         Pandas data frame with coordinates

    max_radius : float
                 Maximum allowable point radius
    bounds : list
             In the form xmin, ymin, xmax, ymax

    min_radius : float
                 The smallest allowable radius size
    xkey : str
           The column name for the x coordinates

    nvalid : int
             The number of valid points after suppression
    ykey : str
           The column name for the y coordinates
    
    k : int
        The number of points to be saved
        The desired number of points after suppression

    error_k : float
              [0,1] the acceptable error in k
              The percentage of allowable error in the domain [0,1]

    domain : tuple
             The (x,y) extent of the input domain
    nsteps : int
             The granularity of the search. This controls the number of
             buckets in the x and y dimension. More granular search adds processing
             time, but can result in a more accurate solution.

    Returns
    -------
    mask : pd.Series
           Boolean suppression mask

    k : int
        The number of unsuppressed observations

    References
    ----------
    [Gauglitz2011]_
    mask : nd.array
           A boolean mask of the valid points

    len(result) : int
                  The numer of valud points
    """
    columns = df.columns
    for i in ['x', 'y', 'strength']:
        if i not in columns:
            raise ValueError('The dataframe is missing a {} column.'.format(i))
    df = df.sort_values(by=['strength'], ascending=False).copy()
    # Compute the bounding area inside of which the suppression will be applied
    minx = min(bounds[0], bounds[2])
    maxx = max(bounds[0], bounds[2])
    miny = min(bounds[1], bounds[3])
    maxy = max(bounds[1], bounds[3])
    domain = (maxx-minx),(maxy-miny)

    min_radius = min(domain) / 20
    max_radius = max(domain)
    mask = pd.Series(False, index=df.index)

    process = True
    if k > len(df):
        warnings.warn('Only {} valid points, but {} points requested'.format(len(df), k))
        k = len(df)
        result = df.index
        process = False
    nsteps = max(domain) * 0.95
    search_space = np.linspace(min_radius, max_radius, nsteps)
    cell_sizes = search_space / math.sqrt(2)
    min_idx = 0
    max_idx = len(search_space) - 1

    # Setup flags to watch for looping
    prev_min = None
    prev_max = None

    while process:
        # Setup to store results
        result = []
    # Sort the dataframe (hard coded to ascending as lower strength (cost) is better)
    df = df.sort_values(by=['strength'], ascending=True).copy()
    df = df.reset_index(drop=True)
    mask = pd.Series(False, index=df.index)

    process = True
    while process:
        # Binary search
        mid_idx = int((min_idx + max_idx) / 2)

        if min_idx == mid_idx or mid_idx == max_idx:
            warnings.warn('Unable to optimally solve.  Returning with {} points'.format(len(result)))
            warnings.warn('Unable to optimally solve.')
            process = False
        else:
            # Setup to store results
            result = []

        # Get the current cell size and grid the domain
        cell_size = cell_sizes[mid_idx]
        n_x_cells = int(domain[0] / cell_size)
        n_y_cells = int(domain[1] / cell_size)
        grid = np.zeros((n_x_cells, n_y_cells), dtype=np.bool)
        n_x_cells = int(round(domain[0] / cell_size, 0)) - 1
        n_y_cells = int(round(domain[1] / cell_size, 0)) - 1

        if n_x_cells <= 0:
            n_x_cells = 1
        if n_y_cells <= 0:
            n_y_cells = 1

        grid = np.zeros((n_y_cells, n_x_cells), dtype=np.bool)
        # Assign all points to bins
        x_edges = np.linspace(0, domain[0], n_x_cells)
        y_edges = np.linspace(0, domain[1], n_y_cells)
        xbins = np.digitize(df['x'], bins=x_edges)
        ybins = np.digitize(df['y'], bins=y_edges)

        # Convert bins to cells
        xbins -= 1
        ybins -= 1
        pts = []
        x_edges = np.linspace(minx, maxx, n_x_cells)
        y_edges = np.linspace(miny, maxy, n_y_cells)
        xbins = np.digitize(df[xkey], bins=x_edges)
        ybins = np.digitize(df[ykey], bins=y_edges)

        # Starting with the best point, start assigning points to grid cells
        for i, (idx, p) in enumerate(df.iterrows()):
            x_center = xbins[i]
            y_center = ybins[i]
            x_center = xbins[i] - 1
            y_center = ybins[i] - 1
            cell = grid[y_center, x_center]

            if cell == False:
                result.append(idx)
                pts.append((p[['x', 'y']]))
                if len(result) > k + k * error_k:
                    # Too many points, break
                    min_idx = mid_idx
                    break

                y_min = y_center - int(round(cell_size, 0))
                if y_min < 0:
                    y_min = 0

                x_min = x_center - int(round(cell_size, 0))
                if x_min < 0:
                    x_min = 0
                # Set the cell to True
                grid[y_center, x_center] = True

                y_max = y_center + int(round(cell_size, 0))
                if y_max > grid.shape[0]:
                    y_max = grid.shape[0]
            # If everything is already 'covered' break from the list
            if grid.all() == False:
                continue

                x_max = x_center + int(round(cell_size, 0))
                if x_max > grid.shape[1]:
                    x_max = grid.shape[1]

                # Cover the necessary cells
                grid[y_min: y_max,
                     x_min: x_max] = True

        #  Check break conditions
        # Check to see if the algorithm is completed, or if the grid size needs to be larger or smaller
        if k - k * error_k <= len(result) <= k + k * error_k:
            # Success, in bounds
            process = False

        elif len(result) < k - k * error_k:
            # The radius is too large
            max_idx = mid_idx
@@ -181,10 +164,14 @@ def spatial_suppression(df, domain, min_radius=1.5, k=250, error_k=0.1):
                process = False
            if min_idx == max_idx:
                process = False
    mask = pd.Series(False, df.index)

        elif len(result) > k + k * error_k:
            # Too many points, break
            min_idx = mid_idx

    mask.loc[list(result)] = True
    
    return mask, k
    return mask, len(result)


def self_neighbors(matches):
+34 −0
Original line number Diff line number Diff line
import numpy as np

from skimage.feature import register_translation

from autocnet.matcher import naive_template
from autocnet.matcher import ciratefi

@@ -55,6 +57,38 @@ def clip_roi(img, center, img_size):
                                             x_stop + 1, y_stop + 1])
    return clipped_img

def subpixel_phase(template, search, **kwargs):
    """
    Apply the spectral domain matcher to a search and template image. To
    shift the images, the x_shift and y_shift, need to be subtracted from
    the center of the search image. It may also be necessary to apply the
    fractional pixel adjustment as well (if for example the center of the
    search is not an integer); this function do not manage shifting.

    Parameters
    ----------
    template : ndarray
               The template used to search

    search : ndarray
             The search image

    Returns
    -------
    x_offset : float
               Shift in the x-dimension

    y_offset : float
               Shift in the y-dimension

    strength : tuple
               With the RMSE error and absolute difference in phase
    """
    if not template.shape == search.shape:
        raise ValueError('Both the template and search images must be the same shape.')

    (y_shift, x_shift), error, diffphase = register_translation(search, template, **kwargs)
    return x_shift, y_shift, (error, diffphase)

def subpixel_offset(template, search, **kwargs):
    """
Loading