文/殷仕淑 武岳 常郝
對多維信號處理的研究工作最早可以追溯到上世紀(jì)60年代。D.P. Petersen 與D. Middleton在文獻(xiàn)[1]中首次開展了對多維頻帶有限信號進(jìn)行抽樣和重建工作的理論研究。此后,在多維數(shù)字濾波器的研究領(lǐng)域里,致力于二維數(shù)字濾波器理論及設(shè)計(jì)方法的研究工作占據(jù)絕大部分分量,發(fā)展很快。相較而言,由于數(shù)據(jù)量大,系統(tǒng)復(fù)雜度增加,處理三維空間信號(例如三維醫(yī)學(xué)影像信號,地球物理空間信號等)的三維數(shù)字濾波器(Three dimensional digital filters: 3D DF)的相關(guān)理論及設(shè)計(jì)工作進(jìn)展不大。設(shè)計(jì)一個(gè)長度為10的FIR一維數(shù)字濾波器,最多有10個(gè)未知變量待求。設(shè)計(jì)一個(gè)三維長度都為10的FIR三維數(shù)字濾波器,待求的未知變量數(shù)為1000。雖然濾波器的設(shè)計(jì)問題描述相同,但龐大的未知變量數(shù)對于算法的有效收斂是個(gè)極大考驗(yàn)。雖然可以利用低維濾波器去處理三維信號,但考慮到數(shù)據(jù)之間的相關(guān)性,尤其是在移動圖像處理、機(jī)器視覺、醫(yī)學(xué)圖像處理、地震信號處理和其他地球物理等三維空間信號處理領(lǐng)域中,發(fā)展三維數(shù)字濾波器的概念和研究三維數(shù)字濾波器的設(shè)計(jì)問題仍舊是不可避免的。
本文主要討論線性相位三維FIR數(shù)字濾波器的分類和設(shè)計(jì)工作。與IIR濾波器相比,F(xiàn)IR系統(tǒng)的單位脈沖長度有限,滿足絕對可和系統(tǒng)穩(wěn)定性判斷條件。其次,F(xiàn)IR數(shù)字濾波器可以實(shí)現(xiàn)線性相位響應(yīng)。最初FIR三維數(shù)字濾波器的設(shè)計(jì)是以一維數(shù)字濾波器為原型,采用McClellan變換(McClellan transform)進(jìn)行轉(zhuǎn)換得到的[2-7]。遺憾的是,這種方法無法精確控制三維數(shù)字濾波器的頻率特性。優(yōu)化迭代算法是另一種比較常見的用來設(shè)計(jì)FIR 三維數(shù)字濾波器方法。濾波器的設(shè)計(jì)問題寫成線性規(guī)劃問題求解,且不用約束其系統(tǒng)穩(wěn)定性[8,9]。在文獻(xiàn)[10,11]中,解析最小二乘法用來求解FIR三維數(shù)字濾波器。為了降低設(shè)計(jì)工作的計(jì)算復(fù)雜度,F(xiàn)IR三維數(shù)字濾波器的設(shè)計(jì)工作充分利用三維空間的對稱特性[4,5,8-13]。即便如此,在上述文獻(xiàn)中所獲得的FIR三維數(shù)字濾波器頻譜特性仍舊不太令人滿意。此外還有一些工作嘗試采用比較新穎的算法思想來解決FIR多維數(shù)字濾波器的設(shè)計(jì)問題,如模擬退火算法[14],凸集投影法[15],但遺憾的是三維數(shù)字濾波器的實(shí)驗(yàn)數(shù)據(jù)并沒有被成功獲取。
研究三維數(shù)字濾波器的對稱特性可以從頻域(幅度響應(yīng),相位響應(yīng))和空域(或稱時(shí)域,即系統(tǒng)單位脈沖響應(yīng)或?yàn)V波器系數(shù))兩個(gè)方面入手,如:球面對稱和立方對稱等屬于頻域概念[16-18]。基于平面對稱概念所定義的多種不同幅度響應(yīng)對稱特性以列表的形式報(bào)告于[19]。在文獻(xiàn)[8]中,作者對濾波器的幅度平方響應(yīng)對稱性進(jìn)行研究,得出三維數(shù)字濾波器的幅度平方響應(yīng)共有47種不同類型對稱性。而在文獻(xiàn)[11]中,基于單位脈沖響應(yīng)三維數(shù)據(jù)的長度奇偶性和對稱奇偶性的不同組合,得出共有64種不同類型的線性相位FIR三維數(shù)字濾波器類型。類似地,在1998年,線性相位一維數(shù)字濾波器單位脈沖響應(yīng)對稱性的四種組合被移植至多維系統(tǒng)[20]。本文中,在文獻(xiàn)[11,20]工作的基礎(chǔ)上,擬對單位脈沖響應(yīng)三維數(shù)據(jù)對稱性的64種組合進(jìn)行分析,剔除冗余數(shù),初步提出20種不同類型的線性相位FIR三維數(shù)字濾波器。結(jié)合零點(diǎn)位置分析,給出設(shè)計(jì)不同種類(低通、高通、帶通和帶阻)線性相位FIR三維數(shù)字濾波器的單位脈沖響應(yīng)對稱特性的選取。為了得到線性相位三維FIR數(shù)字濾波器,本文把此類型濾波器的設(shè)計(jì)問題寫出凸問題,采用約束凸規(guī)劃進(jìn)行求解[21],獲得具有良好頻率特性的線性相位三維FIR數(shù)字濾波器。最后給出設(shè)計(jì)實(shí)例。
對于一個(gè)因果的FIR三維數(shù)字濾波器,其單位脈沖響應(yīng)可以用三維矩陣數(shù)據(jù)來表示:其系統(tǒng)函數(shù)表示為:
系統(tǒng)的頻率響應(yīng)定義為濾波器單位脈沖響應(yīng)的離散時(shí)間傅里葉變換(Discrete-time Fourier transform: DTFT),是式(1)中z變量在z平面單位圓上取值,即,
類似于一維數(shù)字濾波器單位脈沖響應(yīng)的對稱性分析, FIR三維數(shù)字濾波器在三維方向?yàn)V波器長度在每一維度方向上可能是奇數(shù)長度也可能是偶數(shù)長度,且在每一維度方向可能是奇對稱,也可能是偶對稱的。因此,F(xiàn)IR三維數(shù)字濾波器的單位脈沖響應(yīng)每一維度方向有四種可能組合:
表1:單位脈沖響應(yīng)三維方向?qū)ΨQ特性組合種類
(1)長度為奇數(shù)且偶對稱;
(2)長度為偶數(shù)且偶對稱;
(3)長度為奇數(shù)且奇對稱;
(4)長度為偶數(shù)且奇對稱。
也就是文獻(xiàn)[11]中所說的一共有64種組合。我們可以用數(shù)組(111),(112),(121),。。。,(444)來表示這64種組合。
利用線性相位濾波器單位脈沖響應(yīng)的數(shù)據(jù)對稱性,我們可以很方便地寫出表1中三維FIR線性相位濾波器的頻率響應(yīng)。例如,現(xiàn)有某FIR三維數(shù)字濾波器的單位脈沖響應(yīng)三維方向?qū)ΨQ特性為表1中組合“111”,則有:
其中,
表2:例1中三維線性相位FIR濾波器系數(shù)
發(fā)現(xiàn):
圖1:三維FIR數(shù)字濾波器幅頻特性
圖2
再例如:對于組合(444),利用濾波器系數(shù)對稱特性,濾波器的頻率響應(yīng)可以寫成:
其中,
對于一維線性相位FIR濾波器的四種類型,類型1可用于設(shè)計(jì)低通(LP),帶通(BP),高通(HP),和帶阻(BS)濾波器;類型2可用于設(shè)計(jì)LP和BP濾波器;類型3僅能用于設(shè)計(jì)BP濾波器;類型4可用于設(shè)計(jì)HP和BP濾波器。因此,表1中20種類型中適用于設(shè)計(jì)在三個(gè)維度里具有相同通阻特性的線性相位FIR濾波器為:低通(111,112,122,222)、帶通(全部組合)、高通(111,114,144,444)和帶阻(111)。如果濾波器的三個(gè)維度具有不同的通阻特性,則根據(jù)具體要求可選擇不同組合種類。
在接下來的小節(jié),我們把三維線性相位FIR濾波器組的設(shè)計(jì)問題寫出極小化極大值問題,帶入約束凸規(guī)劃求解。
基于CVX工具箱的約束凸規(guī)劃法(DCP),用來求解各類凸問題[21]。與半定規(guī)劃(SDP)法和二階錐規(guī)劃(SOCP)法相比,DCP對問題描述簡單直觀,求解復(fù)雜問題的能力更強(qiáng)。為了采用CVX求解三維線性相位FIR數(shù)字濾波器的設(shè)計(jì)問題,我們先把FIR三維數(shù)字濾波器的設(shè)計(jì)問題表達(dá)成極小化極大值問題:
采用DCP設(shè)計(jì)FIR一維或二維數(shù)字濾波器時(shí),問題描述相同。唯一不同的是在(7)式中待求未知變量數(shù)是三維變量。在濾波器每維度系數(shù)長度一樣時(shí)FIR三維數(shù)字濾波器的待求變量是一維數(shù)字濾波器的倍。式(7)可以用來求任意相位響應(yīng)的FIR濾波器。對于線性相位FIR濾波器,可以利用其濾波器系數(shù)的對稱性,來降低待求變量數(shù)。只是對于表1中的不同組合形式,式5要做不同的改寫。例如組合111,式(7)改寫如下:
設(shè)計(jì)三維線性相位FIR數(shù)字濾波器時(shí),必須預(yù)先了解以下幾點(diǎn):
(1)所期望的濾波器單位脈沖響應(yīng)長度;
(2)所期望的濾波器分段頻率特性(低通、高通等)及截止頻率點(diǎn);
(3)所期望的濾波器通帶/阻帶形狀(球形、立方形等);
(4)所期望的濾波器通帶或阻帶所能接受的最大波紋誤差。
一般來說,濾波器系數(shù)的長度、過渡帶寬度(通帶與阻帶之間的間隔)和通帶、阻帶上的波紋誤差之間是一種相互制衡的關(guān)系。濾波器系數(shù)越長,過渡帶越寬,則通帶和阻帶波紋誤差則越小。而通帶和阻帶波紋誤差之間也有相互制衡的關(guān)系,壓縮通帶的波紋誤差,也就是要獲得更好的通帶平坦性,則在其他參數(shù)不變的前提下,必然以增大阻帶波紋誤差作為代價(jià)。
例1:設(shè)計(jì)一個(gè)長度(5,5,5),球面對稱,線性相位的低通FIR三維數(shù)字濾波器。通帶截止頻率和阻帶截止頻率分別為(0.3π,0.6π),濾波器的阻帶最大波紋誤差設(shè)為-30dB。把濾波器設(shè)計(jì)問題寫成式(8)進(jìn)行求解。待求系數(shù)c(n1,n2,n3)大小為(3×3×3)個(gè)。圖1中顯示其幅度頻率響應(yīng)。
由于所得濾波器在三個(gè)維度方向上具有對稱性,從另外兩個(gè)維度看進(jìn)去的濾波器幅頻特性與圖1中是一致的。其濾波器系數(shù)列于表2中,以供讀者參考。仔細(xì)觀察表2中數(shù)據(jù),可發(fā)現(xiàn)每頁數(shù)據(jù)具有關(guān)于對角線對稱的特性;并且第一頁第二列數(shù)據(jù)與第二頁第一列數(shù)據(jù)相同,以此類推。表2中濾波器系數(shù)的對稱性緣于本例中所要求設(shè)計(jì)的濾波器在三個(gè)維度上通阻帶截止頻率點(diǎn)相同,濾波器長度一致,即濾波器在幅頻上具有對稱性。這種幅頻對稱特性也可以用來降低濾波器設(shè)計(jì)過程的未知變量數(shù)。但這種幅頻對稱特性在單位脈沖響應(yīng)上的影響較為復(fù)雜,尤其是濾波器比較長的時(shí)候。例2:設(shè)計(jì)一個(gè)長度(12,12,12),通帶區(qū)域方形對稱,線性相位的低通FIR三維數(shù)字濾波器。通帶截止頻率和阻帶截止頻率分別為(0.3π,0.6π),濾波器的阻帶最大波紋誤差設(shè)為-30dB。待求系數(shù)c(n1,n2,n3)大小為63個(gè)。圖2顯示的是所得三維FIR數(shù)字濾波器幅頻響應(yīng)圖,其在通帶的最大波紋誤差為2.048e-3。由于篇幅所限,濾波器系數(shù)沒有列出來。
在處理三維數(shù)字信號時(shí),考慮到數(shù)據(jù)之間的相關(guān)性,三維數(shù)字濾波器的設(shè)計(jì)不可避免。三維數(shù)字濾波器設(shè)計(jì)工作的難點(diǎn)在于龐大的待求未知變量數(shù)。本文從一維線性相位FIR濾波器系數(shù)的對稱性出發(fā),首次提出三維線性相位FIR濾波器共有20種不同類型。利用濾波器系數(shù)的對稱性,大幅度減少濾波器設(shè)計(jì)過程中待求未知變量淑,把三維線性相位FIR數(shù)字濾波器的設(shè)計(jì)問題寫成極大值極小化問題,帶入凸規(guī)劃求解。實(shí)驗(yàn)結(jié)果顯示,這種方法可以獲得具有良好頻率特性的三維線性相位FIR濾波器。