張師語,劉章玉,石仁翠
(1.四川農業大學理學院,雅安 625014;2.四川農業大學信息工程學院,雅安 625014;3.四川農業大學理學院,雅安 625014)
有“數學顯微鏡”之稱的小波變換是近年迅速發展的一種新的時頻分析方法,可以用來邊緣檢測和分割,還可以把短時的物理現象作為瞬態過程來進行分析的分析方法是時頻分析;主要區分突變信號和穩定信號還可以定量分析其能量的分析方法是頻域分析[1]。
在傳統的傅里葉分析中,可以發現信號完全是在頻域里展開的,沒有包含任何時域的信息。僅僅分析信號的時域或頻域性質在很多情況下是遠遠不夠的,例如在電力監測系統中對信號分析的要求,既要求監控穩定信號的成分,又要求故障信號的準確定位[1]。后來人們在傳統傅里葉分析的基礎上引入了時域信息的最初嘗試,即短時傅里葉變換[1],但是其對信號進行分析大多時候只能在一個分辨率上,對許多應用場景來說并不夠準確,存在較大的缺陷。
小波變換正是因為這種缺陷引入的,它克服了短時傅里葉變換只能單分辨率分析處理信號的缺點,具有多分辨率分析信號的特點。小波變換不論在時域還是頻域里都能展現信號的局部信息,可以根據信號的具體形態來動態調整其時間窗和頻率窗。通常情況下,在低頻部分,此時信號是較平穩的狀態,不需要太記錄信號隨時間的變化,因此通過降低時間分辨率來提高頻率分辨率;在高頻部分,頻率的變化情況不明顯,此時可以用較低的頻率分辨率來得到時間定位的精確度[1]。因為可以探測到正常信號中的瞬態部分,并且展示其頻率成分的特性,所以小波變換廣泛應用于各個時頻分析領域[2-4]。
人們提出了窗函數,因為傳統傅里葉變換不能表達隨時間變化的頻率,這種函數被稱為“時間-頻率分析”[1]。時間-頻率窗是靈活可變的,這個“窗”不僅能夠體現頻率的細節信息,而且頻率隨時間的變化是由窗與窗之間反映的。綜上所述,短時傅里葉變換屬于一種時間-頻率分析方法,而且固定窗[1]。
把時域和頻域分解成多個大小相等的小窗口是窗口傅里葉變換的實現過程,其中對信號的任何部分采用的都是相同的時間和頻率分辨率。著名的測不準原理[1]是無論采用任何函數作為窗函數,其時間窗和頻率窗寬度的乘積最小值都是2。該定理說明,同時以任意精度在時間和頻率兩個空間去逼近被測信號這是不可能實現的,在信號的分析上一定要對時間或頻率的精度做權衡取舍。
針對上述情況,小波變換提供了可以根據需要來選擇時間或者頻率的精度的一種靈活性很高的方法。可行的是,找到一組在能量有限信號的函數空間上稠密的正交基函數,并且這組函數可以單純地由一個函數的伸縮和平移生成[1]。在不同的分辨率下分解信號由伸縮實現,可以把這組信號作為窗,來觀察關心的部分是平移的結果。另外,滿足上述要求的函數一定是具有緊支集性,可以在有限的區域內迅速衰減到0,這樣的函數稱為母小波,它生成的一組正交基稱為小波函數[1]。
假設Ψ(t)?L2(R)滿足容許性條件:

即Ψ(t)稱作可允許小波或者基小波。并且是Ψ(t)的傅里葉變換。
那么由可允許小波生成的小波系數可表示為:

離散小波函數可表示為:

離散小波變換的系數可表示為:

類似傅里葉級數,其重構公式為:

由小波變換重構原始信號的逆變換為:


從圖1中的(d)可以發現時間窗和頻率窗組合成的面積在小波變換中相同,隨著|a|的增大,頻率窗逐漸減小,時間窗逐漸增大,因此時間窗和頻率窗的比例可以自行調節的途徑是通過選取不同的a。不用像只能獲取固定時間窗內的固定頻率段的短時傅里葉變換那樣,小波變換可以獲取任何時間內任何感興趣的頻譜。

圖1 各種分析方法的時域與頻域分辨率
在小波域中要實現去噪,二維信號的基本思想與一維信號的情況類似。通常來說,在應用中二維信號表現為圖像信號。去噪效果的重要影響因素是閾值選擇,在去噪過程中既可以使用統一固定的全局閾值,也可以從三個方向,水平、豎直和對角,分別作用閾值抑制其成分[1],這樣可以把所有方向的噪聲分離出來。
假設原始的圖像信號為Xj,k(j,k=1,2,…,N),在受到均值為0,噪聲強度為σ的高斯白噪聲τj,k污染后變為Yj,k,這一過程的基本模型如下:

與一維離散小波變換類似,二維小波變換的實現方法是,從原始信號S開始在每個層次只需要分解上一層次的近似系數。在對每個進行分解的時候,需要從兩個維度分別作用兩次濾波器[1],分別是縱向和橫向,這是與一維離散小波變換不同的地方,如圖2。
小波閾值圖像去噪的基本步驟[5]為:
(1)分解過程:根據噪聲特點選擇合適的小波基和分解層數,對含噪圖像作小波變換;
(2)作用閾值過程:尋找一個合適的數值作為閾值,并用該閾值對每一層的高頻系數進行量化;
(3)重建過程:將處理后的小波系數通過小波逆變換,從而恢復原圖像。
以上去噪過程的依據是,表現為能量密集區域的是攜帶信息的原始信號[1],在頻域或小波域中能量相對集中,其信號分解系數的絕對值比較大;而對于能量譜相對分散的噪聲信號,其系數的絕對值較小。概括來說,含噪圖像信號經過小波變換后,信號的小波系數大于噪聲的小波系數,因此通過作用閾值,保留絕對值大于選取閾值的小波系數,過濾絕對值小于選取閾值的小波系數,從而實現去噪。

圖2 二維小波變換的二層分解示意圖
傳統的閾值處理過程中,常采用的是Donoho提出的軟、硬閾值函數[6]。設Wj,k為原始的圖像小波系數,λ為閾值,為閾值量化后的圖像小波系數。則硬閾值函數表達式為:

可以看出,硬閾值是間斷函數,在處不連續。這導致在圖像信號的間斷點處容易產生比較大的震蕩,進行小波重構時容易出現偽吉布斯現象[7],使得重構的圖像失真。
軟閾值函數表達式為:

之后大量研究學者根據上述軟、硬閾值的優點與不足對閾值函數作出改進,如均方根閾值函數,軟硬折中閾值函數[9]等。本文分別對上述具有代表性的兩種改進閾值函數進行改進,并與其進行去噪效果對比分析。
均方根閾值函數表達式為:

該函數雖然使估計小波系數在| |Wj,k逐漸增大時更趨于原始小波系數Wj,k,但在| |Wj,k=λ處間斷。
引入調節參數μ對均方根閾值函數進行改進,提出第一種改進閾值函數:

(1)第一種改進閾值函數的漸近性分析
當Wj,k→+∞時,
當Wj,k→-∞時,

因此,當Wj,k→∞時,

(2)第一種改進閾值函數的連續性分析
當Wj,k→λ時,

折中閾值函數表達式為:

雖然折中閾值函數減少了閾值處理后的小波系數與原始圖像的小波系數之間存在的恒定偏差,恒定偏差為eλ,e小于1。但仍會影響重構原始輸入圖像信號的逼近程度且在| |Wj,k=λ處不連續。
引入調節參數a對折中閾值函數進行改進,提出第二種改進閾值函數:

(1)第二種改進閾值函數的漸近性分析

當Wj,k→-∞時,

因此,當Wj,k→∞時,

(2)第二種改進閾值函數的連續性分析
當Wj,k→λ時,

當Wj,k→-λ時,


(3)可調因子m的影響分析
當m=0時,第二種改進閾值函數變成軟閾值;當m→+∞時,第二種改進閾值函數變成硬閾值。調節因子m使得第二種改進閾值函數在軟、硬閾值函數之間變動,可以根據實際應用進行合理調節,增大或者減小m的值,從而達到最佳的去噪效果。
在小波系數上作用閾值是小波變換去噪的核心步驟。因為去噪的質量直接決定于閾值的選取,如果閾值選取過小,閾值量化后的小波系數中會含有較多的噪聲分量,導致圖像去噪不充分;如果閾值選取過大,會丟失部分有用的圖像信息,導致重構后的圖像失真[10]。
一般是根據原始輸入信號的信號噪聲比來取得對各層小波系數去噪的閾值,常見的小波閾值選取方法有缺省的閾值確定模型、基于Stein無偏似然估計的軟閾值估計、啟發式SURE閾值等。現在有很多種提取噪聲強度σ的方法,假定噪聲為白噪聲,噪聲的數學期望為0,其噪聲強度一般是用原信號的小波分解的各層系數的標準差來衡量[1]。MATLAB提供wnoisest來實現這個功能。缺省的閾值確定模型如下:
噪聲的小波系數與分解尺度成反比,即分解尺度越大,分解的小波系數越小。對圖像進行去噪處理時,不同的分解尺度,閾值的選取應有所不同,閾值的選取準則應是隨著分解尺度的增大而減小[11]。本文對現有閾值選取算法做適當的改進,其表達式為:
缺省閾值確定模型的閾值在不同尺度上是固定的,而在本文實驗中取為可變的,可有效克服上述幾種閾值選取準則的缺點,且計算相對簡單,可以取得較好的去噪效果,實用性也強。
小波去噪中,MATLAB軟件提供了很多既可以通過某種通用的方式直接對信號去噪,也可以對其中的某個環節靈活控制的靈活使用方式。如wdencmp命令可以用于分層閾值去噪,與其他含噪圖像經小波分解后的系數閾值量化不同的地方在于,傳入命令里的不是單一的閾值,而是一個數組,存放了各層閾值,并且不需要指定保留的層數[1]。
由于只關心系數的絕對值,并不關心系數的位置的閾值處理,所以二維小波變換系數的閾值化方法與一維情況類似。MATLAB提供了wthcoef2命令對二維小波變換閾值化,分別從水平、豎直、對角方向對不同層數作用不同的閾值。可見,分層閾值化同全局閾值化相比,因為結合分解層數和分解方向的閾值化方法能夠在進行閾值化處理時利用更精細的細節信息[1],所以在能量損失不是很大的情況下可以獲得更高的去噪比。
(1)根據噪聲特點選擇合適的分解層數和小波基;
(2)對帶噪圖像Yj,k進行多層小波變換,結果得到一組含噪圖像的小波系數Wj,k;
(3)利用MATLAB提供了wnoisest來計算噪聲強度σ,再根據改進閾值選取公式計算λj;
(4)對各層小波系數做閾值量化處理,具體操作是以λj為比較對象,保留絕對值大于λj的小波系數,舍棄絕對值小于λj的小波系數且置為零,得到閾值處理后的小波系數;
(5)根據處理后的小波系數進行小波逆變換,即可得到去噪后的圖像。
本文使用MATLAB軟件編程實現所有程序,選定Lena圖像作為研究的原始圖像,加入均值為0,方差為0.01的高斯噪聲。
為對去噪后的圖像質量進行量化,評價常常有以下3種評定參數:峰值信噪比PSNR、信噪比SNR、均方誤差MSE[12],它們的計算公式如下:

其中f(x,y)和是原來圖像數據和去噪后的圖像數據。M×N是圖像的大小。
假如去噪信號的信噪比SNR越高,原始信號與去噪信號的均方誤差MSE越小,則去噪信號就越接近于原始信號,去噪效果越好[13]。
針對兩種改進函數和分層閾值函數,利用MAT?LAB中的wmaxlev2命令計算在已知信號長度的情況下,選用不同的小波基分解的最大層數,分別以小波基類型和分解層數作為自變量探究每一種改進函數所能實現的最好效果,如表 1、2、3。

表1 小波基對第一種改進閾值去噪影響

表2 小波基對第二種改進閾值去噪影響

表3 小波基對分層閾值去噪影響
通過表中數據對比發現且參考圖像紋理識別中最優小波基的選擇[14-15],選用coif5作為小波基且分解層數為4層,三種改進函數都能達到最優的效果,故選定其來對圖像進行分解與重構。
使用效果最優傳統閾值函數去噪后的結果圖與使用第一種改進閾值函數、第二種改進閾值函數、分層閾值函數去噪后的結果圖對比,如圖3、4、5。

圖3 含噪Lena圖像第一種改進閾值去噪效果對比

圖4 含噪Lena圖像第二種改進閾值去噪效果對比

圖5 含噪Lena圖像分層閾值去噪效果對比

表4 含噪Lena圖像去噪效果評定指標
由表4可知,相比效果最優傳統閾值所得結果,三種改進方法均能使去噪后的圖像更清晰,邊界更清楚。
本文在傳統軟硬閾值基礎上提出了兩個方向的改進措施:改進閾值函數、分層設置閾值,形成了三種改進方法。在原始圖像中加入高斯噪聲,通過分別改變小波基類型和分解層數來探究每一種改進函數所能實現的最好效果,以期去噪效果更好、邊緣信息保留度更高。通過比較去噪效果圖和評價參數PSNR、SNR、MSE和改變比例的數值,發現三種改進方法都比傳統閾值去噪效果更優,并都更好地保留圖像邊緣細節信息,使圖像更接近原始圖像,驗證了本文所提方法的正確性和有效性。