影象增強演算法實現--影象邊緣提取
(1)演算法描述:
Laplacian 運算元是n維歐幾里德空間中的一個二階微分運算元,定義為梯度grad()的散度div()。因此如果f是二階可微的實函式,則f的拉普拉斯運算元定義為:(1) f的拉普拉斯運算元也是笛卡兒座標系xi中的所有非混合二階偏導數求和:(2) 作為一個二階微分運算元,拉普拉斯運算元把C函式對映到C函式,對於k ≥ 2。表示式(1)(或(2))定義了一個運算元Δ : C(R) → C(R),或更一般地,定義了一個運算元Δ : C(Ω) → C(Ω),對於任何開集Ω。
Roberts運算元是一種最簡單的運算元,是一種利用區域性差分運算元尋找邊緣的運算元,他採用對角線方向相鄰兩象素之差近似梯度幅值檢測邊緣。
Prewitt運算元是一種一階微分運算元的邊緣檢測,利用畫素點上下、左右鄰點的灰度差,在邊緣處達到極值檢測邊緣,去掉部分偽邊緣,對噪聲具有平滑作用 。其原理是在影象空間利用兩個方向模板與影象進行鄰域卷積來完成的,這兩個方向模板一個檢測水平邊緣,一個檢測垂直邊緣。對數字影象f(x,y),Prewitt運算元的定義如下:
G(i)=|[f(i-1,j-1)+f(i-1,j)+f(i-1,j+1)]-[f(i+1,j-1)+f(i+1,j)+f(i+1,j+1)]|
G(j)=|[f(i-1,j+1)+f(i,j+1)+f(i+1,j+1)]-[f(i-1,j-1)+f(i,j-1)+f(i+1,j-1)]|
則 P(i,j)=max[G(i),G(j)]或 P(i,j)=G(i)+G(j)
經典Prewitt運算元認為:凡灰度新值大於或等於閾值的畫素點都是邊緣點。即選擇適當的閾值T,若P(i,j)≥T,則(i,j)為邊緣點,P(i,j)為邊緣影象。這種判定是欠合理的,會造成邊緣點的誤判,因為許多噪聲點的灰度值也很大,而且對於幅值較小的邊緣點,其邊緣反而丟失了。
Sobel 運算元有兩個,一個是檢測水平邊緣的 ;另一個是檢測垂直邊緣的 。與Prewitt運算元相比,Sobel運算元對於象素的位置的影響做了加權,可以降低邊緣模糊程度,因此效果更好。該運算元包含兩組3x3的矩陣,分別為橫向及縱向,將之與影象作平面
影象的每一個畫素的橫向及縱向梯度近似值可用以下的公式結合,來計算梯度的大小。
然後可用以下公式計算梯度方向。
在以上例子中,如果上述角度Θ等於零,即代表影象該處擁有縱向邊緣,左方較右方暗。
Canny 運算元的目標是找到一個最優的邊緣檢測演算法,Canny 使用了變分法(calculus of variations),這是一種尋找優化特定功能的函式的方法。最優檢測使用四個指數函式項表示,但是它非常近似於高斯函式的一階導數。Canny邊緣檢測演算法可以分為以下5個步驟:1.應用高斯濾波來平滑影象,目的是去除噪聲。 2.找尋影象的強度梯度(intensity gradients)。3.應用非最大抑制(non-maximum suppression)技術來消除邊誤檢(本來不是但檢測出來是)。4.應用雙閾值的方法來決定可能的(潛在的)邊界。5.利用滯後技術來跟蹤邊界。
(2)實驗結果,介面圖:
整體執行截圖:
左上:原圖 右上:sobel運算元
左中:robert運算元 右中:laplacian運算元
左下:prewitt運算元 右下:差分(對角)
(3)關鍵程式碼:
#include <opencv/cv.h>
#include <opencv/highgui.h>
#include<opencv/cxcore.h>
#include "math.h"
float max(float x, float y);
void prewitt(IplImage *src, IplImage *dst);//Prewitt運算元
void roberts(IplImage *src, IplImage *dst); //roberts運算元
void sobel(IplImage *src, IplImage *dst); //sobel運算元
void canny(IplImage *src, IplImage *dst);//canny運算元
void Laplacian(IplImage *src, IplImage *dst32,IplImage *dst); //Laplacian運算元
int main()
{
//載入影象
IplImage *src = NULL;
src = cvLoadImage("lena.jpg", 0);
cvNamedWindow("【原圖】", 1);
cvShowImage("【原圖】", src);
IplImage *pCannyImg = cvCreateImage(cvGetSize(src), 8, 1);
canny(src, pCannyImg);
IplImage *pLaplaceImg_32 = cvCreateImage(cvGetSize(src), 32, 1);
IplImage *pLaplaceImg = cvCreateImage(cvGetSize(src), 8, 1);
Laplacian(src, pLaplaceImg_32, pLaplaceImg);
IplImage *pSobelImg = cvCreateImage(cvGetSize(src), 8, 1);
sobel(src, pSobelImg);
cvReleaseImage(&pSobelImg);
IplImage *pRobertsImg = cvCreateImage(cvGetSize(src), 8, 1);
roberts(src, pRobertsImg);
cvReleaseImage(&pRobertsImg);
IplImage *pPrewittImg = cvCreateImage(cvGetSize(src), 8, 1);
prewitt(src, pPrewittImg);
cvReleaseImage(&pPrewittImg);
cvWaitKey(0);
//銷燬視窗,釋放影象
cvDestroyWindow("【原圖】");
cvReleaseImage(&src);
return -1;
}
float max(float x, float y)
{
float z;
if (x>y)z = x;
else z = y;
return(z);
}
//Prewitt運算元
void prewitt(IplImage *src, IplImage *dst)
{
//定義prewitt運算元的模板
float prewittx[9] =
{
-1, 0, 1,
-1, 0, 1,
-1, 0, 1
};
float prewitty[9] =
{
1, 1, 1,
0, 0, 0,
-1, -1, -1
};
CvMat px;
px = cvMat(3, 3, CV_32F, prewittx);
CvMat py;
py = cvMat(3, 3, CV_32F, prewitty);
IplImage *dstx = cvCreateImage(cvGetSize(src), 8, 1);
IplImage *dsty = cvCreateImage(cvGetSize(src), 8, 1);
//對影象使用模板,自動填充邊界
cvFilter2D(src, dstx, &px, cvPoint(-1, -1));
cvFilter2D(src, dsty, &py, cvPoint(-1, -1));
//計算梯度,
int i, j, temp;
float tempx, tempy; //定義為浮點型是為了避免sqrt函式引起歧義
uchar* ptrx = (uchar*)dstx->imageData;
uchar* ptry = (uchar*)dsty->imageData;
for (i = 0; i<src->width; i++)
{
for (j = 0; j<src->height; j++)
{
tempx = ptrx[i + j*dstx->widthStep]; //tempx,tempy表示的是指標所指向的畫素
tempy = ptry[i + j*dsty->widthStep];
temp = (int)sqrt(tempx*tempx + tempy*tempy);
dst->imageData[i + j*dstx->widthStep] = temp;
}
}
double min_val = 0, max_val = 0;//取圖並顯示像中的最大最小畫素值
cvMinMaxLoc(dst, &min_val, &max_val);
printf("max_val = %f\nmin_val = %f\n", max_val, min_val);
cvSaveImage("PrewittImg.jpg", dst);//把影象存入檔案
cvReleaseImage(&dstx);
cvReleaseImage(&dsty);
cvNamedWindow("【效果圖】prewitt運算元", 1);
cvShowImage("【效果圖】prewitt運算元", dst);
}
//roberts運算元
void roberts(IplImage *src, IplImage *dst)
{
dst = cvCloneImage(src);
int x, y, i, w, h;
int temp, temp1;
uchar* ptr = (uchar*)(dst->imageData);
int ptr1[4] = { 0 };
int indexx[4] = { 0, 1, 1, 0 };
int indexy[4] = { 0, 0, 1, 1 };
w = dst->width;
h = dst->height;
for (y = 0; y<h - 1; y++)
for (x = 0; x<w - 1; x++)
{
for (i = 0; i<4; i++) //取每個2*2矩陣元素的指標 0 | 1
{ // 3 | 2
ptr1[i] = *(ptr + (y + indexy[i])*dst->widthStep + x + indexx[i]);
}
temp = abs(ptr1[0] - ptr1[2]); //計算2*2矩陣中0和2位置的差,取絕對值temp
temp1 = abs(ptr1[1] - ptr1[3]); //計算2*2矩陣中1和3位置的差,取絕對值temp1
temp = (temp>temp1 ? temp : temp1); //若temp1>temp,則以temp1的值替換temp
temp = (int)sqrt(float(temp*temp) + float(temp1*temp1)); //輸出值
*(ptr + y*dst->widthStep + x) = temp; //將輸出值存放於dst畫素的對應位置
}
double min_val = 0, max_val = 0;//取圖並顯示像中的最大最小畫素值
cvMinMaxLoc(dst, &min_val, &max_val);
printf("max_val = %f\nmin_val = %f\n", max_val, min_val);
cvSaveImage("RobertsImg.jpg", dst);
cvNamedWindow("【效果圖】robert運算元", 1);
cvShowImage("【效果圖】robert運算元", dst);
}
//sobel運算元
void sobel(IplImage *src, IplImage *dst)
{
IplImage *pSobelImg_dx = cvCreateImage(cvGetSize(src), 32, 1);
IplImage *pSobelImg_dy = cvCreateImage(cvGetSize(src), 32, 1);
IplImage *pSobelImg_dxdy = cvCreateImage(cvGetSize(src), 32, 1);
//用sobel運算元計算兩個方向的微分
cvSobel(src, pSobelImg_dx, 1, 0, 3);
cvSobel(src, pSobelImg_dy, 0, 1, 3);
int i, j;
double v1, v2, v;
for (i = 0; i<src->height; i++)
{
for (j = 0; j<src->width; j++)
{
v1 = cvGetReal2D(pSobelImg_dx, i, j);
v2 = cvGetReal2D(pSobelImg_dy, i, j);
v = sqrt(v1*v1 + v2*v2);
/* if(v>100) v = 255;
else v = 0;*/
cvSetReal2D(pSobelImg_dxdy, i, j, v);
}
}
cvConvertScale(pSobelImg_dxdy, dst); //將影象轉化為8位
double min_val = 0, max_val = 0;//取圖並顯示像中的最大最小畫素值
cvMinMaxLoc(pSobelImg_dxdy, &min_val, &max_val);
printf("max_val = %f\nmin_val = %f\n", max_val, min_val);
//歸一化
cvNormalize(dst, dst, 0, 255, CV_MINMAX, 0);
cvReleaseImage(&pSobelImg_dx);
cvReleaseImage(&pSobelImg_dy);
cvReleaseImage(&pSobelImg_dxdy);
cvSaveImage("SobelImg.jpg", dst);//把影象存入檔案
cvNamedWindow("【效果圖】sobel運算元", 1);
cvShowImage("【效果圖】sobel運算元", dst);
}
//canny運算元
void canny(IplImage *src, IplImage *dst)
{
cvCanny(src, dst, 100, 150, 3),
cvNamedWindow("【效果圖】差分(對角)", 1);
cvShowImage("【效果圖】差分(對角)", dst);
cvSaveImage("CannyImg.jpg", dst);//把影象存入檔案
cvReleaseImage(&dst);
}
//Laplacian運算元
void Laplacian(IplImage *src, IplImage *dst32,IplImage *dst)
{
double min_val = 0, max_val = 0;//取圖並顯示像中的最大最小畫素值
cvLaplace(src, dst32, 5);
cvConvertScale(dst32, dst); //將影象轉化為8位
cvMinMaxLoc(dst, &min_val, &max_val);
printf("max_val = %f\nmin_val = %f\n", max_val, min_val);
cvNormalize(dst, dst, 0, 255, CV_MINMAX, 0);
cvNamedWindow("【效果圖】laplacian運算元", 1);
cvShowImage("【效果圖】laplacian運算元", dst);
cvSaveImage("LaplaceImg.jpg", dst);//把影象存入檔案
cvReleaseImage(&dst);
}
(4)完成作業體會,總結:
寫這個程式遇見的最大的難題在於對幾種運算元的原理,不是特別的理解。特別是差分,網上根本沒有相關的演算法介紹,後來看到了canny演算法,發現它是用到了差分演算法。可以說,寫完這個程式後,對每種演算法的具體實現並不能做到百分百的理解,對各種矩陣變換還是一頭霧水。寫這個程式主要的時間都花費在理解各種演算法,實現各種演算法,將各種演算法封裝成函式直接呼叫,這期間也遇到了許多零零碎碎的問題,學到了很多。