一、基本原理
先利用高斯函数对图像进行低通滤波;然后对图像中的每个像素进行处理,寻找边缘的位置及在该位置的边缘法向,并采用一种称之为“非极值抑制”的技术在边缘法向寻找局部最大值;最后对边缘图像做滞后阈值化处理,消除虚假响应。
二、实现流程
2.1 计算步骤
- 先利用高斯平滑滤波器来平滑图像以除去噪声(即用高斯平滑滤波器与图像作卷积);
- 计算梯度的幅值和方向;
- 对梯度幅值进行非极大值抑制;
- 用双阈值检测和连接边缘。
2.2 流程图
三、实现代码
#include "gdal_priv.h"
#include "cpl_conv.h"
#include <iostream>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/opencv.hpp>
#include<cmath>
#include <string>
#include <vector>
//栈与队列
#include<stack>
#include<queue>
using namespace std;
using namespace cv;
#define PI 3.141592
typedef Point Cpoint;
int main(int argc, char** argv)
{
void canny(Mat image);
//将所有需要的库函数等,导入进来并改变位数为X64
//包括VC++目录,C/C++常规、连接器常规等,不同的库路径不同,引入的位置也不同
string filename = "G:/myself/Major/MathPhotograph/Experiment/Picture05.jpg";//斜杠的方向不能反
//Mat image = readImage(argc, argv,filename);
//Mat pInterValue(image.rows, image.cols, CV_8UC1, Scalar::all(0));//创建一个全为0的Mat格式,注意,c:通道数
//argc: 整数, 用来统计你运行程序时送给main函数的命令行参数的个数
//* argv[ ]: 指针数组,用来存放指向字符串参数的指针,每一个元素指向一个参数
if (argc > 1) //如果参数比1大的话,就把第一个参数赋值给filename
{
filename = argv[1]; //将参数赋值给filename
}
Mat image = imread(filename, IMREAD_GRAYSCALE); //灰色样式读取图像到Mat中IMREAD_COLOR//IMREAD_GRAYSCALE//IMREAD_UNCHANGED
if (image.empty()) //如果为空
{
std::cout << "Could not open or find the image." << std::endl;
return -1;
}
std::cout << "图像的类型编号为:" << image.type() << endl; //网上有类型对应的编号
std::cout << "图像的通道数量为:" << image.channels() << endl;
canny(image);
cv::waitKey(0);//等待用户需要多长时间毫秒,零意味着永远等待
return 0;
}
void canny(Mat image)
{
int rows = image.rows - 1, cols = image.cols - 1;
Mat GaussImage(rows + 1, cols + 1, CV_8U, Scalar::all(0));
Mat sobelXImage(rows + 1, cols + 1, CV_32SC1, Scalar::all(0));
Mat sobelYImage(rows + 1, cols + 1, CV_32SC1, Scalar::all(0));
Mat sobelXYImage(rows + 1, cols + 1, CV_32SC1, Scalar::all(0));
Mat sobelRUImage(rows + 1, cols + 1, CV_32SC1, Scalar::all(0));
Mat sobelLUImage(rows + 1, cols + 1, CV_32SC1, Scalar::all(0));
Mat pixDirect(rows + 1, cols + 1, CV_8U, Scalar::all(0));
Mat imageX8UC;
Mat imageY8UC;
Mat imageXY8UC;
int winR, winC, k;
//Step1:高斯滤波
//------------------------------------------------------------------------------
//自己写得高斯滤波(实际用的是系统自带的高斯函数)
//float x, y, sigma;
//vector<float> opera;
float distance;
//int window = 5; sigma = 1.25;
//for (int i = -window / 2; i <= window / 2; i++)
//{
// for (int j = -window / 2; j <= window / 2; j++)
// {
// opera.push_back(Gauss(i, j, sigma));
// cout <<Gauss(i, j, sigma) << endl << endl;
// }
//}
//for (int row = 0; row <= rows; row++)
//{
// for (int col = 0; col <= cols; col++)
// {
// //windowNum = pow(window, 2);
// k = 0;
// for (winR = row - window / 2; winR <= row + window / 2; winR++)
// {
// for (winC = col - window / 2; winC <= col + window / 2; winC++)
// {
// if (winR >= 0 && winR <= rows && winC >= 0 && winC <= cols)
// {
// diffImage.at<int>(row, col) += image.at<uchar>(winR, winC)*opera[k];
// k++;
// //diffImage.at<int>(row, col) = image.at<uchar>(winR, winC);
// }
// else
// {
// k++;
// }//else
// }
// }//for
// }
//}
//用自带的函数进行运算
GaussianBlur(image, GaussImage, cv::Size(5, 5), 0);
//------------------------------------------------------------------------------
//Step2:sobel算子提取边缘
//------------------------------------------------------------------------------
//Mat mtest = (Mat_<float>(4, 1) << -0.055818, -0.734866, -0.675912, 0.506045);
int window = 3;
int sobelX[9] = { -1, 0, 1, -2, 0, 2, -1, 0, 1 };
int sobelY[9] = { 1, 2, 1, 0, 0, 0, -1, -2, -1 };
int sobelRU[9] = { 0, -1, -2, 1, 0, -1, 2, 1, 0 };
int sobelLU[9] = { -2, -1, 0, -1, 0, 1, 0, 1, 2 };
//int pixDirect;
for (int row = 0; row <= rows; row++)
{
for (int col = 0; col <= cols; col++)
{
//windowNum = pow(window, 2);
k = 0;
for (winR = row - window / 2; winR <= row + window / 2; winR++)
{
for (winC = col - window / 2; winC <= col + window / 2; winC++)
{
if (winR >= 0 && winR <= rows && winC >= 0 && winC <= cols)
{
sobelXImage.at<int>(row, col) += GaussImage.at<uchar>(winR, winC) * sobelX[k];
sobelYImage.at<int>(row, col) += GaussImage.at<uchar>(winR, winC) * sobelY[k];
sobelRUImage.at<int>(row, col) += GaussImage.at<uchar>(winR, winC) * sobelRU[k];
sobelLUImage.at<int>(row, col) += GaussImage.at<uchar>(winR, winC) * sobelLU[k];
k++;
//diffImage.at<int>(row, col) = image.at<uchar>(winR, winC);
}
else
{
k++;
}//else
}
}
sobelXImage.at<int>(row, col) = abs(sobelXImage.at<int>(row, col));
sobelYImage.at<int>(row, col) = abs(sobelYImage.at<int>(row, col));
sobelRUImage.at<int>(row, col) = abs(sobelRUImage.at<int>(row, col));
sobelLUImage.at<int>(row, col) = abs(sobelLUImage.at<int>(row, col));
//pixDirect = 0;
if (sobelYImage.at<int>(row, col) > pixDirect.at<uchar>(row, col))
pixDirect.at<uchar>(row, col) = 1;
if (sobelRUImage.at<int>(row, col) > pixDirect.at<uchar>(row, col))
pixDirect.at<uchar>(row, col) = 2;
if (sobelLUImage.at<int>(row, col) > pixDirect.at<uchar>(row, col))
pixDirect.at<uchar>(row, col) = 3;
}
}//for
addWeighted(sobelXImage, 0.5, sobelYImage, 0.5, 0, sobelXYImage); //加权相加
Mat MaxImage(rows + 1, cols + 1, CV_8U, Scalar::all(0));
convertScaleAbs(sobelXYImage, imageXY8UC); //将图像扩展到0--255的范围内
//输出图像
imshow("Canny算子高斯滤波", imageXY8UC);
cv::imwrite("G:/myself/Major/MathPhotograph/Experiment/05_41.png", imageXY8UC);
cv::waitKey(0);
//------------------------------------------------------------------------------
//Step3:非最大值抑制
//------------------------------------------------------------------------------
for (int row = 1; row < rows; row++)
{
for (int col = 1; col < cols; col++)
{
switch (pixDirect.at<uchar>(row, col))
{
case 0:
if (!(sobelXYImage.at<int>(row, col)>sobelXYImage.at<int>(row - 1, col) &&
sobelXYImage.at<int>(row, col)>sobelXYImage.at<int>(row + 1, col)))
imageXY8UC.at<uchar>(row, col) = 0;
break;
case 1:
if (!(sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row, col - 1) &&
sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row, col + 1)))
imageXY8UC.at<uchar>(row, col) = 0;
break;
case 2:
if (!(sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row - 1, col + 1) &&
sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row + 1, col - 1)))
imageXY8UC.at<uchar>(row, col) = 0;
break;
case 3:
if (!(sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row - 1, col + 1) &&
sobelXYImage.at<int>(row, col) > sobelXYImage.at<int>(row + 1, col - 1)))
imageXY8UC.at<uchar>(row, col) = 0;
break;
}
}
}
imshow("Canny算子高斯滤波", imageXY8UC);
cv::imwrite("G:/myself/Major/MathPhotograph/Experiment/05_42.png", imageXY8UC);
cv::waitKey(0);
//------------------------------------------------------------------------------
//Step4:设置双阈值
//------------------------------------------------------------------------------
for (int row = 1; row < rows; row++)
{
for (int col = 1; col < cols; col++)
{
if (imageXY8UC.at<uchar>(row, col) < 80)
{
imageXY8UC.at<uchar>(row, col) = 0;
}
else if (imageXY8UC.at<uchar>(row, col) >= 80 && imageXY8UC.at<uchar>(row, col) < 100)
{
imageXY8UC.at<uchar>(row, col) = 100; //弱边缘
}
else
{
imageXY8UC.at<uchar>(row, col) = 255; //强边缘
}
}
}
imshow("Canny算子高斯滤波", imageXY8UC);
cv::imwrite("G:/myself/Major/MathPhotograph/Experiment/05_43.png", imageXY8UC);
cv::waitKey(0);
//------------------------------------------------------------------------------
//Step5:滞后边缘跟踪
/*//------------------------------------------------------------------------------
//强边缘点可以认为是真的边缘。弱边缘点则可能是真的边缘,
//也可能是噪声或颜色变化引起的。为得到精确的结果,后者引起的弱边缘点应该去掉。
//通常认为真实边缘引起的弱边缘点和强边缘点是连通的,而又噪声引起的弱边缘点则不会。
//所谓的滞后边界跟踪算法检查一个弱边缘点的8连通领域像素,只要有强边缘点存在,
//那么这个弱边缘点被认为是真是边缘保留下来。
//这个算法搜索所有连通的弱边缘,如果一条连通的弱边缘的任何一个点和强边缘点连通,
//则保留这条弱边缘,否则抑制这条弱边缘。搜索时可以用广度优先或者深度优先算法,
//我在这里实现了应该是最容易的深度优先算法。一次连通一条边缘的深度优先算法如下:
*/
vector<Point> queues; //队的集合
bool connected; //联通指示变量
Mat Mark(rows + 1, cols + 1, CV_8U, Scalar::all(0)); //用来标记这个像素是否被访问
Mat borderPoint(rows + 1, cols + 1, CV_8U, Scalar::all(0));
for (int row = 1; row < rows; row++)
{
for (int col = 1; col < cols; col++)
{
vector<Point> queue; //存放一系列边缘点
if (imageXY8UC.at<uchar>(row, col) != 0) //如果是边缘
{
bool highBorder = false;
queue.push_back(Point(col, row)); //把这个点压进去
Search(row, col, imageXY8UC, Mark, queue, highBorder);
if (highBorder == true)
{
for (int i = 0; i < queue.size(); i++)
{
queues.push_back(queue[i]);
}
}//if
}//if
}//for
}
vector<Point>::iterator itr; //vector容器迭代器
for (itr = queues.begin(); itr != queues.end(); itr++) //对待选点循环操作
{
borderPoint.at<uchar>(itr->y, itr->x) = 255;
}
imshow("Canny算子高斯滤波", borderPoint);
cv::imwrite("G:/myself/Major/MathPhotograph/Experiment/05_44.png", borderPoint);
cv::waitKey(0);
}
四、处理结果
总结:canny实现的比较好,很清晰同时也去除了部分噪声,缺点是线段不太连续。代码都是由opencv库编写,本来想用gdal库但是过于麻烦,只表现算法内核的话opencv足够,不用太在图像显示上费心力。