26第二章边缘检测俊[Shen1985，1992]用前两条准则及相同的边缘信号与噪声模型推出最优无限脉冲响应滤波器为对称指数函数ce-ri的一阶微分，并进一步证明了在存在多个阶跃边缘时，上述滤波器也是最优的无限脉冲滤波器；Sarkar与Boyer[Sarkar1991]给出了类似于Can-ny准则的设计最优两阶微分滤波器的三条准则，并证明了最优两阶微分滤波器的近似表达式为一个多项式乘以对称指数函数.计算机视觉界早就发现，高斯函数的两阶微分可以用两个不同尺度高斯函数的差近似表示，而且当两个高斯函数的尺度比，：02为0.628时，近似程度最好，马领德、李炳成等证明了［Ma1995]，虽然两个高斯函数的差并非一个分析函数的两阶微分，但也满足两阶微分滤波器的条件，因此，任何尺度比的两个高斯滤波器的差都是两阶微分滤波器，而且根据Sarkar与Boyer给出的最优准则，证明了当尺度比为0.176时为最优，马颂德、李炳成进一步证明了任意两个不同尺度的平滑滤波器的差在一定条件下都是两阶微分滤波器，并证明了不同尺度的两个对称指数函数的差在尺度比o：o，为0.3时为最优.总之，目前给出了许多“最优滤波器”，但由于最优准则与边缘模型不完全一样，而且有些是针对一阶微分滤波器的，有些是针对两阶微分滤波器的，因此难以比较.下面，我们将所有提出的“最优滤波器”全部转换成两阶微分滤波器（因为在实际边缘抽取中，检测过零点比检测局部最大值方便），给出Sarkar等人的两阶微分滤波器最优准则，并用这些准则对所有这些滤波器进行比较，设边缘模型为阶跃模型，并假设边缘点位于工一0处，噪声为可加性白噪声，方差为e，输出为信号与两阶微分滤波器h（）的卷积，将输出的有用信号看着是过零点处的导数（该导数越大，过零点越易检测），则容易证明信号噪声比为A|h(0)|SNRh'2(r)d由于A，e2为常数，与h（r）的选择无关，于是定义下式为选择h（r）的信号噪声比最优准则：[h(0) M(2.27)h2(x)dr第二条准则为最优边缘点定位准则.在2.3节中我们已推出过零点在有噪声时偏移六为最优定位准则A，并去掉与h（r）选择无关的e，得的方差为8（见式（2.26）），我们取一f(-y)h'(a)da4h2(r)da将阶跃边缘信号f(r）代入上式并去掉与h（r）选择无关的A，便得到选择h(r）的第二条最优准则：[h(0) /(2.28)A-h2(r)dr在有噪声情况下，用h（）对观察信号进行滤波后，在=x=0附近会有多个过零
272.4最优边缘检测滤波器.点，Canny[Canny1983]给出了两个检测到的过零点的距离平均值：3Jh2(r)drE(d) =(2.29)h'2(r)da上述均值越大，过零点间距越大，越容易检测，Sarkar给出了无限脉冲响应滤波器h(r)的有效宽度度量："h(r)drW(2.30)h(α)dr宽度W越小，则h（a）的能量越集中在r=O附近，当||>W时，h（r）迅速降为零，因此，虽然hα)为无限脉冲响应滤波器，窗口为无穷大，但W表示了滤波器的实际宽度，滤波器的宽度越小，能量越集中，即使式（2.29）表示的过零点平均距离较小，过零点也容易检测，Sarkar将式（2.30)的倒数与式（2.29)相乘后得到的描述多过零点分离的难易程度的准则MRC作为选择h（）的第三条最优准则：h2(r)daMRC=(2.31)h2(r)darh(r)dr用变分原理可以求出在MRC给定的约束条件下使ZA最大的滤波器[Sarkar1991]，但由于得不到分析形式，只能用近似式表示，Sarker用一个多项式与对称指数函数的乘积近似表示h（z）.下面我们将近年来提出的各种两阶微分滤波器及它们的ZA与MRC的值列出（见表2.1），由表2.1可见，各种最优滤波器的性能大致上大同小异：（1）高斯两阶微分滤波器（即Marr[Marr1982]与Canny[Canny1983］提出的滤波器）或称高斯-拉普拉斯两阶微分滤波器（简写√G)(2.32)h(r):/2元g3-（2）两个尺度比为t的高斯滤波器之差（简称DoG)121(2.33)h(α)=（1t2)/2元g2元0tA.该滤波器中t的值为0.176时，A·Z·MRC最大[MA1995],表2.1中所列值为当t=0.176时的值.当t=0.625时，DoG与√G很接近（3）两个尺度比为t的指数滤波器之差（简称DoE)[MA1995]1(1-1号(2.34)h(r):at2(1 -t)g2g
28第二章边缘检测-该滤波器当t=0.3时，A·Z·MRC最大，表2.1中所列最优化准则的值均为t=0.30时的值，L(4）Sarkar滤波器1α+1区α(2.35)h(r)=1+2(3α+5)g3202a当α=0.312时，该滤波器为Sarkar根据以上最优化准则并由变分法推出的最优滤波器，以下表2.1中所列值为α=0.312时的值，（5）对称指数滤波器[Shen1985]该滤波器实际上是DoE滤波器t→+0时的极限(e号-8(r)h(r) =(2.36)g22g（6）Deriche滤波器Je-l1h(a)=1-(2.37)4g3当式（2.35)表示的滤波器取α=一1时，即为式（2.37)表示的滤波器将上述所有滤波器表示式代人式（2.27），（2.28），（2.31)三式表示的三个最优准则，得到下表：表2.1各种滤波器性能比较DoESarkerShenDericheDoG性能VG0.890.920.8010.400.684201. 411.681.321. 711.88MRC01. 241. 21.361. 2AE·MRC0. 79从上表可见，各种滤波器大同小异.在实际应用时，还需考虑下一节将讨论的计算复杂性.下一节我们主要介绍快速算法、2.5边缘检测快速算法以上介绍了边缘检测所涉及的滤波器的设计.用滤波器对信号进行卷积的计算量很大，对于定义在n=1～N的离散信号f(n），若假设滤波器h（n）为无限脉冲响应滤波器，假设f(n）当n<1与n>N时为零，则卷积运算f(n）?hn）的计算量为NXN次乘法与LN×N次加法，若信号为二维信号f（m，n）（n=1，，N，m=1，，N），则计算量为N乘nvW法与N加法，对一般图像信号，N=512，因此一次卷积的计算量为236次乘法与加法，当然，在实际使用时，由于即使是无限脉冲响应滤波器，当|nl>w/2时，h(n）～0(w称为滤波器窗口尺寸，例如对高斯滤波器w=6o），故一次卷积的计算量对一维信号是Nw，对二维信号是N2w2，计算量与窗口大小的平方成正比，对512×512的图像信号，若窗口为5!×5，则一次卷积计算量约为600万次为此，许多研究者提出了各种快速算法，本节首先以一维离散信号为例，给出对称指数滤波器的递推算法，算法的优点是计算量与窗口大小无关.对于上一节所介绍的六种滤波器中的四种（除了第一种√G与第二种DoG)都可以
292.5边缘检测快速算法用类似的方法设计出递推算法，对于其他任意形状的滤波器，李炳成、马颂德提出了一种用Lagurre多项式逼近的方法，使它们也可以用指数滤波器的递推算法实现[Li1994]最后我们将介绍在二维情况下的可分离滤波器，用可分离滤波器与上述递推算法结合，计算复杂性可进一步大大减少.2.5.1、滤波器的递推算法沈俊[Shen1985]首先给出了对称指数滤波器的递推算法.式（2.36)为连续形式下的对称指数两阶微分滤波器，以下为它的离散形式：h(n) =ea -)[e -8(n)](2.38)-1+e"a其中有一些系数与连续形式的h（α）不同，这是因为滤波器必须满足离散化后的归一化条件.归一化条件说明如下：对于k阶差分滤波器h()（α），它一定是某平滑滤波器万(n）的f(n)h(m一n），我们希望信号平k阶差分，平滑滤波器五(n）与信号的卷积为F（m）=滑后不改变均值，由此得到平滑滤波器瓦(n)的归一化条件为h(n）=1.容易证明，若h（n)满足归一化条件，则它的k阶差分的i阶矩(i<k）满足以下归一化条件[Ma1995]：i<k0,(2.39)Znhc"(n) :l(-1)i,i=k对两阶差分滤波器，k=2.容易验证，式（2.38）满足式（2.39）表示的归一化条件.对于上一节所给出的连续形式的滤波器离散化时，都要重新调整系数，使它满足归一化条件，将式（2.38)表示的滤波器对信号f(n）卷积得g(n) = f(n)@h(n) = ci(f(n) @ce-l - f(n))其中c1与C2为式（2.38)中的两个系数.由上式可见，对信号用式（2.38)滤波等同于将信号与c2e（它是对称指数平滑滤波器）卷积后与原信号的差.滤波过程可用图2.19表示.因此，我们只需介绍实现c2e-滤波器的快速算法，+输出输人图2.19两阶差分指数滤波器的实现下面，我们将输入信号，输出信号与滤波器c2e-分别记为z(n)，y(n),h(n），将它们的Z变换分别记作X（），Y（）,H(z），于是我们有Y）=H()X（）滤波器h(n)=cze-的双边变换为
30第二章边缘检测e-"=c(..+e+e+1+e-l+e2-2+..)H(z) =zerH+ (z) +H (z)1-e-.-1-e-上式中H(z)与H_(z)分别为Cze(2.40)H+(z) =1-e-zC2(2.41)H_(z) =1e-21因为Y()=H()X(2)=(H(z)+H_()),X(z),将Y(z)记作Y(z)=Y+(z)+Y_(z),其中(2.42)Y(z)=H-(z)X(z)Y+ (2) = H+ (z)X(2),由上式，输出y(t)被分解成两部分，下面分别考察这两部分：由式（2.41)与（2.42)，有c2X(2) =Y. (2)[1 -e--1上式的反变换为y- (n) -e-y- (n - 1) = czr(n)于是，我们得到y-(n)的如下递推公式：(2.43)y- (n) =e-y- (n -1) + czr(n)上式表示，在求出y-(n一1)时，y-(n)可由y-(n—1)与z(n)得到，而且求y-(n)只涉及两次乘法与一次加法，与滤波器的窗口大小无关，类似地，由式（2.40)与（2.42)，有Y+ (2)[1 e-]= Caze-X(2)上式的Z反变换为y+ (n) -e-y+ (n + 1) = ce-a(n + 1)于是得到如下递推公式：(2.44)y+ (n) = e-(y+ (n + 1) + c2a(n + 1))注意，上式的递推方向与式(2.43)相反，即在已知+（n十1)时求y+（n)（自右向左），计算次数为一次乘法(c2与r(n十1)的乘积已由式（2.43)中求出）与一次加法将y+(n)与y_(n)相加，得到输出y(n)=y+(n)+y-(n).在递推计算时，信号两边初始点（n=1，或n=N）处，可以假设y-（0)与y+（N+1)为零.算法可用图2.20表示：上述求（n）的算法中，只包含3N次乘法与2N次加法，与指数滤波器窗口大小无关，这是这种递推算法的最大优点。若p(n)为关于n 的多项式，凡是形式为p(n)e-的滤波器（例如2. 4节中的DoE，Sarker，Shen，Deriche等滤波器）均可用类似以上递推滤波器来实现.这一点可以用以下Z变换的性质看出来：