本發(fā)明涉及聲速測量的技術(shù)領(lǐng)域，尤其是指一種基于三維正交陣的水下探測路徑速度測量方法及裝置。

背景技術(shù)：

海水聲速，即聲波在海洋中的傳播速度，它不僅是水聲學(xué)研究中的一個重要參數(shù)，而且越來越明顯地被列為海洋學(xué)的基本要素之一。水聲探測是獲取海洋信息的基本手段，其中水下目標定位是水聲探測中的一個關(guān)鍵技術(shù)，它能從定位中獲得水下目標的距離和方位等信息。在進行水下目標定位時，探測路徑聲速往往影響著最終的目標定位結(jié)果，對聲速測量的不精準都會導(dǎo)致不精確的定位。而在一些經(jīng)典的水下目標定位方法，如MUSIC算法和ESPRIT算法中，每條探測路徑上的聲速都是利用一片海水區(qū)域內(nèi)的聲速來代替估計得到的，并且是相同的，但海水環(huán)境復(fù)雜，且受水流因素影響較大，每一條探測路徑上的聲速其實是不一樣的，即對探測路徑聲速的估計是不精準的，這使得這些目標定位方法得不到很好的利用。因此，精準的水下探測路徑速度測量對目標定位至關(guān)重要。

目前，海水聲速測量方法按照測量方式主要分為直接法和間接法兩種方法。直接法是用水聲設(shè)備直接在海洋現(xiàn)場測量海水聲速，它可分為時差法、共振聲譜法、駐波干涉法與相位比較法等方法。直接法雖然操作較簡便，但需在目標位置已知的情況下才能對一條路徑上的聲速進行測量，如時差法，它必須在目標與接收源距離已知的情況下進行測量，而在進行水下探測時，目標位置都是未知的，這使得直接法在目標定位中的水下探測路徑聲速測量中失去了實用性。間接法則是根據(jù)海水中的深度、溫度和鹽度三個主要影響海水聲速大小的因素，經(jīng)過成百上千次的分析和實踐，提出一種經(jīng)驗公式，利用此經(jīng)驗公式來測量海水聲速。根據(jù)經(jīng)驗公式的不同，間接法可分為Del Grosso聲速算法、Wilson聲速算法和Chen-Millero-Li聲速算法三種聲速算法。然而間接法測量的是一個聲速分布較均勻的海水環(huán)境中的聲波速度，它測量的并不是一條路徑上的聲速大小，即點對點聲速；同時，間接法也需要在目標位置已知的情況下才能進行測量，這使得間接法在目標定位中的水下探測路徑聲速測量中同樣不具有實用性。

為了克服直接法和間接法存在的局限性，本發(fā)明使用三個兩兩正交的均勻線陣作為接收陣列，利用從目標反射回來的聲波在三個陣列中形成的方向角間的關(guān)系，得到探測路徑聲速與三個旋轉(zhuǎn)算子之間的關(guān)系式，然后計算出在目標和接收源之間的探測路徑上聲速的大小，實現(xiàn)了在目標定位中對水下探測路徑聲速的精準測量。因此，本專利方法在水下探測路徑聲速測量中更具有實用性。

技術(shù)實現(xiàn)要素：

本發(fā)明的目的在于克服現(xiàn)有方法存在的各種不足，提供一種基于三維正交陣的水下探測路徑聲速測量方法，實現(xiàn)在水下目標定位中對探測路徑上聲速的大小進行精確的測量。

另外，本發(fā)明另一目的在于提供一種基于三維正交陣的水下探測路徑聲速測量裝置。

為了實現(xiàn)上述目的，本發(fā)明至少可通過以下技術(shù)方案之一實現(xiàn)。

本發(fā)明方法的原理是基于空間直角坐標系中任一一條直線與三條軸形成的夾角余弦值的平方和等于1的這個關(guān)系式，得出從目標反射回來的聲波在三個兩兩正交的均勻線陣上形成的方向角正弦值的平方和等于1的關(guān)系，然后利用這個關(guān)系推導(dǎo)出含有三個旋轉(zhuǎn)算子的探測路徑聲速表達式，最后通過計算三個均勻線陣上的旋轉(zhuǎn)算子來求得最終探測路徑聲速大小的。

基于三維正交陣的水下探測路徑聲速測量方法，其具體步驟包括以下幾步：

步驟一：建立三個兩兩正交的均勻線陣上的接收信號模型；

步驟二：計算三個兩兩正交的均勻線陣上的旋轉(zhuǎn)算子；

步驟三：建立聲波在三個兩兩正交的均勻線陣上三個方向角之間的關(guān)系；

步驟四：根據(jù)三個方向角之間的關(guān)系求出探測路徑上聲速v的大小。

進一步地，步驟一具體包括：

以三維正交的均勻線陣上的接收信號模型所在的直線上建立空間直角坐標系，空間直角坐標系三個軸分別為x、y、z；每個均勻線陣上都有M個接收陣元，M為正整數(shù)，且陣元間間距為d，窄帶目標聲源為S，中心頻率為f，目標對應(yīng)于x軸方向線陣的方向角為θ_x，對應(yīng)于y軸方向線陣的方向角為θ_y，對應(yīng)于z軸方向線陣的方向角為θ_z，其中x軸方向線陣的M個陣元為x₁,x₂,…,x_M，y軸方向線陣的M個陣元為y₁,y₂,…,y_M，z軸方向線陣為z₁,z₂,…,z_M；

將x軸方向線陣中的M個陣元分為平移矢量為d的兩個子陣列X_h和Y_h；子陣列X_h由x軸方向線陣中的第一到第M-1個陣元組成，即有：

x_h1(t)＝x₁(t),x_h2(t)＝x₂(t),…,x_h(M-1)(t)＝x_M-1(t)

其中，x_h1(t),x_h2(t),…,x_h(M-1)(t)分別是子陣列X_h上第一個陣元到第M-1個陣元接收到的信號；

子陣列Y_h由x軸方向線陣的第二到第M個陣元組成，即有：

y_h1(t)＝x₂(t),y_h2(t)＝x₃(t),…,y_h(M-1)(t)＝x_M(t)

其中，y_h1(t),y_h2(t),…,y_h(M-1)(t)分別是子陣列Y_h上第一個陣元到第M-1個陣元接收到的信號；

x₁(t),x₂(t),…,x_M(t)是x軸方向線陣中第一個陣元到第M個陣元的接收信號；

對于子陣列X_h接收信號，以第一個陣元x_h1為參考點，則第一個陣元接收的信號為：

x_h1(t)＝s(t)+n_hx1(t)，

其中s(t)表示目標信號，n_hx1(t)表示子陣列X_h第一個陣元上的噪聲；

接收信號滿足窄帶條件，即當信號延遲遠小于帶寬倒數(shù)時，延遲作用相當于使基帶信號產(chǎn)生一個相移；那么子陣列X_h上第m個陣元在同一時刻接收到的信號為：

x_hm(t)＝s(t)a_m(θ_x)+n_hxm(t),m＝1,2,…,M-1

其中a_m(θ_x)中v表示聲波在探測路徑上的速度，n_hxm(t)表示子陣列X_h上第m個陣元上的噪聲；

由于子陣列Y_h和子陣列X_h的相對平移矢量為d，那么子陣列Y_h上第m個陣元在同一時刻接收到的信號為：

其中n_hym(t)表示子陣列X_h上第m個陣元上的噪聲；

將子陣列X_h和子陣列Y_h上的各陣元的接收信號排列成列向量形式，則上面兩式可寫成矢量形式：

X_h(t)＝A_xs(t)+N_hx(t) (公式1)

Y_h(t)＝A_xΦ_xs(t)+N_hy(t) (公式2)

其中，為(M-1)×1的子陣列X_h導(dǎo)向矢量矩陣，

X_h(t)＝[x_h1(t),x_h2(t),…,x_h(M-1)(t)]^T為(M-1)×1的子陣列X_h接收信號矩陣，

Y_h(t)＝[y_h1(t),y_h2(t),…,y_h(M-1)(t)]^T為(M-1)×1的子陣列Y_h接收信號矩陣，

N_hx(t)＝[n_hx1(t),n_hx2(t),…,n_hx(M-1)(t)]^T為(M-1)×1的子陣列X_h噪聲矩陣，

N_hy(t)＝[n_hy1(t),n_hy2(t),…,n_hy(M-1)(t)]^T為(M-1)×1的子陣列Y_h噪聲矩陣；

Φ_x是把子陣X_h和Y_h的輸出聯(lián)系起來的一個因子，稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意x軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

對于y軸方向線陣，同樣將其分為兩個平移矢量為d的子陣列X_v和Y_v；子陣列X_v由y軸方向線陣的第一到第M-1個陣元組成，則有：

x_v1(t)＝y(tǒng)₁(t),x_v2(t)＝y(tǒng)₂(t),…,x_v(M-1)(t)＝y(tǒng)_M-1(t)

其中，x_v1(t),x_v2(t),…,x_v(M-1)(t)分別是子陣列X_v上第一個陣元到第M-1個陣元接收到的信號；

子陣列Y_v由y軸方向線陣的第二到第M個陣元組成，則有：

y_v1(t)＝y(tǒng)₂(t),y_v2(t)＝y(tǒng)₃(t),…,y_v(M-1)(t)＝y(tǒng)_M(t)

其中，y_v1(t),y_v2(t),…,y_v(M-1)(t)分別是子陣列Y_v上第一個陣元到第M-1個陣元接收到的信號；

y₁(t),y₂(t),…,y_M(t)是y軸方向線陣中第一個陣元到第M個陣元的接收信號；

以第一個陣元x_v1為參考陣元，那么子陣列X_v和子陣列Y_v中第m個陣元的接收信號分別為：

x_vm(t)＝s(t)a_m(θ_y)+n_vxm(t),m＝1,2,…,M-1

其中n_vxm(t)和n_vym(t)分別為子陣X_v和Y_v上第m個陣元的加性噪聲；將上式寫成矢量形式有：

X_v(t)＝A_ys(t)+N_vx(t) (公式4)

Y_v(t)＝A_yΦ_ys(t)+N_vy(t) (公式5)

其中，為(M-1)×1的子陣列X_v導(dǎo)向矢量矩陣，

X_v(t)＝[x_v1(t),x_v2(t),…,x_v(M-1)(t)]^T為(M-1)×1的子陣列X_v接收信號矩陣，

Y_v(t)＝[y_v1(t),y_v2(t),…,y_v(M-1)(t)]^T為(M-1)×1的子陣列Y_v接收信號矩陣，

N_vx(t)＝[n_vx1(t),n_vx2(t),…,n_vx(M-1)(t)]^T為(M-1)×1的子陣列X_v噪聲矩陣，

N_vy(t)＝[n_vy1(t),n_vy2(t),…,n_vy(M-1)(t)]^T為(M-1)×1的子陣列Y_v噪聲矩陣；

同樣，Φ_y是把子陣X_v和Y_v的輸出聯(lián)系起來的一個因子，也稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意y軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

對于z軸方向線陣，同樣將其分為兩個平移矢量為d的子陣列X_g和Y_g；子陣列X_g由z軸方向線陣的第一到第M-1個陣元組成，則有：

x_g1(t)＝z₁(t),x_g2(t)＝z₂(t),…,x_g(M-1)(t)＝z_M-1(t)

其中，x_g1(t),x_g2(t),…,x_g(M-1)(t)分別是子陣列X_g上第一個陣元到第M-1個陣元接收到的信號；

子陣列Y_g由z軸方向線陣的第二到第M個陣元組成，則有：

y_g1(t)＝z₂(t),y_g2(t)＝z₃(t),…,y_g(M-1)(t)＝z_M(t)

其中，y_g1(t),y_g2(t),…,y_g(M-1)(t)分別是子陣列Y_g上第一個陣元到第M-1個陣元接收到的信號；

z₁(t),z₂(t),…,z_M(t)是z軸方向線陣中第一個陣元到第M個陣元的接收信號；

以第一個陣元x_g1為參考陣元，那么子陣列X_g和子陣列Y_g中第m個陣元的接收信號分別為：

x_gm(t)＝s(t)a_m(θ_z)+n_gxm(t),m＝1,2,…,M-1

其中n_gxm(t)和n_gym(t)分別為子陣X_g和Y_g上第m個陣元的加性噪聲，將上式寫成矢量形式有：

X_g(t)＝A_zs(t)+N_gx(t) (公式7)

Y_g(t)＝A_zΦ_zs(t)+N_gy(t) (公式8)

其中，為(M-1)×1的子陣列X_g導(dǎo)向矢量矩陣，

X_g(t)＝[x_g1(t),x_g2(t),…,x_g(M-1)(t)]^T為(M-1)×1的子陣列X_g接收信號矩陣，

Y_g(t)＝[y_g1(t),y_g2(t),…,y_g(M-1)(t)]^T為(M-1)×1的子陣列Y_g接收信號矩陣，

N_gx(t)＝[n_gx1(t),n_gx2(t),…,n_gx(M-1)(t)]^T為(M-1)×1的子陣列X_g噪聲矩陣，

N_gy(t)＝[n_gy1(t),n_gy2(t),…,n_gy(M-1)(t)]^T為(M-1)×1的子陣列Y_g噪聲矩陣；

同樣，Φ_z是把子陣X_g和Y_g的輸出聯(lián)系起來的一個因子，也稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意z軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

進一步地，步驟二具體包括：

對于x軸方向線陣，計算X_h(t)的協(xié)方差矩陣R_hxx與X_h(t)和Y_h(t)的互協(xié)方差矩陣R_hxy，然后對其進行特征值分解得到最小特征值為σ_x²，利用σ_x²計算C_hxx＝R_hxx-σ_x²I和C_hxy＝R_hxy-σ_x²Z，其中I是(M-1)×(M-1)階的單位矩陣，Z也是(M-1)×(M-1)階的矩陣，并且

最后計算矩陣束{C_hxx,C_hxy}的廣義特征值分解，得到非零特征值λ_x,此特征值就是x軸方向線陣旋轉(zhuǎn)算子的值，即：

對于y軸方向線陣，同樣計算X_v(t)的協(xié)方差矩陣R_vxx與X_v(t)和Y_v(t)的互協(xié)方差矩陣R_vxy，然后對其進行特征值分解得到最小特征值為σ_y²，利用σ_y²計算C_vxx＝R_vxx-σ_y²I和C_vxy＝R_vxy-σ_y²Z；最后計算矩陣束{C_vxx,C_vxy}的廣義特征值分解，得到非零特征值λ_y,此特征值就是y軸方向線陣的旋轉(zhuǎn)算子的值，即：

對于z軸方向線陣，同樣計算X_g(t)的協(xié)方差矩陣R_gxx與X_g(t)和Y_g(t)的互協(xié)方差矩陣R_gxy，然后對其進行特征值分解得到最小特征值為σ_z²，利用σ_z²計算C_gxx＝R_gxx-σ_z²I和C_gxy＝R_gxy-σ_z²Z。最后計算矩陣束{C_gxx,C_gxy}的廣義特征值分解，得到非零特征值λ_z,此特征值就是z軸方向線陣上的旋轉(zhuǎn)算子的值，即：

進一步地，步驟三具體包括：

假設(shè)聲波與x軸方向線陣形成的夾角為α，與y軸方向線陣形成的夾角為β，與z軸方向線陣形成的夾角為γ，空間直角角坐標系的原點為O，三個線陣在坐標軸上分別的設(shè)OD所在直線即為聲波入射方向，將其轉(zhuǎn)換到空間直角坐標系中；

聲波入射方向上的一個點D在各線陣即各坐標軸上的正投影為點A、B、C，則有：

進一步有：

而方向角的大小是聲波與線陣法線之間的夾角大小，并且當聲波到達的第一個陣元不是線陣的參考陣元時，時延參數(shù)τ是小于0，又因為所以此時方向角是負值，于是就有：

根據(jù)以上三個關(guān)系式子，則可得：

cos²α＝sin²θ_x，cos²β＝sin²θ_y，cos²γ＝sin²θ_z

將上式帶入公式(13)得：

sin²θ_x+sin²θ_y+sin²θ_z＝1 (公式14)。

進一步地，步驟四具體包括：

聯(lián)立公式(10)，公式(11)，公式(12)和公式(14)消去θ_x，θ_y和θ_z可得：

整理可得速度v的最終表達為：

利用三個兩兩正交的均勻線陣接收到的聲波信號求出三個旋轉(zhuǎn)算子的值即λ_x、λ_y和λ_z，然后將其帶入公式(16)便可計算出探測路徑上聲速v的大小，即能在目標定位中對水下探測路徑聲速進行精確的測量。

進一步地，當有多個目標即對應(yīng)有多條探測路徑時，每個均勻線陣上的旋轉(zhuǎn)算子是一個對角矩陣，每一個對角元素都對應(yīng)于一個目標，求得的特征值也有多個，它們一一對應(yīng)著對角矩陣上的對角元素，最后利用這些特征值可以求得每條探測路徑上的聲波速度。

實現(xiàn)所述基于三維正交陣的水下探測路徑聲速測量方法的裝置，其包括處理與控制模塊、發(fā)射模塊、接收模塊、A/D和D/A轉(zhuǎn)換模塊、電源模塊和外設(shè)接口；處理與控制模塊由一個處理器組成，并且與A/D和D/A轉(zhuǎn)換模塊、電源模塊和外設(shè)接口相連；處理器將設(shè)置的控制參數(shù)通過D/A轉(zhuǎn)換模塊轉(zhuǎn)換為模擬信號，從而控制發(fā)射模塊，使發(fā)射模塊發(fā)射指定的聲波信號；同時，還能夠?qū)慕邮漳K傳送回來并經(jīng)過A/D轉(zhuǎn)換模塊轉(zhuǎn)換后的數(shù)字信號進行處理，然后計算出精確的水下探測路徑聲速大?。惶幚砥鬟€能夠?qū)⑻幚砗玫臄?shù)據(jù)傳送給外設(shè)接口，再通過有線接口或者無線接口傳送到外部裝置；

發(fā)射模塊由一個阻抗匹配電路和一個超聲波發(fā)送探頭組成，并且與A/D轉(zhuǎn)換模塊和電源模塊相連；發(fā)射模塊能夠根據(jù)處理與控制模塊發(fā)出并經(jīng)過A/D轉(zhuǎn)換模塊轉(zhuǎn)換后的指令，進行超聲發(fā)射；

接收模塊由三維正交的超聲波接收探頭陣列組成，發(fā)射模塊也自身含有阻抗匹配電路，并與A/D轉(zhuǎn)換模塊和電源模塊相連；發(fā)射模塊能夠接收從目標發(fā)射回來的信號，并通過A/D轉(zhuǎn)換模塊傳送給處理與控制模塊；

A/D轉(zhuǎn)換器和D/A轉(zhuǎn)換器與處理與控制模塊、發(fā)射模塊和接收模塊相連；A/D轉(zhuǎn)換器將接收模塊接收到的模擬信號轉(zhuǎn)換為數(shù)字信號并傳遞給處理與控制模塊進行處理；D/A轉(zhuǎn)換器將處理與控制模塊發(fā)出的數(shù)字信號轉(zhuǎn)換為模擬信號并傳送發(fā)射模塊使其發(fā)射指定的聲波；

外設(shè)接口由一個有線接口和無線接口組成，并且與處理與控制模塊和電源模塊相連，夠提供人機交互，將處理與控制模塊中處理好的數(shù)據(jù)通過有線接口或無線接口傳送到外部裝置；

電源模塊與處理與控制模塊、發(fā)射模塊、接收模塊、A/D和D/A轉(zhuǎn)換模塊和外設(shè)接口相連。

本發(fā)明與現(xiàn)有技術(shù)相比，具有如下優(yōu)點：

1、本發(fā)明相比于現(xiàn)有的直接測量法，本發(fā)明可以在目標定位中對水下探測路徑聲速進行測量，克服了直接法需要在目標位置已知情況下才能對一條路徑上的聲速進行測量的局限性，因此本發(fā)明算法在水下探測路徑聲速測量中更具有實用性。同時本發(fā)明的裝置設(shè)備較直接法簡單，且操作方便。

2、本發(fā)明相比于現(xiàn)有的間接測量法，本發(fā)明可以對一條探測路徑上的聲速進行測量，而不是像間接法測量的只是一個聲速分布較均勻的海水環(huán)境中的聲速，因此，在水下探測路徑聲速的測量中本算法更具實用性。同時，本發(fā)明的測量方法不受海水環(huán)境的影響，無需對水的深度、溫度和鹽度進行測量，操作較直接法簡便。

3、本發(fā)明的運算量和復(fù)雜度都相對適中，從而保證了算法的可行性。

4、本發(fā)明裝置可行性強，安裝簡單。除此之外，現(xiàn)代處理器計算處理能力的不斷提高，這使得本發(fā)明所使用的處理器等芯片的集成度高，并且計算能力強，從而保證了本發(fā)明的可行性。

附圖說明

圖1為實例中測量裝置的硬件結(jié)構(gòu)模塊圖。

圖2為實例中三維正交的超聲波接收探頭陣列模型圖。

圖3為本發(fā)明所用的三維正交均勻線陣模型。

圖4為x軸方向線陣中子陣列X_h的接收信號的示意圖。

圖5為空間直角坐標系中聲波與三個線陣所成夾角示意圖。

圖6為本發(fā)明算法的流程圖。

具體實施方式

下面結(jié)合實施例及附圖對本發(fā)明作進一步詳細的描述，但本發(fā)明的實施方式不限于此。

如圖6，基于三維正交陣的水下探測路徑聲速測量方法，具體步驟包括以下幾步：

步驟一：推導(dǎo)三個兩兩正交的均勻線陣上的接收信號模型。

三維正交的均勻線陣模型場景如圖3所示，并以它們所在的直線上建立空間直角坐標系。每個均勻線陣上都有M個接收陣元，且陣元間間距為d，窄帶目標聲源為S，中心頻率為f，目標對應(yīng)于x軸方向線陣的方向角為θ_x，對應(yīng)于y軸方向線陣的方向角為θ_y，對應(yīng)于z軸方向線陣的方向角為θ_z。其中x軸方向線陣的M個陣元為x₁,x₂,…,x_M，y軸方向線陣的M個陣元為y₁,y₂,…,y_M，z軸方向線陣為z₁,z₂,…,z_M。

將x軸方向線陣中的M個陣元分為平移矢量為d的兩個子陣列X_h和Y_h。子陣列X_h由x軸方向線陣中的第一到第M-1個陣元組成，即有：

x_h1(t)＝x₁(t),x_h2(t)＝x₂(t),…,x_h(M-1)(t)＝x_M-1(t)

其中，x_h1(t),x_h2(t),…,x_h(M-1)(t)分別是子陣列X_h上第一個陣元到第M-1個陣元接收到的信號。

子陣列Y_h由x軸方向線陣的第二到第M個陣元組成，即有：

y_h1(t)＝x₂(t),y_h2(t)＝x₃(t),…,y_h(M-1)(t)＝x_M(t)

其中，y_h1(t),y_h2(t),…,y_h(M-1)(t)分別是子陣列Y_h上第一個陣元到第M-1個陣元接收到的信號。

x₁(t),x₂(t),…,x_M(t)是x軸方向線陣中第一個陣元到第M個陣元的接收信號。

子陣列X_h接收信號的示意圖如圖4所示。以第一個陣元x_h1為參考點，則第一個陣元接收的信號為：

x_h1(t)＝s(t)+n_hx1(t)

其中s(t)表示目標信號，n_hx1(t)表示子陣列X_h第一個陣元上的噪聲。

接收信號滿足窄帶條件，即當信號延遲遠小于帶寬倒數(shù)時，延遲作用相當于使基帶信號產(chǎn)生一個相移。那么子陣列X_h上第m個陣元在同一時刻接收到的信號為：

x_hm(t)＝s(t)a_m(θ_x)+n_hxm(t),m＝1,2,…,M-1

其中a_m(θ_x)中v表示聲波在探測路徑上的速度，n_hxm(t)表示子陣列X_h上第m個陣元上的噪聲。

由于子陣列Y_h和子陣列X_h的相對平移矢量為d，那么子陣列Y_h上第m個陣元在同一時刻接收到的信號為：

其中n_hym(t)表示子陣列X_h上第m個陣元上的噪聲。

將子陣列X_h和子陣列Y_h上的各陣元的接收信號排列成列向量形式，則上面兩式可寫成矢量形式：

X_h(t)＝A_xs(t)+N_hx(t) (公式1)

Y_h(t)＝A_xΦ_xs(t)+N_hy(t) (公式2)

其中，為(M-1)×1的子陣列X_h導(dǎo)向矢量矩陣，

X_h(t)＝[x_h1(t),x_h2(t),…,x_h(M-1)(t)]^T為(M-1)×1的子陣列X_h接收信號矩陣，

Y_h(t)＝[y_h1(t),y_h2(t),…,y_h(M-1)(t)]^T為(M-1)×1的子陣列Y_h接收信號矩陣，

N_hx(t)＝[n_hx1(t),n_hx2(t),…,n_hx(M-1)(t)]^T為(M-1)×1的子陣列X_h噪聲矩陣，

N_hy(t)＝[n_hy1(t),n_hy2(t),…,n_hy(M-1)(t)]^T為(M-1)×1的子陣列Y_h噪聲矩陣。

而Φ_x是把子陣X_h和Y_h的輸出聯(lián)系起來的一個因子，稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意x軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

對于y軸方向線陣，同樣將其分為兩個平移矢量為d的子陣列X_v和Y_v。子陣列X_v由y軸方向線陣的第一到第M-1個陣元組成，則有：

x_v1(t)＝y(tǒng)₁(t),x_v2(t)＝y(tǒng)₂(t),…,x_v(M-1)(t)＝y(tǒng)_M-1(t)

其中，x_v1(t),x_v2(t),…,x_v(M-1)(t)分別是子陣列X_v上第一個陣元到第M-1個陣元接收到的信號。

子陣列Y_v由y軸方向線陣的第二到第M個陣元組成，則有：

y_v1(t)＝y(tǒng)₂(t),y_v2(t)＝y(tǒng)₃(t),…,y_v(M-1)(t)＝y(tǒng)_M(t)

其中，y_v1(t),y_v2(t),…,y_v(M-1)(t)分別是子陣列Y_v上第一個陣元到第M-1個陣元接收到的信號。

y₁(t),y₂(t),…,y_M(t)是y軸方向線陣中第一個陣元到第M個陣元的接收信號。

以第一個陣元x_v1為參考陣元，那么子陣列X_v和子陣列Y_v中第m個陣元的接收信號分別為：

x_vm(t)＝s(t)a_m(θ_y)+n_vxm(t),m＝1,2,…,M-1

其中n_vxm(t)和n_vym(t)分別為子陣X_v和Y_v上第m個陣元的加性噪聲。將上式寫成矢量形式有：

X_v(t)＝A_ys(t)+N_vx(t) (公式4)

Y_v(t)＝A_yΦ_ys(t)+N_vy(t) (公式5)

其中，為(M-1)×1的子陣列X_v導(dǎo)向矢量矩陣，

X_v(t)＝[x_v1(t),x_v2(t),…,x_v(M-1)(t)]^T為(M-1)×1的子陣列X_v接收信號矩陣，

Y_v(t)＝[y_v1(t),y_v2(t),…,y_v(M-1)(t)]^T為(M-1)×1的子陣列Y_v接收信號矩陣，

N_vx(t)＝[n_vx1(t),n_vx2(t),…,n_vx(M-1)(t)]^T為(M-1)×1的子陣列X_v噪聲矩陣，

N_vy(t)＝[n_vy1(t),n_vy2(t),…,n_vy(M-1)(t)]^T為(M-1)×1的子陣列Y_v噪聲矩陣。

同樣，Φ_y是把子陣X_v和Y_v的輸出聯(lián)系起來的一個因子，也稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意y軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

對于z軸方向線陣，同樣將其分為兩個平移矢量為d的子陣列X_g和Y_g。子陣列X_g由z軸方向線陣的第一到第M-1個陣元組成，則有：

x_g1(t)＝z₁(t),x_g2(t)＝z₂(t),…,x_g(M-1)(t)＝z_M-1(t)

其中，x_g1(t),x_g2(t),…,x_g(M-1)(t)分別是子陣列X_g上第一個陣元到第M-1個陣元接收到的信號。

子陣列Y_g由z軸方向線陣的第二到第M個陣元組成，則有：

y_g1(t)＝z₂(t),y_g2(t)＝z₃(t),…,y_g(M-1)(t)＝z_M(t)

其中，y_g1(t),y_g2(t),…,y_g(M-1)(t)分別是子陣列Y_g上第一個陣元到第M-1個陣元接收到的信號。

z₁(t),z₂(t),…,z_M(t)是z軸方向線陣中第一個陣元到第M個陣元的接收信號。

以第一個陣元x_g1為參考陣元，那么子陣列X_g和子陣列Y_g中第m個陣元的接收信號分別為：

x_gm(t)＝s(t)a_m(θ_z)+n_gxm(t),m＝1,2,…,M-1

其中n_gxm(t)和n_gym(t)分別為子陣X_g和Y_g上第m個陣元的加性噪聲。將上式寫成矢量形式有：

X_g(t)＝A_zs(t)+N_gx(t) (公式7)

Y_g(t)＝A_zΦ_zs(t)+N_gy(t) (公式8)

其中，為(M-1)×1的子陣列X_g導(dǎo)向矢量矩陣，

X_g(t)＝[x_g1(t),x_g2(t),…,x_g(M-1)(t)]^T為(M-1)×1的子陣列X_g接收信號矩陣，

Y_g(t)＝[y_g1(t),y_g2(t),…,y_g(M-1)(t)]^T為(M-1)×1的子陣列Y_g接收信號矩陣，

N_gx(t)＝[n_gx1(t),n_gx2(t),…,n_gx(M-1)(t)]^T為(M-1)×1的子陣列X_g噪聲矩陣，

N_gy(t)＝[n_gy1(t),n_gy2(t),…,n_gy(M-1)(t)]^T為(M-1)×1的子陣列Y_g噪聲矩陣。

同樣，Φ_z是把子陣X_g和Y_g的輸出聯(lián)系起來的一個因子，也稱旋轉(zhuǎn)算子，它包含了目標反射回來的信號波前在任意z軸方向線陣中一個陣元偶之間的相位延遲信息，表示為：

步驟二：計算三個兩兩正交的均勻線陣上的旋轉(zhuǎn)算子。

對于x軸方向線陣，計算X_h(t)的協(xié)方差矩陣R_hxx與X_h(t)和Y_h(t)的互協(xié)方差矩陣R_hxy，然后對其進行特征值分解得到最小特征值為σ_x²，利用σ_x²計算C_hxx＝R_hxx-σ_x²I和C_hxy＝R_hxy-σ_x²Z。其中I是(M-1)×(M-1)階的單位矩陣，Z也是(M-1)×(M-1)階的矩陣，并且

最后計算矩陣束{C_hxx,C_hxy}的廣義特征值分解，得到非零特征值λ_x,此特征值就是x軸方向線陣旋轉(zhuǎn)算子的值，即：

對于y軸方向線陣，同樣計算X_v(t)的協(xié)方差矩陣R_vxx與X_v(t)和Y_v(t)的互協(xié)方差矩陣R_vxy，然后對其進行特征值分解得到最小特征值為σ_y²，利用σ_y²計算C_vxx＝R_vxx-σ_y²I和C_vxy＝R_vxy-σ_y²Z。最后計算矩陣束{C_vxx,C_vxy}的廣義特征值分解，得到非零特征值λ_y,此特征值就是y軸方向線陣的旋轉(zhuǎn)算子的值，即：

對于z軸方向線陣，同樣計算X_g(t)的協(xié)方差矩陣R_gxx與X_g(t)和Y_g(t)的互協(xié)方差矩陣R_gxy，然后對其進行特征值分解得到最小特征值為σ_z²，利用σ_z²計算C_gxx＝R_gxx-σ_z²I和C_gxy＝R_gxy-σ_z²Z。最后計算矩陣束{C_gxx,C_gxy}的廣義特征值分解，得到非零特征值λ_z,此特征值就是z軸方向線陣上的旋轉(zhuǎn)算子的值，即：

步驟三：建立聲波在三個兩兩正交的均勻線陣上三個方向角之間的關(guān)系。

假設(shè)聲波與x軸方向線陣形成的夾角為α，與y軸方向線陣形成的夾角為β，與z軸方向線陣形成的夾角為γ，將其轉(zhuǎn)換到空間直角坐標系中，如圖5所示。

其中OD所在直線即為聲波入射方向。

由圖可知，于是就有：

而方向角的大小是聲波與線陣法線之間的夾角大小，并且當聲波到達的第一個陣元不是線陣的參考陣元時，時延參數(shù)τ是小于0，又因為所以此時方向角是負值，于是就有：

根據(jù)以上三個關(guān)系式子，則可得：

cos²α＝sin²θ_x，cos²β＝sin²θ_y，cos²γ＝sin²θ_z

將上式帶入公式(13)得：

sin²θ_x+sin²θ_y+sin²θ_z＝1 (公式14)

步驟四：根據(jù)三個方向角之間的關(guān)系求出探測路徑上聲速v的大小。

聯(lián)立公式(10)，公式(11)，公式(12)和公式(14)消去θ_x，θ_y和θ_z可得：

整理可得速度v的最終表達為：

因此，利用三個兩兩正交的均勻線陣接收到的聲波信號求出三個旋轉(zhuǎn)算子的值(即λ_x，λ_y和λ_z)，然后將其帶入公式(16)便可計算出探測路徑上聲速v的大小。

本發(fā)明提出的測量方法可以在目標定位中對水下探測路徑聲速進行精確的測量，克服了現(xiàn)有直接測量法和間接測量法存在的缺陷。

同時，本發(fā)明方法可以推廣到有多個目標的情況中去，這時候的每個均勻線陣上的旋轉(zhuǎn)算子是一個對角矩陣，每一個對角元素都對應(yīng)于一個目標，求得的特征值也有多個，它們一一對應(yīng)著對角矩陣上的對角元素，最后利用這些特征值可以求得每條探測路徑上的聲波速度。

本實例提供的基于三維正交陣的水下探測路徑聲速測量裝置如圖1所示，包括處理與控制模塊、發(fā)射模塊、接收模塊、A/D和D/A轉(zhuǎn)換模塊、電源模塊和外設(shè)接口。

處理與控制模塊由一個處理器組成，并且與A/D和D/A轉(zhuǎn)換模塊、電源模塊和外設(shè)接口相連。它能夠設(shè)置特定參數(shù)，該參數(shù)通過D/A轉(zhuǎn)換模塊轉(zhuǎn)換為模擬信號，從而控制發(fā)射模塊，使發(fā)射模塊發(fā)射指定的聲波信號；同時，它還能夠?qū)慕邮漳K傳送回來，并經(jīng)過A/D轉(zhuǎn)換模塊轉(zhuǎn)換后的數(shù)字信號進行處理，然后運用本發(fā)明的算法計算出精確的水下探測路徑聲速大小。除此之外，它還能夠?qū)⑻幚砗玫臄?shù)據(jù)傳送給外設(shè)接口，使它們通過有線接口或者無線接口傳送到外部裝置。

發(fā)射模塊由一個阻抗匹配電路和一個超聲波發(fā)送探頭組成，并且與A/D轉(zhuǎn)換模塊和電源模塊相連。它能夠根據(jù)處理與控制模塊發(fā)出，并經(jīng)過A/D轉(zhuǎn)換模塊轉(zhuǎn)換后的指令，獲取指定的信號并對其進行超聲發(fā)射；

接收模塊由三維正交的超聲波接收探頭陣列組成，具體模型場景如圖2所示，與發(fā)射模塊一樣，它同樣含有阻抗匹配電路，并與A/D轉(zhuǎn)換模塊和電源模塊相連。它能夠接收從目標發(fā)射回來的信號，并通過A/D轉(zhuǎn)換模塊將其傳送給處理與控制模塊。

A/D轉(zhuǎn)換器和D/A轉(zhuǎn)換器與處理與控制模塊、發(fā)射模塊和接收模塊相連。A/D轉(zhuǎn)換器將接收模塊接受到的模擬信號轉(zhuǎn)換為數(shù)字信號并傳遞給處理與控制模塊進行處理；D/A轉(zhuǎn)換器將處理與控制模塊發(fā)出的數(shù)字信號轉(zhuǎn)換為模擬信號并傳送發(fā)射模塊使其發(fā)射指定的聲波。

外設(shè)接口由一個有線接口和無線接口組成，并且與處理與控制模塊和電源模塊相連。它能夠提供人機交互，將處理與控制模塊中處理好的數(shù)據(jù)通過有線接口或無線接口傳送到外部裝置。

電源模塊由一個電源組成，并且與處理與控制模塊、發(fā)射模塊、接收模塊、A/D和D/A轉(zhuǎn)換模塊和外設(shè)接口相連。它能夠為這些模塊供電。

本發(fā)明裝置的主要工作流程如下：操作者通過處理與控制模塊輸入對應(yīng)的參數(shù)，使該模塊產(chǎn)生所需要的數(shù)字信號，該數(shù)字信號經(jīng)過D/A轉(zhuǎn)換模塊將其轉(zhuǎn)換成模擬信號后，傳送給發(fā)發(fā)射模塊，超聲波發(fā)射探頭將對其進行發(fā)射。接收模塊中的三維正交超聲波接收探頭陣列接收到從目標源發(fā)射回來的模擬信號后將其傳送給A/D轉(zhuǎn)換模塊，信號轉(zhuǎn)換為數(shù)字信號后傳送給處理與控制模塊，然后模塊中的處理器利用本發(fā)明算法對其進行計算，得出水下探測路徑聲速的大小。最后處理與控制模塊將計算得出的速度值傳給外設(shè)接口，使其通過有線接口或者無線接口傳送給外部裝置。工作的全程中，電源模塊為所有模塊供電。

本發(fā)明裝置包括處理與控制模塊、發(fā)射模塊、接收模塊、A/D和D/A轉(zhuǎn)換模塊、電源模塊和外設(shè)接口。處理與控制模塊和A/D轉(zhuǎn)換模塊可以用DSP芯片實現(xiàn)(如：TI公司的TMS320VC5509A型號的DSP芯片)，此DSP芯片可實現(xiàn)A/D轉(zhuǎn)換和D/A轉(zhuǎn)換的功能，并能夠?qū)崿F(xiàn)正交均勻線陣的旋轉(zhuǎn)算子和最終聲波傳播路徑速度的計算；發(fā)射模塊使用一個超聲波發(fā)射探頭；接收模塊使用三個兩兩正交的均勻線陣，其中每個陣列包括多個超聲接收探頭，并且數(shù)量相同；輸出模塊使用一個USB有線接口和W-CDMA無線接口。圖1即為本發(fā)明所述裝置的硬件結(jié)構(gòu)模塊圖。

本實例的主要工作步驟具體如下：

步驟1：在一個淡水庫中放置一個目標S。設(shè)定陣列接收模塊和發(fā)射模塊中的各項參數(shù)：將接收模塊中的每個均勻線陣中的陣元個數(shù)M定為9，則每個均勻線陣中前8個陣元為子陣X_h，X_v和X_g，后8個陣元為子陣Y_h，Y_v和Y_g；兩個子陣列之間的距離d設(shè)為5mm；通過處理與控制模塊控制發(fā)射模塊使待發(fā)射的超聲頻率為100kHz，脈沖長度5ms。

步驟2：對超聲接收探頭接收到的目標聲源信號進行采樣；x軸方向均勻線陣接收到的信號為x₁(t),x₂(t),…,x₉(t)，y軸方向均勻線陣接收的信號為y₁(t),y₂(t),…,y₉(t)，z軸方向均勻線陣接收的信號為z₁(t),z₂(t),…,z₉(t)共采樣接收200次，并將接收到的信號通過AD轉(zhuǎn)換模塊傳遞給處理與控制模塊模塊做運算處理。

步驟3：信號在數(shù)據(jù)采集處理與控制模塊中的處理步驟具體如下：

1)將處在x軸方向上的均勻線陣接收到的信號排成矢量形式X_h(t)和Y_h(t)，計算X_h(t)的協(xié)方差矩陣R_hxx與X_h(t)和Y_h(t)之間的互協(xié)方差矩陣R_hxy。同時對y軸和z軸方向上的均勻線陣接收到的信號也進行相同處理，得到R_vxx和R_vxy與R_gxx和R_gxy。

2)對x軸方向線陣中的兩個協(xié)方差矩陣R_hxx和R_hxy進行特征值分解，得到最小的特征值從而有和同時對y軸方向線陣和z軸方向線陣中的兩個協(xié)方差矩陣進行相同的處理，得到C_vxx和C_vxy與C_gxx和C_gxy。

3)分別計算矩陣束{C_hxx,C_hxy}，{C_vxx,C_vxy}和{C_gxx,C_gxy}的廣義特征值分解，得到特征值λ_x，λ_y和λ_z。

4)根據(jù)3)中的特征值，聯(lián)立式子對v進行求解，最終求得

步驟4：將計算出的水下探測路徑上的聲速信息存儲下來，并傳送給外設(shè)接口，使其通過USB有線接口或者W-CDMA無線接口傳送給外部裝置。根據(jù)本發(fā)明算法測量出的水下探測路徑聲速大小為1452.8m/s。