增强RANSAC, Adaptive Parameter Calculation

# -*- coding: utf-8 -*-
"""
Created on Fri May  4 00:31:43 2018

@author: Sikasjc
"""
import numpy as np
import random, sys

class Model(object):
    def fit(self, data):
        raise NotImplementedError
        
    def distance(self, data):
        raise NotImplementedError

class LineModel(Model):
    """
    A 2D line model.
    """
    def __init__(self):
        self.params = None
        self.dists = None
        
    def fit(self, data):
        """
        Fits the model to the data, minimizing the sum of absolute errors.
        """
        X = data[:,0]
        Y = data[:,1]
#        denom = (X[-1] - X[0])
#        if denom == 0:
#            raise ZeroDivisionError
#        k = (Y[-1] - Y[0]) / denom
#        m = Y[0] - k * X[0]
        A = np.vstack([X, np.ones(len(X))]).T
        k, m = np.linalg.lstsq(A, Y, None)[0]
        X_ = np.mean(X)
        Y_ = np.mean(Y)
        self.params = [k, m]
        self.residual = abs(k * X_ + m - Y_)
    
    def distance(self, samples):
        """
        Calculates the vertical distances from the samples to the model.
        """
        X = samples[:,0]
        Y = samples[:,1]
        k = self.params[0]
        m = self.params[1]
        dists = abs(k * X + m - Y)
#        dists = abs(-k * X + Y - m) / math.sqrt(k**2 + 1)
        
        return dists

def ransac_APC(data, model, min_samples, min_inliers, eps=1e-10, P = 0.99, random_seed=42):
    """
    Fits a model to observed data.
    
    Uses the RANSC iterative method of fitting a model to observed data.
    """
    random.seed(random_seed)
    
    if len(data) <= min_samples:
        raise ValueError("Not enough input data to fit the model.")
        
    iterations = sys.maxsize
    count = 0
    min_inliers = 0
    
    best_params = None
    best_inliers = None
    best_residual = np.inf
    
    while iterations > count: 
        if 0 < min_inliers and min_inliers < 1:
            min_inliers = int(min_inliers * len(data))
        
        indices = list(range(len(data)))
        random.shuffle(indices)
        inliers = np.asarray([data[i] for i in indices[:min_samples]])
        shuffled_data = np.asarray([data[i] for i in indices[min_samples:]])
        
        try:
            model.fit(inliers)
            dists = model.distance(shuffled_data)
            more_inliers = shuffled_data[np.where(dists <= eps)[0]]
            inliers = np.concatenate((inliers, more_inliers))
            
            if len(inliers) >= min_inliers:
                model.fit(inliers)
                if model.residual < best_residual:
                    best_params = model.params
                    best_inliers = inliers
                    best_residual = model.residual
                                        
        except ZeroDivisionError as e:
            print(e)
            
        ratio = len(inliers) / len(data)
        if min_inliers/len(data) < ratio:
            best_residual = np.inf
            min_inliers = ratio
            iterations = int(np.log(1 - P) / np.log(1 - min_inliers**min_samples)) + 1 
        count += 1 
        
    if best_params is None:
        raise ValueError("RANSAC failed to find a sufficiently good fit for the data.")
    else:
        return (best_params, best_inliers, best_residual, count)

# -*- coding: utf-8 -*-
"""
Created on Fri May  4 00:17:50 2018

@author: Sikasjc
"""
import numpy as np
import matplotlib.pyplot as plt

import random
import ransac_APC
import time

num_samples = 1000
noise_ratio = 0.8
num_noise = int(noise_ratio * num_samples)

def setup():
    global ax1
    X = np.asarray(range(num_samples))
    Y = 0.5 * X
    noise = [random.randint(0, 2 * (num_samples - 1)) for i in range(num_noise)]
    Y[random.sample(range(len(Y)), num_noise)] = noise
    data = np.asarray([X, Y]).T
    model = ransac_APC.LineModel()
    plt.figure(figsize=(10,5))
    ax1 = plt.subplot(1,2,1)
    plt.plot(X, Y, 'bx')
    
    return data, model

def run(data, model):
    random_seed = random.randint(0,100)
    start_time = time.time()
    (params, inliers, residual, iterations) = ransac_APC.ransac_APC(data, model, 2, min_inliers = 0.1, eps=1e-10, P=0.99, random_seed=random_seed)
    end_time = time.time()
    mean_time = (end_time - start_time) / iterations
    return params, residual, mean_time, iterations
    
    
def summary(params, residual, mean_time, iterations):
    print(" Paramters ".center(40, '='))
    print(params)
    print(" Residual ".center(40, '='))
    print(residual)
    print(" Iteration ".center(40, '='))
    print(iterations, "the number of iterations in the running")
    print(" Time ".center(40, '='))
    print("%.1f msecs mean time spent per iteration" % (1000 * mean_time))
    X = np.asanyarray([0, num_samples - 1])
    Y = params[0] * X + params[1]
    plt.subplot(1,2,2, sharey = ax1)
    plt.plot(X, Y, 'g-', linewidth=2)
    
    plt.show()
    
if __name__ == '__main__':
    data, model = setup()
    try:
        params, residual, mean_time, iterations = run(data, model)
        summary(params, residual, mean_time, iterations)
    except ValueError as e:
       print(e)