利用壓縮感知理論的聲速剖面反演方法研究

張東玉，邢傳璽，吳耀文

(云南民族大學 電氣信息工程學院，云南 昆明 650500)

海洋擁有著豐富的資源，如礦物資源、生物資源、水資源，海洋資源的開發、利用和保護是水聲學研究的最終目的，水下目標探測、導航、定位以及聲吶系統性能的預測等研究成為近幾年的熱點，這些研究的開展都必須首先獲得海洋環境參數[1].聲速作為水聲研究中重要的環境參數之一，是了解和研究水聲信道和水聲傳播的基礎，是所有水聲學研究領域都不可缺少的重要環境參數，海洋聲速剖面在聲場建模和水聲設備性能的評估方面占據著十分重要的地位[2].

采用直接測量的方法測量水聲環境中的聲速，需要多次實地測量，需要耗費大量的物力人力，相較而言，反演的方法只需一次測量就可多次重復使用，能降低成本，節省資源[3].隨著水聲理論日漸成熟，水聲設備技術不斷完善，為了更好的研究、利用和保護海洋，大量學者提出利用己有的歷史聲速剖面，通過反演的方法獲得海洋的垂直聲速剖面，最常用的方法是使用經驗正交函數模型反演聲速剖面，但經驗正交函數模型反演聲速剖面的方法容易遺漏SSP中微小擾動的特征點，同時計算效率低下，反演誤差較大[5].SSP反演問題傳統上是通過最小化最小二乘代價函數的解的能量來正則化，需要對復雜的SSP結構進行低采樣或使用少量基函數來解釋SSP的變化[6].這種分辨率的降低會導致SSP的變化不確定，當內部波或水流產生強烈的、隨時間變化時，SSP的變化具有不確定性，這種不確定性會嚴重影響其他參數反演的準確性.

近年來，有學者開始將壓縮感知理論應用于海洋聲學中，如應用于波束形成，地聲參數反演.壓縮感知理論能夠將可壓縮信號在降采樣觀測后以很高的精度重建，打破了傳統采樣定理的局限.目前，壓縮感知理論已經廣泛應用于醫學、導航、遙感等領域.文獻[5]證明了CS在海洋聲學中的應用，可以將重構問題可以轉化為求解最小l0范數問題.在與距離無關的淺海環境中，SSP反演可以用CS來解決.在稀疏假設下，海洋中聲學SSP反演可以表述為一個不確定的線性問題，在稀疏域中的SSP擾動可以用形狀基函數來表示.

本文采用Kraken模型對淺海水聽器的垂直線陣壓力觀測進行了正演模擬.用一階泰勒展開方法對正演模型的非線性響應進行線性化[7].在稀疏域中使用形狀基函數表示SSP的擾動.并用帶有l1稀疏約束的最小二乘代價函數的凸優化來求解字典系數向量，最后使用CS方法反演得到SSP.實驗結果表明，使用CS方法可以很好地估計SSP，并驗證了CS方法對于聲速剖面反演的優點.

1 基本原理

1.1 傳播模型

在反演海洋環境參數的研究中，選取適當的傳播模型對水聲傳播實驗中的聲場進行建模是反演研究的重要前提，反演結果的精度由傳播模型的計算精度決定.從水聽器的垂直線陣列(vertical line array，VLA)觀察到海洋環境中的聲壓，然后利用簡正波傳播模型對海洋聲場進行正向建模.簡正波傳播模型的基礎是求解式(1)的波動方程[3].

(1)

(2)

Kraken方法使用有限差分來對簡正波方程進行求解，波動方程的解可以寫為Φ(r)和Ψ(z)的乘積：

P(r,z)=Φ(r)Ψ(z).

(3)

將(3)式代入聲壓方程，可得：

(4)

(5)

把式(5)代入式(4)可得到深度方向的模式方程如式(6)：

(6)

假設海面是絕對軟邊界，即，ψ(0)=0，海底為剛性邊界，水深為D，滿足條件：

(7)

式(7)的解有無窮多個，每一階簡正波都有對應的模態函數Ψm(z)和相應的一個水平波數krm.Ψm(z)也可叫做為本征函數，krm也可稱做本征值.各階簡正波的波數都是正交的，因此有：

(8)

聲壓函數可由所有簡正波疊加表示：

(9)

利用貝賽爾函數推導出簡正波的時域聲壓表達式為:

(10)

水聲傳播信道用沖擊響應h(t)來描述，信號x(t)=e-iω0t作為激勵信號，用來激勵信道h(t)，則相應的輸出為：

y(t)=x(t)*h(t)=p(r,z)e-iω0t.

(11)

對上式進行傅里葉變換得到：

Y(ω)=X(ω)*H(ω)=2πδ(ω+ω0)p(r,z).

(12)

(13)

即求得信道沖擊響應也即聲壓值P.

1.2 參數稀疏化和壓縮感知

由于CS是一種尋找非定值線性系統稀疏解的技術，因此利用一階泰勒展開對SSP的非線性響應進行線性化[7].對于k點離散化海洋SSP，c(x)∈RK建模為：

c(x)=c0+Qx,

(14)

其中，c0∈RK是使用CTD探測儀測量經轉換得出的離散參考聲速.Q=[q1,…,qn]∈RK×N是形狀基函數qn的字典，x∈RN是字典系數向量.SSP由形狀基函數建模，形狀基函數描述了在x中具有少數的非零系數的SSP擾動[9].在M元VLA處接收的聲壓為pobs.

pobs=g(x)+n,

(15)

其中，pobs∈CM，g(x)∈CM是模型傳播函數，n∈CM是高斯白噪聲，聲源和其他海洋環境已知.

模型傳播函數的一階泰勒展開的線性化響應為[8]：

(16)

矩陣D=[d1,…,dn]∈CK×N包含了M個壓力觀測值相對于字典Q中的N個形狀基函數的導數.dn通過對參考剖面進行qn的一部分擾動,用雙邊有限差分計算.

假設dn充分非相干，利用l1范數凸優化方法找到x的稀疏估計.式(15)的稀疏解可以表示為

(17)

(18)

(19)

(20)

其中，λ控制解能量和度量擬合的相對重要性.

2 方法實現

使用Kraken模型對簡正波聲場進行仿真計算，首先編寫程序提前設置好式(10)中的密度、聲源距離、聲源深度、接收陣元深度、衰減系數等參數，調用軟件計算聲場可以得到模態函數和水平波數，再把得到的相關值和參數代入求聲壓的式中即可得到某一頻率對應的聲壓值.

采用Kraken模型對水深為 30 m 的海洋聲場進行了數值模擬.聲場由 100 Hz 的聲源產生，且聲源位于水下深度 10 m 處，在 1 km 范圍內采樣，采用M列8個均勻間隔的VLA單元，跨度6～20 m.實驗參數如下：海底聲速 1 700 m/s，密度 1.5 g/cm3，衰減 0.1 dB/λ.環境圖如圖1所示.

使用2013年在中國某海域實測的聲速剖面數據作為SSP訓練數據集，通過對SSP訓練數據集進行主成分分析，計算出一組EOF，前3個EOF如圖2所示.通過傳播模型計算獲得聲壓場的聲壓值，使用CS方法結合壓縮估算得到的稀疏系數x對SSP進行反演.

圖1 實驗環境圖 圖2 前3個EOF

圖3 SSP反演流程圖

聲速剖面反演步驟的過程流程圖如圖3所示，主要用到的聲速剖面反演具體算法步驟如下：

1) 采用Kraken模型對水深 30 m 的海洋聲場進行了數值模擬，得到水聽器垂直陣列(VLA)中的復聲壓，并使用一階泰勒展開對其進行線性化.

pobs=g(x)+n,

(21)

(22)

2) 從SSP訓練集中使用數據壓縮算法主成分分析(principal component analysis,PCA)方法計算得到訓練集的字典系數向量X和特征向量V.

3) 離散的K個點從CTD轉換得出的線性壓力值離散的K個點的參考聲速結合傳播模型得N個離散為K個點的形狀基函數qn的字典Q和M個壓力觀測值相對于Q中的N個形狀基函數的字典D.

Q=[q1,…,qn]∈RK×N,

(23)

D=[d1,...,dn]∈CK×N.

(24)

(25)

3 結果分析

圖4 發射信號 圖5 水聽器接收的聲壓

發射信號為頻率 120 Hz，持續時間 2 s 的加Blackman窗進行幅度調制后的線性調頻脈沖，時域信號及頻譜如圖4所示.圖5是水聽器陣列接收到的聲壓值P，圖6是8個陣元接收到的信號，陣元接收到的信號來自同一信號，大致輪廓和信號形式是統一的，但由于陣元擺放位置不同，信號傳播的路徑不同，最終接收到的信號有誤差.由接收信號圖可以看出接收信號的時間約為 6.8 s，結合水下的平均聲速 1 500 m/s，可計算得到聲源水平距離在 10 000 m 附近，與 10 000 m 的聲源位置一致，計算結果證明，仿真得到的水聽器陣列接收到的信號的聲壓值可進行下一步的反演實驗.

圖6 8個水聽器的接收信號

圖7 聲速剖面反演結果圖

為了驗證實驗結果的可靠性和有效性，進一步驗證CS方法對于聲速剖面反演的優點，分析CS的性能.對圖7中的實測聲速剖面進行30次聲速剖面計算反演過程，對最優估計值統計分析，求解獲得反演聲速剖面的置信區間，置信區間定義為均值±標準差的結果，近似代表置信度為70%[10-11].反演得到的聲速剖面均值與測量聲速剖面、置信區間與測量聲速剖面對比如圖8所示.

從圖8看出，反演得到的聲速剖面均值與實際測量的聲速剖面非常接近，并且實際測量的反演聲速剖面置信區間內.在30組反演計算中，得到聲速剖面最大誤差為 0.642 8 m/s，最小誤差為 0.442 2 m/s，平均誤差 0.509 5 m/s.

圖8 反演結果對比分析

使用2013年在黃海實驗的真實海試數據作為待反演聲速，分別使用傳統的EOF方法和本文方法進行聲速剖面反演實驗，反演結果圖如圖9所示，(a)圖是使用EOF方法選前3階經驗正交函數作為基函數得到的反演剖面.(b)圖使用EOF結合壓縮感知方法得到的反演剖面.從圖給出的2種方法的仿真結果來看，壓縮感知方法的重構效果明顯優于傳統的EOF方法.使用壓縮感知方法反演得到的聲速剖面和實測的聲速剖面，對 30 m 深度的 120 Hz 的單頻信號的傳播進行波形預報，波形預報對比圖如圖10所示，反演聲速剖面得到的預報波形與實測聲速剖面得到的預報波形基本相同，研究結果表明使用文中方法對聲速進行反演是一種行之有效的方法，在水聲學研究領域具有重要作用.

圖9 反演結果對比分析

圖10 波形預報對比圖

4 結語

本文使用一階泰勒展開對SSP的非線性響應進行線性化，使用海洋聲速統計數據的先驗知識以稀疏表示SSP的形狀基函數字典，將PCA方法應用于海洋SSP數據以計算EOF，結合海洋聲傳播模型通過適當更新的EOF訓練字典，使用l1范數凸優化方法找到稀疏估計值，結合訓練字典使用CS算法對聲速剖面進行反演.計算了SSP估計誤差，誤差的70%置信區間在 ±0.6 m/s 之內.實驗結果表明，本文算法對SSP的反演是一種行之有效的方法.