图像处理－线性滤波－2 图像微分（1、2阶导数和拉普拉斯算子）

逃离我推掉我的手 2022-07-24 09:26 177阅读 0赞

From:http://www.cnblogs.com/pegasus/archive/2011/05/20/2051780.html

更复杂些的滤波算子一般是先利用高斯滤波来平滑，然后计算其1阶和2阶微分。由于它们滤除高频和低频，因此称为带通滤波器（band-pass filters）。

在介绍具体的带通滤波器前，先介绍必备的图像微分知识。

## 1 一阶导数 ##

连续函数，其微分可表达为[![image][]][image 1] ，或[![image][image 2]][image_image 2]                         （1.1）

**对于离散情况（图像），其导数必须用差分方差来近似**，有

[ ![image][image 3]][_image_image 3]，前向差分 forward differencing                  （1.2）

[ ![image][image 4]][_image_image 4] ，中心差分 central differencing                     （1.3）

**1）前向差分的Matlab实现**

> <table> 
>  <tbody> 
>   <tr> 
>    <td> 
>     <div>
>       1 
>     </div> 
>     <div>
>       2 
>     </div> 
>     <div>
>       3 
>     </div> 
>     <div>
>       4 
>     </div> 
>     <div>
>       5 
>     </div> 
>     <div>
>       6 
>     </div> 
>     <div>
>       7 
>     </div> 
>     <div>
>       8 
>     </div> 
>     <div>
>       9 
>     </div> 
>     <div>
>       10 
>     </div> 
>     <div>
>       11 
>     </div> 
>     <div>
>       12 
>     </div> 
>     <div>
>       13 
>     </div> 
>     <div>
>       14 
>     </div> 
>     <div>
>       15 
>     </div> 
>     <div>
>       16 
>     </div> 
>     <div>
>       17 
>     </div> 
>     <div>
>       18 
>     </div> 
>     <div>
>       19 
>     </div> 
>     <div>
>       20 
>     </div> 
>     <div>
>       21 
>     </div> 
>     <div>
>       22 
>     </div> 
>     <div>
>       23 
>     </div> 
>     <div>
>       24 
>     </div> 
>     <div>
>       25 
>     </div> 
>     <div>
>       26 
>     </div> 
>     <div>
>       27 
>     </div> </td> 
>    <td> 
>     <div> 
>      <div> 
>       <code>function dimg = mipforwarddiff(img,direction)</code> 
>      </div> 
>      <div> 
>       <code>% MIPFORWARDDIFF&nbsp;&nbsp;&nbsp;&nbsp; Finite difference calculations </code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; DIMG = MIPFORWARDDIFF(IMG,DIRECTION)</code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; Calculates the forward-difference </code> 
>       <code>for</code> 
>       <code>a given direction</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; IMG&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; : input image</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; DIRECTION : </code> 
>       <code>'dx'</code> 
>       <code>or </code> 
>       <code>'dy'</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; DIMG&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; : resultant image</code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; See also MIPCENTRALDIFF MIPBACKWARDDIFF MIPSECONDDERIV</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; MIPSECONDPARTIALDERIV</code> 
>      </div> 
>      <div>
>        &nbsp; 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; Omer Demirkaya, Musa Asyali, Prasana Shaoo, ... 9/1/06</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; Medical Image Processing Toolbox</code> 
>      </div> 
>      <div>
>        &nbsp; 
>      </div> 
>      <div> 
>       <code>imgPad = padarray(img,[1 1],</code> 
>       <code>'symmetric'</code> 
>       <code>,</code> 
>       <code>'both'</code> 
>       <code>);%将原图像的边界扩展</code> 
>      </div> 
>      <div> 
>       <code>[row,col] = size(imgPad);</code> 
>      </div> 
>      <div> 
>       <code>dimg = zeros(row,col);</code> 
>      </div> 
>      <div> 
>       <code>switch</code> 
>       <code>(direction)&nbsp;&nbsp; </code> 
>      </div> 
>      <div> 
>       <code>case</code> 
>       <code> 'dx'</code> 
>       <code>,</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;</code> 
>       <code>dimg(:,1:col-1) = imgPad(:,2:col)-imgPad(:,1:col-1);%x方向差分计算，</code> 
>      </div> 
>      <div> 
>       <code>case</code> 
>       <code> 'dy'</code> 
>       <code>,</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;</code> 
>       <code>dimg(1:row-1,:) = imgPad(2:row,:)-imgPad(1:row-1,:); </code> 
>      </div> 
>      <div> 
>       <code>otherwise, disp(</code> 
>       <code>'Direction is unknown'</code> 
>       <code>);</code> 
>      </div> 
>      <div> 
>       <code>end;</code> 
>      </div> 
>      <div> 
>       <code>dimg = dimg(2:end-1,2:end-1);</code> 
>      </div> 
>     </div> </td> 
>   </tr> 
>  </tbody> 
> </table>

**2）中心差分的Matlab实现**

> <table> 
>  <tbody> 
>   <tr> 
>    <td> 
>     <div>
>       1 
>     </div> 
>     <div>
>       2 
>     </div> 
>     <div>
>       3 
>     </div> 
>     <div>
>       4 
>     </div> 
>     <div>
>       5 
>     </div> 
>     <div>
>       6 
>     </div> 
>     <div>
>       7 
>     </div> 
>     <div>
>       8 
>     </div> 
>     <div>
>       9 
>     </div> 
>     <div>
>       10 
>     </div> 
>     <div>
>       11 
>     </div> 
>     <div>
>       12 
>     </div> 
>     <div>
>       13 
>     </div> 
>     <div>
>       14 
>     </div> 
>     <div>
>       15 
>     </div> 
>     <div>
>       16 
>     </div> 
>     <div>
>       17 
>     </div> 
>     <div>
>       18 
>     </div> 
>     <div>
>       19 
>     </div> 
>     <div>
>       20 
>     </div> 
>     <div>
>       21 
>     </div> 
>     <div>
>       22 
>     </div> 
>     <div>
>       23 
>     </div> 
>     <div>
>       24 
>     </div> 
>     <div>
>       25 
>     </div> 
>     <div>
>       26 
>     </div> 
>     <div>
>       27 
>     </div> 
>     <div>
>       28 
>     </div> </td> 
>    <td> 
>     <div> 
>      <div> 
>       <code>function dimg = mipcentraldiff(img,direction)</code> 
>      </div> 
>      <div> 
>       <code>% MIPCENTRALDIFF&nbsp;&nbsp;&nbsp;&nbsp; Finite difference calculations </code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; DIMG = MIPCENTRALDIFF(IMG,DIRECTION)</code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; Calculates the central-difference </code> 
>       <code>for</code> 
>       <code>a given direction</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; IMG&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; : input image</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; DIRECTION : </code> 
>       <code>'dx'</code> 
>       <code>or </code> 
>       <code>'dy'</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp; DIMG&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; : resultant image</code> 
>      </div> 
>      <div> 
>       <code>%</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; See also MIPFORWARDDIFF MIPBACKWARDDIFF MIPSECONDDERIV</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; MIPSECONDPARTIALDERIV</code> 
>      </div> 
>      <div>
>        &nbsp; 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; Omer Demirkaya, Musa Asyali, Prasana Shaoo, ... 9/1/06</code> 
>      </div> 
>      <div> 
>       <code>%&nbsp;&nbsp; Medical Image Processing Toolbox</code> 
>      </div> 
>      <div>
>        &nbsp; 
>      </div> 
>      <div> 
>       <code>img = padarray(img,[1 1],</code> 
>       <code>'symmetric'</code> 
>       <code>,</code> 
>       <code>'both'</code> 
>       <code>);</code> 
>      </div> 
>      <div> 
>       <code>[row,col] = size(img);</code> 
>      </div> 
>      <div> 
>       <code>dimg = zeros(row,col);</code> 
>      </div> 
>      <div> 
>       <code>switch</code> 
>       <code>(direction)</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>case</code> 
>       <code>'dx'</code> 
>       <code>,</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>dimg(:,2:col-1) = (img(:,3:col)-img(:,1:col-2))/2;</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>case</code> 
>       <code>'dy'</code> 
>       <code>,</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>dimg(2:row-1,:) = (img(3:row,:)-img(1:row-2,:))/2;</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>otherwise,</code> 
>      </div> 
>      <div> 
>       <code>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;</code> 
>       <code>disp(</code> 
>       <code>'Direction is unknown'</code> 
>       <code>);</code> 
>      </div> 
>      <div> 
>       <code>end</code> 
>      </div> 
>      <div> 
>       <code>dimg = dimg(2:end-1,2:end-1);</code> 
>      </div> 
>     </div> </td> 
>   </tr> 
>  </tbody> 
> </table>
> 
> <table> 
>  <tbody> 
>   <tr> 
>    <td> 
>     <div>
>       1 
>     </div> </td> 
>    <td> 
>     <div> 
>      <div> 
>       <code>&nbsp;</code>&nbsp; 
>      </div> 
>     </div> </td> 
>   </tr> 
>  </tbody> 
> </table>

**实例：技术图像x方向导数**

> <table> 
>  <tbody> 
>   <tr> 
>    <td> 
>     <div>
>       1 
>     </div> 
>     <div>
>       2 
>     </div> </td> 
>    <td> 
>     <div> 
>      <div> 
>       <code>I = imread(</code> 
>       <code>'coins.png'</code> 
>       <code>); figure; imshow(I);</code> 
>      </div> 
>      <div> 
>       <code>Id = mipforwarddiff(I,</code> 
>       <code>'dx'</code> 
>       <code>); figure, imshow(Id);</code> 
>      </div> 
>     </div> </td> 
>   </tr> 
>  </tbody> 
> </table>

[ ![image][image 5]][_image_image 5] [ ![image][image 6]][_image_image 6]

>     原图像                                                   x方向1阶导数

## 2 图像梯度（Image Gradient） ##

图像I的梯度定义为[![image][image 7]][image_image 7]  ，其幅值为[![image][image 8]][image_image 8] 。出于计算性能考虑，幅值也可用[![image][image 9]][image_image 9] 来近似。

Matlab函数

> 1）gradient：梯度计算
> 
> 2）quiver：以箭头形状绘制梯度。注意放大下面最右侧图可看到箭头，由于这里计算横竖两个方向的梯度，因此箭头方向都是水平或垂直的。

**实例：仍采用上面的原始图像**

> <table> 
>  <tbody> 
>   <tr> 
>    <td> 
>     <div>
>       1 
>     </div> 
>     <div>
>       2 
>     </div> 
>     <div>
>       3 
>     </div> 
>     <div>
>       4 
>     </div> 
>     <div>
>       5 
>     </div> </td> 
>    <td> 
>     <div> 
>      <div> 
>       <code>I = </code> 
>       <code>double</code> 
>       <code>(imread(</code> 
>       <code>'coins.png'</code> 
>       <code>));</code> 
>      </div> 
>      <div> 
>       <code>[dx,dy]=gradient(I);</code> 
>      </div> 
>      <div> 
>       <code>magnitudeI=sqrt(dx.^2+dy.^2);</code> 
>      </div> 
>      <div> 
>       <code>figure;imagesc(magnitudeI);colormap(gray);%梯度幅值</code> 
>      </div> 
>      <div> 
>       <code>hold </code> 
>       <code>on</code> 
>       <code>;quiver(dx,dy);%叠加梯度方向</code> 
>      </div> 
>     </div> </td> 
>   </tr> 
>  </tbody> 
> </table>

[ ![image][image 10]][_image_image 10] [ ![image][image 11]][_image_image 11]

梯度幅值                                   梯度幅值+梯度方向

## 3 二阶导数 ##

对于一维函数，其二阶导数[![image][image 12]][image_image 12] ，即[![image][image 13]][image_image 13] 。它的差分函数为

[ ![image][image 14]][_image_image 14]                  （3.1）

## 3.1 普拉斯算子（laplacian operator） ##

### 3.1.2 概念 ###

**拉普拉斯算子是n维欧式空间的一个二阶微分算子**。它定义为两个梯度向量算子的内积

[ ![image][image 15]][_image_image 15]       （3.2）

其在二维空间上的公式为：    ![image][image 16]                （3.3)

对于1维离散情况，其二阶导数变为二阶差分

> 1）首先，其一阶差分为[![image][image 17]][image_image 17]

> 2）因此，二阶差分为

>            [ ![image][image 18]][_image_image 18]
> 
> 3）因此，**1维拉普拉斯运算可以通过1维卷积核**[**![image][image 19]**][image_image 19] **实现**
> 
>

对于2维离散情况（图像），拉普拉斯算子是2个维上二阶差分的和（见式3.3），其公式为：

[![image][image 20]][image_image 20]   （3.4）

上式对应的卷积核为

[ ![image][image 21]][_image_image 21]

常用的拉普拉斯核有：

[ ![image][image 22]][_image_image 22]

### 3.1.2 应用 ###

拉普拉斯算子会突出像素值快速变化的区域，因此常用于边缘检测。

Matlab里有两个函数

1）del2

> 计算公式：[![image][image 23]][image_image 23] ，[![image][image 24]][image_image 24]

2）fspecial：图像处理中一般利用Matlab函数fspecial

h = fspecial('laplacian', alpha) returns a 3-by-3 filter approximating the shape of the two-dimensional Laplacian operator.  
The parameter alpha controls the shape of the Laplacian and must be in the range 0.0 to 1.0. The default value for alpha is 0.2.

### 3.1.3 资源 ###

[http://fourier.eng.hmc.edu/e161/lectures/gradient/node8.html][http_fourier.eng.hmc.edu_e161_lectures_gradient_node8.html] （非常清晰的Laplacian Operator介绍，本文的主要参考）

[http://homepages.inf.ed.ac.uk/rbf/HIPR2/log.htm][http_homepages.inf.ed.ac.uk_rbf_HIPR2_log.htm]

[image]: /images/20220722/6da6e8ba8cf140f4a6dd22302cc20354.png
[image 1]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201051588495.png
[image 2]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201051582922.png
[image_image 2]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201051589018.png
[image 3]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201051593968.png
[_image_image 3]: /images/20220722/226abba04ff54644bb050f81eb4e1309.png
[image 4]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201051598919.png
[_image_image 4]: /images/20220722/cb4a51ac5c8a408e9307b59371e12875.png
[image 5]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/20110520105201706.png
[_image_image 5]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052006279.png
[image 6]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261626257058.png
[_image_image 6]: /images/20220722/3c404fc876384776b0c7815ff327f2c9.png
[image 7]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052023638.png
[image_image 7]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052016495.png
[image 8]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/20110520105202574.png
[image_image 8]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052029494.png
[image 9]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052024968.png
[image_image 9]: /images/20220722/cb097dd5d72046708bf6624e91d35c62.png
[image 10]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052036853.png
[_image_image 10]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052033822.png
[image 11]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261626301513.png
[_image_image 11]: /images/20220722/5abc476047cc4aedb5d6ed51d1d30db0.png
[image 12]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052052752.png
[image_image 12]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052051672.png
[image 13]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201052068192.png
[image_image 13]: /images/20220722/bc08762bc40a466fa521bdba7d96dcb8.png
[image 14]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201053147439.png
[_image_image 14]: /images/20220722/b1a8b644fd4d4d5bb6189e4673d49aed.png
[image 15]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617359687.png
[_image_image 15]: /images/20220722/52fd4bfbed1e41fabce6fa2efb1cfc22.png
[image 16]: /images/20220722/a2e7b827c16a4360aebccff491b428dd.png
[image 17]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617367080.png
[image_image 17]: /images/20220722/c30826fb0c6c4c568485484d6b0dab61.png
[image 18]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617386010.png
[_image_image 18]: /images/20220722/5f0c512950ec49d88aefae112f20d4d0.png
[image 19]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617399291.png
[image_image 19]: /images/20220722/499527a45d2c40a787d5d426391f6beb.png
[image 20]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617442384.png
[image_image 20]: /images/20220722/0bedf851cd514c8b85ad18d204f2f1af.png
[image 21]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617441238.png
[_image_image 21]: /images/20220722/3a8d02a3e9044e02b5075a841421e982.png
[image 22]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105261617457791.png
[_image_image 22]: /images/20220722/47cc818c7509412fa10fb842f8421d6b.png
[image 23]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201053189518.png
[image_image 23]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201053173663.png
[image 24]: http://images.cnblogs.com/cnblogs_com/pegasus/201105/201105201053187467.png
[image_image 24]: /images/20220722/533a6672f14d42709033ae7a46911ed7.png
[http_fourier.eng.hmc.edu_e161_lectures_gradient_node8.html]: http://fourier.eng.hmc.edu/e161/lectures/gradient/node8.html
[http_homepages.inf.ed.ac.uk_rbf_HIPR2_log.htm]: http://homepages.inf.ed.ac.uk/rbf/HIPR2/log.htm