关于canny算子边缘检测的原理，此处不再做详细说明，这里直接将实现代码附在下面。读者可对照代码，进行学习。
//canny边缘检测算子
/*
#include <highgui.h>
#include <math.h>
#include<iostream>
using namespace std;
//函数声明
void ImageGray(IplImage *sourceImage, IplImage *grayImage);
void GaussFliter(IplImage *grayImage,IplImage *gaussImage);
void CalGradient(IplImage *gaussImage,IplImage *edgeImage);
void TraceEdge(int y,int x,int ThrLow, unsigned char *pResult, int *pMag,CvSize temp);
//定义存放高斯滤波像素的数组
char *pCanny;
int main()
{
//确定要读入的AVI视频
CvCapture *capture = cvCreateFileCapture("E:\\新建文件夹\\行车视频.avi");
//将视频文件的下一帧加载到内存
IplImage *tempFrame= cvQueryFrame(capture);
IplImage *grayImage,*gaussImage,*edgeImage;
//创建头并分配数据,取与这一帧图大小一样的尺寸，数据类型为无符号8位整型，单通道
grayImage = cvCreateImage(cvSize(tempFrame->width, tempFrame->height),8,1);
gaussImage=cvCreateImage(cvSize(tempFrame->width,tempFrame->height),8,1);
edgeImage=cvCreateImage(cvSize(tempFrame->width,tempFrame->height),8,1);
while (IplImage *tempFrame = cvQueryFrame(capture))
{
cvShowImage("原始视频",tempFrame);//播放原视频
ImageGray(tempFrame, grayImage);//得到灰度图
cvShowImage("灰度视频",grayImage);
GaussFliter(grayImage,gaussImage);
cvShowImage("高斯视频",gaussImage);
CalGradient(gaussImage,edgeImage );
cvShowImage("边缘视频",edgeImage);
char c = cvWaitKey(33);
if (c == 27)
{
return 0;
}
}
cvReleaseCapture(&capture);
cvReleaseImage(&grayImage);
cvReleaseImage(&gaussImage);
cvReleaseImage(&edgeImage);
cvDestroyWindow("原始视频");
cvDestroyWindow("灰度视频");
cvDestroyWindow("高斯视频");
cvDestroyWindow("边缘视频");
}
void ImageGray(IplImage *sourceImage, IplImage *grayImage)
{
for(int i = 0; i < sourceImage->height; i++)
{
for(int j = 0; j < sourceImage->width; j++)
{
unsigned char b = (unsigned char)sourceImage->imageData[i * sourceImage->widthStep + j * 3 + 0];
unsigned char g = (unsigned char)sourceImage->imageData[i * sourceImage->widthStep + j * 3 + 1];
unsigned char r = (unsigned char)sourceImage->imageData[i * sourceImage->widthStep + j * 3 + 2];
unsigned char gray = 0.299 * r + 0.587 * g + 0.114 * b;
grayImage->imageData[i * sourceImage->width + j] = gray;
}
}
}
void GaussFliter(IplImage *grayImage,IplImage *gaussImage)
{
//定义一个数组存放高斯模板
int temp[9];
temp[0] = 1;
temp[1] = 2;
temp[2] = 1;
temp[3] = 2;
temp[4] = 4;
temp[5] = 2;
temp[6] = 1;
temp[7] = 2;
temp[8] = 1;
//申请图像的高度乘以宽度个char类型的内存单元，
//并把这块内存单元的首地址赋给一个char类型的指针pCanny
pCanny = new char[grayImage->height *grayImage->width];
uchar *data0=(uchar*)(grayImage->imageData );
int step0=grayImage->width;
for(int i= 1; i < grayImage->height-1; i++)
{
for(int j= 1; j <grayImage->width-1; j++)
{
double gauss=
( temp[0]*data0[(i-1)*step0+j-1]+ temp[1]*data0[(i-1)*step0+j]+ temp[2]*data0[(i-1)*step0+j+1]
+temp[3]*data0[i*step0+j-1]+temp[4]*data0[i*step0+j]+temp[5]*data0[i*step0+j+1]
+temp[6]*data0[(i+1)*step0+j-1]+temp[7]*data0[(i+1)*step0+j]+temp[8]*data0[(j+1)*step0+i+1])/16;
pCanny[i*grayImage->width + j]=gauss;
}
}
gaussImage->imageData = pCanny;
}
//图像增强——计算图像梯度及其方向
/***************************************************************
代码功能：边缘检测
*****************************************************************/
/*
void CalGradient(IplImage *gaussImage,IplImage *edgeImage)
{
int i,j;
double *P = new double[gaussImage->height * gaussImage->width];//存放x向偏导数
double *Q = new double[gaussImage->height * gaussImage->width];//存放y向偏导数
int *M = new int[gaussImage->height * gaussImage->width ]; //存放梯度幅值
double *Theta = new double[gaussImage->height * gaussImage->width]; //存放梯度方向
//计算x、y向偏导数
for( i = 0; i < (gaussImage->height-1); i++)
{
for( j = 0; j<(gaussImage->width-1); j++)
{
P[i * gaussImage->width + j] = (double)( pCanny[i * gaussImage->width + min(j+1,gaussImage->width-1)]
-pCanny[i * gaussImage->width + ( j )]
+ pCanny[min(i+1, gaussImage->height-1) * gaussImage->width + min(j+1,gaussImage->width-1)] - pCanny[min(i+1, gaussImage->height-1) * gaussImage->width + ( j )] )/2;
Q[i * gaussImage->width + j] = (double)( pCanny[i * gaussImage->width + ( j )]-pCanny[min(i+1, gaussImage->height-1) * gaussImage->width + ( j )]
+pCanny[( i ) * gaussImage->width + min(j+1, gaussImage->width-1)]-pCanny[min(i+1, gaussImage->height-1)*gaussImage->width+min(j+1, gaussImage->width-1)] )/2;
}
}
//计算梯度幅值和梯度的方向
for( i = 0; i < gaussImage->height; i++)
{
for( j = 0; j < gaussImage->width; j++)
{
M[i * gaussImage->width + j]=(int)(sqrt(P[i * gaussImage->width + j] * P[i * gaussImage->width + j]
+Q[i * gaussImage->width + j] * Q[i * gaussImage->width + j]) + 0.5);
Theta[i * gaussImage->width + j] = atan( Q [i * gaussImage->width + j]/P[i * gaussImage->width + j] ) * 57.3;
if( Theta[i * gaussImage->width + j] < 0 )
Theta[i * gaussImage->width + j] += 360;//将这个角度转换到0~360范围
}
}
//非极大值抑制
unsigned char *N = new unsigned char[gaussImage->height * gaussImage->width];//存放非极大值抑制后的结果
double g1,g2,g3,g4;
double Tmp1 = 0.0,Tmp2 = 0.0;//保存求得的插值后的像素值
double weight;//进行插值时的权重
//对边界进行初始化
//对第一行和最后一行进行初始化
for( i = 0; i < gaussImage->width; i++)
{
N[i] = 0;
N[(gaussImage->height - 1) * gaussImage->width + i] = 0;
}
//对第一列和最后一列进行初始化
for( j = 0; j < gaussImage->height; j++)
{
N[j * gaussImage->width] = 0;
N[j * gaussImage->width + (gaussImage->width - 1)] = 0;
}
for( i = 1; i < gaussImage->height - 1 ; i++)
{
for( j = 1; j < gaussImage->width - 1 ; j++)
{
int nowPoint = i * gaussImage->width + j;
if(M[nowPoint] == 0)
N[nowPoint] = 0;
else
{
////////////首先判断属于那种情况，然后根据情况插值///////
////////////////////第一种情况///////////////////////
///////// g1 g2 /////////////
///////// C /////////////
///////// g3 g4 /////////////
/////////////////////////////////////////////////////
if( ((Theta[nowPoint] >= 90) && (Theta[nowPoint] < 135)) ||
((Theta[nowPoint] >= 270) && (Theta[nowPoint] < 315)) )
{
g1 = M[(i-1) * gaussImage->width + (j-1)];
g2 = M[(i-1) * gaussImage->width + ( j )];
g3 = M[(i+1) * gaussImage->width + (j+1)];
g4 = M[(i+1) * gaussImage->width + ( j )];
weight = fabs(P[nowPoint])/fabs(Q[nowPoint]);//反正切
Tmp1 = g1 * weight + g2 * (1-weight);
Tmp2 = g4 * weight + g3 * (1-weight);
}
////////////////////第二种情况///////////////////////
///////// g1 /////////////
///////// g2 C g3 /////////////
///////// g4 /////////////
/////////////////////////////////////////////////////
else if( ((Theta[nowPoint] >= 135) && (Theta[nowPoint] < 180)) ||
((Theta[nowPoint] >= 315) && (Theta[nowPoint] < 360)) )
{
g1 = M[(i-1) * gaussImage->width + (j-1)];
g2 = M[( i ) * gaussImage->width + (j-1)];
g3 = M[( i ) * gaussImage->width + (j+1)];
g4 = M[(i+1) * gaussImage->width + (j+1)];
weight = fabs(Q[nowPoint])/fabs(P[nowPoint]);//正切
Tmp1 = g1 * weight + g2 * (1-weight);
Tmp2 = g4 * weight + g3 * (1-weight);
}
////////////////////第三种情况///////////////////////
///////// g1 g2 /////////////
///////// C /////////////
///////// g4 g3 /////////////
/////////////////////////////////////////////////////
else if( ((Theta[nowPoint] >= 45) && (Theta[nowPoint] < 90)) ||
((Theta[nowPoint] >= 225) && (Theta[nowPoint] < 270)) )
{
g1 = M[(i-1) * gaussImage->width + ( j )];
g2 = M[(i-1) * gaussImage->width + (j+1)];
g3 = M[(i+1) * gaussImage->width + (j-1)];
g4 = M[(i+1) * gaussImage->width + ( j )];
weight = fabs(P[nowPoint])/fabs(Q[nowPoint]);//反正切
Tmp1 = g2 * weight + g1 * (1-weight);
Tmp2 = g3 * weight + g4 * (1-weight);
}
////////////////////第四种情况///////////////////////
///////// g1 /////////////
///////// g4 C g2 /////////////
///////// g3 /////////////
/////////////////////////////////////////////////////
else if ( ((Theta[nowPoint] >= 0) && (Theta[nowPoint] < 45)) ||
((Theta[nowPoint] >= 180) && (Theta[nowPoint] < 225)) )
{
g1 = M[(i-1) * gaussImage->width + (j+1)];
g2 = M[( i ) * gaussImage->width + (j+1)];
g3 = M[(i+1) * gaussImage->width + (j-1)];
g4 = M[( i ) * gaussImage->width + (j-1)];
weight = fabs(Q[nowPoint])/fabs(P[nowPoint]);//正切
Tmp1 = g1 * weight + g2 * (1-weight);
Tmp2 = g3 * weight + g4 * (1-weight);
}
}
if ( (M[nowPoint] > Tmp1) && (M[nowPoint] > Tmp2) )
N[nowPoint] = 128;
else
N[nowPoint] = 0;
}
}
//*************************************************************
//双阈值检测实现
int Hist[1024];
int EdgeNum;//可能边界数
int MaxMag = 0;//最大梯度值
int HighCount;
//*********************************************************
//构造灰度图的统计直方图
for( i = 0; i < 1024; i++)
Hist[i] = 0;
for( i = 0; i < gaussImage->height; i++)
{
for( j = 0; j < gaussImage->width; j++)
{
if( N[i * gaussImage->width + j]==128 )
Hist[ M[i * gaussImage->width + j] ]++;
}
}
//******************************************8888****************
//获取最大梯度幅值及潜在边缘点个数
EdgeNum = Hist[0];
MaxMag = 0;//获取最大的梯度幅值
for( i = 1; i < 1024; i++)//统计经过 "非最大值抑制" 后有多少像素
{
if( Hist[i] != 0 )//梯度为零的点不可能为边界点
{
MaxMag = i;
}
EdgeNum += Hist[i];//经过non-maximum suppression后有多少像素
}
//计算两个阈值
double gRatHigh = 0.79;
double ThrHigh;
double ThrLow;
double gRatLow = 0.5;
//按照灰度值从低到高的顺序，选取前79%个灰度值中的最大
//的灰度值为高阈值，低阈值大约为高阈值的一半。
HighCount = (int)(gRatHigh * EdgeNum + 0.5);
j=1;
EdgeNum = Hist[1];
while ( (j < (MaxMag-1)) && (EdgeNum < HighCount) )
{
j++;
EdgeNum += Hist[j];
}
ThrHigh = j;
ThrLow = (int)((ThrHigh) * gRatLow + 0.5);
//***********************************************************************
//进行边缘检测
//以上代码在非极大值抑制产生的二值灰度矩阵的潜在点中按照高阈值寻找边缘，
//并以所找到的点为中心寻找邻域内满足低阈值的点，从而形成一个闭合的轮廓。
//然后对于不满足条件的点，可用如下代码直接删除掉。
//********************************************************************
CvSize temp = cvSize(gaussImage->width, gaussImage->height);
for( i = 0; i < gaussImage->height; i++)
{
for( j = 0; j < gaussImage->width; j++)
{
if( (N[i * gaussImage->width + j]==128) && (M[i * gaussImage->width + j] >= ThrHigh) )
{
N[i * gaussImage->width + j] = 255;
//TraceEdge(i,j,ThrLow,N,M,temp);
}
if( (N[i * gaussImage->width + j]==128) && (M[i * gaussImage->width + j] >= ThrLow) )
{
N[i * gaussImage->width + j] = 255;
}
}
}
//*****************************************************************
//将还没有设置为边界的点设置为非边界点
for( i = 0; i < gaussImage->height; i++)
{
for( j = 0; j < gaussImage->width; j++)
{
if( N[i * gaussImage->width + j] != 255)
{
N[i * gaussImage->width + j] = 0;
}
}
}
edgeImage->imageData = (char *)N;
//释放申请的动态内存
if(gaussImage==NULL)
{
delete pCanny;
delete [] N;
}
delete [] Theta;
delete [] P;
delete [] Q;
delete [] M;
}
/*
void TraceEdge(int y,int x,int ThrLow, unsigned char *pResult, int *pMag,CvSize temp)
{
//对8邻域像素进行查询
int xNum[8] = {1,1,0,-1,-1,-1,0,1};
int yNum[8] = {0,1,1,1,0,-1,-1,-1};
long yy,xx,k;
for( k = 0; k < 8; k++)
{
yy = y + yNum[k];
xx = x + xNum[k];
if(pResult[yy * temp.width + xx] == 128 && pMag[yy * temp.width + xx] >= ThrLow )
{
//该点设为边界点
pResult[yy * temp.width + xx] = 255;
//以该点为中心再进行跟踪
TraceEdge(yy,xx,ThrLow,pResult,pMag,temp);
}
}
}