基于Field II的非衍射陣列波聲場(chǎng)的仿真研究*

劉冬雪，彭虎，韓志會(huì)

(合肥工業(yè)大學(xué)，合肥 230000)

1 引 言

常用的非衍射聲場(chǎng)有Bessel聲場(chǎng)，X_wave聲場(chǎng)，陣列波聲場(chǎng)等[1]。Bessel聲場(chǎng)是連續(xù)波聲場(chǎng)，X_wave聲場(chǎng)激勵(lì)信號(hào)是復(fù)數(shù)且能量集中在z軸附近很難實(shí)現(xiàn)波束偏轉(zhuǎn)，難以應(yīng)用到實(shí)際中。盧建宇教授提出了陣列波，這種波束可以用于高幀率超聲成像。其中分層陣列波的聲場(chǎng)分布與縱坐標(biāo)無(wú)關(guān)，因?yàn)椴ㄊ螤钕穹謱拥臇艡诙妹鸞2]。聲場(chǎng)在傳播過(guò)程中無(wú)法抑制噪聲的干擾導(dǎo)致成像效果不夠理想，因此，研究如何使聲場(chǎng)保持非衍射特性傳播更遠(yuǎn)為非衍射波應(yīng)用于實(shí)際提供了理論支持。

波動(dòng)場(chǎng)常用積分衍射公式法進(jìn)行仿真，但是這種方法計(jì)算量大，精度低[3]。為了提高效率，本研究使用超聲領(lǐng)域中常用的仿真軟件Field II進(jìn)行仿真研究。為了使聲場(chǎng)傳播更遠(yuǎn)，一般采用幅度變跡技術(shù)進(jìn)行加權(quán)，通常使用三角函數(shù)作為加權(quán)函數(shù)。常見(jiàn)的有Hanning窗、Hamming窗、布萊克曼窗、余弦窗等，加權(quán)值是從陣元中間向陣元兩側(cè)緩慢變小的。通過(guò)幅度變跡技術(shù)可以增加成像的分辨率，這種方法雖然抑制了旁瓣，但是卻加寬了主瓣，空間分辨率也下降了。究其原因是幅度變跡都在發(fā)射聲場(chǎng)以前進(jìn)行設(shè)置，這種設(shè)置是固定的，與聲場(chǎng)本身是沒(méi)有關(guān)系的[4-6]。因此，人們提出了自適應(yīng)波束形成的方式進(jìn)行加權(quán)，該方法是根據(jù)接收聲場(chǎng)信號(hào)計(jì)算最優(yōu)權(quán)并對(duì)接收信號(hào)進(jìn)行加權(quán)，從而抑制噪聲獲取更高的成像分辨率[7]。

本研究通過(guò)對(duì)Capon波束形成算法的研究提出一種適用于發(fā)射聲場(chǎng)的Capon窗加權(quán)方式，并結(jié)合幅度變跡技術(shù)，使陣列波在更大的范圍內(nèi)保持非衍射特性。

2 理論介紹

2.1 陣列波理論

陣列波是非衍射波的一種，基于一般標(biāo)量波動(dòng)方程推導(dǎo)而來(lái)[8]。一般N階齊次標(biāo)量波動(dòng)方程見(jiàn)式(1)，其中xi(j=1，2，…，N)代表N維空間直角坐標(biāo)，t代表時(shí)間，N是大于0的整數(shù)，c是聲速，Φ(x1，x2，……，xN；t)=f(s)是N維聲場(chǎng)函數(shù)。

(1)

令N=2，x1=x，x2=z解得，

ΦL(x，z-c1t)=(cosKxx)eiKz(z-c1t)

(2)

2.2 Capon算法應(yīng)用到發(fā)射聲場(chǎng)

標(biāo)準(zhǔn)Capon算法(standard capon beamformer，SCB)也稱(chēng)之為最小方差無(wú)畸變響應(yīng)，其基本目的是尋找到最佳加權(quán)向量W=[W1，W2，…，WN]，在保證期望方向θ增益最大的條件下，使干擾信號(hào)和噪聲信號(hào)的輸出功率最小[9]。Capon算法一般用于接收聲場(chǎng)的處理，本研究用于發(fā)射聲場(chǎng)，因此假設(shè)Capon算法處理的聲場(chǎng)是某深度下傳感器接收到的信號(hào)。見(jiàn)圖1，仿照產(chǎn)生分層陣列波的發(fā)射傳感器，假設(shè)某深度zt有N個(gè)線(xiàn)性傳感器陣列組成的二維傳感器進(jìn)行發(fā)射聲場(chǎng)接收。每個(gè)線(xiàn)性傳感器有M個(gè)陣元。因此本研究某時(shí)刻下發(fā)射聲場(chǎng)系統(tǒng)的功率可表達(dá)為：

(3)

圖1 陣列波傳播空間示意圖Fig 1 The diagram of layered array beams propagation

式(3)中，

(4)

設(shè)Φxx為深度zt下矩形傳感器上接收信號(hào)X(n)的自協(xié)方差。

Capon算法的代價(jià)函數(shù)[7]為minW*ΦxxWsubjecttoW*a(θ)=1，a(θ)表示指向方向θ的方向向量，W*表示加權(quán)向量共軛轉(zhuǎn)置。上式可以通過(guò)拉格朗日因子的方法求解，將這個(gè)問(wèn)題轉(zhuǎn)化為無(wú)約束最優(yōu)化問(wèn)題：

L(w)=w*Φxxw+λ(w*a(θ)-1)

(5)

式中λ是拉格朗日因子。至此就可以求出加權(quán)向量W的值。

(6)

利用延時(shí)疊加法中的延時(shí)處理，將方向向量表示為a=[1，1，…，1]T，T表示矩陣轉(zhuǎn)置[11-12]。

2.3 復(fù)合權(quán)加權(quán)方式

圖2 傳感器復(fù)合權(quán)加權(quán)方式Fig 2 The composite weighting method of transducer

3 Field II真

3.1 Field II仿真函數(shù)

對(duì)聲場(chǎng)仿真主要使用以下函數(shù)：

(1)xdc_2d_array()，該函數(shù)功能是建模2d傳感器[13]。xdc_2d_array()函數(shù)中特殊參數(shù)enabled表示陣元是否采用，1代表有效陣元，0代表無(wú)效未被選用陣元[14]?？捎迷摵瘮?shù)設(shè)計(jì)同心圓環(huán)傳感器。

(2)xdc_get()，該函數(shù)功能是獲取傳感器上每個(gè)陣元的物理信息。xdc_get()函數(shù)會(huì)獲得26行、物理陣元總數(shù)列大小的矩陣，每一行表示的是物理陣元的實(shí)際信息。

(3)ele_waveform()，該函數(shù)的功能是對(duì)每個(gè)陣元逐一激勵(lì)。由于仿真分層陣列波不同橫坐標(biāo)激勵(lì)信號(hào)不同，可以使用該函數(shù)對(duì)每個(gè)陣元逐一激勵(lì)[15]。

(4)calc_hp()，該函數(shù)功能是計(jì)算空間中各位置聲場(chǎng)的強(qiáng)度。calc_hp()函數(shù)逐點(diǎn)對(duì)空間點(diǎn)進(jìn)行聲場(chǎng)計(jì)算，每點(diǎn)返回值分別為“hp”、“start_time”，其中start_time為起始采樣時(shí)刻，hp為從采樣時(shí)刻開(kāi)始以采樣間隔1/fs采樣到的所有聲場(chǎng)信號(hào)。因?yàn)楦鼽c(diǎn)采樣時(shí)刻不同，需做時(shí)間對(duì)齊才可得到該時(shí)刻下聲場(chǎng)分布。

本研究使用150×150個(gè)陣元的傳感器，陣元高、寬均為0.2 mm，陣元x、y方向間隔均為0.1 mm，聲速1 540 m/s，采樣頻率100 MHz，傳感器脈沖響應(yīng)為中心頻率f0=2.5 MHz、脈沖寬度2個(gè)周期的正弦信號(hào)。

3.2 加權(quán)仿真研究

由上圖可見(jiàn)分層陣列波聲場(chǎng)有清晰的6條明亮豎條紋波束，分層陣列波類(lèi)似于平面波的傳播方式，波束以圖中形狀在x-y平面沿z軸向前傳播[2]。圖中加了噪聲的聲場(chǎng)尾影明顯，中間聲場(chǎng)部分也受到了噪聲影響。在實(shí)際應(yīng)用中聲場(chǎng)傳播過(guò)程無(wú)法避免噪聲的干擾，為了使陣列波保持非衍射特性傳播更遠(yuǎn)是將陣列波聲場(chǎng)應(yīng)用于實(shí)際非常重要的一個(gè)待解決問(wèn)題。

本研究通過(guò)Field II分別進(jìn)行了幅度變跡、單一Capon窗加權(quán)、復(fù)合窗加權(quán)的聲場(chǎng)仿真，討論加各種權(quán)值對(duì)含高斯白噪聲的聲場(chǎng)產(chǎn)生的效果見(jiàn)圖4。

圖3分層陣列波近場(chǎng)區(qū)聲場(chǎng)分布

Fig3Thenearfielddistributionofthelayeredarraybeams

圖4 含噪聲場(chǎng)不同加權(quán)方式仿真結(jié)果Fig 4 The simulation results of different weighted methods of noise field

通過(guò)觀察上圖的仿真結(jié)果可知，含有噪聲的聲場(chǎng)在z=150 mm時(shí)波束已經(jīng)完全發(fā)散，波束受噪聲干擾嚴(yán)重已經(jīng)看不到明顯的分層形狀。使用余弦權(quán)進(jìn)行加權(quán)可見(jiàn)在該深度下聲場(chǎng)有分層形狀，該方式對(duì)噪聲有一定的抑制作用。z=150 mm時(shí)加Capon窗的聲場(chǎng)和加復(fù)合權(quán)的聲場(chǎng)效果相近都有很好的非衍射特性，因此，在近場(chǎng)區(qū)域可以使用Capon算法代替復(fù)合權(quán)進(jìn)行激勵(lì)。對(duì)于聲束沿y軸方向尾影的抑制，加復(fù)合權(quán)方式比加單一的Capon窗效果好。復(fù)合權(quán)方式有效的抑制噪聲的同時(shí)減少了聲場(chǎng)邊緣旁瓣幅度，使聲場(chǎng)能量更集中從而可以保持良好的非衍射特性，因此可以傳播的更遠(yuǎn)。在一些工業(yè)檢測(cè)以及醫(yī)學(xué)中較厚組織器官的檢查要求更深的聲場(chǎng)傳播范圍，此時(shí)可以使用復(fù)合權(quán)對(duì)陣列波進(jìn)行加權(quán)處理。

4 結(jié)論

本研究基于Capon算法提出一種適用于發(fā)射聲場(chǎng)的Capon窗加權(quán)方法，并使用FieldII仿真分析各種加權(quán)方式對(duì)含噪聲場(chǎng)的作用效果。通過(guò)仿真結(jié)果得知，使用本研究提出的Capon算法比余弦窗更好的保持了非衍射波的非衍射特性。單一Capon窗仿真的聲場(chǎng)邊緣開(kāi)始發(fā)散，復(fù)合權(quán)方法抑制噪聲的同時(shí)也對(duì)邊緣發(fā)散聲場(chǎng)進(jìn)行了有效抑制。本研究進(jìn)行了理論探索研究，通過(guò)仿真驗(yàn)證了該方法的可行性，在實(shí)際應(yīng)用中可以通過(guò)計(jì)算機(jī)模擬實(shí)際介質(zhì)參數(shù)，并通過(guò)本研究仿真方法獲得發(fā)射聲場(chǎng)信號(hào)計(jì)算Cpaon權(quán)值。也可以使用水聽(tīng)器直接獲取任意空間位置發(fā)射聲場(chǎng)信號(hào)。

該方法具有如下特點(diǎn)：(1)使用超聲領(lǐng)域內(nèi)常用的仿真工具Field II進(jìn)行仿真研究，仿真結(jié)果更可靠，相比傳統(tǒng)的計(jì)算仿真可以避免計(jì)算過(guò)程產(chǎn)生的誤差。軟件仿真能夠更快速便捷的進(jìn)行參數(shù)調(diào)試，節(jié)省了成本，為應(yīng)用研發(fā)提供了理論依據(jù)。(2)本研究基于Field II開(kāi)發(fā)出一套聲場(chǎng)仿真系統(tǒng)，該系統(tǒng)可以實(shí)現(xiàn)任意聲場(chǎng)仿真，并顯示聲場(chǎng)任意時(shí)刻任意平面波形灰度圖。(3)在實(shí)際應(yīng)用中只有有限孔徑的傳感器，激勵(lì)在傳感器邊界會(huì)驟減為零，一些頻率不連續(xù)點(diǎn)附近會(huì)引起較大誤差，為了避免這種誤差的產(chǎn)生，本研究基于Capon波束形成提出復(fù)合窗算法，傳感器中間范圍采用本研究提出的Capon窗加權(quán)激勵(lì)，邊緣采用余弦窗加權(quán)激勵(lì)，使其譜函數(shù)的主瓣包含更多能量，相應(yīng)旁瓣幅度更小。通過(guò)仿真結(jié)果可知加復(fù)合權(quán)激勵(lì)能夠使非衍射波保持非衍射特性傳播更遠(yuǎn)。