SIFT匹配中使用RANSAC算法剔除外點 MATLAB實現

　　該算法轉載自：https://blog.csdn.net/czl389/article/details/60580574

　　RANSAC為Random Sample Consensus的縮寫，它是根據一組包含異常數據的樣本數據集，計算出數據的數學模型參數，得到有效樣本數據的算法。它於1981年由Fischler和Bolles最先提出。RANSAC算法經常用於計算機視覺中。例如，在立體視覺領域中同時解決一對相機的匹配點問題及基本矩陣的計算。

　　RANSAC算法的基本假設是樣本中包含正確數據(inliers，可以被模型描述的數據)，也包含異常數據(outliers，偏離正常范圍很遠、無法適應數學模型的數據)，即數據集中含有噪聲。這些異常數據可能是由於錯誤的測量、錯誤的假設、錯誤的計算等產生的。同時RANSAC也假設，給定一組正確的數據，存在可以計算出符合這些數據的模型參數的方法。

 

RANSAC基本思想描述如下：

　　①考慮一個最小抽樣集的勢為n的模型(n為初始化模型參數所需的最小樣本數)和一個樣本集P，集合P的樣本數#(P)>n，從P中隨機抽取包含n個樣本的P的子集S初始化模型M；

　　②余集SC=P\S中與模型M的誤差小於某一設定閾值t的樣本集以及S構成S*。S*認為是內點集，它們構成S的一致集(Consensus Set)；

　　③若#(S*)≥N，認為得到正確的模型參數，並利用集S*(內點inliers)采用最小二乘等方法重新計算新的模型M*；重新隨機抽取新的S，重復以上過程。

　　④在完成一定的抽樣次數后，若未找到一致集則算法失敗，否則選取抽樣后得到的最大一致集判斷內外點，算法結束。

 

　　假設已知一組特征匹配點，使用RANSAC算法剔除其中外點。建立RANSACFit函數如下：

　　RANSACFit(p1, p2, match, maxIter, seedSetSize, maxInlierError, goodFitThresh )

　　其中輸入參數：p1，N1×2矩陣，每行都是一個特征點坐標；p2，N2×2矩陣，每行都是一個特征點坐標；match，M×2矩陣，每行代表一組匹配關系[index of p1, index of p2]；maxIter，RANSAC算法將執行的迭代次數；seedSetSize，隨機抽取的樣本點數；maxInlierError，最大內點偏差，如果這個匹配不在內點集中，將其看作內點最大可以容忍的誤差；goodFitThresh，用於評判得到的模型是否足夠優異的閾值。

　　輸出：H，p1對p2進行魯棒性仿射變換的單應估計。

　　MATLAB代碼：

function H = RANSACFit(p1, p2, match, maxIter, seedSetSize, maxInlierError, goodFitThresh )

    N = size(match, 1);
    if N<3
        error('not enough matches to produce a transformation matrix')
    end
    if ~exist('maxIter', 'var'),
        maxIter = 200;
    end
    if ~exist('seedSetSize', 'var'),
        seedSetSize = ceil(0.2 * N);
    end
    seedSetSize = max(seedSetSize,3);
    if ~exist('maxInlierError', 'var'),
        maxInlierError = 30;
    end
    if ~exist('goodFitThresh', 'var'),
        goodFitThresh = floor(0.7 * N);
    end
    H = eye(3);
    
    iota = Inf;
    kappa = 0;
    lambda = iota;
    alpha = seedSetSize;
    for i = 1 : maxIter,
        [beta, gamma] = part(match, alpha);
        eta = ComputeAffineMatrix(p1(beta(:, 1), :), p2(beta(:, 2), :));
        delta = ComputeError(eta, p1, p2, gamma);
        epsilon = (delta <= maxInlierError);
        if sum(epsilon(:)) + alpha >= goodFitThresh,
            zeta = [beta; gamma(epsilon, :)];
            eta = ComputeAffineMatrix(p1(zeta(:, 1), :), p2(zeta(:, 2), :));
            theta = sum(ComputeError(eta, p1, p2, zeta));
            if theta < iota,
                H = eta;
                kappa = lambda;
                iota = theta;
            end
        end
    end

    if sum(sum((H - eye(3)).^2)) == 0,
        disp('No RANSAC fit was found.')
    end
end

function dists = ComputeError(H, pt1, pt2, match)

    dists = zeros(size(match,1),1);
    transform_pt1=H*[pt1(match(:,1),:)';ones(1,size(match,1))];
    subtract=pt2(match(:,2),:)-transform_pt1(1:2,:)';
    dists=sqrt(subtract(:,1).^2+subtract(:,2).^2);


    if size(dists,1) ~= size(match,1) || size(dists,2) ~= 1
        error('wrong format');
    end
end

function [D1, D2] = part(D, splitSize)
    idx = randperm(size(D, 1));
    D1 = D(idx(1:splitSize), :);
    D2 = D(idx(splitSize+1:end), :);
end

2021-05-16 18:10:23 哲