奚 暢,蔡志明,袁 駿
(海軍工程大學 電子工程學院, 湖北 武漢 430033)
定常直航是拖線陣聲吶的常規工況,但在實際探測過程中拖船常需要進行必要的戰術機動,這勢必造成拖線陣陣形畸變,從而偏離假設的直線形態[1-2]。大多數波束形成算法都是以精確已知陣形為前提,高分辨空間譜估計算法對陣形誤差尤其敏感[3],實際陣形與理論假設之間的失配不可避免地導致信號處理增益的下降和方位估計性能的惡化。因此,有效地估計陣形對于提高轉向機動條件下拖線陣探測性能具有重要意義。
water-pulley模型[4]是一種簡單有效的陣形估計方法[5-6],模型認為拖船機動較為平緩時拖船的橫向位移沿陣無衰減地傳遞到尾端,陣上每個點都沿著前一個點的運動軌跡運動,只需利用拖船當前及歷史時刻的位置擬合航跡即可得到拖線陣陣形。其原理較為簡單,運算負擔小,易于工程實現。但缺陷同樣明顯,即適用范圍有限,拖船機動較劇烈時模型失效,成為制約其工程應用前景的最顯著因素。
利用water-pulley模型估計陣形的過程中,對拖船位置坐標進行平滑是擬合航跡前的必要步驟。Gerstoft等[7]根據經驗將平滑窗寬度設置為陣長的1.5倍,在海試中利用water-pulley模型校正陣形,得到較好的陣形估計效果,但此平滑經驗缺乏理論支持,無法將其應用于其他物理屬性的拖線陣。目前尚無公開發表的文獻對航跡平滑方法進行研究。
針對water-pulley模型天然存在的適用范圍有限問題以及航跡平滑研究領域的真空狀態,本文假設在拖船常規機動情況下,聲陣段上各點總是沿同一軌跡運動,提出將拖點航跡平滑為聲陣段航跡的方法,實現陣形實時估計。
Paidoussis[8]在慣性坐標系下分析柔性細長圓柱體微元段的受力平衡情況,建立的流體中零浮力纜的運動方程(Paidoussis方程)為:
(1)
式中:m是單位長度的纜質量;M是單位長度纜等體積的流體質量;dc是纜直徑;L是纜長;U是纜軸向水流速度;Ct和Cn分別是纜的切向和法向阻力系數;C′t是纜尾部的形阻系數,當尾部處于自由狀態時該系數為零。
拖船做簡諧運動且纜達到穩態時,纜上各點均做相同頻率、不同振幅的簡諧運動。將Paidoussis方程無因次化,再代入拖船位移方程,可得零浮力纜的穩態振蕩響應公式[9]:
(2)

利用式(2)計算不同振蕩頻率情況下纜上各位置的歸一化振幅可知,流體中的纜可看作低通濾波器,無因次頻率較大時,拖船處的振蕩在向纜尾傳導的過程中顯著衰減;無因次頻率較小時,纜上各位置的振幅與拖船振幅相差較小,可以近似認為纜上各點均沿拖船航跡運動,此即water-pulley模型[4],又稱worm-in-a-hole運動[10]。并且,振幅從拖點到纜尾呈類指數下降,即越靠近纜尾,振幅變化越平緩,纜尾附近振幅不變,下面針對此問題進行分析。
由式(2)整理可得:

(3)
式中前兩項與纜上位置無關,令第三項等于g(z),纜物理屬性和振蕩頻率一定時,振幅|v(z)|可表示為|g(z)|與常數相乘的形式,因此可利用|g(z)|分析振幅隨纜上位置的變化情況。由于b0是復數,為討論方便,不失一般性,假設b0=2reiφ,r和φ可為任意實數,將b0=2reiφ代入g(z)可得:

(4)
式中:
對式(4)取模并求導可得:
(5)
(6)

(7)



(8)

(9)
(10)


由此得出結論,拖點簡諧運動時陣上各點的振幅從拖點到纜尾逐漸減小,且減小速度越來越慢,纜尾附近振幅不變,理論分析結論與計算結果一致。
拖線陣包括拖纜段和聲陣段兩部分,拖纜段兩端分別連接拖船和聲陣段頭部,聲陣段長度通常小于拖纜段。由上述分析可知,拖點簡諧運動時,靠近拖線陣尾端的聲陣段振幅隨纜上位置變化較平緩,即使拖點振幅與尾端振幅相差較大,聲陣段上各點振幅依然較為一致。因此,拖點劇烈機動時water-pulley模型失效,但可認為聲陣段上各點沿同一軌跡運動,若已知拖纜段尾端航跡,即可得到聲陣段運動軌跡。
式(2)可計算零浮力纜上各點的穩態振蕩響應,但拖纜段密度大于流體密度且尾端非自由狀態,應采用文獻[9]方法調整法向阻力系數Cn并代入形阻系數C′t,即可計算得到拖點簡諧運動時拖纜段尾端運動特性。

如圖1所示,對回轉運動情況下拖船航跡進行平滑可以得到航跡的同心圓,平滑結果的半徑與平滑窗寬度相關。假設聲陣段上各點均沿聲陣段頭部的航跡運動,且聲陣段頭部回轉半徑Ra可通過穩態振蕩響應公式計算得到,則可利用合適寬度的平滑窗對拖船航跡進行平滑,得到聲陣段運動軌跡。

圖1 拖船回轉運動時拖線陣示意Fig.1 Array during circular maneuvering
確定平滑窗寬度時,可利用不同寬度的平滑窗分別對半徑為R的圓進行平滑,找到平滑結果半徑最接近Ra的平滑窗即可。當采用矩形平滑窗時,寬度l可由式(11)近似得到:
(11)
實際偵測過程中,常用的戰術機動模式是在直行的基礎上調整操舵角度,令拖船以一固定轉彎半徑進行轉向,調整到指定航向后將舵角歸零繼續直航,完整航跡由圓弧和圓弧兩端的切線組成,本文將此稱為轉向機動,并針對此機動模式進行分析。
考慮轉彎角度較大的情況,在拖船轉向過程中,拖線陣上各點依次脫離直行穩態并逐漸進入回轉穩態;拖船結束轉向開始直行后,陣上各點依次脫離回轉穩態,拖線陣逐漸被拖直。因此,可將拖船轉向機動過程中拖線陣上各點的運動狀態變化情況歸納為“直行穩態—過渡態—回轉穩態—過渡態—直行穩態”。
用回轉運動中得到的聲陣段運動軌跡的平滑窗對圓弧運動拖船航跡進行處理,可以得到“直線—弧線—圓弧—弧線—直線”,其中圓弧部分的半徑是拖線陣達到回轉穩態時聲陣段運動半徑。
由此可發現,平滑結果的分段線形屬性與聲陣段分段運動狀態屬性一致,且平滑結果中的直線和圓弧可以較好描述聲陣段在直行和圓周穩態下的特性。雖然平滑得到的各分段長度和銜接點位置與真實情況并非完全吻合,但可以將此平滑結果近似認為是聲陣段的運動軌跡,示意圖如圖2所示。

圖2 拖船轉向機動時拖線陣示意Fig.2 Array during course-change maneuvering
拖船機動過程中需要實時估計拖線陣流形,利用上節所述方法對拖船當前及歷史時刻的航跡點進行平滑可以得到聲陣段運動軌跡,但無法確定當前聲陣段的具體位置。因此需先估計拖纜段的陣形,由拖纜段尾端位置確定聲陣段頭部位置,進而沿聲陣段軌跡確定當前聲陣段流形。
如圖2所示,如果有一段軸向長度等于拖纜段長度的弧線光滑地連接拖船和聲陣段運動軌跡上某點,則可認為這段弧線是拖纜段陣形。連接拖船與聲陣段運動軌跡的曲線,可以用完整平滑窗的一部分對距離拖船不足平滑窗寬度的航跡部分進行平滑得到。
為使完整平滑窗和不完整平滑窗的平滑結果銜接部分光滑,應采用三角平滑窗。設窗寬度是W,待平滑航跡的一端為拖船,另一端與拖船沿航跡距離為D,對于所有滿足D≤W的航段內坐標點分別進行加權求和,即可得到拖船附近航跡的平滑結果。當D≤(W+1)/2時,平滑窗可表示為:
(12)
當(W+1)/2 (13) 綜上所述,拖船轉向機動過程中的陣形估計方法可簡述如下:第一步,根據拖船轉彎半徑及拖線陣物理屬性,利用穩態振蕩響應公式計算拖纜段尾端振幅;第二步,確定三角窗寬度,利用其對拖船轉彎半徑對應的圓進行平滑,得到以拖纜段尾端振幅為半徑的圓;第三步,利用此三角窗對拖船轉向機動軌跡進行平滑,并利用平滑窗的一部分對拖船附近航跡進行平滑,得到拖線陣流形。 討論兩種特殊情況。第一種情況:當拖船轉彎角度較小時,拖線陣尚未進入回轉穩態就受到拖船直行的影響逐漸拉直進入直行穩態,拖線陣上各點狀態變化情況為“直行穩態—過渡態—直行穩態”。在此情況下,拖船航跡的圓弧部分距離較短,若小于平滑窗寬度則無法平滑得到圓弧,平滑結果為“直線—弧線—直線”,陣形估計結果與實際變化趨勢相符。 第二種情況:第1節證明越靠近纜尾,振幅變化越平緩,進而認為位于纜尾的聲陣段上各點沿同一軌跡運動。然而當拖纜段較短時,聲陣段頭部可能靠近全陣中部甚至位于全陣前部,此時聲陣段頭部和尾端振幅存在較大差異,不能認為聲陣段上各點沿同一軌跡運動。 陣形估計時,用平滑窗的一部分和完整平滑窗對拖船航跡進行平滑得到拖線陣流形,不同時刻完整平滑窗的平滑結果均在同一軌跡上(圖3中點線部分),而不完整平滑窗并非如此。如圖3所示,當拖纜段較長時,估計的聲陣段首端位于完整平滑窗的平滑結果上,聲陣段上各點沿同一軌跡運動;當拖纜段長度較短時,估計的聲陣段首端位于不完整平滑窗的平滑結果上,聲陣段首端和尾部并非沿同一軌跡運動。拖纜段長度變化時的陣形估計結果與實際趨勢相符。 圖3 不同拖纜段長度的聲陣段首端示意Fig.3 Top of sonar array for different cable length 采用表1所示拖線陣參數,假設拖纜段和聲陣段均為光滑圓柱體,法向及切向阻力系數分別為1.275 8及0.004 6。利用穩態振蕩響應公式計算不同無因次頻率纜上各位置的振幅以及相對于拖船的相位差,結果如圖4、圖5所示。 表1 拖線陣參數 圖4 纜上各位置振幅Fig.4 Normalized amplitude of each point on array 圖5 纜上各位置與拖船之間的相位差Fig.5 Phase difference between array and tow point 由圖4可知,拖點簡諧運動時陣上位置越靠近纜尾振幅變化越平緩,即使拖點振幅與尾端振幅相差較大,聲陣段上各點振幅依然較為一致,可知第1節的假設合理。 由圖5可知,當無因次頻率小于5時,纜尾與拖船相位差小于π,滿足拖纜坐標系條件,穩態振蕩響應模型適用于拖船回轉運動的情況。 Ablow模型[11]經過海試數據驗證,具有一定的可靠性,可以將其計算結果作為真實值。利用表1所示拖線陣參數,令拖船以不同無因次頻率對應的轉彎半徑轉過150°,用Ablow模型計算轉向機動過程中的真實陣形,并用本文方法和不對航跡進行平滑的樸素water-pulley模型估計陣形。假設400 m長的聲陣段上均勻嵌入81個陣元,將估計陣形、直線陣對應的陣形與真實陣形對比,認為估計陣形與真實陣形各自首尾陣元連線法向之差是陣形估計誤差導致的方位估計偏差,分別計算轉向機動過程中的平均陣增益損失和平均方位估計偏差,結果如圖6、圖7所示。 圖6 不同無因次頻率的陣增益損失Fig.6 Gain loss for different dimensionless frequency 圖7 不同無因次頻率的方位估計偏差Fig.7 Direction of arrival error for different dimensionless frequency 由圖6~7可知,本文方法與樸素water-pulley模型相比可以使轉向過程中的平均陣增益提高約0.8 dB,平均方位估計偏差減小約4.7°,與不做陣形估計相比可使陣增益提高約1.9 dB,方位估計偏差減小約39.2°,間接證明了算法的有效性。 進一步分析可知,無因次頻率小于4時,本文方法估計陣形依然會造成一定的陣增益損失和方位估計偏差,這是由于本文假設聲陣段上各點沿同一軌跡運動,而真實情況是聲陣段首尾存在一定振幅差(如圖4所示)。無因次頻率為5時,轉彎半徑過小導致真實陣形未進入回轉穩態就被重新拖直,本文方法的估計效果未出現明顯惡化,依然較為穩健有效,這與第3節的討論結果一致。 本文陣形估計方法的關鍵是確定三角平滑窗寬度,它由拖船轉彎半徑和拖線陣物理屬性共同決定。拖線陣在實際使用的過程中,拖纜段長度可隨應用場景變化,其余物理屬性均固定不變。因此,算法的輸入變量是拖船轉彎半徑和拖纜段長度,下面分析平滑窗寬度對輸入的敏感性。 令轉彎半徑等于500 m,拖纜段長度從400 m到800 m變化,其余拖線陣參數如表1所示,計算三角窗寬度,結果如圖8(a)所示;利用表1所示拖線陣參數,令轉彎半徑從300 m到700 m變化,計算三角窗寬度,結果如圖8(b)所示。 (a) 拖纜段長度變化(a) Length of cable varies (b) 轉彎半徑變化(b) Turning radius varies圖8 平滑窗寬度變化情況Fig.8 Change of smooth window width 由圖8分析可知,三角平滑窗寬度對拖纜段長度敏感,而轉彎半徑對窗寬度影響較小。拖線陣實際使用過程中,拖纜段長度可以由絞車記錄的全陣長度減去聲陣段長度得到,而轉彎半徑很難測量。因此,可以將實際的拖纜段長度和拖船常規操舵角度對應的轉彎半徑用于三角平滑窗寬度的計算,從而將本文方法簡化為一種更易于工程實現的方法。 海試數據來源于2018年1月在我國東海海域實施的一次綜合性水聲試驗,試驗中拖船航速大約5.5 kn,在60 s時間內發生一次約20°的轉向機動,用此船常規舵角對應的轉彎半徑計算三角平滑窗寬度。不做陣形校正以及用本文方法校正陣形后的方位歷程如圖9、圖10所示。 圖9中不做陣形校正的左右舷方位歷程圖互為鏡像,從中選取兩個目標進行分析,目標1的舷角從57°變化到73°,目標2的舷角從95°變化到63°。由圖10可知,用本文方法進行陣形估計后,目標1在右舷歷程圖上的軌線以及目標2在左舷歷程圖上的軌線均變得聚集且清晰,可以判斷目標1位于右舷且目標2位于左舷,波束形成效果間接證明了簡化算法的有效性和可行性。 需要說明的是,本節在波束形成時采用基于遠場平面波假設的波束形成算法,圖9中目標3是一個位于近場的目標,因此陣形校正反而造成對于此目標波束形成效果的惡化。 圖9 不做陣形估計的方位歷程Fig.9 Bearing time record without array shape estimation 圖10 用本文方法校正陣形的方位歷程Fig.10 Bearing time record using proposed method 通過分析纜的穩態振蕩響應特性,近似認為聲陣段上各點沿同一軌跡運動,設計平滑窗將拖船航跡平滑為聲陣段運動軌跡,利用平滑窗的一部分對拖船附近航跡進行平滑得到當前拖纜段陣形,從而實現轉向機動過程中陣形的實時估計。仿真結果表明:對于陣元間隔為5 m的81個陣元,所提出的方法與樸素water-pulley模型相比可使轉向機動過程中的平均陣增益提高約0.8 dB,平均方位估計偏差減小約4.7°,與不做陣形估計相比可使陣增益提高約1.9 dB,方位估計偏差減小約39.2°,間接證明了算法的有效性。通過輸入敏感性分析,對陣形估計方法進行簡化使其更易于工程實現,海試數據驗證表明簡化的方法可行、有效。
4 實驗結果與分析
4.1 算法仿真驗證





4.2 輸入敏感性分析


4.3 海試數據運用


5 結論