OpenCV 真圆度测量
最近一个项目需要在图像上测量一些小孔的真圆度。因此专门研究了一下真圆度计算问题。
对于一个轮廓,我们可以求出这个轮廓的外接圆和内切圆。这两个圆的半径差定义为真圆度。这个数值越小,表示这个圆越标准。
外接圆在 OpenCV 中有现成的函数来计算,但是内切圆是没有的。去算内切圆难度还是蛮大的。
为此,我采用了个变通的方法。对这个轮廓先计算最小二乘拟合圆。之后计算轮廓上的各个点到圆心的距离。最大距离减去最小距离的差可以作为真圆度的一个近似的计算。
真圆度的计算代码封装到了一个类里。这个类的头文件如下:
#ifndef TRUEROUNDNESS_H
#define TRUEROUNDNESS_H
#include <opencv2/opencv.hpp>
/**
* @brief The TrueRoundness class 真圆度计算。只实现了一种最基本的算法。
*/
class TrueRoundness
{
public:
/**
* @brief TrueRoundness 真圆度计算的类的构造函数。考虑到除去噪声影响,可以排除部分数据点。
* @param top 计算真圆度时,最大直径排除前 top 个点。
* @param bottom 计算真圆度时,最小直径排除最后 bottom 个点。
*/
TrueRoundness(unsigned int top = 0, unsigned int bottom = 0);//可以去除两头的一部分数据
/**
* @brief TrueRoundness 真圆度计算的类的构造函数。考虑到除去噪声影响,可以排除部分数据点。
* @param top_percent 计算真圆度时,最大直径排除前 top_percent % 的点。
* @param bottom_percent 计算真圆度时,最小直径排除最后 bottom_percent % 的点。
*/
TrueRoundness(double top_percent, double bottom_percent); //去除两头一定比例的数据
/**
* @brief operator ()计算真圆度。首先计算出拟合圆,然后计算轮廓上到点到圆心的最远距离和最近距离。之差就是真圆度。
* @param [in] contour 轮廓上的点。
* @param [out] r_min
* @param [out] r_max
* @return 返回真圆度。如果结果小于 0 ,则表示计算失败。通常是拟合圆失败。
*/
double operator() (const std::vector<cv::Point2i> &contour, double &r_min, double &r_max);
double operator() (const std::vector<cv::Point2f> &contour, double &r_min, double &r_max);
/**
* @brief circleLeastFit 最小二乘法拟合圆。这个函数是辅助函数,通常不需要使用。
* @param points 圆的轮廓上的点
* @param center_x 计算出的圆心坐标 x
* @param center_y 计算出的圆心坐标 y
* @param radius 计算出的半径值
* @return true 表示拟合成功。false 表示失败。
*/
bool circleLeastFit(const std::vector<cv::Point2i> &points, double ¢er_x, double ¢er_y, double &radius);
bool circleLeastFit(const std::vector<cv::Point2f> &points, double ¢er_x, double ¢er_y, double &radius);
/**
* @brief sort 辅助函数,对边缘数据进行排序,离圆心远的数据排在前面。
* @param contour 轮廓上的点
* @param center_x 圆心坐标 x
* @param center_y 圆心坐标 y
*/
void sort( std::vector<cv::Point2i> &contour, double center_x, double center_y );
void sort( std::vector<cv::Point2f> &contour, double center_x, double center_y );
/**
* @brief distance 辅助函数,计算一个点到圆心的距离
* @param p
* @param center_x
* @param center_y
* @return
*/
double distance(cv::Point2i p, double center_x, double center_y);
double distance(cv::Point2f p, double center_x, double center_y);
private:
int m_type = 0;
int m_top = 0;
int m_bottom = 0;
double m_top_percent = 0;
double m_bottom_percent = 0;
};
#endif // TRUEROUNDNESS_H
实现代码如下:
#include "TrueRoundness.h"
TrueRoundness::TrueRoundness(unsigned int top, unsigned int bottom)
{
m_type = 0;
m_top = top;
m_bottom = bottom;
}
TrueRoundness::TrueRoundness(double top_percent, double bottom_percent)
{
m_type = 1;
m_top_percent = top_percent;
m_bottom_percent = bottom_percent;
}
double TrueRoundness::operator()(const std::vector<cv::Point2i> &contour, double &r_min, double &r_max)
{
double center_x;
double center_y;
double radius;
bool ret = circleLeastFit(contour, center_x, center_y, radius);
if(!ret) return -1;
std::vector<cv::Point2i> newcontour = contour;
sort(newcontour, center_x, center_y);
int N = newcontour.size();
if(m_type == 0)
{
if(m_top > N || m_bottom > N)
{
return -2;
}
int N1 = m_top;
int N2 = N - 1 - m_bottom;
if(distance(newcontour[N1], center_x, center_y) < radius || distance(newcontour[N2], center_x, center_y) > radius)
{
return -3;
}
r_max = distance(newcontour[N1], center_x, center_y);
r_min = distance(newcontour[N2], center_x, center_y);
return r_max - r_min ;
}
else if(m_type == 1)
{
if(m_top_percent >= 1 || m_bottom_percent >= 1) return -4;
int N1 = 0, N2 = 0;
for(int i = 0; i < N; i++)
{
double r = distance(newcontour[i], center_x, center_y);
if(r > radius) N1 ++;
if(r < radius) N2 ++;
}
N1 = N1 * (m_top_percent);
N2 = N - 1 - N2 * (m_bottom_percent);
r_max = distance(newcontour[N1], center_x, center_y);
r_min = distance(newcontour[N2], center_x, center_y);
return r_max - r_min ;
}
return -5;
}
double TrueRoundness::operator()(const std::vector<cv::Point2f> &contour, double &r_min, double &r_max)
{
double center_x;
double center_y;
double radius;
bool ret = circleLeastFit(contour, center_x, center_y, radius);
if(!ret) return -1;
std::vector<cv::Point2f> newcontour = contour;
sort(newcontour, center_x, center_y);
int N = newcontour.size();
if(m_type == 0)
{
if(m_top > N || m_bottom > N)
{
return -2;
}
int N1 = m_top;
int N2 = N - 1 - m_bottom;
if(distance(newcontour[N1], center_x, center_y) < radius || distance(newcontour[N2], center_x, center_y) > radius)
{
return -3;
}
r_max = distance(newcontour[N1], center_x, center_y);
r_min = distance(newcontour[N2], center_x, center_y);
return r_max - r_min ;
}
else if(m_type == 1)
{
if(m_top_percent >= 1 || m_bottom_percent >= 1) return -4;
int N1 = 0, N2 = 0;
for(int i = 0; i < N; i++)
{
double r = distance(newcontour[i], center_x, center_y);
if(r > radius) N1 ++;
if(r < radius) N2 ++;
}
N1 = N1 * (m_top_percent);
N2 = N - 1 - N2 * (m_bottom_percent);
r_max = distance(newcontour[N1], center_x, center_y);
r_min = distance(newcontour[N2], center_x, center_y);
return r_max - r_min ;
}
return -5;
}
bool TrueRoundness::circleLeastFit(const std::vector<cv::Point2i> &points, double ¢er_x, double ¢er_y, double &radius)
{
center_x = 0.0;
center_y = 0.0;
radius = 0.0;
if (points.size() < 3)
{
return false;
}
double sum_x = 0.0, sum_y = 0.0;
double sum_x2 = 0.0, sum_y2 = 0.0;
double sum_x3 = 0.0, sum_y3 = 0.0;
double sum_xy = 0.0, sum_x1y2 = 0.0, sum_x2y1 = 0.0;
int N = points.size();
for (int i = 0; i < N; i++)
{
double x = points[i].x;
double y = points[i].y;
double x2 = x * x;
double y2 = y * y;
sum_x += x;
sum_y += y;
sum_x2 += x2;
sum_y2 += y2;
sum_x3 += x2 * x;
sum_y3 += y2 * y;
sum_xy += x * y;
sum_x1y2 += x * y2;
sum_x2y1 += x2 * y;
}
double C, D, E, G, H;
double a, b, c;
C = N * sum_x2 - sum_x * sum_x;
D = N * sum_xy - sum_x * sum_y;
E = N * sum_x3 + N * sum_x1y2 - (sum_x2 + sum_y2) * sum_x;
G = N * sum_y2 - sum_y * sum_y;
H = N * sum_x2y1 + N * sum_y3 - (sum_x2 + sum_y2) * sum_y;
a = (H * D - E * G) / (C * G - D * D);
b = (H * C - E * D) / (D * D - G * C);
c = -(a * sum_x + b * sum_y + sum_x2 + sum_y2) / N;
center_x = a / (-2);
center_y = b / (-2);
radius = sqrt(a * a + b * b - 4 * c) / 2;
return true;
}
bool TrueRoundness::circleLeastFit(const std::vector<cv::Point2f> &points, double ¢er_x, double ¢er_y, double &radius)
{
center_x = 0.0;
center_y = 0.0;
radius = 0.0;
if (points.size() < 3)
{
return false;
}
double sum_x = 0.0, sum_y = 0.0;
double sum_x2 = 0.0, sum_y2 = 0.0;
double sum_x3 = 0.0, sum_y3 = 0.0;
double sum_xy = 0.0, sum_x1y2 = 0.0, sum_x2y1 = 0.0;
int N = points.size();
for (int i = 0; i < N; i++)
{
double x = points[i].x;
double y = points[i].y;
double x2 = x * x;
double y2 = y * y;
sum_x += x;
sum_y += y;
sum_x2 += x2;
sum_y2 += y2;
sum_x3 += x2 * x;
sum_y3 += y2 * y;
sum_xy += x * y;
sum_x1y2 += x * y2;
sum_x2y1 += x2 * y;
}
double C, D, E, G, H;
double a, b, c;
C = N * sum_x2 - sum_x * sum_x;
D = N * sum_xy - sum_x * sum_y;
E = N * sum_x3 + N * sum_x1y2 - (sum_x2 + sum_y2) * sum_x;
G = N * sum_y2 - sum_y * sum_y;
H = N * sum_x2y1 + N * sum_y3 - (sum_x2 + sum_y2) * sum_y;
a = (H * D - E * G) / (C * G - D * D);
b = (H * C - E * D) / (D * D - G * C);
c = -(a * sum_x + b * sum_y + sum_x2 + sum_y2) / N;
center_x = a / (-2);
center_y = b / (-2);
radius = sqrt(a * a + b * b - 4 * c) / 2;
return true;
}
void TrueRoundness::sort(std::vector<cv::Point2i> &contour, double center_x, double center_y)
{
auto cmp = [center_x, center_y](cv::Point2i p1, cv::Point2i p2)
{double r1 = hypot(p1.x - center_x, p1.y - center_y);
double r2 = hypot(p2.x - center_x, p2.y - center_y);
return r1 > r2;};
std::sort( contour.begin(), contour.end(), cmp );
}
void TrueRoundness::sort(std::vector<cv::Point2f> &contour, double center_x, double center_y)
{
auto cmp = [center_x, center_y](cv::Point2i p1, cv::Point2i p2)
{double r1 = hypot(p1.x - center_x, p1.y - center_y);
double r2 = hypot(p2.x - center_x, p2.y - center_y);
return r1 > r2;};
std::sort( contour.begin(), contour.end(), cmp );
}
double TrueRoundness::distance(cv::Point2i p, double center_x, double center_y)
{
return hypot(p.x - center_x, p.y - center_y);
}
double TrueRoundness::distance(cv::Point2f p, double center_x, double center_y)
{
return hypot(p.x - center_x, p.y - center_y);
}
下面是个简单的测试用例。首先我们需要个测试图片:
下面是测试代码:
#include <QCoreApplication>
#include <QDebug>
#include <opencv2/opencv.hpp>
#include "TrueRoundness.h"
bool circleLeastFit(const std::vector<cv::Point2i> &points, double ¢er_x, double ¢er_y, double &radius);
int main(int argc, char *argv[])
{
QCoreApplication a(argc, argv);
cv::Mat image;
image = cv::imread("E:\\testimage\\round.png", cv::IMREAD_COLOR);
cv::Mat gray;
cv::cvtColor(image, gray, cv::COLOR_BGRA2GRAY);
cv::Canny(gray ,gray, 40, 120);
std::vector<std::vector<cv::Point2i>> contours;
std::vector<cv::Vec4i> hierarchy;
cv::findContours(gray, contours, hierarchy, cv::RETR_EXTERNAL, cv::CHAIN_APPROX_NONE);
double center_x, center_y;
double radius, r_min, r_max;
TrueRoundness tr(1U, 1U);
tr.circleLeastFit(contours[0], center_x, center_y, radius);
double t = tr(contours[0], r_min, r_max);
cv::circle(image, cv::Point(center_x, center_y), r_min, cv::Scalar(0, 255, 0), 1);
cv::circle(image, cv::Point(center_x, center_y), r_max, cv::Scalar(0, 0, 255), 1);
qDebug() << t;
cv::imshow("circle", image);
return a.exec();
}
输出的结果如下: