引导线的提取预处理(3)——颜色阈值、霍夫直线及其聚类

主要流程

分为三个步骤

1、颜色通道分离,引导线(红色)颜色提取(详情参考前面的博文);
2、霍夫直线检测;
3、霍夫直线聚类;

霍夫直线

检测函数

  • 第一个检测函数
void HoughLines(InputArray image, OutputArray lines, 
				double rho, double theta, int threshold, 
				double srn=0, double stn=0, 
				double min_theta = 0, double max_theta = CV_PI )

参数解释:
第一个参数:InputArray类型的image,输入图像,即源图像,需为8位的单通道二进制图像,可以将任意的源图载入进来后由函数修改成此格式后,再填在这里。
第二个参数:InputArray类型的lines,经过调用HoughLines函数后储存了霍夫线变换检测到线条的输出矢量。每一条线由具有两个元素的矢量表示,其中,是离坐标原点((0,0)(也就是图像的左上角)的距离。 是弧度线条旋转角度(0垂直线,π/2水平线)。
第三个参数:double类型的rho,以像素为单位的距离精度。另一种形容方式是直线搜索时的进步尺寸的单位半径。PS:Latex中/rho就表示 。
第四个参数:double类型的theta,以弧度为单位的角度精度。另一种形容方式是直线搜索时的进步尺寸的单位角度。
第五个参数:int类型的threshold,累加平面的阈值参数,即识别某部分为图中的一条直线时它在累加平面中必须达到的值。大于阈值threshold的线段才可以被检测通过并返回到结果中。
第六个参数:double类型的srn,有默认值0。对于多尺度的霍夫变换,这是第三个参数进步尺寸rho的除数距离。粗略的累加器进步尺寸直接是第三个参数rho,而精确的累加器进步尺寸为rho/srn。
第七个参数:double类型的stn,有默认值0,对于多尺度霍夫变换,srn表示第四个参数进步尺寸的单位角度theta的除数距离。且如果srn和stn同时为0,就表示使用经典的霍夫变换。否则,这两个参数应该都为正数。
第八个参数:double类型的 min_theta,对于标准和多尺度Hough变换,检查线条的最小角度。必须介于0和max_theta之间。
第九个参数:double类型的 max_theta, 对于标准和多尺度Hough变换,检查线条的最大角度。必须介于min_theta和CV_PI之间.

  • 第二个检测函数

直线描述

霍夫直线参数含义

霍夫直线参数位rho,theta。
其中,
rho代表图像原点至直线的距离;
theta代表直线的角度,垂直为0度,顺时针旋转逐渐增大,水平位90度。
theta的取值范围:大于等于0度, 小于180度。

数学方程表示直线

已知: rho,theta;
直线系数:a = cos(theta), b = sin(theta);
垂点坐标:x0 = a * rho, y0 = b * rho;
直线方程:x = x0 - bt, y = y0 + at;

直线反侧变换

标准直线参数的theta范围位[0,PI),在直线之间进行运算时会用到theta在此范围外的值,例如将1度的直线变换为181度,或者将179度的直线变换为-1度,这种变换成为反侧变换,反侧变换前后的直线等效,为同一直线。

  • 反侧变换一:
new_theta = old_theta - PI;
new_rho = - old_rho;
  • 反侧变换二:
new_theta = old_theta + PI;
new_rho = - old_rho;
直线规范化操作

将直线参数的角度值控制在[0,PI)范围内。

while (center_theta > PI)
{
	center_theta -= PI;
	center_rho = - center_rho;
}

while (center_theta < 0)
{
	center_theta += PI;
	center_rho = - center_rho;
}
霍夫直线间求差

在霍夫直线聚类时,需要将直线参数看作二维点,在此计算二维点之间的距离。
设两直线参数为(rho1, theta1)、(rho2, theta2),求取绝对误差:

error_theta = abs(theta1-theta2);
if error_theta <= PI/2 时:
	error_rho = abs(rho1 - rho2);
if error_theta > PI/2 时:
	error_theta = PI - error_theta;
	error_rho = abs(rho1 + rho2);

计算归一化直线误差(距离):

d0 = drho/Rho_max;
d1 = dtheta/(CV_PI*0.5);
error_cj = sqrtf(d0*d0 + d1*d1);   
霍夫直线求和与均值

在霍夫直线聚类求取中心点时,需要计算多直线参数累加和、均值。

  • 计算累加和:
    当样本直线距离该聚类中心直线的角度差大于90度时,需要做反侧变换,伪代码描述过程如下。
sum_theta = 0;
sum_rho = 0;
this_cluster_num = 0
for every lines in this cluster center:
	this_cluster_num ++;
	# 直线角度距离大于90度,则做反侧变换
	if (this_theta - center_theta) > PI/2:
		this_theta = this_theta - PI;
		this_rho = - this_rho;
	else if (this_theta - center_theta) < - PI/2:
		this_theta = this_theta + PI;
		this_rho = - this_rho;
	sum_theta += this_theta;
	sum_rho += this_rho;
  • 求取均值(中心点)
center_theta = sum_theta/this_cluster_num;
center_rho = sum_rho/this_cluster_num;
  • 结果规范化操作
    为了保持中心直线参数在多伦的聚类迭代运算时,始终保持在一定的范围内,需要对其进行常规化操作,将其范围控制在[0,PI)之间。
霍夫直线滤波

在直线跟踪时,为了保证目标直线的平滑行,需要加入滤波处理。
普通一阶低通滤波的方法: y = y_1 + k*(x - y_1);
霍夫直线滤波时需要考虑直线反侧变换,过程用伪代码表示:

# 偏差检测与反侧变换
error_theta = x_theta - y_theta;	
if (error_theta > PI/2):
	# 直线反侧变换
	x_u_theta = x_theta - PI;
	x_u_rho = - x_rho;
else if (error_theta < -PI/2):
	# 直线反侧变换
	x_u_theta = x_theta + PI;
	x_u_rho = - x_rho;
else:
	x_u_theta = x_theta;
	x_u_rho = -x_rho;

# 滤波处理
y_line = y_line + k*(x_u_line - y_line);

霍夫直线聚类

直线聚类的算法流程

聚类中心数量cluster_center_num = 1;
while 聚类未完成:	
	for 聚类N次数:
		聚类初始化工作;
		随机分布聚类中心点:在数据样本集中随机选取cluster_center_num 个点,作为数据中心点;
		while 本次聚类未达到收敛状态:
			根据当前聚类中心值,计算每一个点归属于哪一个聚类中心;
			根据归属关系,计算新的聚类中心值;
			当新的聚类中心与上次聚类中心一致时,本次聚类达到收敛状态;
		计算本次聚类完成后的样本点均方差,并记录;
	选取N次聚类结果中,均方差最小的那一次数据,作为聚类结果。
	若聚类均方差小于设定阈值,则聚类完成;否则,cluster_center_num  ++,开启新的N轮聚类。

随机分布聚类中心点的方法:
(1)在取值空间中随机取值坐标点,产生cluster_center_num 个数据点;
(2)在数据样本集中随机选取cluster_center_num 个点,作为数据中心点;
如果样本点在空间中分布非常不均匀的情况下,方法(1)产生的数据中心点有些偏离数据较远,则很难达到较好的收敛状态,即某些聚类中心可能没有隶属样本点。

需要考虑的问题:

  • 聚类的总类别个数不确定
    因为事先不知道图形中存在多少真实的直线,所以要自适应寻找聚类的类别数,判断依据是聚类后的均方差是否足够小;
  • 霍夫直线参数误差计算时的特殊处理
    (1) 角度处理
    由霍夫直线参数物理含义可知,1度的直线与179度的直线,在斜率上仅仅差别2度。所以,
    当两直线的误差error = | theta1 - theta2 |大于90度时,
    其真实误差是| PI - error|。
    (2) 距离处理
    当两个直线的角度偏差计算经过折算处理时,其斜率计算也需要进行折算处理。

效果展示

图1

原图

在这里插入图片描述

霍夫直线检测与聚类

在这里插入图片描述

原图添加效果

在这里插入图片描述

图2

原图

在这里插入图片描述

霍夫直线检测与聚类

在这里插入图片描述

原图添加效果

在这里插入图片描述

图3

原图

在这里插入图片描述

霍夫直线检测与聚类

在这里插入图片描述

原图添加效果

在这里插入图片描述

程序示例

程序主体

#include <opencv2/opencv.hpp>
#include <opencv2/xfeatures2d.hpp>
#include <iostream>
 
using namespace cv;
using namespace cv::xfeatures2d;         // 不要忘了导入扩展模块
using namespace std;
 
Mat src_img, gray_img;
const string output_name = "二值化";
int minHessian = 100;                // 定义SURF中的hessian阈值特征点检测算子
int max_value = 250;

Mat bool_img;

// =================================================
//      颜色通道分离 + 霍夫直线检测 + 直线绘制
// =================================================
int Rho_max = 500;
void cluster_lines(vector<Vec2f> lines, vector<Vec2f> *lines_center);
string Name = "6";
int main() 
{
    // Mat src = imread("/home/parobot/Pictures/rgb.jpeg");
    // resize(src, src_img, Size(round(1.0*src.cols), round(1.0*src.rows)));
    Mat src = imread("/home/parobot/Pictures/line/" + Name + ".jpg");
	resize(src, src_img, Size(round(0.15*src.cols), round(0.15*src.rows)));
    // resize(src, src_img, Size(round(1.0*src.cols), round(1.0*src.rows)));
	
	if (src_img.empty())
	{
		printf("could not load the image...\n");
		return -1;
	}
	namedWindow("原图", WINDOW_AUTOSIZE);
	imshow("原图", src_img);
	cvtColor(src_img, gray_img, COLOR_BGR2GRAY);

    vector<Mat> channels;
    split(src_img, channels);
    // imshow("red",   channels.at(2));
    // imshow("green", channels.at(1));
    // imshow("blue",  channels.at(0));
    cout << channels.at(0).size << endl;

    Rho_max = sqrtf(src_img.rows*src_img.rows + src_img.cols*src_img.cols);
    
    for (int i = 0; i < channels.at(0).rows; i ++)
    {
        for (int j = 0; j < channels.at(0).cols; j ++)
        {
            // 提取红色
            // 红色通道值,应该大于其它通道值           
            int red = channels.at(2).at<uchar>(i,j) - 30;
            if ((red < channels.at(0).at<uchar>(i,j))
                ||(red < channels.at(1).at<uchar>(i,j)))
            {
                channels.at(0).at<uchar>(i,j) = 0;
                channels.at(1).at<uchar>(i,j) = 0;  
                channels.at(2).at<uchar>(i,j) = 0;
            }  
        }
    }

    Mat red_line;
    merge(channels, red_line);

    // 霍夫直线检测
    vector<Vec2f> lines;
    HoughLines(channels.at(2), lines, 1, CV_PI/180, 200, 0, 0);
    cout << lines.size() << endl;
    for (int i = 0; i < lines.size(); i ++)
    {
        float rho = lines[i][0], theta = lines[i][1];
        Point p1, p2;
        float a = cos(theta), b = sin(theta);
        float x0 = a*rho, y0 = b*rho;
        p1.x = cvRound(x0 + 1000*(-b));
        p1.y = cvRound(y0 + 1000*a);
        p2.x = cvRound(x0 - 1000*(-b));
        p2.y = cvRound(y0 - 1000*a);
        line(red_line, p1, p2, Scalar(255, 0, 0), 1, LINE_AA);
    }
    
    // 概率霍夫直线检测
    // vector<Vec4i> lines;
    // HoughLinesP(channels.at(2), lines, 1, CV_PI/180, 80, 50, 10);
    // cout << lines.size() << endl;
    // for (int i = 0; i < lines.size(); i ++)
    // {
    //    Vec4i l = lines[i];
    //     line(red_line, Point(l[0],l[1]), Point(l[2], l[3]), Scalar(255, 0, 0), 1, LINE_AA);
    // }

    imwrite("/home/parobot/Pictures/line/orgin" + Name + ".jpg", src_img);

    // 聚类
    vector<Vec2f> lines_center;
    cluster_lines(lines, &lines_center);
    cout << "lines_center.size : " << lines_center.size() << endl;
    for (int i = 0; i < lines_center.size(); i ++)
    {
        float rho = lines_center[i][0], theta = lines_center[i][1];
        Point p1, p2;
        float a = cos(theta), b = sin(theta);
        float x0 = a*rho, y0 = b*rho;
        p1.x = cvRound(x0 + 1000*(-b));
        p1.y = cvRound(y0 + 1000*a);
        p2.x = cvRound(x0 - 1000*(-b));
        p2.y = cvRound(y0 - 1000*a);
        line(red_line, p1, p2, Scalar(0, 255, 0), 2, LINE_AA);
        line(src_img, p1, p2, Scalar(0, 255, 0), 2, LINE_AA);        
    }
    imshow("cluster",  red_line);
    imshow("result",  src_img);
    
    imwrite("/home/parobot/Pictures/line/cluster" + Name + ".jpg", red_line);
    imwrite("/home/parobot/Pictures/line/result" + Name + ".jpg", src_img);
    
	waitKey(0);
	return 0;
}

聚类算法

void cluster_lines(vector<Vec2f> lines, vector<Vec2f> *lines_center)
{
    int cluster_num = 1;
    bool cluster_done = false;
    vector<Vec2f> center_this;
    vector<Vec2f> center_best;

    vector<vector<Vec2f>> center_episode;
    center_episode.resize(3);

    while (!cluster_done)
    {
        cout << " ================================================ " << endl;
        cout << "新的一次聚类中心点数量, cluster_num = " << cluster_num << endl;
        center_this.resize(cluster_num);
        center_best.resize(cluster_num);

        float error_min_rec = 0.0;
        float error_points[lines.size()];
        float error_center[cluster_num];
        float error_center_sum_min;

        // 多次聚类取最优
        int episode_best_index = 0;
        for (int episode = 0; episode < 3; episode ++)
        {
            cout << " ------------------------------------------ " << endl;
            cout << "cluster begin, episode = " << episode << endl;
            // 初始化:给中心点赋随机值 
            // for (int i = 0; i < cluster_num; i ++)
            // {
            //     unsigned seed = time(0);
            //     srand(seed);
            //     int rho = rand()%(Rho_max + 1);
            //     int angle = rand()%180;

            //     if (angle > 90)
            //     {
            //         center_this[i][0] =  -(float)rho;
            //     }
            //     else
            //     {
            //         center_this[i][0] =  (float)rho;
            //     }
            //     center_this[i][1] =  ((float)angle)*CV_PI/180.0;                
            // }
            unsigned seed = time(0);
            srand(seed);
            for (int i = 0; i < cluster_num; i ++)
            {                
                int index = rand()%lines.size();
                center_this[i] = lines[index];
                cout << "center_this[" << i << "]" << center_this[i] << endl;
            }

            // 循环多次,直到聚类中心收敛到稳态
            bool center_stable = false;
            int Point_Belong[lines.size()];
            int Center_Points_Num[cluster_num];
            //while (!center_stable)
            for (int time = 0; time < 3; time ++)
            {         
                cout << "update the cluster center ... ... ..." << endl;
                // 新一轮聚类的初始化
                for (int cj = 0; cj < cluster_num; cj ++)
                {
                    Center_Points_Num[cj] = 0.0;
                }

                // 分类:遍历所有点,判断点归属于哪个中心                
                Vec2f Center_Points_Vaule_Sum[cluster_num];

                for (int pi = 0; pi < lines.size(); pi ++)
                {
                    float moni_d0_min , moni_d1_min;
                    float error_min = 1.0e10;
                    // 计算中心距离:遍历所有中心点
                    for (int cj = 0; cj < cluster_num; cj ++)
                    {
                        // 折算最近方向的距离
                        float dtheta = fabsf(lines[pi][1] - center_this[cj][1]);
                        float drho;
                        if (dtheta > CV_PI*0.5)
                        {
                            dtheta = fabsf(CV_PI - dtheta);
                            drho = fabsf(lines[pi][0] + center_this[cj][0]);
                        }
                        else
                        {
                            drho = fabsf(lines[pi][0] - center_this[cj][0]);
                        }

                        // 计算归一化距离
                        float d0 = drho/Rho_max;
                        float d1 = dtheta/(CV_PI*0.5);

                        float error_cj = sqrtf(d0*d0 + d1*d1);   
                        if (error_min > error_cj)
                        {
                            error_min = error_cj;
                            Point_Belong[pi] = cj;

                            moni_d0_min = d0;
                            moni_d1_min = d1;
                        }                    
                    }

                    error_points[pi] = error_min;
                    
                    cout << "d0 = " << moni_d0_min     << ", d1 = " << moni_d1_min << endl;

                    // 计算累加和、个数
                    Center_Points_Num[Point_Belong[pi]] ++;
                    // 计算累加和:需要折算
                    Vec2f re_line = lines[pi];
                    float theta_err = lines[pi][1] - center_this[Point_Belong[pi]][1];
                    if (theta_err > CV_PI*0.5)
                    {
                        re_line[1] -= CV_PI;
                        re_line[0] = - lines[pi][0];
                    }
                    else if (theta_err < -CV_PI*0.5)
                    {
                        re_line[1] += CV_PI;
                        re_line[0] = - lines[pi][0];
                    }
                    
                    Center_Points_Vaule_Sum[Point_Belong[pi]] += re_line;

                    // cout << "lines[" << pi << "] " << lines[pi] << endl;
                    // cout << "Center_Points_Vaule_Sum[" << Point_Belong[pi]<< "] " << Center_Points_Vaule_Sum[Point_Belong[pi]] << endl;
                }

                // 更新中心点坐标
                bool center_change = false;
                for (int cj = 0; cj < cluster_num; cj ++)
                {
                    cout << "            center_this[" << cj << "] " << center_this[cj] << endl;
                    cout << "Center_Points_Vaule_Sum[" << cj << "] " << Center_Points_Vaule_Sum[cj] << endl;
                    cout << "      Center_Points_Num[" << cj << "] " << Center_Points_Num[cj] << endl;
                    Vec2f center_temp = Center_Points_Vaule_Sum[cj]/Center_Points_Num[cj];
                    if (center_this[cj] != center_temp)
                    {
                        center_change = true;
                        center_this[cj] = center_temp;
                    }
                    cout << "               center_temp " << center_temp << endl;
                    
                }
                center_stable = !center_change;
            }

            cout << "this elipse cluster is done " << center_stable << endl;
            // 记录本次
            center_episode[episode] = center_this;
            

            // 更新方差            
            for (int cj = 0; cj < cluster_num; cj ++)
            {
                error_center[cj] = 0;
            }
            for (int pi = 0; pi < lines.size(); pi ++)
            {
                error_center[Point_Belong[pi]] += error_points[pi]*error_points[pi]/Center_Points_Num[Point_Belong[pi]];
            }

            // 选择最小方差
            float error_center_sum = 0;
            for (int cj = 0; cj < cluster_num; cj ++)
            {
                error_center_sum += error_center[cj];
                cout << "error_center[" << cj << "] " << error_center[cj] << endl;
            }
            if (episode == 0)
            {
                error_center_sum_min = error_center_sum;
                episode_best_index = 0;
            }
            else
            {
                if (error_center_sum_min < error_center_sum)
                {
                    error_center_sum_min = error_center_sum;
                    episode_best_index = episode;
                }             
            }            
        }

        center_best = center_episode[episode_best_index];

        if((error_center_sum_min/cluster_num) < 0.05)
        {
            cluster_done = true;           
        }
        else
        {
            cluster_num ++;
        }

        cout << "error_center_sum_min = " << error_center_sum_min << ", avg = " << error_center_sum_min/cluster_num << endl;

    }
    
    lines_center->resize(center_best.size());
    (*lines_center) = center_best;

    for (int i = 0; i < cluster_num; i ++)
    {
        cout << "center_best[i] " << center_best[i] << endl;
    }
    
}
import cv2 as cv import numpy as np #直线检测 #使用霍夫直线变换做直线检测,前提条件:边缘检测已经完成 #标准霍夫线变换 def line_detection(image): gray = cv.cvtColor(image, cv.COLOR_RGB2GRAY) edges = cv.Canny(gray, 50, 150) #apertureSize参数默认其实就是3 cv.imshow("edges", edges) #cv.HoughLines参数设置:参数1,灰度图像;参数二,以像素为单位的距离精度(一般都是1,进度高,但是速度会慢一点) #参数三,以弧度为单位的角度精度(一般是1rad);参数四,阈值,大于阈值threshold的线段才可以被检测通过并返回到结果中 #该函数返回值为rho与theta lines = cv.HoughLines(edges, 1, np.pi/180, 200) for line in lines: rho, theta = line[0] #line[0]存储的是点到直线的极径和极角,其中极角是弧度表示的。 a = np.cos(theta) #theta是弧度 b = np.sin(theta) x0 = a * rho #代表x = r * cos(theta) y0 = b * rho #代表y = r * sin(theta) x1 = int(x0 + 1000 * (-b)) #计算直线起点横坐标 y1 = int(y0 + 1000 * a) #计算起始起点纵坐标 x2 = int(x0 - 1000 * (-b)) #计算直线终点横坐标 y2 = int(y0 - 1000 * a) #计算直线终点纵坐标 注:这里的数值1000给出了画出的线段长度范围大小,数值越小,画出的线段越短,数值越大,画出的线段越长 cv.line(image, (x1, y1), (x2, y2), (0, 0, 255), 2) #点的坐标必须是元组,不能是列表。 cv.imshow("image-lines", image) #统计概率霍夫线变换 def line_detect_possible_demo(image): gray = cv.cvtColor(image, cv.COLOR_RGB2GRAY) edges = cv.Canny(gray, 50, 150, apertureSize=3) # apertureSize参数默认其实就是3 lines = cv.HoughLinesP(edges, 1, np.pi / 180, 60, minLineLength=60, maxLineGap=5) for line in lines: x1, y1, x2, y2 = line[0] cv.line(image, (x1, y1), (x2, y2), (0, 0, 255), 2) cv.imshow("line_detect_possible_demo",image) src = cv.imread("E:/opencv/picture/track.jpg") print(src.shape) cv.namedWindow('input_image', cv.WINDOW_AUTOSIZE) cv.imshow('input_image', src) line_detection(src) src = cv.imread("E:/opencv/picture/track.jpg") #调用上一个函数后,会把传入的src数组改变,所以调用下一个函数时,要重新读取图片 line_detect_possible_demo(src) cv.waitKey(0) cv.destroyAllWindows() 霍夫检测直线原理: 关于hough变换,核心以及难点就是关于就是有原始空间到参数空间的变换上。以直线检测为例,假设有一条直线L,原点到该直线的垂直距离为p,垂线与x轴夹角为θθ,那么这条直线是唯一的,且直线的方程为 ρ=xcosθ+ysinθρ=xcosθ+ysinθ, 如下图所
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值