王 顏
濾波是數(shù)字信號(hào)處理中的常見(jiàn)過(guò)程。通過(guò)選頻濾波可以提取時(shí)域信號(hào)有用的頻率分量,壓制干擾信號(hào)。選頻濾波通常在頻域進(jìn)行,其具體步驟是:
(1)對(duì)長(zhǎng)度為N 的時(shí)域信號(hào)x(n)作離散傅里葉變換(DFT),得到對(duì)應(yīng)的頻域信號(hào)X(k)。通常使用快速傅里葉變換算法(FFT)進(jìn)行。
(2)將頻域信號(hào)X(k)各點(diǎn)處的值與頻域?yàn)V波器F(k)相同頻率處的值相乘,得到濾波后的頻域信號(hào)Y(k):
(3)對(duì) Y(k)做離散傅里葉反變換(IDFT),即得到頻域?yàn)V波后的時(shí)域信號(hào) y(n)。通常使用快速傅里葉反變換算法(IFFT)進(jìn)行。
對(duì)時(shí)域信號(hào)進(jìn)行頻域?yàn)V波后所得到的信號(hào)可能在時(shí)域的一端被截?cái)?,被截去的部分在另一端出現(xiàn),造成時(shí)域信號(hào)失真。本文稱這種時(shí)域信號(hào)截?cái)嗪鸵莆辉斐傻男盘?hào)失真為頻域?yàn)V波的時(shí)域端點(diǎn)效應(yīng)。下文通過(guò)編寫MATLAB 程序,應(yīng)用式(2)所表示的帶余弦鑲邊的帶通濾波器對(duì)地震勘探數(shù)據(jù)處理中廣泛應(yīng)用的雷克子波[1]進(jìn)行頻域?yàn)V波,說(shuō)明頻域?yàn)V波的流程,并指出這一過(guò)程產(chǎn)生的時(shí)域截?cái)嗪鸵莆滑F(xiàn)象,即頻域?yàn)V波的時(shí)域端點(diǎn)效應(yīng)。
式(2)中,l1和h1分別是低頻端和高頻端的截止頻率,l1-l2和h1-h2分別是低頻端和高頻端的過(guò)渡帶,在過(guò)渡帶內(nèi)濾波器的幅值以余弦函數(shù)的形式從1 衰減到0. 顯然,對(duì)于連續(xù)的頻率,該濾波器的頻譜是光滑的。
現(xiàn)將濾波器低頻端和高頻端的截止頻率分別設(shè)為20Hz,60Hz,過(guò)渡帶寬度設(shè)為 3Hz,即 l2=16,l1=20,h1=60,h2=64,則濾波器的頻譜如圖1 所示。
圖2 和圖 3 分別是主頻 40Hz 的雷克子波信號(hào)及其頻譜。
圖2 40Hz 的雷克子波信號(hào)
圖4 和圖5 分別是應(yīng)用圖1 所示濾波器進(jìn)行步驟(2)所述濾波后的信號(hào)的頻譜和經(jīng)過(guò)離散傅里葉反變換后所得的時(shí)域信號(hào)。
圖3 40Hz 對(duì)應(yīng)的頻譜
圖4 濾波后的信號(hào)的頻譜
圖5 經(jīng)過(guò)離散傅里葉反變換后所得的時(shí)域信號(hào)
從圖5 可以看出,濾波后的信號(hào)在t=0 處被截?cái)啵唤厝サ牟糠殖霈F(xiàn)在了時(shí)間軸的末端,即出現(xiàn)了前文所述的端點(diǎn)效應(yīng)。
根據(jù)離散傅里葉變換和反變換的性質(zhì)[2],對(duì)頻域信號(hào)進(jìn)行點(diǎn)數(shù)為N 的等間隔采樣,其離散傅里葉反變換對(duì)應(yīng)的時(shí)域信號(hào)為
由(3)式可知,如果原來(lái)的時(shí)域信號(hào)x(n)經(jīng)過(guò)頻域?yàn)V波后對(duì)應(yīng)的延伸至a ≤n <0 范圍內(nèi),則這段信號(hào)會(huì)被截除,同時(shí)因?yàn)楸唤爻男盘?hào)會(huì)出現(xiàn)在時(shí)間軸末端的a+N ≤n ≤N-1 位置。類似的,如果延伸至N ≤n ≤b 范圍內(nèi),這段信號(hào)同樣會(huì)被截除,并出現(xiàn)在0≤n ≤b-N 位置。
本文現(xiàn)從時(shí)域和頻域運(yùn)算間的關(guān)系角度分析為避免時(shí)域端點(diǎn)效應(yīng),頻域?yàn)V波過(guò)程應(yīng)滿足的條件。
由傅里葉變換的性質(zhì)可知,頻域的乘積運(yùn)算等價(jià)于時(shí)域的卷積運(yùn)算,所以。
由上式可知經(jīng)過(guò)頻域?yàn)V波后的信號(hào)相當(dāng)于原信號(hào)在時(shí)域與濾波器的時(shí)域響應(yīng)做卷積運(yùn)算的結(jié)果。假設(shè)x(n)和 f(n)在時(shí)間軸上的第一個(gè)和最后一個(gè)非零值出現(xiàn)的位置分別為 N1、N2和 N3、N4,濾波后的信號(hào)為 xf(n),則:
若xf(x)≠0,則x(m)≠0 且f(n-m)≠0,所以m ≤N2,nm ≤ N4,兩式相加得。
由(7)式可知,原信號(hào)的采樣總時(shí)長(zhǎng)N 滿足(8)式時(shí),可以確保濾波后的信號(hào)在時(shí)間軸末端不被截?cái)唷?/p>
(8)式可以變形為
因?yàn)閷?duì)于傅里葉反變換,時(shí)間軸的兩個(gè)方向是等價(jià)的,所以由(9)可以得到確保濾波后的信號(hào)在時(shí)間軸起始端不被截?cái)嗟臈l件。
(9)式和(11)式指出了原信號(hào)和濾波器的參數(shù)滿足什么條件時(shí),可以保證經(jīng)過(guò)頻域?yàn)V波的信號(hào)在時(shí)域不出現(xiàn)端點(diǎn)效應(yīng)。從推導(dǎo)過(guò)程可以看出,這兩個(gè)不等式表達(dá)的是避免時(shí)域端點(diǎn)效應(yīng)的充分條件,而非必要條件。所以若頻域?yàn)V波過(guò)程中出現(xiàn)了時(shí)域端點(diǎn)效應(yīng),可以按滿足(8)(11)兩式的方向調(diào)整相關(guān)參數(shù),但相關(guān)參數(shù)并非必須滿足兩式才能避免端點(diǎn)效應(yīng)。
本節(jié)介紹了避免頻域?yàn)V波的時(shí)域端點(diǎn)效應(yīng)的兩種具體的方法,并以濾波后的信號(hào)在時(shí)間軸起始端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)為例,通過(guò)MATLAB 程序進(jìn)行了實(shí)驗(yàn)驗(yàn)證。
通過(guò)在原信號(hào)的末端增加M 個(gè)零值點(diǎn),可使式(8)中的N增大M 而N2、N4不變,各參數(shù)整體向滿足該不等式變化,因此當(dāng)M 大于某個(gè)值時(shí),信號(hào)在時(shí)間軸末端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)會(huì)得到避免。
類似的,通過(guò)在原信號(hào)的起始端增加M 個(gè)零值點(diǎn),可使式(11)中的N、N1各增大M。因?yàn)轭l域采樣點(diǎn)數(shù)增加了M,所以濾波器的離散傅里葉反變換時(shí)長(zhǎng)增加了M,由離散傅里葉反變換結(jié)果的對(duì)稱性可知,N3增大了M/2。因此(11)式左端增大了M,右端增大了3M/2,各參數(shù)整體向滿足該不等式變化,因此當(dāng)M大于某個(gè)值時(shí),信號(hào)在時(shí)間軸起始端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)會(huì)得到避免。
在圖2 所示信號(hào)起始端增加100 個(gè)零值點(diǎn),得到圖6 所示延時(shí)0.25s 后的信號(hào)。
圖6 延時(shí)0.25s 后的信號(hào)
根據(jù)離散傅里葉變換的性質(zhì),信號(hào)的頻譜(振幅譜)不發(fā)生變化,且仍為均勻采樣,但采樣數(shù)比原來(lái)增加100,采樣間隔變小。因此要通過(guò)頻域?yàn)V波獲取與第一節(jié)中相同的頻率成分,所使用的濾波器應(yīng)滿足連續(xù)頻譜與第一節(jié)實(shí)驗(yàn)中的濾波器相同,但采樣位置與圖6 信號(hào)頻譜的采樣位置相同。該濾波器的頻譜如圖 7 所示。
圖7 與圖6 信號(hào)頻譜的采樣位置相同下的濾波器頻譜
圖8 紅線和藍(lán)線分別是第一節(jié)的實(shí)驗(yàn)和本實(shí)驗(yàn)所用濾波器的時(shí)域響應(yīng),可以看出藍(lán)線第一個(gè)非零值的位置比紅線增大了50,與理論分析相符。
圖8 不同實(shí)驗(yàn)所用濾波器的時(shí)域響應(yīng)
圖9 是經(jīng)過(guò)頻域?yàn)V波得到的時(shí)域信號(hào),可見(jiàn)圖5 中的端點(diǎn)效應(yīng)消失,證明了時(shí)域補(bǔ)零方法的有效性。
圖9 經(jīng)過(guò)頻域?yàn)V波得到的時(shí)域信號(hào)
通過(guò)拓寬濾波器頻譜,使其時(shí)域響應(yīng)變窄,(8)式中的N4變小,(11)式中的N3 變大,兩式中的參數(shù)都向使不等式成立的方向變化。因此信號(hào)在時(shí)間軸起始端和末端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)都可以通過(guò)拓寬濾波器頻譜來(lái)解決。
將第一節(jié)實(shí)驗(yàn)所用濾波器的過(guò)渡帶寬度增加到9Hz,其頻譜如圖10 所示。
圖1 0 濾波器過(guò)渡帶寬度增加到9Hz 時(shí)的頻譜
圖11 紅線和藍(lán)線分別是第一節(jié)的實(shí)驗(yàn)和本實(shí)驗(yàn)所用濾波器的時(shí)域響應(yīng),可以看出本實(shí)驗(yàn)濾波器的時(shí)域響應(yīng)變窄,第一個(gè)零值點(diǎn)位置后移,最后一個(gè)零值點(diǎn)位置前移。
圖12 是用過(guò)渡帶寬度增加的濾波器對(duì)圖2 所示信號(hào)進(jìn)行頻域?yàn)V波的結(jié)果,可見(jiàn)圖5 中的端點(diǎn)效應(yīng)消失,證明了拓寬濾波器頻譜方法的有效性。同時(shí)因?yàn)闉V波器的過(guò)渡帶變得更加平緩,吉布斯效應(yīng)也得以減弱。
圖11 不同實(shí)驗(yàn)所用濾波器的時(shí)域響應(yīng)
圖12 經(jīng)過(guò)頻域?yàn)V波得到的時(shí)域信號(hào)
對(duì)數(shù)字信號(hào)進(jìn)行頻域?yàn)V波所得到的時(shí)域信號(hào)可能在一端被截?cái)?,被截去的部分出現(xiàn)在另一端,造成時(shí)域信號(hào)失真。本文稱這一現(xiàn)象為頻域?yàn)V波的時(shí)域端點(diǎn)效應(yīng)。通過(guò)理論分析和數(shù)值實(shí)驗(yàn),得出以下結(jié)論:
(1)頻域?yàn)V波的時(shí)域端點(diǎn)效應(yīng)產(chǎn)生的原因是原信號(hào)經(jīng)過(guò)頻域?yàn)V波后時(shí)寬變長(zhǎng),因采樣總時(shí)長(zhǎng)有限而在一端或兩端被截?cái)唷M瑫r(shí)因?yàn)轭l域信號(hào)對(duì)應(yīng)的時(shí)域信號(hào)的周期性,在一端被截除的信號(hào)在另一端重現(xiàn)。
(2)假設(shè)原信號(hào)采樣總時(shí)長(zhǎng)為N,原信號(hào)和濾波器的時(shí)域響應(yīng)在時(shí)間軸上的第一個(gè)和最后一個(gè)非零值出現(xiàn)的位置分別為 N1、N2和 N3、N4,則上述參數(shù)滿足 N-1≤N1+N3可以確保不出現(xiàn)信號(hào)在起始端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng);滿足N-1≥N2+N4可以確保不出現(xiàn)信號(hào)在末端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)。
(3)若頻域?yàn)V波過(guò)程中出現(xiàn)了時(shí)域端點(diǎn)效應(yīng),可以按滿足結(jié)論(2)中的兩個(gè)不等式的方向調(diào)整相關(guān)參數(shù)。但因?yàn)檫@兩個(gè)不等式表達(dá)的是避免時(shí)域端點(diǎn)效應(yīng)的充分條件而非必要條件,所以調(diào)整后參數(shù)并非必須滿足兩式。
(4)通過(guò)在原信號(hào)的起始端補(bǔ)零,可以避免濾波后的信號(hào)在起始端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng);在原信號(hào)的末端補(bǔ)零,可以避免濾波后的信號(hào)在末端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)。
(5)通過(guò)拓寬濾波器頻譜,信號(hào)在時(shí)間軸起始端和末端被截?cái)嘣斐傻亩它c(diǎn)效應(yīng)都可以得到解決。