RANSAC估计——以直线拟合为例
RANSAC(RANdom SAmple Consensus),即随机采样一致性。该方法最早是由Fischler和Bolles提出的一种鲁棒估计方法,最早用于计算机视觉中位姿估计问题,现在已广泛应用于已知模型的参数估计问题中。对于数据中存在大比例外点(错误数据、野值点)时,RANSAC方法十分有效。下面就以直线估计的例子来说明RANSAC的基本思想。
直线拟合RANSAC估计基本原理
RANSAC的思想比较简单,主要有以下几步:
- 随机选择两点(确定一条直线所需要的最小点集);由这两个点确定一条线ll;
- 使用几何距离,求最大一致集的最佳拟合直线,作为数据点的最佳匹配直线。
如下图所示,首先在点集中随机选择两个点,求解它们构成的直线参数,再计算点集中剩余点到该直线的距离,距离小于设置的距离阈值的点为内点,统计内点个数;接下来再随机选取两个点,同样统计内点个数,以此类推;其中内点个数最多的点集即为最大一致集,最后将该最大一致集里面的点利用最小二乘拟合出一条直线。
参数设置
根据上面RANSAC基本原理的介绍,在这算法流程中存在两个重要的参数需要设置,采样次数和距离阈值。
为了保证算法的效率,当数据量比较大时,我们不可能采集所有的样本,将所有的可能情况都计算一遍。如何在抽样当中保证至少有一个好样本(不包含外点)是我们需要考虑的问题。Fischler和Bolles在原文中通过概率统计的方法得出了采样次数与数据中外点比例和得到一个好样本概率之间的数学关系,在这里我直接给出结论,有兴趣推导该关系的博友可以自己阅读原文或者由吴福朝编写的《计算机视觉中的数学方法》。
K=log(1−z)log(1−wn)K=log(1−z)log(1−wn)
其中,KK为模型参数估计需要的最小点个数,直线拟合最少需要2个点。
距离阈值一般需要根据经验来设定,在直线拟合实验中,我设置的参数为2.99。当观测误差符合0均值和sigmasigma。
RANSAC直线拟合代码
结合OpenCV库(2.0版本以上)实现了RANSAC直线拟合,在主函数中给出了一个例子,仅供大家参考。
#include <opencv2\opencv.hpp>
#include <iostream>
#include <ctime>using namespace std;
using namespace cv;//生成[0,1]之间符合均匀分布的数
double uniformRandom(void)
{return (double)rand() / (double)RAND_MAX;
}//生成[0,1]之间符合高斯分布的数
double gaussianRandom(void)
{/* This Gaussian routine is stolen from Numerical Recipes and is theircopyright. */static int next_gaussian = 0;static double saved_gaussian_value;double fac, rsq, v1, v2;if (next_gaussian == 0) {do {v1 = 2 * uniformRandom() - 1;v2 = 2 * uniformRandom() - 1;rsq = v1*v1 + v2*v2;} while (rsq >= 1.0 || rsq == 0.0);fac = sqrt(-2 * log(rsq) / rsq);saved_gaussian_value = v1*fac;next_gaussian = 1;return v2*fac;}else {next_gaussian = 0;return saved_gaussian_value;}
}//根据点集拟合直线ax+by+c=0,res为残差
void calcLinePara(vector<Point2d> pts, double &a, double &b, double &c, double &res)
{res = 0;Vec4f line;vector<Point2f> ptsF;for (unsigned int i = 0; i < pts.size(); i++)ptsF.push_back(pts[i]);fitLine(ptsF, line, CV_DIST_L2, 0, 1e-2, 1e-2);a = line[1];b = -line[0];c = line[0] * line[3] - line[1] * line[2];for (unsigned int i = 0; i < pts.size(); i++){double resid_ = fabs(pts[i].x * a + pts[i].y * b + c);res += resid_;}res /= pts.size();
}//得到直线拟合样本,即在直线采样点集上随机选2个点
bool getSample(vector<int> set, vector<int> &sset)
{int i[2];if (set.size() > 2){do{for (int n = 0; n < 2; n++)i[n] = int(uniformRandom() * (set.size() - 1));} while (!(i[1] != i[0]));for (int n = 0; n < 2; n++){sset.push_back(i[n]);}}else{return false;}return true;
}//直线样本中两随机点位置不能太近
bool verifyComposition(const vector<Point2d> pts)
{cv::Point2d pt1 = pts[0];cv::Point2d pt2 = pts[1];if (abs(pt1.x - pt2.x) < 5 && abs(pt1.y - pt2.y) < 5)return false;return true;
}//RANSAC直线拟合
void fitLineRANSAC(vector<Point2d> ptSet, double &a, double &b, double &c, vector<bool> &inlierFlag)
{double residual_error = 2.99; //内点阈值bool stop_loop = false;int maximum = 0; //最大内点数//最终内点标识及其残差inlierFlag = vector<bool>(ptSet.size(), false);vector<double> resids_(ptSet.size(), 3);int sample_count = 0;int N = 500;double res = 0;// RANSACsrand((unsigned int)time(NULL)); //设置随机数种子vector<int> ptsID;for (unsigned int i = 0; i < ptSet.size(); i++)ptsID.push_back(i);while (N > sample_count && !stop_loop){vector<bool> inlierstemp;vector<double> residualstemp;vector<int> ptss;int inlier_count = 0;if (!getSample(ptsID, ptss)){stop_loop = true;continue;}vector<Point2d> pt_sam;pt_sam.push_back(ptSet[ptss[0]]);pt_sam.push_back(ptSet[ptss[1]]);if (!verifyComposition(pt_sam)){++sample_count;continue;}// 计算直线方程calcLinePara(pt_sam, a, b, c, res);//内点检验for (unsigned int i = 0; i < ptSet.size(); i++){Point2d pt = ptSet[i];double resid_ = fabs(pt.x * a + pt.y * b + c);residualstemp.push_back(resid_);inlierstemp.push_back(false);if (resid_ < residual_error){++inlier_count;inlierstemp[i] = true;}}// 找到最佳拟合直线if (inlier_count >= maximum){maximum = inlier_count;resids_ = residualstemp;inlierFlag = inlierstemp;}// 更新RANSAC迭代次数,以及内点概率if (inlier_count == 0){N = 500;}else{double epsilon = 1.0 - double(inlier_count) / (double)ptSet.size(); //野值点比例double p = 0.99; //所有样本中存在1个好样本的概率double s = 2.0;N = int(log(1.0 - p) / log(1.0 - pow((1.0 - epsilon), s)));}++sample_count;}//利用所有内点重新拟合直线vector<Point2d> pset;for (unsigned int i = 0; i < ptSet.size(); i++){if (inlierFlag[i])pset.push_back(ptSet[i]);}calcLinePara(pset, a, b, c, res);
}void main()
{//demoint width = 640;int height = 320;//直线参数double a = 1, b = 2, c = -640;//随机获取直线上20个点int ninliers = 0;vector<Point2d> ptSet;srand((unsigned int)time(NULL)); //设置随机数种子while (true){double x = uniformRandom()*(width - 1);double y = -(a*x + c) / b;//加0.5高斯噪声x += gaussianRandom()*0.5;y += gaussianRandom()*0.5;if (x >= 640 && y >= 320)continue;Point2d pt(x, y);ptSet.push_back(pt);ninliers++;if (ninliers == 20)break;}int nOutliers = 0;//随机获取10个野值点while (true){double x = uniformRandom() * (width - 1);double y = uniformRandom() * (height - 1);if (fabs(a*x + b*y + c) < 10) //野值点到直线距离不小于10个像素continue;Point2d pt(x, y);ptSet.push_back(pt);nOutliers++;if (nOutliers == 10)break;}//绘制点集中所有点Mat img(321, 641, CV_8UC3, Scalar(255, 255, 255));for (unsigned int i = 0; i < ptSet.size(); i++)circle(img, ptSet[i], 3, Scalar(255, 0, 0), 3, 8);double A, B, C;vector<bool> inliers;fitLineRANSAC(ptSet, A, B, C, inliers);B = B / A;C = C / A;A = A / A;//绘制直线Point2d ptStart, ptEnd;ptStart.x = 0;ptStart.y = -(A*ptStart.x + C) / B;ptEnd.x = -(B*ptEnd.y + C) / A;ptEnd.y = 0;line(img, ptStart, ptEnd, Scalar(0, 255, 255), 2, 8);cout << "A:" << A << " " << "B:" << B << " " << "C:" << C << " " << endl;imshow("line fitting", img);waitKey();
}
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
- 24
- 25
- 26
- 27
- 28
- 29
- 30
- 31
- 32
- 33
- 34
- 35
- 36
- 37
- 38
- 39
- 40
- 41
- 42
- 43
- 44
- 45
- 46
- 47
- 48
- 49
- 50
- 51
- 52
- 53
- 54
- 55
- 56
- 57
- 58
- 59
- 60
- 61
- 62
- 63
- 64
- 65
- 66
- 67
- 68
- 69
- 70
- 71
- 72
- 73
- 74
- 75
- 76
- 77
- 78
- 79
- 80
- 81
- 82
- 83
- 84
- 85
- 86
- 87
- 88
- 89
- 90
- 91
- 92
- 93
- 94
- 95
- 96
- 97
- 98
- 99
- 100
- 101
- 102
- 103
- 104
- 105
- 106
- 107
- 108
- 109
- 110
- 111
- 112
- 113
- 114
- 115
- 116
- 117
- 118
- 119
- 120
- 121
- 122
- 123
- 124
- 125
- 126
- 127
- 128
- 129
- 130
- 131
- 132
- 133
- 134
- 135
- 136
- 137
- 138
- 139
- 140
- 141
- 142
- 143
- 144
- 145
- 146
- 147
- 148
- 149
- 150
- 151
- 152
- 153
- 154
- 155
- 156
- 157
- 158
- 159
- 160
- 161
- 162
- 163
- 164
- 165
- 166
- 167
- 168
- 169
- 170
- 171
- 172
- 173
- 174
- 175
- 176
- 177
- 178
- 179
- 180
- 181
- 182
- 183
- 184
- 185
- 186
- 187
- 188
- 189
- 190
- 191
- 192
- 193
- 194
- 195
- 196
- 197
- 198
- 199
- 200
- 201
- 202
- 203
- 204
- 205
- 206
- 207
- 208
- 209
- 210
- 211
- 212
- 213
- 214
- 215
- 216
- 217
- 218
- 219
- 220
- 221
- 222
- 223
- 224
- 225
- 226
- 227
- 228
- 229
- 230
- 231
- 232
- 233
- 234
- 235
- 236
- 237
- 238
- 239
- 240
- 241
- 242
- 243
- 244
- 245
- 246
- 247
- 248
- 249
- 250
- 251
- 252
- 253
- 254
- 255
- 256
- 257
- 258
实验结果如下:
RANSAC估计——以直线拟合为例相关推荐
- RANSAC算法做直线拟合
RANSAC算法之前了解过相关的原理,这两天利用晚上闲暇的时间,看了一下RANSAC算法的Python代码实现,这方面的资料很多了,这里就不在重复.在分析该RANSAC.py代码之前,想用自己的对RA ...
- RANSAC算法实现 + 直线拟合
一.RANSAC算法 1.参考资料 [1]题目来源与解析:商汤科技SLAM算法岗的RANSAC编程题 [2]牛客网题目:[编程题]线性回归 [3]牛客网解答参考:商汤科技某算法岗的编程题有点过分了啊 ...
- RANSAC原理及直线拟合(python动态图解)
一.简介 随机采样一致性(Random Sample Consensus,RANSAC)由斯坦福国际研究院的Fischler和Bolles于1981年首次提出[1].RANSAC算法是一种随机参数估计 ...
- RANSAC算法(附RANSAC直线拟合C++与Python版本)
文章目录 RANSAC算法简介 RANSAC算法基本思想和流程 迭代次数推导 RANSAC与最小二乘区别 RANSAC直线拟合代码(C++及Python版本) C++版本代码 Python版本代码如下 ...
- RANSAC 直线拟合算法
1. 参考文献 RANSAC 直线拟合算法 2. 算法实现 #include <iostream> #include <random> #include <vector& ...
- RANSAC直线拟合和最小二乘直线拟合的简单介绍
RANSAC是"RANdom SAmple Consensus(随机抽样一致)"的缩写.它可以从一组包含"局外点"的观测数据集中,通过迭代方式估计 ...
- 用RANSAC算法实现干扰严重的直线拟合~
1.说到直线拟合,一般是用最小二乘啦,在opencv里面就是用cv.fitLine来完成,首先简单介绍一下该函数: cv.fitLine(points, distType, param, reps, ...
- RANSAC 直线拟合
result: code: #include "stdafx.h". #include <opencv2\opencv.hpp> #include <iostre ...
- PCL:多直线拟合(RANSAC)
文章目录 1 RANSAC 空间直线拟合 2 RANSAC 多直线拟合 1 RANSAC 空间直线拟合 具体可参考以下博客. https://blog.csdn.net/weixin_46098577 ...
最新文章
- python 自带虚拟环境
- 大自然创作的分形艺术
- 在Sql2005中,向表中插入数据时遇到uniqueidentifier列,如何插入数据?
- 广域网设备:远程路由器/交换机/服务器—Vecloud
- python字典get计数_Python内部是如何存储GC引用变量的计数的?
- [机器学习]梯度提升决策树--GBDT
- nothing comes free
- 文本怎么整体向左缩进_如何设置Word文档的左缩进或者右缩进呢?并在缩进区域输入文字...
- Android 实现ListView圆角效果
- 136. 只出现一次的数字 golang
- 17秋 SDN课程 第一次作业
- 西门子pcs7 V9.1用Sim_EKB_Install安装长密匙失败
- idea配置tomcat热部署
- java apktool if_apktool反编译详细使用教程
- matlab 植被指数,基于matlab的MOD13A2-NDVI的植被指数重建-SG滤波与质量控制文件
- MATLAB创建三对角线矩阵-211103
- 特种部队2全面反击bt电影种子下载
- 汉洛塔问题(c解决)
- Docer实现Django Uwsgi部署
- office快捷键设置
热门文章
- Nginx的11个phases
- AOP Aspect Oriented Programming 面向切面编程 Spring
- 【python】Macbook的Anaconda查看、创建和管理python环境
- Upload-Labs(1-5)
- 计费系统设计_Web设计人员的按小时计费与基于价值的定价
- 跨域 cors 请求两次_请求两次的故事-CORS
- azure多功能成像好用吗_了解Azure持久功能
- 连续两天高烧_连续工作两天,可以看电视11秒
- 杭电oj 1001 c++版本
- 雨滴桌面时间插件_只需要五步,让你的桌面变得干净漂亮