一種用于超聲成像的零相位濾波后處理方法?

韓曉麗 吳文燾 曹 政

(中國(guó)科學(xué)院聲學(xué)研究所 北京 100190)

0 引言

超聲相控陣系統(tǒng)中，聲波信號(hào)在經(jīng)過(guò)隔離、限幅、放大、濾波等模擬調(diào)理電路時(shí)必然會(huì)引入一定的噪聲干擾[1]，因此，由A/D 采樣后的數(shù)字信號(hào)在成像之前，必須要經(jīng)過(guò)一定的濾波處理，以提高信噪比和成像質(zhì)量[2?3]，也就是說(shuō)，濾波是一個(gè)必備的處理過(guò)程。在對(duì)接收的回波信號(hào)進(jìn)行分析的過(guò)程中，通常需要將信號(hào)不同成分進(jìn)行分離，然后對(duì)不同成分做不同的處理或調(diào)整，從而最終提高對(duì)于缺陷或目標(biāo)的檢測(cè)能力[4]。大部分具有頻率依賴性的檢測(cè)識(shí)別算法的實(shí)現(xiàn)都需要輔以濾波處理的過(guò)程，包括超聲檢測(cè)系統(tǒng)廣泛采用裂譜處理(Split spectrum processing, SSP)方法，從本質(zhì)上講，也是采用若干帶通有限脈沖響應(yīng)(Finite impulse response, FIR)濾波器進(jìn)行濾波的過(guò)程[4?5]。

在超聲無(wú)損檢測(cè)領(lǐng)域，超聲波遇到缺陷時(shí)會(huì)發(fā)生散射，在接收端會(huì)接收到經(jīng)缺陷散射的回波信號(hào)，與缺陷有關(guān)的信號(hào)出現(xiàn)的時(shí)間即到達(dá)時(shí)間(Time of arrival, TOA)或傳播時(shí)間(Time of flight, TOF)[4]恰恰能夠反映缺陷在檢測(cè)對(duì)象中所處的位置，以此來(lái)對(duì)缺陷進(jìn)行定位。同時(shí)，為了保證系統(tǒng)檢測(cè)的靈敏度和能夠通過(guò)缺陷回波幅度在一定程度上估計(jì)缺陷的大小，要求缺陷回波波幅在濾波前后保持不變。基于上述兩點(diǎn)，在對(duì)接收的回波信號(hào)進(jìn)行濾波處理時(shí)，要保證濾波過(guò)程不影響主要頻率成分的到達(dá)時(shí)間，這就需要濾波器具有零相位的相位響應(yīng)特性。雖然FIR濾波器在單位脈沖響應(yīng)序列h(n)滿足一定的對(duì)稱性時(shí)，具有線性相位特性[6?7]，但是濾波前后必然會(huì)出現(xiàn)一定的時(shí)延，難以滿足零相位的要求。而無(wú)限脈沖響應(yīng)(Infinite impulse response,IIR)濾波器不能實(shí)現(xiàn)嚴(yán)格意義上的線性相位特性。

除了直接設(shè)計(jì)具有零相位響應(yīng)的數(shù)字濾波器以外，零相位濾波的過(guò)程還可以通過(guò)一種比較巧妙的方式實(shí)現(xiàn)，即采用前向-后向?yàn)V波-反向輸出(Forward-reverse filtering, reverse output, FRR)[8]或反向?yàn)V波-反向輸出(Reverse filtering,reverse filtering,RRF)的方法[9?10]。這種方法借助于兩次濾波和信號(hào)翻轉(zhuǎn)的組合，實(shí)現(xiàn)對(duì)輸入信號(hào)的零相位濾波。鑒于FIR濾波器具有便于實(shí)現(xiàn)線性相位特性的優(yōu)越性、系統(tǒng)穩(wěn)定性好和使用上的靈活性[11]，通常大家習(xí)慣性地選用FIR 濾波器實(shí)現(xiàn)濾波過(guò)程，包括零相位濾波，也主要以FIR 的方式實(shí)現(xiàn)[10]。然而，由于FIR 濾波器傳遞函數(shù)的極點(diǎn)固定在原點(diǎn)，只能靠改變零點(diǎn)位置來(lái)改變它的性能。所以要達(dá)到高的選擇性，必須用較高的階數(shù)；對(duì)于同樣的濾波器設(shè)計(jì)指標(biāo)，F(xiàn)IR 濾波器所要求的階數(shù)可能比IIR 濾波器高5～10倍，信號(hào)延時(shí)也較大，處理效率低[11]。

近年來(lái)，在超聲檢測(cè)領(lǐng)域，實(shí)時(shí)缺陷檢測(cè)的需求越來(lái)越高，要求系統(tǒng)能夠盡快地實(shí)時(shí)采集、分析和處理超聲波回波信號(hào)，這對(duì)于超聲成像系統(tǒng)的實(shí)時(shí)性提出了更高的要求[5]。許多學(xué)者致力于探索如何高效地實(shí)現(xiàn)超聲回波信號(hào)采集之后的濾波過(guò)程[12]，也包括開(kāi)辟濾波器的高效硬件實(shí)現(xiàn)的新思路[5,10]。快速高效地實(shí)現(xiàn)零相位濾波，將有利于超聲檢測(cè)系統(tǒng)信號(hào)處理效率的提升，確保檢測(cè)系統(tǒng)具有更高的實(shí)時(shí)性能。

本文針對(duì)在超聲相控陣成像領(lǐng)域，在采用基于FIR濾波器的FRR零相位濾波時(shí)，由于FIR濾波器階數(shù)高導(dǎo)致運(yùn)算量大、處理效率低、難以滿足實(shí)時(shí)性要求高的系統(tǒng)需求的問(wèn)題，探討在實(shí)際超聲相控陣系統(tǒng)中使用IIR 濾波器實(shí)現(xiàn)FRR 零相位濾波的可行性，比較基于FIR和IIR 濾波器實(shí)現(xiàn)FRR零相位濾波的性能，探索IIR 濾波器實(shí)現(xiàn)FRR 零相位濾波時(shí)的效率提升程度，并且把IIR 零相位濾波應(yīng)用于實(shí)時(shí)超聲成像后處理系統(tǒng)中，通過(guò)實(shí)際的成像實(shí)驗(yàn)，驗(yàn)證使用IIR 濾波器實(shí)現(xiàn)FRR 零相位濾波的效果。

1 基于IIR的零相位濾波

根據(jù)線性系統(tǒng)理論，對(duì)信號(hào)進(jìn)行經(jīng)典線性濾波的過(guò)程對(duì)原始輸入信號(hào)主要產(chǎn)生兩大方面的影響[9]：一是對(duì)信號(hào)不同頻率成分的幅值進(jìn)行調(diào)整，通帶內(nèi)的頻率成分一般不衰減甚至得到放大，而阻帶內(nèi)的成分通常被衰減，以減少這些成分在原始信號(hào)中所占的比重或影響，這是通過(guò)濾波器的幅頻響應(yīng)來(lái)實(shí)現(xiàn)的；二是對(duì)信號(hào)的不同頻率成分產(chǎn)生相應(yīng)的相移，這是通過(guò)濾波器的相頻響應(yīng)Arg[H(ejw)]來(lái)實(shí)現(xiàn)的。濾波器的線性相位指的是不同頻率成分的時(shí)間延時(shí)是相同的，則總體上輸出信號(hào)相對(duì)于輸入信號(hào)不發(fā)生相位失真，僅發(fā)生一定的延遲。在超聲相控陣成像系統(tǒng)中借助于零相位濾波正是巧妙地實(shí)現(xiàn)了對(duì)超聲回波信號(hào)的濾波處理，同時(shí)保證群延時(shí)為零，即處理后的回波信號(hào)仍具有正確的到達(dá)時(shí)間。

借助于兩次濾波和信號(hào)翻轉(zhuǎn)的組合以實(shí)現(xiàn)零相位濾波的具體方法有兩種：FRR 和RRF。FRR零相位濾波過(guò)程如圖1所示。設(shè)濾波器單位脈沖響應(yīng)序列h(n)，對(duì)應(yīng)的Z變換(Ztransformation)為H(Z)。FRR先將輸入序列x(n)按順序輸入到濾波器H(Z)中進(jìn)行濾波，即前向?yàn)V波；然后將所得的結(jié)果序列y1(n)，翻轉(zhuǎn)后得到序列y1(?n)，y1(?n)再次輸入到濾波器H(Z)中，即反向?yàn)V波；再將所得結(jié)果y2(n)逆轉(zhuǎn)后輸出，即反向輸出；最后得到精確的零相位失真的輸出序列y(n)。RRF 先將輸入信號(hào)序列反轉(zhuǎn)后通過(guò)濾波器，然后將所得結(jié)果逆轉(zhuǎn)后再次通過(guò)濾波器，這樣所得結(jié)果即為精確零相位失真的輸出序列[9]。不管是哪種具體方法都需要預(yù)先設(shè)計(jì)具有特定幅頻特性的濾波器H(ejw)。

圖1 FRR 零相位濾波過(guò)程Fig.1 Zero phase filtering based on FRR

FRR簡(jiǎn)單的推導(dǎo)過(guò)程[5]如下：

聯(lián)立以上4 式，可以推得整個(gè)零相位過(guò)程的等效傳輸函數(shù)Heq(Z)：

即等效系統(tǒng)頻率響應(yīng)為

也就是說(shuō)，F(xiàn)RR 濾波過(guò)程真正實(shí)現(xiàn)了對(duì)輸入信號(hào)x(n)的群延時(shí)為零的濾波處理，僅對(duì)其各頻率成分進(jìn)行了幅度調(diào)整。RRF 的推導(dǎo)過(guò)程與FRR 類似，這里不再贅述。因此，設(shè)計(jì)零相位濾波時(shí)，主要是設(shè)計(jì)滿足特定需求的具有特定幅頻響應(yīng)特性的基礎(chǔ)濾波器H(Z)。在上述推導(dǎo)過(guò)程中，并未指定濾波器H(Z)為有限脈沖響應(yīng)濾波器(FIR)還是無(wú)限脈沖響應(yīng)濾波器(IIR)。本文在后續(xù)章節(jié)中，將基于FIR 的FRR零相位濾波稱為FIR FRR，將基于IIR的FRR零相位濾波稱為IIR FRR。

2 IIR FRR零相位濾波的設(shè)計(jì)與分析

FRR 或RRF 的過(guò)程對(duì)輸入信號(hào)的相位最終不產(chǎn)生影響，因此我們?cè)谠O(shè)計(jì)基礎(chǔ)濾波器H(Z)時(shí)主要關(guān)注其幅頻特性，使其幅頻特性滿足我們的具體應(yīng)用需求即可。下面以帶通FRR 零相位濾波為例，設(shè)計(jì)IIR FRR和FIR FRR的基礎(chǔ)濾波器H(Z)。FIR FRR主要用于與IIR FRR進(jìn)行對(duì)比。

2.1 基礎(chǔ)濾波器的設(shè)計(jì)

分別設(shè)計(jì)幅頻特性參數(shù)如表1所示的FIR 和IIR 濾波器。前者采用等波紋(Equiripple)法，后者采用契比雪夫II (Chebyshev II)法。

設(shè)計(jì)的FIR 和IIR 濾波器的幅頻響應(yīng)曲線如圖2所示。前者為57階，即N=57，后者為12 階，即N=12，其中Nb=12。在表1所示的幅頻響應(yīng)特性參數(shù)的約束下，兩種情況設(shè)計(jì)的幅頻響應(yīng)曲線在通帶內(nèi)的形狀差不多，主要差異在于低頻范圍內(nèi)的通帶到阻帶的過(guò)渡帶上，IIR 濾波器更窄。其中IIR 濾波器的零極點(diǎn)分布如圖3所示。距離單位圓最近的極點(diǎn)是0.9484±0.2190 i。

表1 帶通濾波器幅頻特性描述參數(shù)Table1 Bandpass filter amplitude-frequency characteristic parameters

圖2 IIR 和FIR 濾波器幅頻特性曲線Fig.2 Amplitude-frequency curves of IIR and FIR filters

圖3 IIR 濾波器零極點(diǎn)圖Fig.3 Pole-zero diagram of IIR filter

2.2 FRR運(yùn)算量分析

設(shè)x(n)為長(zhǎng)度為Nx的因果序列， 即n ∈[0,Nx ?1]。對(duì)于一個(gè)N階的直接型FIR 濾波器，根據(jù)公式(7)可以推知，完成N階的FIR 濾波需要的運(yùn)算量為(N+ 1)(Nx ?N+1)次乘法，N(Nx ?N+1)次加法。

對(duì)于一個(gè)N階的直接I 型IIR 濾波器，根據(jù)遞推公式(8)可以推知，完成N階的IIR 濾波需要的運(yùn)算量[13]為(N+Nb+1)(Nx ?N+1)次乘法，(N+Nb)(Nx ?N+1)次加法。其中，Nb為IIR濾波器傳輸函數(shù)H(Z) 的分子多項(xiàng)式的階數(shù)，需要注意的是，這里假定兩種情況下輸出序列的長(zhǎng)度都僅取決于已知的輸入序列x(n)和濾波器階數(shù)，即輸出序列y(n)的序列號(hào)范圍均為n ∈[N,Nx]。

因此，上述的57階FIR和12階IIR濾波器直接I 型實(shí)現(xiàn)時(shí)的運(yùn)算量如表2所示，其中，輸出序列長(zhǎng)度Ny=16384，輸入序列長(zhǎng)度Nx=16384+N ?1。由此看出，理論上，前者的運(yùn)算量約是后者的2倍多。

表2 FIR 和IIR 運(yùn)算量理論值Table2 Amount of calculation comparison between FIR and IIR

圖4 FIR FRR 和IIR FRR 零相位濾波結(jié)果比較Fig.4 Comparison of zero-phase filtering results of FIR FRR and IIR FRR

在Matlab 中，利用設(shè)計(jì)的FIR 和IIR 濾波器對(duì)射頻回波信號(hào)進(jìn)行FRR 零相位濾波處理，對(duì)原始回波信號(hào)進(jìn)行FIR FRR 和IIR FRR 的濾波結(jié)果如圖4(a)和圖4(b)所示。由圖4可知，IIR FRR同F(xiàn)IR FRR 一樣，也可以實(shí)現(xiàn)零相位濾波的功能。濾波前后，F(xiàn)IR FRR 和IIR FRR 的細(xì)節(jié)比較如圖5所示，兩種情況輸出信號(hào)的時(shí)間位置相同，在幅度上只有較小的差異，造成差異的原因是兩種濾波器的幅頻特性曲線并不是完全相同的。因此，在幾乎一致的幅頻響應(yīng)特性情況下，F(xiàn)IR 濾波器的階數(shù)約為IIR濾波器階數(shù)的5 倍，理論上的運(yùn)算量約是IIR 濾波器的2倍多。

圖5 FIR FRR 和IIR FRR 濾波結(jié)果的細(xì)節(jié)比較Fig.5 Detailed comparison of FIR FRR and IIR FRR results

3 IIR FRR用于超聲成像后處理

為了比較IIR FRR 和FIR FRR 的超聲成像效果，設(shè)計(jì)兩組成像實(shí)驗(yàn)，分別利用FIR FRR 和IIR FRR對(duì)原始信號(hào)進(jìn)行零相位濾波，然后經(jīng)過(guò)相同的其他后處理后獲得超聲成像結(jié)果，并對(duì)結(jié)果進(jìn)行比較。從原始數(shù)據(jù)到成像的整個(gè)后處理過(guò)程如圖6所示。兩組成像實(shí)驗(yàn)分別為：(1)單孔成像，(2)多孔成像。采用由課題組自主研發(fā)的UA special I 超聲相控陣系統(tǒng)采集一幀原始回波數(shù)據(jù)，采集時(shí)使用的試塊和換能器的放置形式如圖7和圖8所示。圖7中，試塊在深度25 mm 處有一個(gè)直徑為1 mm 的通孔，相控陣換能器通過(guò)耦合劑與試塊表面直接耦合，聚焦深度為25 mm。圖8中的試塊為標(biāo)準(zhǔn)的B 型便攜式相控陣試塊，由濟(jì)寧模具廠生產(chǎn)，檢測(cè)時(shí)換能器中心位于18??2 排孔最左側(cè)兩個(gè)孔的中間連線上，以18??2排孔為檢測(cè)時(shí)的主要目標(biāo)，聚焦深度設(shè)為50 mm。數(shù)據(jù)采樣率為50 MSPS，超聲相控陣換能器中心頻率為5 MHz。

圖6 超聲成像后處理過(guò)程Fig.6 Post-processing of ultrasonic imaging involving zero-phase filtering

圖7 數(shù)據(jù)采集時(shí)試塊和換能器位置示意圖——單孔Fig.7 Test block and transducer position during data acquisition (single hole)

圖8 數(shù)據(jù)采集時(shí)試塊和換能器位置示意圖——多孔Fig.8 Test block and transducer position during data acquisition (a row of holes)

3.1 基于IIR FRR 的超聲成像

利用圖7情況下采集的回波數(shù)據(jù)集，經(jīng)過(guò)圖6所示的處理過(guò)程，其中的零相位濾波的過(guò)程分別利用FIR FRR 和IIR FRR 完成，最終形成如圖9所示的單孔成像結(jié)果。由圖9可知，經(jīng)FIR FRR 和IIR FRR 處理后的數(shù)據(jù)的成像效果幾乎完全相同，視覺(jué)上難以分辨差異。

圖9 FIR FRR 和IIR FRR 成像結(jié)果——單孔Fig.9 Imaging based on FIR FRR and IIR FRR (single hole)

圖10 FIR FRR 和IIR FRR 成像結(jié)果——多孔Fig.10 Imaging based on FIR FRR and IIR FRR (a row of holes)

另外，利用圖8情況下采集的回波數(shù)據(jù)集，經(jīng)過(guò)圖6所示的處理過(guò)程，其中的零相位濾波分別利用FIR FRR 和IIR FRR 完成，最終對(duì)18??2 排孔形成如圖10所示的成像結(jié)果。由圖10可知，在復(fù)雜的多孔成像情況下，經(jīng)FIR FRR 和IIR FRR 處理后的數(shù)據(jù)的成像效果也幾乎完全相同，僅在18??2后方的偽影處存在極其細(xì)微的差別，這是由于兩種濾波器的幅頻特性曲線并不是完全相同導(dǎo)致的。

由成像結(jié)果可以看出，IIR FRR 可以實(shí)現(xiàn)與FIR FRR幾乎相同的零相位濾波性能，得到幾乎一致的超聲成像效果。

3.2 算法實(shí)時(shí)性的驗(yàn)證

系統(tǒng)后端采用8 核處理器TMS320C6678 進(jìn)行處理，利用定時(shí)器統(tǒng)計(jì)連續(xù)執(zhí)行1000 次FIR FRR和IIR FRR 零相位濾波的時(shí)間，然后計(jì)算平均一次FRR 零相位濾波的時(shí)間。57 階FIR 和12 階IIR 采用直接I 型網(wǎng)絡(luò)結(jié)構(gòu)實(shí)現(xiàn)，數(shù)據(jù)格式采用浮點(diǎn)格式，其中，輸出序列y(n)的長(zhǎng)度Ny= 16384，輸入序列x(n)長(zhǎng)度Nx=16384+N ?1。統(tǒng)計(jì)結(jié)果如表3所示。需要注意的是，在統(tǒng)計(jì)時(shí)間時(shí)，保持兩種濾波過(guò)程采用相同的程序結(jié)構(gòu)來(lái)實(shí)現(xiàn)，保證調(diào)用相同的函數(shù)，以及在調(diào)用函數(shù)時(shí)保證采用相同的調(diào)用格式，以保證在兩種情況下，濾波算法的結(jié)構(gòu)、實(shí)現(xiàn)的具體方式是相同的，從而盡可能減少其他因素對(duì)處理時(shí)間的影響。

從表3中可以看出：基于12 階IIR 的FRR (直接I) 比FIR FRR (直接I) 耗時(shí)少，后者是前者的1.71 倍。

表3 FIR FRR 和IIR FRR 耗費(fèi)的時(shí)間對(duì)比Table3 Time consumed during FIR FRR and IIR FRR

4 討論

(1)處理時(shí)間

從理論角度分析，實(shí)現(xiàn)幅頻響應(yīng)類似的濾波，IIR 確實(shí)相比FIR運(yùn)算量低很多，如表2所示。保守來(lái)講，IIR甚至僅為FIR的1/2。

由表3可知，采用直接型I 結(jié)構(gòu)實(shí)現(xiàn)高階IIR 濾波器時(shí)，相對(duì)于采用直接型結(jié)構(gòu)的FIR，IIR FRR的處理效率優(yōu)于FIR FRR。對(duì)比表2和表3，相比基于57 階FIR 的FRR，基于12 階的IIR(直接I 型)的FRR 處理時(shí)間確實(shí)更短，但并沒(méi)有達(dá)到約2 倍多的處理效率的提升。這是因?yàn)镕RR 的過(guò)程包括兩次濾波和兩次序列的翻轉(zhuǎn)，濾波耗時(shí)只是FRR處理總時(shí)間的一部分，在硬件處理器上執(zhí)行FRR算法的時(shí)間，不僅包括完成兩次FIR 或IIR濾波的時(shí)間，還包括序列翻轉(zhuǎn)的時(shí)間。這些時(shí)間導(dǎo)致最終統(tǒng)計(jì)的時(shí)間效率比值不可能達(dá)到實(shí)現(xiàn)IIR 和FIR濾波器本身的時(shí)間效率比。

(2)穩(wěn)定性

由于IIR 濾波器存在極點(diǎn)，如圖3所示。當(dāng)極點(diǎn)取值在單位圓內(nèi)時(shí)，理論上IIR 濾波器是穩(wěn)定的。但是在硬件處理器上實(shí)現(xiàn)時(shí)，濾波參數(shù)的量化或截?cái)嘈?yīng)可能會(huì)引起實(shí)際濾波器極點(diǎn)位置的偏移，從而會(huì)引起濾波器出現(xiàn)不穩(wěn)定的情況。在定點(diǎn)處理器中，這種問(wèn)題尤為突出。但是FIR 濾波器由于只存在零點(diǎn)，不會(huì)出現(xiàn)不穩(wěn)定的問(wèn)題。從這個(gè)角度上說(shuō)，F(xiàn)IR 濾波器實(shí)現(xiàn)起來(lái)更容易，更穩(wěn)定，對(duì)參數(shù)數(shù)據(jù)占位寬度或精度的要求更低。因此基于IIR 的FRR更適用于浮點(diǎn)處理器。

5 結(jié)論

基于IIR FRR 零相位濾波可以實(shí)現(xiàn)與FIR FRR 相當(dāng)?shù)奶幚硇Ч?。?dāng)IIR 濾波器采用直接I型結(jié)構(gòu)直接實(shí)現(xiàn)時(shí)，相對(duì)于FIR FRR，IIR FRR 可以在保證成像效果的前提下，將處理效率提升近1倍。在實(shí)際實(shí)時(shí)超聲成像系統(tǒng)中，為了獲得更高的處理效率，實(shí)現(xiàn)更高的幀率，可以選擇采用IIR FRR零相位濾波取代FIR FRR。尤其在浮點(diǎn)型硬件處理器中，采用直接型IIR 實(shí)現(xiàn)FRR 零相位濾波，可以一定程度上提高處理效率。在實(shí)際應(yīng)用中，需要仔細(xì)權(quán)衡FIR FRR 和IIR FRR 零相位濾波器在性能、時(shí)效、精度、穩(wěn)定性、可控性、占用內(nèi)存等多方面的要求，合理選擇采用何種方式實(shí)現(xiàn)零相位濾波處理。