RANSAC算法应用

问题描述

在对数据进行建模时，我们往往会遇到数据中存在异常数据（outlier，一般称为外点）的情况，如图，用一条直线拟合图中的点，其中左下和右上的点明显是外点。

此时使用最小二乘法，则会由于外点的存在偏离正确模型较多，如下图，因此我们需要找出并剔除这些外点，从而得到一个较好的模型。

RANSAC算法原理

RANSAC（Random Sample Consensus）基本思想是通过多次随机抽取部分数据进行模型估计，然后使用全部数据评估这些模型的正确性并进行比较，得出最好的模型。具体步骤如下。

随机抽取数据集中的一些样本 Si ，估计一个模型 Mi 。通过一定标准检验能较好的符合模型 Mi 的数据集 S′i ，称为内点（inliner）。若 size(Si)>N ，N为人为设定的阈值，则算法结束，得到模型M否则重新进行步骤1～3.若一定重复次数后算法依然没结束，则选取内点数目最多的模型。使用得到模型的内点再次计算更优的模型 M′

使用RANSAC算法得到的拟合直线如下

具体实现

RANSAC进行直线拟合的代码如下：

---

ransac.hpp

#include <opencv2\opencv.hpp> #include <vector> class RansacModelEstimator { public: RansacModelEstimator() {}; cv::Mat run(std::vector<cv::Point3d> point_set); int iterate(std::vector<int>& mask); void estimateParm(std::vector<cv::Point3d>& initial_point); void estimateInliners(std::vector<int>& mask); cv::Mat estimateParm1(std::vector<cv::Point3d>& initial_point); double estimateDeparture(cv::Point3d& p); private: std::vector<cv::Point3d> point_set; cv::Point3d theta; std::vector<int> rand_index; };

---

ransac.cpp

#include "ransac.hpp" #include <random> #include <vector> using namespace std; using namespace cv; double RansacModelEstimator::estimateDeparture(Point3d& p) { return fabs(p.dot(theta)/ norm(Point2d(theta.x,theta.y))); } void RansacModelEstimator::estimateInliners(vector<int>& mask) { mask.clear(); for (size_t i = 0; i < point_set.size(); i++) { if (estimateDeparture(point_set[i]) < 0.3) mask.push_back(i); } } void RansacModelEstimator::estimateParm(vector<cv::Point3d>& initial_point) { theta = initial_point[0].cross(initial_point[1]); } Mat RansacModelEstimator::estimateParm1(vector<cv::Point3d>& inliner_point) { Mat constrains(inliner_point.size(),3,CV_64F); Mat b(inliner_point.size(), 1, CV_64F); b.setTo(0); for (int i = 0;i<inliner_point.size();i++) { constrains.at<double>(i, 0) = inliner_point[i].x; constrains.at<double>(i, 1) = inliner_point[i].y; constrains.at<double>(i, 2) = inliner_point[i].z; } Mat tconstrains; Mat result; cv::SVD::solveZ(constrains, result); return result; } int RansacModelEstimator::iterate(vector<int>& mask) { static std::random_device rd; static std::default_random_engine e(rd()); static std::uniform_int_distribution<> u(0, point_set.size()-1); vector<cv::Point3d> initial_point; int last=-1; for (size_t i = 0; i < 2; i++) { int index = u(e); while (index == last) index = u(e); initial_point.push_back(point_set[index]); last = index; } estimateParm(initial_point); estimateInliners(mask); return mask.size(); } Mat RansacModelEstimator::run(vector<cv::Point3d> point_set_) { point_set = point_set_; vector<int> optimal_mask; vector<int> current_mask; int optimal_inliner_num = 0; int inliner_num; for (size_t i = 0; i < 300; i++) { inliner_num = iterate(current_mask); if (optimal_inliner_num < inliner_num) { optimal_mask = current_mask; optimal_inliner_num = inliner_num; } } vector<Point3d> inliners; for (size_t i = 0; i < optimal_inliner_num; i++) { inliners.push_back(point_set[optimal_mask[i]]); } return estimateParm1(inliners); } int main() { RansacModelEstimator est; Mat result = est.run({ {1,1,1},{2,2,1},{3,4,1},{5,4,1},{7,6.9,1},{6,6.1,1},{1,111,1},{23,2,1},{11.8,12.1,1} }); std::cout << result << endl; return 0; }