Canny算子是John Canny在1986年提出的,那年老大爺才28歲,該文章發表在PAMI頂級期刊上的(1986. A computational approach to edge detection. IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 8, 1986, pp. 679-698)。老大爺目前在加州伯克利作machine learning,80-90年代視覺都是圖像處理,如今作視覺都是機器學習的天下,大爺的主頁(http://www.cs.berkeley.edu/~jfc/)。html
Canny算子與Marr(LoG)邊緣檢測方法相似(Marr大爺號稱計算機視覺之父),也屬因而先平滑後求導數的方法。John Canny研究了最優邊緣檢測方法所需的特性,給出了評價邊緣檢測性能優劣的三個指標:算法
1 好的信噪比,即將非邊緣點斷定爲邊緣點的機率要低,將邊緣點判爲非邊緣點的機率要低;網絡
2 高的定位性能,即檢測出的邊緣點要儘量在實際邊緣的中心;app
3 對單一邊緣僅有惟一響應,即單個邊緣產生多個響應的機率要低,而且虛假響應邊緣應該獲得最大抑制。機器學習
用一句話說,就是但願在提升對景物邊緣的敏感性的同時,能夠抑制噪聲的方法纔是好的邊緣提取方法。ide
Canny算子求邊緣點具體算法步驟以下:函數
1. 用高斯濾波器平滑圖像.性能
2. 用一階偏導有限差分計算梯度幅值和方向學習
3. 對梯度幅值進行非極大值抑制spa
4. 用雙閾值算法檢測和鏈接邊緣.
具體的步驟是能容易理解,如今就是用C語言怎麼實現了,在參考了網上諸多教程的基礎下,寫了個代碼給你們參考,確定有很多問題,但願能獲得你們的指點。
首先用白話敘述下Canny算子的原理:看做者寫的paper題目就是邊緣檢測,何爲邊緣?圖象局部區域亮度變化顯著的部分,對於灰度圖像來講,也就是灰度值有一個明顯變化,既從一個灰度值在很小的緩衝區域內急劇變化到另外一個灰度相差較大的灰度值。依據我僅有的數學知識,怎麼表徵這種灰度值的變化呢?導數就是表徵變化率的,可是數字圖像都是離散的,也就是導數確定會用差分來代替。也就是具體算法中的步驟2。可是在真實的圖像中,通常會有噪聲,噪聲會影響梯度(換不嚴格說法 偏導數)的計算,因此步驟1上來先濾波。理論上將圖像梯度幅值的元素值越大,說明圖像中該點的梯度值越大,但這不能說明該點就是邊緣。在Canny算法中,非極大值抑制(步驟3)是進行邊緣檢測的重要步驟,通俗意義上是指尋找像素點局部最大值,沿着梯度方向,比較它前面和後面的梯度值進行了。步驟4,是一個典型算法,有時候咱們並不像一刀切,也就是超過閾值的都是邊緣點,而是設爲兩個閾值,但願在高閾值和低閾值之間的點也多是邊緣點,並且這些點最好在高閾值的附近,也就是說這些中間閾值的點是高閾值邊緣點的一種延伸。因此步驟4用了雙閾值來檢測和鏈接邊緣。
基本原理簡單說完,上代碼,代碼按照下面6步驟敘述。
第一步:灰度化
第二步:高斯濾波
第三步:計算梯度值和方向
第四步:非極大值抑制
第五步:雙閾值的選取
第六步:邊緣檢測
1 把彩色圖像變成灰度圖像。該部分主要是爲像我這樣的小菜鳥準備的。
該部分是按照Canny算法一般處理的圖像爲灰度圖,若是獲取的彩色圖像,那首先就得進行灰度化。以RGB格式的彩圖爲例,一般灰度化採用的公式是:
Gray=0.299R+0.587G+0.114B;
說個我常常出問題的代碼:OpenCvGrayImage->imageData[i*OpenCvGrayImage->widthStep+j] 這是opencv iplimage格式經過直接訪問內存讀取像素值的方式,我一直搞不清楚,i*widthStep仍是j*widthStep。
記住一點,是高*widthStep就行。並且是*widthStep,而不是乘以width.若是圖像的寬度不是4的倍數,opencv貌似還有補齊這一說法,因此mark一下。
代碼以下:
1 ////////第一步:灰度化 2 IplImage * ColorImage=cvLoadImage("c:\\photo.bmp",1); 3 if (ColorImage==NULL) 4 { 5 printf("image read error"); 6 return 0; 7 } 8 cvNamedWindow("Sourceimg",0); 9 cvShowImage("Sourceimg",ColorImage); // 10 IplImage * OpenCvGrayImage; 11 OpenCvGrayImage=cvCreateImage(cvGetSize(ColorImage),ColorImage->depth,1); 12 float data1,data2,data3; 13 for (int i=0;i<ColorImage->height;i++) 14 { 15 for (int j=0;j<ColorImage->width;j++) 16 { 17 data1=(uchar)(ColorImage->imageData[i*ColorImage->widthStep+j*3+0]); 18 data2=(uchar)(ColorImage->imageData[i*ColorImage->widthStep+j*3+1]); 19 data3=(uchar)(ColorImage->imageData[i*ColorImage->widthStep+j*3+2]); 20 OpenCvGrayImage->imageData[i*OpenCvGrayImage->widthStep+j]=(uchar)(0.07*data1 + 0.71*data2 + 0.21*data3); 21 } 22 } 23 cvNamedWindow("GrayImage",0); 24 cvShowImage("GrayImage",OpenCvGrayImage); //顯示灰度圖 25 cvWaitKey(0); 26 cvDestroyWindow("GrayImage");
2 對圖像高斯濾波,圖像高斯濾波的實現能夠用兩個一維高斯核分別兩次加權實現,也就是先一維X方向卷積,獲得的結果再一維Y方向卷積。固然也能夠直接經過一個二維高斯核一次卷積實現。也就是二維卷積模板,因爲水平有限,只說二維卷積模板怎麼算。
首先,一維高斯函數:
二維高斯函數:
是否是和一維很像,其實就是兩個一維乘積就是二維,有的書和文章中將r2=x2+y2 來表示距離的平方,也就是當前要卷積的點離核心的距離的平方,若是是3*3的卷積模板,核心就是中間的那個元素,那左上角的點到核心的距離是多少呢,就是sqrt(1+1)=sqrt(2),距離的平方 r2=2。基於這個理論,那麼模板中每個點的高斯係數能夠由上面的公式計算,這樣獲得的是否是最終的模板呢?答案不是,須要歸一化,也便是每個點的係數要除以全部係數之和,這樣纔是最終的二維高斯模板。
這個裏面有個小知識點,要想計算上面的係數,須要知道高斯函數的標準差σ (sigma),還須要知道選3*3仍是5*5的模板,也就是模板要多大,實際應用的時候,這二者是有關係的,根據數理統計的知識,高斯分佈的特色就是數值分佈在(μ—3σ,μ+3σ)中的機率爲0.9974,也就是模板的大小其實就是6σ這麼大就OK了,可是6σ可能不是奇數,由於咱們必定要保證有核心。因此模板窗口的大小通常採用1+2*ceil(3*nSigma) ceil是向上取整函數,例如ceil(0.6)=1。
計算獲得模板,那就是直接卷積就OK,卷積的意思就是圖像中的點附近的模板大小區域乘以高斯模板區域,獲得的結果就是該點卷積後的結果。卷積的核心意義就是獲取原始圖像中像模板特徵的性質。插一句話,目前深度學習中很火的一個卷積神經網絡CNN,就是利用卷積的原理,經過學習出這些卷積模板來識別檢測。
代碼以下:
////////第二步:高斯濾波 /////// double nSigma=0.2; int nWindowSize=1+2*ceil(3*nSigma);//經過sigma獲得窗口大小 int nCenter=nWindowSize/2; int nWidth=OpenCvGrayImage->width; int nHeight=OpenCvGrayImage->height; IplImage * pCanny; pCanny=cvCreateImage(cvGetSize(ColorImage),ColorImage->depth,1); //生成二維濾波核 double *pKernel_2=new double[nWindowSize*nWindowSize]; double d_sum=0.0; for(int i=0;i<nWindowSize;i++) { for (int j=0;j<nWindowSize;j++) { double n_Disx=i-nCenter;//水平方向距離中心像素距離 double n_Disy=j-nCenter; pKernel_2[j*nWindowSize+i]=exp(-0.5*(n_Disx*n_Disx+n_Disy*n_Disy)/(nSigma*nSigma))/(2.0*3.1415926)*nSigma*nSigma; d_sum=d_sum+pKernel_2[j*nWindowSize+i]; } } for(int i=0;i<nWindowSize;i++)//歸一化處理 { for (int j=0;j<nWindowSize;j++) { pKernel_2[j*nWindowSize+i]=pKernel_2[j*nWindowSize+i]/d_sum; } } //輸出模板 for (int i=0;i<nWindowSize*nWindowSize;i++) { if (i%(nWindowSize)==0) { printf("\n"); } printf("%.10f ",pKernel_2[i]); } //濾波處理 for (int s=0;s<nWidth;s++) { for (int t=0;t<nHeight;t++) { double dFilter=0; double dSum=0; //當前座標(s,t) //獲取8鄰域 for (int x=-nCenter;x<=nCenter;x++) { for (int y=-nCenter;y<=nCenter;y++) { if ((x+s>=0)&&(x+s<nWidth)&&(y+t>=0)&&(y+t<nHeight))//判斷是否越界 { double currentvalue=(double)OpenCvGrayImage->imageData[(y+t)*OpenCvGrayImage->widthStep+x+s]; dFilter+=currentvalue*pKernel_2[x+nCenter+(y+nCenter)*nCenter]; dSum+=pKernel_2[x+nCenter+(y+nCenter)*nCenter]; } } } pCanny->imageData[t*pCanny->widthStep+s]=(uchar)(dFilter/dSum); } } cvNamedWindow("GaussImage",0); cvShowImage("GaussImage",pCanny); //顯示高斯圖 cvWaitKey(0); cvDestroyWindow("GaussImage");
3 計算梯度值和方向,圖像灰度值的梯度通常使用一階有限差分來進行近似,這樣就能夠得圖像在x和y方向上偏導數的兩個矩陣。本文采用很是簡單的算子,固然你也能夠選擇高大上的算子或者本身寫的:
其中f爲圖像灰度值,P表明X方向梯度幅值,Q表明Y方向 梯度幅值,M是該點幅值,Θ是梯度方向,也就是角度。
代碼以下:
1 ////////////////////////////////////////////////////////// 2 ////////第三步:計算梯度值和方向 3 //////////////////一樣能夠用不一樣的檢測器////////////////加上把圖像放到0-255之間////// 4 ///// P[i,j]=(S[i+1,j]-S[i,j]+S[i+1,j+1]-S[i,j+1])/2 ///// 5 ///// Q[i,j]=(S[i,j]-S[i,j+1]+S[i+1,j]-S[i+1,j+1])/2 ///// 6 ///////////////////////////////////////////////////////////////// 7 double *P=new double[nWidth*nHeight]; 8 double *Q=new double[nWidth*nHeight]; 9 int *M=new int[nWidth*nHeight]; 10 //IplImage * M;//梯度結果 11 //M=cvCreateImage(cvGetSize(ColorImage),ColorImage->depth,1); 12 double *Theta=new double[nWidth*nHeight]; 13 int nwidthstep=pCanny->widthStep; 14 for(int iw=0;iw<nWidth-1;iw++) 15 { 16 for (int jh=0;jh<nHeight-1;jh++) 17 { 18 P[jh*nWidth+iw]=(double)(pCanny->imageData[min(iw+1,nWidth-1)+jh*nwidthstep]-pCanny->imageData[iw+jh*nwidthstep]+ 19 pCanny->imageData[min(iw+1,nWidth-1)+min(jh+1,nHeight-1)*nwidthstep]-pCanny->imageData[iw+min(jh+1,nHeight-1)*nwidthstep])/2; 20 Q[jh*nWidth+iw]=(double)(pCanny->imageData[iw+jh*nwidthstep]-pCanny->imageData[iw+min(jh+1,nHeight-1)*nwidthstep]+ 21 pCanny->imageData[min(iw+1,nWidth-1)+jh*nwidthstep]-pCanny->imageData[min(iw+1,nWidth-1)+min(jh+1,nHeight-1)*nwidthstep])/2; 22 } 23 } 24 //計算幅值和方向 25 for(int iw=0;iw<nWidth-1;iw++) 26 { 27 for (int jh=0;jh<nHeight-1;jh++) 28 { 29 M[jh*nWidth+iw]=(int)sqrt(P[jh*nWidth+iw]*P[jh*nWidth+iw]+Q[jh*nWidth+iw]*Q[jh*nWidth+iw]+0.5); 30 Theta[jh*nWidth+iw]=atan2(Q[jh*nWidth+iw],P[jh*nWidth+iw])*180/3.1415; 31 if (Theta[jh*nWidth+iw]<0) 32 { 33 Theta[jh*nWidth+iw]=360+Theta[jh*nWidth+iw]; 34 } 35 } 36 }
4 非極大值抑制是進行邊緣檢測的一個重要步驟,通俗意義上是指尋找像素點局部最大值。沿着梯度方向,比較它前面和後面的梯度值進行了。見下圖。
上圖中左右圖:g一、g二、g三、g4都表明像素點,很明顯它們是c的八領域中的4個,左圖中c點是咱們須要判斷的點,藍色的直線是它的梯度方向,也就是說c要是局部極大值,它的梯度幅值M須要大於直線與g1g2和g2g3的交點,dtmp1和dtmp2處的梯度幅值。可是dtmp1和dtmp2不是整像素,而是亞像素,也就是座標是浮點的,那怎麼求它們的梯度幅值呢?線性插值,例如dtmp1在g一、g2之間,g一、g2的幅值都知道,咱們只要知道dtmp1在g一、g2之間的比例,就能獲得它的梯度幅值,而比例是能夠靠夾角計算出來的,夾角又是梯度的方向。
寫個線性插值的公式:設g1的幅值M(g1),g2的幅值M(g2),則dtmp1能夠很獲得:
M(dtmp1)=w*M(g2)+(1-w)*M(g1)
其中w=distance(dtmp1,g2)/distance(g1,g2)
distance(g1,g2) 表示兩點之間的距離。實際上w是一個比例係數,這個比例係數能夠經過梯度方向(幅角的正切和餘切)獲得。
右邊圖中的4個直線就是4個不一樣的狀況,狀況不一樣,g一、g二、g三、g4表明c的八領域中的4個座標會有所差別,可是線性插值的原理都是一致的。
代碼以下:
1 ////////第四步:非極大值抑制 2 //注意事項 權重的選取,也就是離得近權重大 3 ///////////////////////////////////////////////////////////////// 4 IplImage * N;//非極大值抑制結果 5 N=cvCreateImage(cvGetSize(ColorImage),ColorImage->depth,1); 6 IplImage * OpencvCannyimg;//非極大值抑制結果 7 OpencvCannyimg=cvCreateImage(cvGetSize(ColorImage),ColorImage->depth,1); 8 int g1=0, g2=0, g3=0, g4=0; //用於進行插值,獲得亞像素點座標值 9 double dTmp1=0.0, dTmp2=0.0; //保存兩個亞像素點插值獲得的灰度數據 10 double dWeight=0.0; //插值的權重 11 12 for(int i=1;i<nWidth-1;i++) 13 { 14 for(int j=1;j<nHeight-1;j++) 15 { 16 //若是當前點梯度爲0,該點就不是邊緣點 17 if (M[i+j*nWidth]==0) 18 { 19 N->imageData[i+j*nwidthstep]=0; 20 }else 21 { 22 ////////首先判斷屬於那種狀況,而後根據狀況插值/////// 23 ////////////////////第一種狀況/////////////////////// 24 ///////// g1 g2 ///////////// 25 ///////// C ///////////// 26 ///////// g3 g4 ///////////// 27 ///////////////////////////////////////////////////// 28 if((Theta[i+j*nWidth]>=90&&Theta[i+j*nWidth]<135)||(Theta[i+j*nWidth]>=270&&Theta[i+j*nWidth]<315)) 29 { 30 g1=M[i-1+(j-1)*nWidth]; 31 g2=M[i+(j-1)*nWidth]; 32 g3=M[i+(j+1)*nWidth]; 33 g4=M[i+1+(j+1)*nWidth]; 34 dWeight=fabs(P[i+j*nWidth])/fabs(Q[i+j*nWidth]); 35 dTmp1=g1*dWeight+(1-dWeight)*g2; 36 dTmp2=g4*dWeight+(1-dWeight)*g3; 37 ////////////////////第二種狀況/////////////////////// 38 ///////// g1 ///////////// 39 ///////// g2 C g3 ///////////// 40 ///////// g4 ///////////// 41 ///////////////////////////////////////////////////// 42 }else if((Theta[i+j*nWidth]>=135&&Theta[i+j*nWidth]<180)||(Theta[i+j*nWidth]>=315&&Theta[i+j*nWidth]<360)) 43 { 44 g1=M[i-1+(j-1)*nWidth]; 45 g2=M[i-1+(j)*nWidth]; 46 g3=M[i+1+(j)*nWidth]; 47 g4=M[i+1+(j+1)*nWidth]; 48 dWeight=fabs(Q[i+j*nWidth])/fabs(P[i+j*nWidth]); 49 dTmp1=g1*dWeight+(1-dWeight)*g2; 50 dTmp2=g4*dWeight+(1-dWeight)*g3; 51 ////////////////////第三種狀況/////////////////////// 52 ///////// g1 g2 ///////////// 53 ///////// C ///////////// 54 ///////// g4 g3 ///////////// 55 ///////////////////////////////////////////////////// 56 }else if((Theta[i+j*nWidth]>=45&&Theta[i+j*nWidth]<90)||(Theta[i+j*nWidth]>=225&&Theta[i+j*nWidth]<270)) 57 { 58 g1=M[i+1+(j-1)*nWidth]; 59 g2=M[i+(j-1)*nWidth]; 60 g3=M[i+(j+1)*nWidth]; 61 g4=M[i-1+(j+1)*nWidth]; 62 dWeight=fabs(P[i+j*nWidth])/fabs(Q[i+j*nWidth]); 63 dTmp1=g1*dWeight+(1-dWeight)*g2; 64 dTmp2=g4*dWeight+(1-dWeight)*g3; 65 ////////////////////第四種狀況/////////////////////// 66 ///////// g1 ///////////// 67 ///////// g4 C g2 ///////////// 68 ///////// g3 ///////////// 69 ///////////////////////////////////////////////////// 70 }else if((Theta[i+j*nWidth]>=0&&Theta[i+j*nWidth]<45)||(Theta[i+j*nWidth]>=180&&Theta[i+j*nWidth]<225)) 71 { 72 g1=M[i+1+(j-1)*nWidth]; 73 g2=M[i+1+(j)*nWidth]; 74 g3=M[i-1+(j)*nWidth]; 75 g4=M[i-1+(j+1)*nWidth]; 76 dWeight=fabs(Q[i+j*nWidth])/fabs(P[i+j*nWidth]); 77 dTmp1=g1*dWeight+(1-dWeight)*g2; 78 dTmp2=g4*dWeight+(1-dWeight)*g3; 79 80 } 81 82 } 83 84 if ((M[i+j*nWidth]>=dTmp1)&&(M[i+j*nWidth]>=dTmp2)) 85 { 86 N->imageData[i+j*nwidthstep]=200; 87 88 }else N->imageData[i+j*nwidthstep]=0; 89 90 } 91 } 92 93 94 //cvNamedWindow("Limteimg",0); 95 //cvShowImage("Limteimg",N); //顯示非抑制 96 //cvWaitKey(0); 97 //cvDestroyWindow("Limteimg");
5 雙閾值的選取。
雙閾值的選取是按照直方圖來選擇的,首先把梯度幅值的直方圖(扯點題外話:梯度的幅值直方圖和角度直方圖也是SIFT算子中的一個環節)求出來,選取佔直方圖總數%多少(本身定,代碼中定義70%)所對應的梯度幅值爲高閾值,高閾值的一半爲低閾值,這只是一種簡單策略。也能夠採用其餘的。
代碼以下:
1 ///////第五步:雙閾值的選取 2 //注意事項 注意是梯度幅值的閾值 3 ///////////////////////////////////////////////////////////////// 4 int nHist[1024];//直方圖 5 int nEdgeNum;//全部邊緣點的數目 6 int nMaxMag=0;//最大梯度的幅值 7 for(int k=0;k<1024;k++) 8 { 9 nHist[k]=0; 10 } 11 for (int wx=0;wx<nWidth;wx++) 12 { 13 for (int hy=0;hy<nHeight;hy++) 14 { 15 if((uchar)N->imageData[wx+hy*N->widthStep]==200) 16 { 17 int Mindex=M[wx+hy*nWidth]; 18 nHist[M[wx+hy*nWidth]]++;//獲取了梯度直方圖 19 20 } 21 } 22 } 23 nEdgeNum=0; 24 for (int index=1;index<1024;index++) 25 { 26 if (nHist[index]!=0) 27 { 28 nMaxMag=index; 29 } 30 nEdgeNum+=nHist[index];//通過non-maximum suppression後有多少邊緣點像素 31 } 32 //計算兩個閾值 注意是梯度的閾值 33 int nThrHigh; 34 int nThrLow; 35 double dRateHigh=0.7; 36 double dRateLow=0.5; 37 int nHightcount=(int)(dRateHigh*nEdgeNum+0.5); 38 int count=1; 39 nEdgeNum=nHist[1]; 40 while((nEdgeNum<=nHightcount)&&(count<nMaxMag-1)) 41 { 42 count++; 43 nEdgeNum+=nHist[count]; 44 } 45 nThrHigh=count; 46 nThrLow= (int)(nThrHigh*dRateLow+0.5); 47 printf("\n直方圖的長度 %d \n ",nMaxMag); 48 printf("\n梯度的閾值幅值大 %d 小閾值 %d \n ",nThrHigh,nThrLow);
6 邊緣檢測,也便是根據步驟5獲得的雙閾值進行鏈接。
首先判斷該點是否超太高閾值,而後判斷該點的8鄰域點中尋找知足超太低閾值的點,再根據此點收集新的邊緣,直到整個圖像邊緣閉合。整個圖像找完後,將非邊緣點剔除,即灰度值置0.
代碼以下:
1 ////////////////////////////////////////////////////////// 2 ////////第六步:邊緣檢測 3 ///////////////////////////////////////////////////////////////// 4 5 for(int is=1;is<nWidth-1;is++) 6 { 7 for (int jt=1;jt<nHeight-1;jt++) 8 { 9 //CvScalar s=cvGet2D(N,jt,is); 10 //int currentvalue=s.val[0]; 11 int currentvalue=(uchar)(N->imageData[is+jt*N->widthStep]); 12 if ((currentvalue==200)&&(M[is+jt*nWidth]>=nThrHigh)) 13 //是非最大抑制後的點且 梯度幅值要大於高閾值 14 { 15 N->imageData[is+jt*nwidthstep]=255; 16 //鄰域點判斷 17 TraceEdge(is, jt, nThrLow, N, M); 18 } 19 } 20 } 21 for (int si=1;si<nWidth-1;si++) 22 { 23 for (int tj=1;tj<nHeight-1;tj++) 24 { 25 if ((uchar)N->imageData[si+tj*nwidthstep]!=255) 26 { 27 N->imageData[si+tj*nwidthstep]=0; 28 } 29 } 30 31 } 32 33 cvNamedWindow("Cannyimg",0); 34 cvShowImage("Cannyimg",N);
其中,鄰域跟蹤算法給出了兩個,一種是遞歸,一種是非遞歸。
遞歸算法。解決了堆棧溢出問題。以前找了不少canny代碼參考,其中有一個版本流傳很廣,可是對於大圖像堆棧溢出。
1 int TraceEdge(int w, int h, int nThrLow, IplImage* pResult, int *pMag) 2 { 3 int n,m; 4 char flag = 0; 5 int currentvalue=(uchar)pResult->imageData[w+h*pResult->widthStep]; 6 if ( currentvalue== 0) 7 { 8 pResult->imageData[w+h*pResult->widthStep]= 255; 9 flag=0; 10 for (n= -1; n<=1; n++) 11 { 12 for(m= -1; m<=1; m++) 13 { 14 if (n==0 && m==0) continue; 15 int curgrayvalue=(uchar)pResult->imageData[w+n+(h+m)*pResult->widthStep]; 16 int curgrdvalue=pMag[w+n+(h+m)*pResult->width]; 17 if (curgrayvalue==200&&curgrdvalue>nThrLow) 18 if (TraceEdge(w+n, h+m, nThrLow, pResult, pMag)) 19 { 20 flag=1; 21 break; 22 } 23 } 24 if (flag) break; 25 } 26 return(1); 27 } 28 return(0); 29 }
非遞歸算法。以下:
1 void TraceEdge(int w, int h, int nThrLow, IplImage* pResult, int *pMag) 2 { 3 //對8鄰域像素進行查詢 4 int xNum[8] = {1,1,0,-1,-1,-1,0,1}; 5 int yNum[8] = {0,1,1,1,0,-1,-1,-1}; 6 int xx=0; 7 int yy=0; 8 bool change=true; 9 while(change) 10 { 11 change=false; 12 for(int k=0;k<8;k++) 13 { 14 xx=w+xNum[k]; 15 yy=h+yNum[k]; 16 // 若是該象素爲可能的邊界點,又沒有處理過 17 // 而且梯度大於閾值 18 int curgrayvalue=(uchar)pResult->imageData[xx+yy*pResult->widthStep]; 19 int curgrdvalue=pMag[xx+yy*pResult->width]; 20 if(curgrayvalue==200&&curgrdvalue>nThrLow) 21 { 22 change=true; 23 // 把該點設置成爲邊界點 24 pResult->imageData[xx+yy*pResult->widthStep]=255; 25 h=yy; 26 w=xx; 27 break; 28 } 29 } 30 } 31 }
到此,整個算法寫完了。打擊下信心,整個算法跑起來沒問題,可是沒有opencv 的cvCanny 一個函數效果好。分析了下緣由,一個是梯度算子選的太簡單,opencv通常選用的是3*3 sobel。二是邊緣鏈接性仍是不夠好,出現了不少斷的,也就是鄰域跟蹤算法不夠好。但願有高手能改進。
附上代碼:http://www.pudn.com/downloads726/sourcecode/graph/detail2903773.html
整篇博客參考過網上不少canny算法的總結,在此謝過!