一种自适应的全变分图像去噪算法

一种自适应的全变分图像去噪算法

1非线性图像去噪方法图像噪声去除是数据处理的重要部分，在实际的图像采集、传输和存储过程中不可避免地会引入各种噪声。噪声有很多类型和原因。在大多数情况下，需要处理图像以消除噪声，以便将处理后的图像适合分析和提取信息。有很多方法可以去除噪声。在频域处理中，有以高滤波为代表的线性滤波和以中值滤波为代表的非线性滤波。由于其简单、复杂、低噪声计算，高滤波是一种广泛应用的噪声去除方法，但由于噪声和信号的处理，它可以去除噪声，削弱图像的亮度，并在图像边缘移动。因此，在这个问题上，提出了许多非线性噪声去除方法。其中，基于微分法的图像消噪是一种具有代表性的方法。由于这一方法因为偏离差分方程是数学的一个重要分支，因此形成了一个完整的理论体系和数值方法。同时，该方法从新的角度（如能量扩散、横向发展等）来解释图像噪声的去除过程。其中，这种方法的代表性是pm方程和全变分算法（图像遗传，模型）。本文主要讨论的是TV算法,该方法自1992年由Rudin、Osher和Fatime首次提出后,因其去噪效果明显而引起广泛关注,但这种方法也有自身的缺点,如阶梯效应、需要已知图像的噪声方差等.针对这些问题,本文提出一种自适应全变分图像去噪算法,该算法主要是修改传统TV算法中的逼近项,用受噪声影响较小的模糊图像代替原始的含噪声图像,使得在整个图像处理的计算过程中不需要已知图像的噪声方差,便可较好地对图像进行去噪处理,达到显著的去噪效果,同时依据边缘信息对模型中拉格朗日因子进行加权处理,使得处理后不仅不会模糊图像的边缘,而且较好地抑制传统TV算法所产生的阶梯效应.2算法的设计和实现2.1tv模型中的抽象函数传统TV算法就是最小化能量泛函:E=∫Ω(|∇u|+12λ(u-u0)2)dxdy,(1)E=∫Ω(|∇u|+12λ(u−u0)2)dxdy,(1)其中,∫Ω|∇u|dxdy∫Ω|∇u|dxdy和∫Ω12λ(u-u0)2dxdy∫Ω12λ(u−u0)2dxdy分别为TV模型的正则项(RegularizationTerm)和逼近项(FidelityTerm),Ω为图像区域,λ为拉格朗日因子,u0为初始含有噪声的图像,u为经处理后的图像,式(1)的欧拉-拉格朗日(Euler-Lagrange,E-L)方程如下:div(∇u|∇u|)+λ(u0-u)=0.(∇u|∇u|)+λ(u0−u)=0.由最速下降法(SteepestDescentMethod)可解得图像以时间为演化参数的演化公式如下:{ut=div(∇u|∇u|)+λ(u0-u)u|t=0=u0λ=1σ2|Ω|∫Ωdiv(∇u|∇u|)(u-u0)dxdy⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪ut=div(∇u|∇u|)+λ(u0−u)u|t=0=u0λ=1σ2|Ω|∫Ωdiv(∇u|∇u|)(u−u0)dxdy其中,|Ω||Ω|为图像面积,λ也被称为是全局的尺度因子.从上式可看出,TV模型中关于逼近项的演化公式存在着2个问题:1)λ的求解需要已知图像噪声方差σ2,并且λ是一个全局变量,但图像中不同部分(如边缘、背景等)的信噪比是不同的,用相同的λ作为逼近项解前的系数显然不是一个理想的选择;2)当原图像噪声较大(σ2较大)时,u0及其相关的逼近项就会存在较大误差,同时σ2与λ成反比,λ值相对较小,使得演化公式主要由正则项的解div(∇u|∇u|)(∇u|∇u|)决定,所以在处理边缘时很容易出现阶梯效应.实际情况也是当图像噪声越大时,阶梯效应越明显.2.2拉格朗日因子求解本文解决上述问题的出发点是逼近项中的u0,因为计算λ需要的噪声方差是由u0决定,当噪声较大时,u0还会引起较大的误差,所以本文对u0进行一次高斯滤波得到一个较模糊的图像ug=u0*G(G为高斯核),用ug代替u0,这样逼近项就变成∫Ω12λ(u-ug)2dxdy∫Ω12λ(u−ug)2dxdy,新的能量泛函如下:Enew=∫Ω(|∇u|+12λ(u-ug)2)dxdyEnew=∫Ω(|∇u|+12λ(u−ug)2)dxdy,其欧拉-拉格朗日方程如下:div(∇u|∇u|)+λ(ug-u)=0,这样得到新的演化公式:ut=div(∇u|∇u|)+λ(ug-u).(2)如2.1节小结所述,将λ作为一个全局变量不是一个理想的选择,λ应该根据图像局部信息的不同而不同,所以将上式两边同时乘上(ug-u)并在局部图像区域Ω0(Ω0是以当前像素点为中心点大小为(2m+1)×(2m+1)的窗口)进行积分,当达到稳定解时,ut将趋于0,所以得0=∫Ω0div(∇u|∇u|)(ug-u)dxdy+λ∫Ω0(ug-u)2dxdy.(3)下面的重点是如何求解拉格朗日因子λ,而图像局部信息主要是边缘信息,所以接下来本文分边缘区域和非边缘区域进行讨论.1)边缘区域.不失一般性,假设u0组成形式如下:u0=u*+n,其中,u*为无噪声有效信号,n为均值为0的噪声信号,由于高斯滤波的去噪作用,且在边缘区域,图像的有效信号功率通常要远大于噪声信号功率,即P(u*)≫P(n),所以可得ug=u0*G≈u**G.(4)同时对图像信号的高斯滤波可用能量扩散的角度来解释,其过程同样可用以时间为参数的演化公式来表示:u*t=u*xx+u*yy,(5)其中,u*t表示u*对时间的导数,u*xx和u*yy表示u*对图像坐标方向的二阶导数,由式(4)和式(5)可得到一个ug的近似表达式:ug≈u**G=u*+u*t=u*+u*xx+u*yy,将上式代入式(3)等号右边的第二项,则第二项变为λ∫Ω0(u*+u*xx+u*yy-u)2dxdy,将上式分解为注意到u*为无噪声有效信号,显然u演化的目标便是不断地趋近于u*,所以当达到稳定解的时候,新的模型中存在一个约束条件:∫Ω0(u-u*)2dxdy=0,上式的条件表示当ut将趋于0时,u将趋于u*,等价对于任一常数C,∫Ω0C(u-u*)dxdy=0,此时,式(6)的第一项和第三项被消除,同时当达到稳定解的时候,u可近似代替u*,所以式(6)可近似等价为λ∫Ω0(uxx+uyy)2dxdy,将上式替换式(3)等号右边的第二项,可得到求解λ表达式为λ=∫Ω0div(∇u|∇u|)(u-ug)dxdy∫Ω0(uxx+uyy)2dxdy+ε‚(7)其中,ε是一个大于0的常数,引入该常数是为了避免出现∫Ω0(uxx+uyy)2dxdy取值接近0的情况.2)非边缘区域.非边缘区域是图像的背景和平坦区域,其信号主要是由噪声组成,而且不需要考虑边缘区域存在的阶梯效应.在非边缘区域,并不能按式(7)计算拉格朗日因子λ,因为式(7)的成立依赖于式(4)的近似u0*G≈u**G,而这一近似表达成立的前提条件是在边缘区域图像的有效信号功率通常要远大于噪声信号功率,即P(u*)≫P(n),而在非边缘区域这一条件很可能是不成立的,尤其是在噪声较大的时候,所以在非边缘区域图像的演化公式应该有所不同.本文从演化公式中正则项与逼近项的作用来考虑非边缘区域的图像演化.图像噪声的去除是通过正则项来实现,它现实对噪声信号的平滑.虽然在平滑过程中正则项容许边缘的不连续,但逼近项的作用则是要求经过演化后的图像近似地趋于指定图像ug(针对演化式(2)而言),从某种程度来说,正则项与逼近项作用是相反的,两项经过演化后以拉格朗日因子λ为比例权重到达平衡.在非边缘区域中逼近项的作用实际上是多余的,因为在这些区域中并不存在图像的边缘信息和纹理信息,所以使这些区域中的像素点与其周围像素点共同趋于某一相同的灰度值是合理的,完全没有必要趋于指定图像ug对应的灰度值.因此在图像的非边缘区域,演化公式可完全由正则项部分(div(∇u|∇u|))决定,这样不仅可以避免拉格朗日因子λ不能由式(7)计算的问题,而且消除逼近项在演化公式中的作用,从而加速图像的演化,有利于噪声的去除,同时也不必考虑边缘区域才会出现的阶梯效应.2.3保证边缘信号响应的区域因为前面是分边缘区域与非边缘区域讨论演化公式,所以如何区分边缘区域与非边缘区域就成了一个需要解决的问题.但是否需要准确区分两个区域,答案是否定的.首先前面所讨论的演化公式并不是基于准确区分两个区域的,所以准确区分没有必要;其次,如果准确区分两个区域会直接导致两个问题:两个区域的演化公式不能统一;两个区域的相邻区域很可能存在相邻且灰度值相近的像素点却按照差异很大的演化公式进行演化,从而使演化的结果差异很大,显然这也是不合理的.所以本文需要一种度量来表征图像中某一像素点属于边缘区域或非边缘区域的程度,在无噪声情况下,这种度量必须满足2个性质:1)单调性,它需要保证其对边缘信号的响应在边缘区域大于(或小于)在非边缘区域;2)连续性,即当从边缘区域过渡到非边缘区域(或相反方向),该度量在数值是连续变化的.在有噪声的情况下,这种度量即使不能严格满足这两个性质,也需要近似满足这两个性质.从本文后面的结果可看出,采用满足这两个性质的度量对演化公式进行修正,可最好地统一边缘区域与非边缘区域的演化公式,同时也可避免相邻且灰度值相近的像素点因演化公式的巨大差异而产生差异很大的演化结果.本文用局部方差作为上述的度量,显然局部方差在边缘区域对边缘信号的响应要大于非边缘区域,且其也满足上述所说的连续性,当然为了最大程度避免噪声的影响,同时基于我们不需要准确区分边缘区域与非边缘区域这一前提,本文并不对u0计算局部方差,而对模糊图像ug计算局部方差,计算结果为V,计算局部方差的区域是以当前像素点为中心点大小为(2l+1)×(2l+1)的窗口,如果当前像素点的图像坐标为(x,y),则V(x,y)计算表达式如下:V(x,y)=1(2l+1)2i=l∑i=-lj=l∑j=-l[ug(x+i,y+j)-1(2l+1)2i=l∑i=-lj=l∑j=-lug(x+i,y+j)]2.(8)但不同图像的局部方差在数值上存在明显的差异,所以本文对上述的局部方差V作了归一化处理,使V线性映射到[01]的区间内,其最终结果为—V.显然当归一化后的局部方差越趋于0,说明当前位置的像素点越可能属于非边缘区域,反之当该值越趋于1,说明当前位置的像素点越可能属于边缘区域.所以用—V乘以演化公式中的逼近项部分,便得到了统一的演化公式:{ut=div(∇u|∇u|)+λ—V(ug-u)u|t=0=u0λ=∫Ω0div(∇u|∇u|)(u-ug)dxdy∫Ω0(uxx+uyy)2dxdy+ε(9)从上式可看出,在边缘区域时,—V→1,上式由式(2)和式(7)决定,符合边缘区域的演化过程,在非边缘区域时,—V→0,上式由演化公式的正则项部分决定,这也同样符合非边缘区域的演化过程.3边缘区域犯罪性检验为验证本文算法,本文用两幅8位的灰度图像(图1(b))进行了实验((b)分别是在(a)上加入均值为0,方差为50的高斯噪声).并将本文算法与传统TV算法在整幅图像与相同位置的局部放大图像的结果进行对比,如图1所示.需要说明的是,在本文的实验中式(9)中局部图像区域Ω0与式(8)计算局部方差的区域均选择相同大小为5×5的窗口,式(9)中常量ε取值为1,并且实验的计算过程中本文算法与传统TV算法均选取相同的迭代规则和终止条件,在同一个λ下迭代10次,迭代步长为0.2,当图像平均每一个像素的灰度值改变量小于0.08时终止迭代,各结果图λ的迭代次数如下所示.图1(c)左图λ的迭代次数为76,右图为74;图1(d)左图λ的迭代次数为77,右图为70.假设图像的大小为M×N,则传统TV算法与本文算法每迭代一次λ所对应的计算复杂度为O(MN)和O((2m+1)2MN).从实验结果图可很显然地看出本文算法在不需要已知图像噪声方差的情况,与传统TV算法(如果不知道噪声方差,传统TV算法是无法达到如图1(c)所示的去噪效果)一样可较好去除图像的噪声,不会像线性滤波那样模糊图像边缘,较好地保持边缘的对比度,同时还可看出传统TV算法结果的阶梯效应要比本文算法明显的多,这一点从图1中两者的局部放大图(e)和(f)可很明显的看出.所以本文算法从已知条件约束和最终的结果图来看要明显地优于传统TV算法.接下来本文对图像处理的结果进行评价,显然信噪比(SignaltoNoiseRatio,SNR)是评价图像去噪效果的一个重要指标,但仅用这一个指标是不够的,因为SNR并不能反应出图像边缘区域的图像质量,如果去噪算法产生较明显的阶梯效应,就会产生不规则的人造边缘,这样图像质量明显是不好的,但由于阶梯效应的产生,使得这些边缘处会产生更大的灰度值变化,因此它的信噪比也会更大.所以本文需要一个辅助指标,在图像噪声得到明显的改善,SNR有了大幅度提升,图像边缘没有被模糊的前提下,来评价图像边缘的平滑情况,因为阶梯效应越明显,边缘平滑性越差.以∫Ω|∇u|2dxdy作为一幅图像平滑性的度量是一种合理的选择,因为它不充许边缘有跳变,实际如果以∫Ω|∇u|2dxdy作为图像的能量泛函并使其最小化,那么图像的演化公式就是线性扩散方程,如式(5),即高斯滤波,所以∫Ω|∇u|2dxdy的值越小说明图像越平滑.但本文需要考虑的是边缘区域的平滑情况,所以需要对∫Ω|∇u|2dxdy进行加权平均处理,边缘区域权重高,非边缘区域权重低,2.3节所提到的归一化局部方差ˉV满足这权重关系.结合上面的叙述,本文提出边缘区域的平滑性(SmoothnessofEd