應用于聲納的高性能扇形插值方法研究

席偉光，童大為，費玉杰，徐天柏

（海鷹企業(yè)集團有限責任公司，江蘇無錫，214115）

0 引言

高頻多波束圖像聲納因其卓越的圖像質(zhì)量，在海洋工程、環(huán)境監(jiān)測等民用領(lǐng)域以及反水雷、反潛等軍用領(lǐng)域都得到了廣泛的應用[1]。

高頻多波束圖像聲納的探測頁面是一個以聲基陣[2]開角為圓心角的扇形區(qū)域，只有采用扇形顯示方式，才能最大限度地還原水下環(huán)境的真實面貌，因此需要將信號處理輸出的二維矩陣圖像轉(zhuǎn)換成扇形圖像。聲納作為一種專用設備，集成度較高，配備的軟硬件資源有限，很難適配通用扇形插值軟件[3]，且大部分通用軟件的扇形插值功能非專為聲納設計，在圖像質(zhì)量和實時性[4]方面無法滿足使用要求。

RTheta數(shù)字掃描變換算法[5]的插值效果逼真，能滿足聲納圖像質(zhì)量要求，但其復雜的計算過程無法滿足聲納的實時性要求，故優(yōu)化算法是本文研究的重要內(nèi)容。

本文以高頻多波束聲納圖像為研究對象，以R-Theta數(shù)字掃描變換算法為基礎(chǔ)建立扇形插值的數(shù)學模型；使用“查表法”[6]“圖形對稱法”“通用化設計”優(yōu)化數(shù)學模型；使用面向?qū)ο蟮某绦蛟O計方法完成了軟件設計；最后通過聲納設備的試驗結(jié)果，說明方法的有效性。

1 扇形插值原理

本文討論的扇形插值原理基于R-Theta數(shù)字掃描變換算法，它是一種二維線性插值算法[7]，該算法最先應用于B超儀，其原理是利用極坐標與直角坐標的轉(zhuǎn)換關(guān)系，將B超探頭掃描獲得的波束數(shù)據(jù)[8]映射到矩形的顯示區(qū)域中。

高頻多波束圖像聲納與B超儀的工作原理類似。本文以某型前視聲納為例進行說明，該聲納的探測開角為100度，探測扇面內(nèi)均勻分布了200個波束，波束寬度[9]為0.5°，每個波束上等間隔分布600個樣本點，假設扇形的右側(cè)半徑為角度參考線，即θ=0°，由此可求得任意波束上的任意樣本點的徑向距離R和波束偏移角（波束與參考線的夾角）θ。

算法的目的是把如圖1所示的扇形波束圖映射到存貯器矩陣中。為了說明上的方便，將圖1中的局部存貯空間進行放大，局部放大圖見圖2，每個存貯單位（+號所在位置）中應填入合適的像素值。

圖1 扇形波束圖

圖2 扇形圖像的局部放大圖

對于處在相鄰兩個波束之間的所有像素值，使用下列3個方程來插值。

(x,y)是扇形顯示區(qū)域內(nèi)的任意一點，s(*,*)是(x,y)相鄰的四個樣本點的值。而z(i)和z(i+1)是計算f(x,y)的過渡值，Rerr是徑向距離誤差，θerr是方位角相對誤差，由下式確定：

其中R'為點(x,y)的徑向距離，θ'為(x,y)與角度參考線的夾角，ΔR為同一波束上相鄰采樣點(遠端)的徑向距離差，Δθ為波束寬度。

綜上，將(x,y)相鄰兩個波束上的4個樣本點s(i,j)、s(i,j+1)、s(i+1,j)、s(i+1,j+1)代入公式⑴⑵⑶，可得：

由公式⑹可求得顯示區(qū)域扇形范圍內(nèi)任意一點的像素值。

2 方法設計

2.1 建立數(shù)學模型

結(jié)合某型前視聲納，以扇形的圓心為原點建立直角坐標系，如圖3所示。

圖3 建立直角坐標系的扇形波束圖

在直角坐標系中，矩形顯示區(qū)中的任意一點都可以用(x,y)表示。根據(jù)扇形插值原理：

直角坐標系內(nèi)的任意一點(x,y)可以轉(zhuǎn)化為極坐標形式(R',θ')。

θ'的取值范圍為根據(jù)圖2設s(i+1,j+1)對應點的極坐標為(R,θ)，則s(i,j+1)對應的極坐標為(R,θ-Δθ)，(i,j)對應的極坐標為(R-ΔR,θ-Δθ)，(i+1,j)對應的極坐標為(R-ΔR,θ)，根據(jù)圖2中(x,y)和s(i+1,j+1)的位置關(guān)系易推出：

將⑺⑻⑼⑽代入⑾⑿可知s(i+1,j+1)的極坐標可由(x,y)確定，同理f(x,y)周圍的其余三個樣本點的極坐標均可由(x,y)確定。

將⑺⑻⑼⑽⑾⑿代入⑷⑸可知Rerr和θerr由(x,y)確定。由⑹式可知，f(x,y)的值可由周圍四點s(*,*)的值、Rerr和θerr確定，因此只需求得s(*,*)的值，即可得到f(x,y)的值。下文介紹了通過s(*,*)的極坐標獲取s(*,*)的值的方法。

如圖4所示，將待進行扇形插值的200×600（W×H）個樣本點的聲納波束數(shù)據(jù)按順序分布到內(nèi)存中。

設圖2中點s(i+1,j+1)，s(i,j+1)，s(i,j)，s(i+1,j)對應圖4中點A、B、C、D；波束樣本點j的取值范圍為（0,599），波束號i的取值范圍為（0,199）。

設內(nèi)存中波束數(shù)據(jù)的起始地址為pBmembuf[0],則圖4中點A、B、C、D對應的內(nèi)存值如下：

圖4 聲納波束數(shù)據(jù)內(nèi)存分布圖

A：pBmembuf[200(j + 1)+ 200 - i - 1];

B：pBmembuf[200(j + 1)+ 200 - i ];

C：pBmembuf[200j + 200 - i];

D：pBmembuf[200j + 200 - i - 1];

上文設s(i+1,j+1)對應的極坐標為(R,θ)，且由⑾⑿可知R是ΔR的整數(shù)倍，θ是Δθ的整數(shù)倍,結(jié)合⑺⑻⑼⑽可進一步推出：

由⒀⒁可知i和j可通過(x,y)計算獲得，則點A、B、C、D的內(nèi)存序號也可通過(x,y)計算獲得。

由此可實現(xiàn)將直角坐標系的扇形波束圖中的任意一點(x,y)周圍的四點映射為波束域數(shù)據(jù)中對應的內(nèi)存值，再根據(jù)公式⑹計算得到f(x,y)的值。

將圖3扇形顯示區(qū)內(nèi)的像素點按此方法逐一賦值，即可完成扇形插值。

2.2 優(yōu)化數(shù)學模型

在用軟件實現(xiàn)R-Theta算法中，需要多次調(diào)用三角函數(shù)、浮點運算等大開銷算子，為了節(jié)省運算資源，提高圖像刷新率和系統(tǒng)實時性，需優(yōu)化數(shù)學模型。由⒀⒁可知，計算任意一個f(x,y)的值，至少需要進行一次和一次計算，以圖3顯示區(qū)尺寸為800×650（W×H）為例，需分別進行52萬次計算，采用查表法和圖形對稱法可大幅降低計算量。通過上述優(yōu)化，在奔4處理器上測試扇形插值算法，150米量程的單圖刷新時間從140ms降為80ms，性能明顯提高。

本例列舉的數(shù)學模型是根據(jù)某型前視聲納的專用需求設計的，不同聲納對扇形顯示的需求不同，具體體現(xiàn)在扇形圓心角、扇形半徑、顯示區(qū)尺寸、顯示量程等參數(shù)有差異。為了使數(shù)學模型更具通用性，應進行通用化設計。

2.2.1 查表法

查表法的原理是通過預處理[10]的方式提前計算復雜算法在一定范圍的運算結(jié)果，并按規(guī)律儲存到指定內(nèi)存中，當程序需要調(diào)用算法時，查詢并訪問內(nèi)存中的運算結(jié)果即可，因計算機訪問內(nèi)存地址的速度遠高于執(zhí)行算法本身[11]，故可達到降低計算量的目的。使用查表法的約束條件是：算法輸入值的范圍有限且可確定。

圖3顯示區(qū)尺寸為800×650（W×H），則x的取值范圍為[-400,400]，y的取值范圍為[0,650]，可通過以下方法將取值范圍內(nèi)的進行預處理。

初始化預處理結(jié)果存儲數(shù)組：float m_sqrtTable [401][651]; float m_atanTable [401][651];

預處理偽代碼如下：

2.2.2 圖形對稱法

圖形對稱法的原理是利用顯示區(qū)左右兩側(cè)以y軸對稱的特性，首先使用數(shù)學模型計算得到顯示區(qū)一側(cè)像素點相鄰四點的內(nèi)存序號以及Rerr和θerr的值，再利用軸對稱性質(zhì)，推算出另一側(cè)對稱點相鄰四點的內(nèi)存序號及Rerr和θerr的值，該方法可將復雜算法和的運算次數(shù)降低50%，大幅提升運算效率。

舉例說明：圖4中點A、B、C、D為(x,y)的四個相鄰點，由軸對稱性質(zhì)可知，一定存在對稱點A'、B'、C'、D'、(x',y')，已知點A對應的內(nèi)存序號為200(j+1)+200-i-1，則很容易推出點A'的內(nèi)存序號為200(j+1)+i，同理可推出點B'、C'、D'的內(nèi)存序號，且R'err=Rerr；θ'err=θerr。然后根據(jù)公式⑹計算得到，f(x',y')的值。

2.2.3 通用化設計

高頻多波束圖像聲納的探測開角、波束數(shù)、采樣點數(shù)均為聲納固有參數(shù)，不同型號的聲納設備或有差異；而顯示量程為聲納工作參數(shù)，可在聲納使用過程中根據(jù)按需設置。以上參數(shù)均為扇形插值的輸入值，參與了扇形插值的計算過程，只有將這些參數(shù)對用戶開放，才能使扇形插值算法更具通用性。

下文介紹了聲納固有參數(shù)通用化和顯示量程參數(shù)通用化的詳細設計方法。

2.2.3.1 聲納固有參數(shù)通用化

本節(jié)只討論顯示起始量程為0米的情況，非0情況在“顯示量程參數(shù)通用化”中介紹。

探測開角決定了扇形圓心角，扇形圓心角結(jié)合波束數(shù)決定了Δθ，采樣點數(shù)決定了扇形半徑和ΔR，聲納圖像顯示區(qū)尺寸可通過以上參數(shù)估算獲得。

在程序中，聲納固有參數(shù)的用戶接口，可通過宏定義設定：

2.2.3.2 顯示量程參數(shù)通用化

圖像聲納的距離分辨率[12]較高，可達厘米級，假設某聲納的距離分辨率為2.5cm，則在200m量程下每個波束包含的采樣點數(shù)為8000，遠高于普通顯示器的分辨率，為了實現(xiàn)滿量程下的圖像顯示，一般會對樣本點進行抽點處理，以輸出符合圖像顯示區(qū)尺寸的圖像數(shù)據(jù)，即犧牲了距離分辨率來滿足全量程顯示的需求；但也可通過縮小顯示量程范圍的方法，提升聲納圖像的距離分辨率。

如圖5中(1)所示，聲納全量程為0～200m，當目標處于50～100m范圍內(nèi)時，為了提高目標的距離分辨率，可將50～100m量程區(qū)間（陰影部分）的圖像進行放大顯示，50m量程定義為顯示起始量程，100m量程定義為顯示終止量程。

圖5 顯示量程放大圖

一般有兩種方法實現(xiàn)顯示量程范圍內(nèi)的圖像放大：

方法1：如圖5中(2)所示，首先將顯示量程內(nèi)的樣本點重新抽點處理，使之與扇形半徑的像素點一致，然后將顯示量程內(nèi)的數(shù)據(jù)按照扇形的半徑、圓心角進行插值，從圖上可以發(fā)現(xiàn)，顯示起始量程50m處的圓弧已縮小成了一個點，丟失了大量有效數(shù)據(jù)，越接近顯示起始量程，丟失的有效數(shù)據(jù)元素越多，同時顯示終止量程100m處的圓弧被拉伸了一倍，降低了方位分辨率。可見本方法與聲納實際應用場景有較大差異，容易造成目標變形。

方法2：如圖5中(3)所示，同樣，首先將顯示量程內(nèi)的樣本點重新抽點處理，使之與原扇形半徑的像素點一致，然后保持顯示起始量程50m處的弧長和顯示終止量程100m處的弧長不變，并將兩圓弧間的環(huán)寬拉伸至與原扇形半徑一致，如陰影部分所示，陰影部分是新生成扇形中的一部分，根據(jù)上文數(shù)學模型，只需知道扇形的圓心角、半徑、波束數(shù)即可完成扇形插值。新扇形的波束數(shù)不變，圓心角β和半徑R可通過以下方法計算獲得：

r為原扇形半徑，endrange為顯示終止量程，startrange為顯示起始量程，α為原扇形圓心角。由公式（15）可知endrange與startrange越接近新扇形的半徑R越大，應用中可確定endrange與startrange間的最小差值，即聲納在最大距離分辨率時，顯示區(qū)能夠顯示的最大量程范圍。

方法2最大限度還原了聲納探測扇面，聲納目標變形小。

2.3 軟件設計

軟件設計基于C++編程語言，建立了扇形圖像插值類：PolarImage。共由三個文件組成：PolarImage.h、PolarImage.cpp、UsrIf.h，PolarImage.h聲明了PolarImage類，PolarImage.cpp實現(xiàn)了類成員函數(shù)，UsrIf.h為用戶提供了聲納固有參數(shù)和顯示尺寸的配置接口。

如圖6所示為PolarImage類的結(jié)構(gòu)組成。

從圖6可知transform()函數(shù)是PolarImage類提供的唯一用戶接口，其完整聲明應為：static bool transform (unsigned char *InBuf,unsigned char *OutBuf, float startrange, float endrange)；InBuf為待插值的波束數(shù)據(jù)所在的內(nèi)存地址，OutBuf為插值結(jié)果存放地址，startrange為顯示起始量程，endrange為顯示終止量程。

圖6 PolarImage類結(jié)構(gòu)圖

如圖7所示為聲納應用程序及PolarImage類成員之間調(diào)用關(guān)系。聲納應用程序利用UsrIf.h完成對靜態(tài)成員變量的初始化。

圖7 函數(shù)調(diào)用關(guān)系圖

聲納應用程序需要進行扇形插值時調(diào)用transform()函數(shù)，該函數(shù)先對輸入?yún)?shù)進行有效性檢查，通過后調(diào)用instance()函數(shù)獲取PolarImage類對象（若對象不存在，則instance()函數(shù)創(chuàng)建對象再返回，創(chuàng)建對象時調(diào)用構(gòu)造函數(shù)PolarImage()對m_sqrtTable和m_atanTable表空間進行初始化），再利用獲取的PolarImage類對象調(diào)用_transform()函數(shù)，該函數(shù)完成了扇形插值數(shù)學模型的代碼實現(xiàn)，是核心函數(shù)。

聲納應用程序刪除PolarImage類對象時調(diào)用析構(gòu)函數(shù)～PolarImage()，該函數(shù)釋放了m_sqrtTable和m_atan Table表空間。

3 試驗效果

本文研究的扇形插值方法已成功應用于多個型號的聲納設備。

3.1 某型前視聲納應用效果圖

通過水池試驗[13]驗證了扇形顯示效果，試驗水池為矩形，長寬約為20×6m，在聲基礎(chǔ)陣正前方約15m處放置了一個直徑約60cm的環(huán)型鐵塊，正前方1m處放置了一個直徑約2cm的鐵球，分別采用矩形顯示和扇形顯示兩種方式進行了對比。

從圖8矩形顯示圖可見，原本矩形的水池壁，被顯示成塔形，且越接近基陣目標變形越明顯，1m處的小球尺寸已嚴重失真，遠大于實際尺寸，15m處的環(huán)型鐵塊失真不明顯。說明在矩形顯示模式下，目標越近失真越明顯。

圖8 前視聲納矩形顯示圖

從圖9扇形顯示圖可見，矩形的水池壁無變形，1m處的小球呈現(xiàn)為1個亮點，15m處的環(huán)型鐵塊與矩形圖差別不大。可見扇形模式下，近距離目標無失真。

圖9 前視聲納扇形顯示圖

經(jīng)測試，在奔四處理器上完成圖像扇形插值時間約為80ms；在i7 2.4G處理器上完成圖像扇形插值時間約為30ms。均能滿足前視聲納的實時性要求。

3.2 某型分類識別聲納應用效果

該型聲納的固定參數(shù)如下：

探測開角：30°；

波束數(shù)：512；

波束上樣本點數(shù)：1024；

扇形圖像顯示區(qū)大小：540×1038（W×H）；

矩形圖像顯示區(qū)大小：512×1024（W×H）；

通過湖試[14]驗證了扇形顯示效果，試驗目標為一根φ0.066m×64m的陣纜[15]，將陣纜懸浮（利用浮球）于水中，兩端用繩子固定，使其姿態(tài)保持直線，左端距離聲基陣40余米，右端距離聲基陣60余米。

圖10為全量程[16]（0～350m）下的矩形顯示圖和扇形顯示圖，圖中白圈標注部分為目標陣纜，可見在全量程模式下聲納圖像顯示的范圍較大，但距離分辨率較低；矩形顯示模式下的陣纜圖像被拉伸發(fā)生了變形，扇形顯示模式下的陣纜目標圖像較小，真實還原了聲基陣的水中姿態(tài)。

圖10 識別聲納全量程顯示圖

圖11 為顯示量程 (0～80m)下的矩形顯示圖和扇形顯示圖。通過減小顯示量程范圍提升了目標的距離分辨率；矩形顯示模式下的陣纜圖像被拉伸發(fā)生了一定程度的變形，呈現(xiàn)為弧形；扇形顯示模式下能夠還原基陣的真實姿態(tài)。

圖11 識別聲納0～80m顯示量程顯示圖

圖12 為顯示量程(30～80m)下的矩形顯示圖和扇形顯示圖。顯示量程范圍進一步縮小，距離分辨率得到更大提高；矩形顯示模式下的陣纜圖像被拉伸發(fā)生了一定程度的變形，呈現(xiàn)為弧形，扇形顯示模式下基本能夠還原基陣的真實姿態(tài)，但較之圖11的扇形圖像有了些許彎曲。

圖12 識別聲納30～80顯示量程圖

經(jīng)測試在TMS320C6678處理器[17]上完成圖像扇形插值的時間約為50ms。

4 結(jié)語

本文以高頻多波束聲納圖像為研究對象，采用R-Theta數(shù)字掃描變換算法建立了扇形插值的數(shù)學模型；根據(jù)裝備實際使用要求，提出三種優(yōu)化模型的方法：（1）查表法；（2）圖形對稱法；（3）通用化設計；并針對數(shù)學模型完成了軟件設計。通過兩型聲納設備的試驗驗證結(jié)果表明，該扇形插值方法平臺適應性強，且圖像質(zhì)量和實時性均能滿足聲納使用要求，可在行業(yè)內(nèi)進一步推廣應用。