[转] 图像处理中的拉普拉斯算子

http://book.51cto.com/art/200808/84592.htm

5.5.2  拉普拉斯掩模锐化（1）

1．基本理论

拉普拉斯算子是最简单的各向同性微分算子，具有旋转不变性。一个二维图像函数 的拉普拉斯变换是各向同性的二阶导数，定义为：

  

（5-11）

为了更适合于数字图像处理，将该方程表示为离散形式： 

   

（5-12）

另外，拉普拉斯算子还可以表示成模板的形式，如图5-9所示。图5-9（a）表示离散拉普拉斯算子的模板，图5-9（b）表示其扩展模板，图5-9（c）则分别表示其他两种拉普拉斯的实现模板。从模板形式容易看出，如果在图像中一个较暗的区域中出现了一个亮点，那么用拉普拉斯运算就会使这个亮点变得更亮。因为图像中的边缘就是那些灰度发生跳变的区域，所以拉普拉斯锐化模板在边缘检测中很有用。一般增强技术对于陡峭的边缘和缓慢变化的边缘很难确定其边缘线的位置。但此算子却可用二次微分正峰和负峰之间的过零点来确定，对孤立点或端点更为敏感，因此特别适用于以突出图像中的孤立点、孤立线或线端点为目的的场合。同梯度算子一样，拉普拉斯算子也会增强图像中的噪声，有时用拉普拉斯算子进行边缘检测时，可将图像先进行平滑处理。

   

图5-9  拉普拉斯的4种模板

图像锐化处理的作用是使灰度反差增强，从而使模糊图像变得更加清晰。图像模糊的实质就是图像受到平均运算或积分运算，因此可以对图像进行逆运算，如微分运算能够突出图像细节，使图像变得更为清晰。由于拉普拉斯是一种微分算子，它的应用可增强图像中灰度突变的区域，减弱灰度的缓慢变化区域。因此，锐化处理可选择拉普拉斯算子对原图像进行处理，产生描述灰度突变的图像，再将拉普拉斯图像与原始图像叠加而产生锐化图像。拉普拉斯锐化的基本方法可以由下式表示：

  

这种简单的锐化方法既可以产生拉普拉斯锐化处理的效果，同时又能保留背景信息，将原始图像叠加到拉普拉斯变换的处理结果中去，可以使图像中的各灰度值得到保留，使灰度突变处的对比度得到增强，最终结果是在保留图像背景的前提下，突现出图像中小的细节信息。

 

2．算法实现

CImgEnhance类中的成员函数LapTemplate()实现拉普拉斯锐化操作，在对灰度图像进行梯度锐化时直接调用LapTemplate()函数即可。以下是LapTemplate()函数的代码实现。

 /**************************************
********************************
*
* 函数名称：
*     LapTemplate(int inputH, int 
inputW, int inputX0, int inputY0, float 
*     pTemplate, float fCoef)
*
* 参数：
*     int inputH-模板的高度
*     int inputW-模板的宽度
*     int inputX0-模板的中心元素X坐标
*     int inputY0-模板的中心元素Y坐标
*     float *pTemplate-指向模板数组的指针
*     float f-模板系数
*
* 返回值：
*     void
*
* 说明：
*     该函数用指定的模板（任意大小）来对
图像进行操作，参数inputH指定模板
*     的高度，参数inputW指定模板的宽度，
参数inputX0和iinputY0指定模板的中心
*     元素坐标，参数pTemplate指定模板元素，f指定系数
*     拉普拉斯锐化
***********************************
***********************************/
void CImgEnhance::LapTemplate(int 
inputH, int inputW, int inputX0, int 
inputY0, float *pTemplate, float f)
{
unsigned char *pSrc, *pDst; 
int  i,j,k,l;
float value;

if(m_pImgDataOut != NULL)
{
delete []m_pImgDataOut;
m_pImgDataOut = NULL;
}
 int lineByte = (m_imgWidth * 
m_nBitCount / 8 + 3) / 4 * 4;

if(m_nBitCount != 8)
{
AfxMessageBox("只能处理8位灰度图像!");
return ;
}
//创建要复制的图像区域
m_nBitCountOut = m_nBitCount;
int lineByteOut = (m_imgWidth * 
m_nBitCountOut / 8 + 3) / 4 * 4;
if (!m_pImgDataOut)
{
m_pImgDataOut = new unsigned char
[lineByteOut * m_imgHeight];
}

int pixelByte = m_nBitCountOut / 8;
for(i = 0; i < m_imgHeight; i++){
for(j = 0; j < m_imgWidth * pixelByte; j++)
*(m_pImgDataOut + i * lineByteOut + j)
= *(m_pImgData + i * lineByteOut + j);
}

//行处理
for (i = inputY0; i < m_imgHeight - 
inputH + inputY0 + 1; i++)
{
//列处理
for (j = inputX0; j 
inputW + inputX0 + 1; j++)
{
//指向新DIB第i行第j列像素的指针
pDst = m_pImgDataOut + lineByte * 
(m_imgHeight -1 - i) + j;

value=0;
//计算
for (k = 0; k < inputH; k++)
{
for (l = 0; l < inputW; l++)
{
pSrc = m_pImgData + lineByte * 
(m_imgHeight - 1 - i + inputY0 - k)+ j - inputX0 + l;
//计算加权平均
value += (*pSrc) * pTemplate[k * inputW + l];
}
}
//乘以系数
value *= f;
//取结果的绝对值
value = (float)fabs(value);
if (value > 255)
{
*pDst = 255; 
}
else
{
*pDst = (unsigned char)(value+0.5);
}
}
}
}

   5.5.2  拉普拉斯掩模锐化（2）

 

3．函数调用

在视图类CDemoView中映射“拉普拉斯锐化”事件的处理函数:OnLaplacienSharp()，以下是CDemoView:: :OnLaplacienSharp()函数的代码实现。

void CDemoView::OnLaplacienSharp() 
{
CDemoDoc *pDoc=GetDocument();
ImgCenterDib *pDib=pDoc->GetPDib();

if(pDib->m_nBitCount!=8&{
::MessageBox(0,"只处理灰度图像",MB_OK,0);
return ;
}
CImgEnhance imgnoise(pDib->GetDimensions(),pDib->m_nBitCount,
pDib->m_lpColorTable, pDib->m_pImgData);
int TempH, TempW, TempCX, TempCY;
float fCoef;
float Temp[9];
//设模板为平均模板
Temp[0]=-1.0;
Temp[1]=-1.0;
Temp[2]=-1.0;
Temp[3]=-1.0;
Temp[4]=9.0;
T[5]=-1.0;
T[6]=-1.0;
T[7]=-1.0;
T[8]=-1.0;
//初始化对话框变量
TH=3;
TW=3;
TCX=1;
TCY=1;
fCoef=1.0;
imgnoise.LapTlate(TH,TW,TCX,TCY,T,fCoef);

CMainFrame* pFrame = (CMainFrame *)(AfxGetApp()->m_pMainWnd);
pFrame->SendMessage(WM_COMMAND, ID_FILE_NEW);

CDemoView* pView=(CDemoView*)pFrame->MDIGetActive()->GetActiveView();
CDemoDoc* pDocNew=pView->GetDocument();
ImgCenterDib *dibNew=pDocNew->GetPDib();

dibNew->ReplaceDib(imgnoise.GetDimensions(),imgnoise.m_nBitCountOut,imgnoise.m_lpColorTable, imgnoise.m_pImgDataOut);
pDocNew->SetModifiedFlag(TRUE);
pDocNew->UpdateAllViews(pView);

Invalidate();
}

4．实验结果与分析

图5-10（a）显示了一幅花朵的图片，图5-10（b）显示了用图5-9（a）所示的拉普拉斯模板对该图像滤波后的结果。由图可以看出，将原始图像通过拉普拉斯变换后增强了图像中灰度突变处的对比度，使图像中小的细节部分得到增强并保留了图像的背景色调，使图像的细节比原始图像更加清晰。基于拉普拉斯变换的图像增强已成为图像锐化处理的基本工具。

 

 

图5-10  拉普拉斯方法的锐化结果

 

5.6  本章小结

本章主要介绍了几种图像的空域增强方法及其每种方法的基本原理，在程序中用类CImgEnhance实现了图像添加噪声和常用的图像增强的算法。本章从不同角度介绍了图像中噪声的来源和模型，以及如何计算信噪比及编程实现随机噪声和椒盐噪声的添加，为研究后续的图像平滑和滤波方法奠定基础。本章还简要介绍了图像的灰度修正方法，其中包括灰度校正、灰度变换及直方图修正；详细介绍了图像的空间域平滑方法，包括邻域平均法、加权平均法和选择式掩模平滑法；介绍了传统中值滤波方法和一种适合并行实现的快速中值滤波方法；在“图像的锐化”一节中重点介绍了两种微分锐化方法：梯度锐化方法和拉普拉斯掩模锐化方法。上述方法都给出了相应的实验结果，便于读者根据不同的应用场合和应用需求选择合适的算法。

 

Matlab拉普拉斯算子

http://blog.csdn.net/superdont/archive/2007/04/08/1556934.aspx

%function PicOut=Lap_edge(PicInput,thresh) 
PicInput=imread('coins.png');
% PicInput=imread('rice.png');
thresh=48;
% 本程序能够将BMP格式的黑白灰度图像用拉普拉斯算子进行边缘检测
% 生物图像处理作业2
% 格式为 a=Lap_edge(PicInput,thresh) 或者嵌套为 Lap_edge(imread('rice.tif'),15)
% thresh参数可自选，对于rice.tif这张图来说最合适的值大约为14到18
% 使用例子：PicInput=imread('rice.tif');
% a=Lap_edge(PicInput,15);
%－－－－－－－－－－－－－－－－ BEGIN CODE －－－－－－－－－－－－－－－－

% 一，原图像预处理，读入黑白图片并确定长和宽
[m,n]=size(PicInput); %确定图片的长和宽

% 二，拉普拉斯变换预处理，定义镜框矩阵和输出矩阵
r=m+2; %把图片的长和宽各加2
c=n+2;
PicFrame=zeros(r,c); %定义二维数组“PicFrame”，长、宽比“Input”各多2，成为镜框的尺寸
b=zeros(m,n); %定义滤波后的数组

% 三，拉普拉斯运算的三个矩阵
Temp=zeros(3); %定义三阶方阵“Temp”，为临时矩阵
op=[0-10;-14-1;0-10]; %定义拉普拉斯算子
Result=zeros(3); %定义三阶方阵“Result”，为运算结果矩阵

% 四，原图像矩阵处理，做一个“像框” 
%PicFrame(2:m+1,2:n+1)=PicInput; %把原图的矩阵放到新的矩阵“PicFrame”中心，它的第一行、最后一行、第一列、
PicFrame(2:m+1,2:n+1)=PicInput; %最后一列都是“0”，即原图矩阵周围有一圈“0”的边缘，好像给图像加一个像框
figure,
imshow(PicFrame,[]);
PicFrame(1,:)=PicFrame(2,:); %把第二行的值赋给第一行
PicFrame(r,:)=PicFrame(r-1,:); %把倒数第二行的值赋给最后一行
PicFrame(:,1)=PicFrame(:,2); %把第二列的值赋给第一列
PicFrame(:,c)=PicFrame(:,c-1); %把倒数第二列的值赋给最后一列
figure,imshow(PicFrame,[]);

% 五，用拉普拉斯算子进行滤波
for i=1:m
for j=1:n
Temp=PicFrame(i:i+2,j:j+2); %从“PicFrame”矩阵中依次取出三阶方阵，赋值给临时矩阵“Temp”
Result=Temp.*op; %临时矩阵与拉普拉斯算子“点乘”，赋值给结果矩阵“Result”
b(i,j)=sum(sum(Result)); 
%结果矩阵中“十”字线上元素相加，赋值给输出矩阵中相应的位置，
%即临时矩阵中心元素所对应的位
end
end
figure,imshow(b);

% 六，设定阈值，将图像二值化
% thresh=1.618*mean2(abs(b)) 可用黄金分割的比例选阈值优点是边缘清晰 
e=repmat(logical(uint8(0)),m,n); %创建数组
e(find(b>thresh))=1; %阈值判断二值化
PicOut=e; %函数输出
figure,subplot(1,2,1),imshow(PicInput); %显示原图片
title('原图像');
subplot(1,2,2),imshow(e); %显示拉普拉斯边缘检测后的图片
title('自编函数边缘检测后的图像');
%－－－－－－－－－－－－－－－－END OF CODE －－－－－－－－－－－－－－－－