ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

RANSAC圆拟合算法详解:C++实现与工业视觉应用

RANSAC圆拟合算法详解:C++实现与工业视觉应用 简介这是一份基于C实现的RANSAC算法示例工程用于从给定的二维点集中寻找n个最佳拟合圆主要面向计算机视觉、图像处理与三维重建方向的开发者也适合希望掌握鲁棒参数估计方法的中高级C学习者。压缩包内共12个文件以5个cpp源码和4个hpp头文件为主体分别实现点数据管理、线性方程求解、圆模型拟合、距离度量及RANSAC迭代主循环另附1个Markdown说明文档和1个文本格式的测试点集便于对照理解与直接运行。算法实现覆盖了RANSAC的关键流程随机选取最小样本集、通过代数方法拟合圆参数、按距离阈值筛选内点并构建共识集、迭代更新最优模型同时展示了针对噪声数据的异常值剔除与多圆查找策略。整个资源包仅26KB结构简洁清晰适合快速移植到实际项目中也可作为学习C算法工程组织的参考。目前已有1042人学习下载对想从代码层面吃透RANSAC原理的开发者而言是一份轻量且具有实用价值的参考资料。1. 先用 RANSAC 把「噪声里的圆」捞出来做视觉定位或点云处理的人大概率遇到过这种场景工件上的圆形轮廓被反光、遮挡和随机噪点污染直接最小二乘拟合出的圆心偏出去好几个像素换了几种滤波都不收敛。项目里这套RANSAC-Algorithm-master解决的正是这类问题——从带有大量离群点的 2D 点集中稳健地找出一个或多个最佳拟合圆。它不依赖人工剔除坏点而是通过反复随机采样三个点构成候选圆再用全局评估选出内点最多的模型。对 C 开发者来说这份代码的价值不只是某个main.cpp能跑通而是把随机采样、线性方程求解、模型评估和迭代控制完整地串了起来适合直接移植到工业视觉或三维重建的预处理链路里。2. 拟合圆的几何基础三点定圆与线性化求解2.1 为什么最小二乘对离群点无能为力先看一个直观现象。假设点集里有 100 个真实圆上的点带轻微高斯噪声和 30 个随机离群点。普通最小二乘拟合的目标是让所有点到圆的代数距离平方和最小离群点误差大、权重高会把圆心往错误方向拖拽。RANSAC 的策略完全不同它先假设一小撮点完全正确推导出模型再去验证这个模型能覆盖多少其他点。只要采样到的三个点都属于真实圆这个假设就能获得大量支持。三点确定一个圆的几何事实是 RANSAC 在这里的立身之本不共线的三个点 $P_1(x_1,y_1)$、$P_2(x_2,y_2)$、$P_3(x_3,y_3)$ 唯一确定一个圆。圆心 $C(x_c,y_c)$ 是三条垂直平分线的交点但直接求垂直平分线涉及浮点除法和斜率无穷大的边界处理。更稳妥的做法是解线性方程组$$ \begin{cases} (x_1-x_2)x_c (y_1-y_2)y_c \frac{x_1^2 - x_2^2 y_1^2 - y_2^2}{2} \ (x_1-x_3)x_c (y_1-y_3)y_c \frac{x_1^2 - x_3^2 y_1^2 - y_3^2}{2} \end{cases} $$这个形式的好处是直接处理坐标差规避了垂直于 x 轴的圆切线这类特殊情况。2.2 Kåsa 方法与代数距离的工程权衡上面的线性方程组能精确求圆但 RANSAC 做完后通常还要对内点做一次精拟合。最常用的精拟合是 Kåsa 方法它的思想是把圆方程 $x^2 y^2 Dx Ey F 0$ 视为关于 $D,E,F$ 的线性问题。把 $x^2y^2$ 移到等号右边就成了标准的最小二乘形式求解一个 $3 \times 3$ 线性系统即可得到圆心 $(-D/2, -E/2)$ 和半径 $\sqrt{D^2E^2-4F}/2$。这个方法的代价是引入了平方项会放大远点噪声但对工业场景中半径在几十到几千像素范围的圆误差完全可接受而且实现简洁、不需要迭代初值。下面的代码对应项目里的LinearEqu.cpp本质就是解多元一次方程组#include cmath struct Circle3 { double cx, cy, r; }; // 用 Kåsa 线性化方法从内点集合中精拟合圆 // 输入 pts: 内点集n 3 Circle3 fitCircleKasa(const std::vectorPoint2f pts) { double sumX 0, sumY 0; double sumX2 0, sumY2 0; double sumXY 0, sumX3 0, sumY3 0; for (const auto p : pts) { double x p.x, y p.y; double x2 x * x, y2 y * y; sumX x; sumY y; sumX2 x2; sumY2 y2; sumXY x * y; sumX3 x * x2; sumY3 y * y2; } int n (int)pts.size(); // 正规方程组的三个等式 // D*sumX2 E*sumXY F*sumX -sumX3 - sumXY2 // D*sumXY E*sumY2 F*sumY -sumY3 - sumX2Y // D*sumX E*sumY F*n -sumX2 - sumY2 // 实际工程一般用高斯消元或 Cramer 法则求解 Circle3 c; // 省略 3x3 求解细节见 LinearEqu.cpp 中的 Solve3x3 return c; }这里把大段矩阵求解细节留给了LinearEqu.cpp文章关注的是方法选择理由为什么不直接用三点外接圆做精拟合因为三点外接圆吃满了采样噪声内点多了反而浪费Kåsa 方法能利用全部内点做统计平均把单个点的噪声抹平。2.3 距离度量的单位问题评估点是否属于某个圆需要算点到圆心的距离与半径之差。这个差的量纲是长度如果半径是 200 像素、阈值给了 0.5 像素那标准非常严苛反之如果数据来自世界坐标系单位米噪声可能有 0.01 米阈值就得相应放大。项目里的pts_to_ransac.txt是二维笛卡尔坐标通常已经归一化到像素单位所以阈值直接按像素设。// 评估点到圆的几何误差返回带符号的距离差 // 正数表示点在圆外负数表示点在圆内 inline double pointCircleError(const CPoint2f p, const Circle3 c) { double dist std::hypot(p.x - c.cx, p.y - c.cy); return std::fabs(dist - c.r); }std::hypot比手动sqrt(x*x y*y)更稳能避免中间结果溢出值得在 C 实现里养成习惯。3. 数据结构与模块划分这套代码是怎么组织的3.1 点与圆的类型设计项目里的mypoint2f.hpp和Circle.hpp承担了最基础的数据封装。点类型不必花哨x、y两个浮点成员加一个默认构造就够圆类型需要保存cx、cy、r三个量并提供一个计算点到圆误差的成员函数方便后续评估循环直接调用。// mypoint2f.hpp struct CPoint2f { double x, y; CPoint2f() : x(0.0), y(0.0) {} CPoint2f(double px, double py) : x(px), y(py) {} };// Circle.hpp #include mypoint2f.hpp struct CCircle { double cx, cy, r; CCircle() : cx(0.0), cy(0.0), r(1.0) {} CCircle(double x, double y, double radius) : cx(x), cy(y), r(radius) {} // 点到圆轮廓的几何误差 double error(const CPoint2f p) const { return std::fabs(std::hypot(p.x - cx, p.y - cy) - r); } // 判断点是否为内点 bool isInlier(const CPoint2f p, double threshold) const { return error(p) threshold; } };这里把error和isInlier都放进CCircle调用方便语义也贴近数学定义。我自己在重构这类算法代码时会避免把点集和圆模型混在一个大类里——它们各自独立演进混在一起会让后续改造成多模型检测时非常痛苦。3.2 线性方程求解模块的作用看到项目里有LinearEqu.cpp时不要以为它只是为了 RANSAC 服务。三点定圆要用它精拟合也要用它。工程上这个模块解的是如下形式的 2x2 系统// LinearEqu.cpp 核心Cramer 法则解二元一次方程组 // a1*x b1*y c1 // a2*x b2*y c2 bool Solve2x2(double a1, double b1, double c1, double a2, double b2, double c2, double x, double y) { double det a1 * b2 - a2 * b1; if (std::fabs(det) 1e-12) return false; // 共线或退化 x (c1 * b2 - c2 * b1) / det; y (a1 * c2 - a2 * c1) / det; return true; }行列式判据很关键。三点近共线时det趋近于零求出的圆心接近无穷远半径极大——这类候选圆本身没有意义应该在拟合后直接丢弃。RANSAC 主循环里最容易踩的坑就是忘了检查这个返回值导致 NaN 污染后续比较。3.3 从文件读入点集的格式约定项目提供了pts_to_ransac.txt这是一份标准的测试输入。读文件时注意点坐标分隔符可能是空格、逗号或 Tab实现里最好两类都处理。下面这段读取函数考虑了空行和注释行的干扰#include fstream #include sstream #include string #include vector bool loadPointsFromFile(const std::string path, std::vectorCPoint2f pts) { std::ifstream in(path); if (!in.is_open()) return false; std::string line; while (std::getline(in, line)) { if (line.empty() || line[0] #) continue; std::stringstream ss(line); double x, y; // 逗号和空格都作为分隔符 while (ss x) { if (ss.peek() , || ss.peek() ) ss.ignore(); if (!(ss y)) break; pts.emplace_back(x, y); } } return !pts.empty(); }这个函数把容错放在了解析层保证后续所有算法代码只面对干净的std::vectorCPoint2f关注点分离后续做边缘检测提取到的轮廓点也能直接复用。4. RANSAC 主循环与多圆检测从三万行代码里提炼出的骨架4.1 单模型 RANSAC 的完整实现主干逻辑不复杂但每个分支都要想到最大迭代次数、内点阈值、最小内点数、模型更新策略。下面是精简过的核心循环对应main.cpp的主要流程#include random #include vector #include algorithm // 从点集中随机采样三个不同点 bool sampleTriplet(const std::vectorCPoint2f pts, std::mt19937 rng, CPoint2f a, CPoint2f b, CPoint2f c) { int n (int)pts.size(); if (n 3) return false; std::uniform_int_distributionint dist(0, n - 1); int i dist(rng), j dist(rng), k dist(rng); // 简单去重最多重试 10 次 int guard 0; while ((i j || j k || i k) guard 10) { j dist(rng); k dist(rng); } if (i j || j k || i k) return false; a pts[i]; b pts[j]; c pts[k]; return true; } // 核心 RANSAC找到最佳拟合圆 CCircle ransacCircle(const std::vectorCPoint2f pts, double threshold, int maxIterations, int minInliers, std::vectorint bestInliers, unsigned int seed 42) { std::mt19937 rng(seed); CCircle bestCircle; bestInliers.clear(); int bestCnt 0; for (int iter 0; iter maxIterations; iter) { CPoint2f a, b, c; if (!sampleTriplet(pts, rng, a, b, c)) break; // 三点拟合圆用第 3 章提到的 Solve2x2 double x1 a.x, y1 a.y; double x2 b.x, y2 b.y; double x3 c.x, y3 c.y; double A1 2 * (x2 - x1), B1 2 * (y2 - y1); double C1 x2*x2 - x1*x1 y2*y2 - y1*y1; double A2 2 * (x3 - x1), B2 2 * (y3 - y1); double C2 x3*x3 - x1*x1 y3*y3 - y1*y1; double cx, cy; if (!Solve2x2(A1, B1, C1, A2, B2, C2, cx, cy)) continue; double r std::hypot(a.x - cx, a.y - cy); // 半径合理性检查防止退化的超大半圆 if (!std::isfinite(cx) || !std::isfinite(cy) || !std::isfinite(r)) continue; if (r 1e6) continue; CCircle cand(cx, cy, r); // 评估一致性 std::vectorint inliers; inliers.reserve(pts.size() / 2); for (int idx 0; idx (int)pts.size(); idx) { if (cand.isInlier(pts[idx], threshold)) inliers.push_back(idx); } if ((int)inliers.size() bestCnt) { bestCnt (int)inliers.size(); bestCircle cand; bestInliers.swap(inliers); // 早期终止内点比例已经足够高 if (bestCnt minInliers) break; } } return bestCircle; }这段代码有几个细节值得解释。第一随机数引擎用std::mt19937而不是rand()因为rand()的周期短且低位随机性差在高迭代次数下可能产生重复三元组。项目中若看到srand(time(NULL))的痕迹建议改成带种子的std::mt19937在调试阶段固定种子可以让算法结果可复现。第二半径上限检查被我加上了。真实工业场景中圆的半径通常在 1 到 5000 像素之间如果三点近似共线拟合出的圆心在百万像素距离外这种候选模型不可能代表真实几何。工程代码里保留这个检查能省掉大量无效的内点计数。第三bestInliers通过swap赋值而不是拷贝减少了大向量复制开销。当点集有几十万个点、迭代几千次时这个优化是实打实的。换用std::move也可以但swap在循环里语义更清晰。4.2 找到 n 个最佳拟合圆的内点移除策略标题里「n 个最佳拟合圆」有两种理解一是拟合一个圆最少需要 n 个采样点n3二是从数据中分离出 n 个不同的圆。从工程角度多圆场景更常见——电路板上有多个圆形焊盘或者工件上有多个定位孔。策略很直接跑一次 RANSAC 找到最佳圆记录其内点然后把这部分内点从数据集中剔除对剩余点再次执行 RANSAC重复直到找到指定数量的圆或剩余点数不足。// 多圆提取反复执行单模型 RANSAC 并移除内点 std::vectorCCircle extractMultipleCircles( std::vectorCPoint2f pts, double threshold, int maxIterations, int minInliers, int targetCircleCount, unsigned int seedBase 100) { std::vectorCCircle circles; std::vectorint inliers; for (int idx 0; idx targetCircleCount; idx) { if ((int)pts.size() 3) break; // 每轮用不同种子避免重复采样路径 CCircle c ransacCircle(pts, threshold, maxIterations, minInliers, inliers, seedBase idx); // 内点太少说明已经提不出有效圆 if ((int)inliers.size() minInliers) break; // 用 Kåsa 方法基于内点精修 std::vectorCPoint2f inlierPts; inlierPts.reserve(inliers.size()); for (int i : inliers) inlierPts.push_back(pts[i]); CCircle refined fitCircleKasa(inlierPts); circles.push_back(refined); // 移除已分配的内点保留离群点参与下一轮 std::vectorCPoint2f remaining; remaining.reserve(pts.size() - inliers.size()); std::vectorbool removeFlag(pts.size(), false); for (int i : inliers) removeFlag[i] true; for (size_t i 0; i pts.size(); i) { if (!removeFlag[i]) remaining.push_back(pts[i]); } pts.swap(remaining); } return circles; }这里有个关键选择移除内点时是否把离群点保留下来我的做法是只移除内点、保留离群点因为本轮未被某个圆接受的点有可能是另一个圆的真实成员尤其当两个圆靠得较近、阈值设置偏严时。如果开头就把所有非内点当垃圾删掉两个交叉圆中靠后的那个基本无法被识别。4.3 自适应迭代次数估计固定maxIterations不够优雅因为内点比例 $\varepsilon$ 未知。每次采样 3 个点若想以概率 $p$ 保证至少有一次采样全部落在内点上迭代次数 $K$ 满足$$ K \ge \frac{\log(1-p)}{\log(1-\varepsilon^3)} $$把这个公式落地成代码后当数据质量好内点比例高时算法会自动提前收敛当数据质量差时迭代次数也不会拍脑袋乱设。int computeAdaptiveIterations(double inlierRatio, double confidence, int sampleSize 3) { double eps std::max(inlierRatio, 1e-6); double denominator std::log(1.0 - std::pow(eps, sampleSize)); if (std::fabs(denominator) 1e-12) return 1; int k (int)std::ceil(std::log(1.0 - confidence) / denominator); return std::max(k, 10); // 保底下限 }实践中可以先用 50 次快速采样估算内点比例再代入这个函数求出总迭代次数。比如内点比例 0.6、置信度 0.99 时大约需要 $\log(0.01)/\log(1-0.216) \approx 19$ 次内点比例降到 0.3同样置信度需要约 148 次。参数含义和推荐范围汇总在下表。参数含义推荐范围设定依据threshold点到圆轮廓的最大允许误差0.5 ~ 5.0像素由噪声标准差 $\sigma$ 决定通常取 $3\sigma$maxIterations采样次数上限500 ~ 5000内点比例越低需要越大minInliers接受一个模型的最小内点数点集总数的 $10% \sim 30%$低于 10% 的模型基本是随机凑出来的confidence至少找到一次全内点采样的概率0.99固定值工程上不需要再高5. 用合成数据检验阈值、噪声与多圆分离的调参实战5.1 构造可控的测试基准写一个生成器合成三个互不重叠的圆加高斯噪声和随机离群点这是验证 RANSAC 实现正确性最直接的办法。有了合成数据就可以把「期望找 3 个圆」作为 ground truth量化检测成功率和圆心偏差。#include random // 生成一个圆心 (cx, cy)、半径 r、噪声 sigma 的圆环点集 // 外加 ratio 比例的全局随机离群点 void generateSyntheticData(std::vectorCPoint2f pts, double cx, double cy, double r, int pointsOnCircle, double sigma, int outliers, unsigned int seed) { std::mt19937 rng(seed); std::normal_distributiondouble noise(0.0, sigma); std::uniform_real_distributiondouble angleDist(0.0, 2.0 * 3.141592653589793); for (int i 0; i pointsOnCircle; i) { double theta angleDist(rng); double x cx (r noise(rng)) * std::cos(theta); double y cy (r noise(rng)) * std::sin(theta); pts.emplace_back(x, y); } std::uniform_real_distributiondouble box(-200.0, 200.0); for (int i 0; i outliers; i) { pts.emplace_back(box(rng), box(rng)); } }建议把含离群点的数据写入文件后再用项目里的main读取这样能验证loadPointsFromFile的解析逻辑。我在复现时习惯性会看一眼pts_to_ransac.txt的点数规模和数据范围确认源文件的量纲后再设阈值而不是盲目套用代码里的默认参数。5.2 阈值设置的统计学视角阈值直接决定内点集合的大小。取太小真实圆上的点会被大量标记为外点最佳模型的内点数不足取太大远距离的离群点混入精拟合结果被污染。常用的经验法则是把阈值设为点到圆距离标准差的 3 倍但如果想更精准可以从卡方分布的角度推导对于一元高斯噪声距离差的平方服从自由度为 1 的卡方分布95% 置信水平对应约 $1.96 \sigma$99% 对应约 $2.58 \sigma$。下表给出了不同噪声水平下建议的阈值。噪声标准差 σ理论 99% 阈值 (2.58σ)工程推荐值0.2 像素0.520.60.5 像素1.291.51.0 像素2.583.02.0 像素5.165.5RANSAC 对阈值的敏感度是「先平后陡」阈值从一个较小值逐步放大时内点数增长会经历一个平台期平台期之后继续增大离群点开始大量混入。实际操作中我会以 0.2 像素为步长扫描 0.5 到 4.0 的阈值区间观察内点数曲线的拐点拐点对应的阈值就是当前数据的最优设置。5.3 多圆场景的分离技巧与失败模式多圆提取最怕的是两个圆靠得过近内点区域交叠。这时的处理技巧是第一轮 RANSAC 找到主圆后不要急着移除全部内点而是把距离该圆轮廓 2 倍阈值以内的点打上「暂缓」标记第二轮优先从未标记点中采样。具体到代码就是给extractMultipleCircles增加一个softRemove模式把硬删除改为权重降低。另一个常见失败模式是半径退化。离群点分布均匀时三个离群点偶尔能凑出一个看似合理的圆虽然内点数不多但在内点比例极低的数据上这一小撮假内点可能超过真实圆的得分。排查方法是看模型得分分布真实圆的得分通常远高于随机组合如果两个候选圆的内点数接近多半是阈值过松或者数据本身不存在明显圆形结构。最后提一个通用排查技巧把 RANSAC 的每一次迭代中间结果采样点、拟合圆、内点数写入日志然后用 Python 脚本或矢量工具可视化。你会发现失败的迭代往往集中在初始几次采样时因半径上限检查被跳过、或是三点中混入了离群点导致模型偏移——这些问题一旦可视化几分钟内就能定位比盯着源码猜要快得多。本文还有配套的精品资源点击获取
返回列表