于東康，楊功流，謝祖輝

（北京航空航天大學(xué)儀器科學(xué)與光電工程學(xué)院，北京100191）

0 引言

現(xiàn)代化大型艦船通常裝備有雷達(dá)、光學(xué)觀瞄設(shè)備、跟蹤測量設(shè)備、武器發(fā)射設(shè)備等［1］。這些設(shè)備正常工作的前提，是其擁有統(tǒng)一的時空基準(zhǔn)。一般而言，時間基準(zhǔn)通常為GPS時鐘或船上當(dāng)?shù)貢r鐘信號，而空間基準(zhǔn)則是船舶在陸地建造期間將各設(shè)備的空間姿態(tài)坐標(biāo)進(jìn)行高精度匹配而建立起來的［2］。

由于艦船并非絕對剛體，在實(shí)際航行過程中難免受日曬、海浪沖擊等外力作用而產(chǎn)生角形變，這使得艦上各設(shè)備失去了統(tǒng)一的空間基準(zhǔn)，影響了其信息共享與協(xié)同作戰(zhàn)的能力。學(xué)者M(jìn)ochalov等人所做的研究表明，艦船的靜態(tài)變形角最大可達(dá)1°～1.5°，其動態(tài)撓曲變形角可達(dá)數(shù)角分。為保證設(shè)備的正常工作，需對船體形變角進(jìn)行精確測量并進(jìn)行補(bǔ)償［3?4］。

近年來，各國不斷加深在船體變形測量領(lǐng)域的研究，并相繼出現(xiàn)了一系列的測量方法。其中，美國先后提出激光跟蹤法、應(yīng)變測量法等方法；俄羅斯Mochalov團(tuán)隊(duì)開展了基于激光陀螺組合體（LGU）的形變角測量方法研究。國內(nèi)方面，秦石喬教授提出了自準(zhǔn)直平行光管法；長春光機(jī)所采用大鋼管法研制了一套光學(xué)形變監(jiān)測系統(tǒng)，并進(jìn)行了實(shí)船驗(yàn)證；學(xué)者鄭佳興針對基于姿態(tài)匹配的變形角解算方法，進(jìn)行了細(xì)致的研究。

相較于傳統(tǒng)的光學(xué)測量法、應(yīng)變測量法等，利用多套捷聯(lián)慣導(dǎo)系統(tǒng)、基于慣性量匹配的船體變形角測量方法擁有測量精度高、易于實(shí)現(xiàn)等優(yōu)點(diǎn)，是時下研究的熱點(diǎn)，而基于角速度匹配［5］、慣性姿態(tài)匹配［6］等的多種匹配算法也陸續(xù)出現(xiàn)。要得到高精度的測量結(jié)果，必須解決的一項(xiàng)難點(diǎn)就是動態(tài)撓曲變形模型參數(shù)的辨識。對此，傳統(tǒng)的做法是通過處理一段變形測量的歷史數(shù)據(jù)來解出相應(yīng)的參數(shù)，即所謂的 “先驗(yàn)數(shù)據(jù)分析法”。但實(shí)際上，受環(huán)境影響，動態(tài)模型參數(shù)是實(shí)時變化的，由先驗(yàn)數(shù)據(jù)解出的既定參數(shù)無法滿足需要。文獻(xiàn)［7］提出了基于交互式多模型濾波的變形測量方案，通過自適應(yīng)估計(jì)的方法解決了參數(shù)變化的問題。文獻(xiàn)［8］提出，在2套INS的角增量輸出差值中包含動態(tài)撓曲角信息，并推導(dǎo)得到了角增量差值的自相關(guān)函數(shù)值，其為指數(shù)衰減正弦信號（EDS）形式，最終利用基于后向預(yù)測及奇異值截?cái)嗟姆椒ǎ═?K法）解出了動態(tài)模型的參數(shù)，具有較高的精度。

受文獻(xiàn)［8］啟發(fā)，本文利用2套INS的角增量差值自相關(guān)函數(shù)為EDS形式建立了模型，采用粒子群算法進(jìn)行動態(tài)撓曲變形模型參數(shù)的在線估計(jì)，并將所估參數(shù)應(yīng)用于基于角速度匹配的Kalman濾波器中，進(jìn)行變形角的估計(jì)。上述方法在仿真實(shí)驗(yàn)中取得了良好的求參效果，且在不同程度噪聲的干擾下均獲得了較高的求參精度。結(jié)合該參數(shù)辨識方法的船體變形角估計(jì)誤差能夠控制在10″左右，滿足了工程應(yīng)用的需要。

1 角速度匹配法及動態(tài)撓曲角模型

1.1 角速度匹配方程

光纖陀螺捷聯(lián)慣導(dǎo)系統(tǒng)（FOG SINS）主要由3只光纖陀螺與3只加速度計(jì)組成，它們分別安裝在由精密機(jī)械結(jié)構(gòu)支撐的3個正交面上。光纖陀螺能夠?qū)崟r敏感載體繞3個正交軸方向的旋轉(zhuǎn)角速度。

假定2套INS分別安裝在艦艇的相應(yīng)位置，如圖1所示。其中，INS1坐標(biāo)系oxyz與船體坐標(biāo)系（b 系）取齊，INS2 坐標(biāo)系 ox′y′z′與艦載設(shè)備坐標(biāo)系（s系）取齊。

圖1 INS安裝坐標(biāo)示意圖Fig.1 Installation schematic diagram of INS

總形變角φ由2部分組成，分別為靜態(tài)形變角Φ和動態(tài)形變角θ（t），即有

在認(rèn)為變形角為小角度的前提下，有

在式（4）基礎(chǔ)上建立量測方程，以角速度差Δω作為濾波觀測量，合理設(shè)計(jì)線性Kalman濾波器，可以實(shí)現(xiàn)對變形角的精確估計(jì)。

1.2 動態(tài)撓曲角模型

對于短時間的觀測而言，艦艇的靜態(tài)變形角可被視作常值，而動態(tài)撓曲變形角來自于外界的動態(tài)擾動。更進(jìn)一步說，其產(chǎn)生的主要原因?yàn)榕烍w在航行過程中所受到的隨機(jī)風(fēng)力和流體力的作用，因此它是一種由外界隨機(jī)力驅(qū)動的隨機(jī)變量，與由白噪聲驅(qū)動的Markov過程類似。故可以用3個方向獨(dú)立的二階Markov過程來建立模型，這種簡單的二階模型適用于大多數(shù)應(yīng)用場合［9］。

動態(tài)撓曲變形的自相關(guān)函數(shù)為

式（6）中，i＝x，y，z分別代表3 個軸向，μi為衰減因子，λi為支配頻率，σi為幅度標(biāo)準(zhǔn)差。自相關(guān)函數(shù)對應(yīng)的濾波器方程為

變形角的估計(jì)精度與動態(tài)撓曲模型的3個參數(shù)（σ2、μ、λ）的準(zhǔn)確性有著直接的關(guān)系。由于艦船的航行環(huán)境具有多變性和不可預(yù)測性的特點(diǎn)，受其影響的動態(tài)模型參數(shù)也是實(shí)時變化的。如何精確辨識出模型參數(shù)并將其應(yīng)用于Kalman濾波，是一項(xiàng)極具挑戰(zhàn)性的工作。

1.3 Kalman濾波設(shè)計(jì)

考慮光纖陀螺漂移的影響，式（4）變化為

引入陀螺漂移誤差模型

式（10）中，κi為隨機(jī)漂移不規(guī)則系數(shù)，σ′i為表征漂移不規(guī)則程度的均方差。由式（7）、式（10）可給出系統(tǒng)狀態(tài)方程的矩陣形式

2 基于粒子群算法的參數(shù)辨識

2.1 角增量差值自相關(guān)函數(shù)

將角速度輸出在采樣時間內(nèi)進(jìn)行積分，推導(dǎo)2套INS的角增量匹配關(guān)系。記采樣周期為t0，考慮第k次采樣，對式（4）兩邊同時積分得到

整理后，得到2套INS在第k個采樣周期內(nèi)的角增量匹配方程

式（14）中，δΘ1k和 δΘ2k分別為INS1和 INS2在單采樣周期內(nèi)的角增量。由φk＝Φ ＋θk，可將其進(jìn)一步寫作

通過所估靜態(tài)變形角 Φ，可對式（15）中的-δΘ1k×Φ項(xiàng)進(jìn)行在線補(bǔ)償，從而將其消掉，即化為

在式（17）的基礎(chǔ)上，文獻(xiàn)［8］詳細(xì)推導(dǎo)了 ΔΘ的自相關(guān)函數(shù)形式，這里直接給出

2.2 粒子群算法的參數(shù)辨識原理

現(xiàn)在面臨的工作是，如何利用式（18）辨識出其中的參數(shù)。學(xué)者吳偉提出，式（18）的形式為指數(shù)衰減正弦波（EDS）疊加的形式，即

式（23）中，信號y（n）由M個EDS信號構(gòu)成，as為復(fù)振幅，w（n）為Gauss白噪聲。

基于這點(diǎn)，吳偉選用了一種解EDS信號的方法——T?K法，來進(jìn)行參數(shù)的求解。其做法是利用RΔΘ（n）的復(fù)共軛數(shù)據(jù)建立后向預(yù)測方程，通過求解預(yù)測誤差濾波多項(xiàng)式的根得到未知參數(shù)。在確定預(yù)測誤差濾波多項(xiàng)式系數(shù)的過程中，還需用到奇異值截?cái)嗟姆椒?。作為一種純數(shù)學(xué)手段，T?K法能夠得到較高的求參精度，但它的短板是計(jì)算量大，對信噪比的門限要求高［10?11］。

作為一種智能算法，粒子群算法（PSO）因其全局搜索能力強(qiáng)，被廣泛用于參數(shù)估計(jì)，并表現(xiàn)出了良好的估計(jì)性能［12］。PSO算法的基本思路為：首先在n維空間內(nèi)初始化s個隨機(jī)粒子（隨機(jī)解），定義為X＝（X1，X2，…，Xs），并定義其中第i個粒子的位置為Xi，速度為Vi，其個體極值為pi，粒子的群體極值為BestSi。在每次迭代過程中，粒子通過式（24）、式（25）完成位置和速度的更新

式（24）中，w（t）為時變慣性權(quán)重，kg＝1，2，…，G，為當(dāng)前進(jìn)化代數(shù)，r1、r2為分布于區(qū)間［0，1］的隨機(jī)數(shù)，c1為局部學(xué)習(xí)因子，c2為全局學(xué)習(xí)因子。

利用PSO算法進(jìn)行動態(tài)變形角模型參數(shù)的辨識，首先需要明確算法模型。回顧式（18），已知它的形式為指數(shù)衰減正弦信號的疊加，借助式（23），可將其進(jìn)一步簡化為

這里，i＝1，2，…，N為數(shù)據(jù)樣本的個數(shù)，a、b、c為待辨識參數(shù)。以式（27）作為算法模型進(jìn)行編碼，并將適應(yīng)度函數(shù)設(shè)為如下形式

式（28）中，y（i）為數(shù)據(jù)樣本真值，J 值越小，越接近真實(shí)值，模型參數(shù)越準(zhǔn)確。PSO算法的具體實(shí)現(xiàn)流程如圖2所示。

圖2 PSO算法流程圖Fig.2 Flow chart of PSO

2.3 動態(tài)變形模型參數(shù)辨識的實(shí)現(xiàn)流程

對利用PSO算法進(jìn)行動態(tài)變形模型參數(shù)辨識的思路進(jìn)行整理：首先推導(dǎo)了2套INS輸出角增量差值的自相關(guān)函數(shù)形式，即式（18）。其中，在簡化角增量差值ΔΘk的過程中，設(shè)定了對 -δΘ1k×Φ項(xiàng)進(jìn)行在線補(bǔ)償。接著，通過對比式（18）和式（21），得出可以利用角增量差值的自相關(guān)函數(shù)對動態(tài)變形模型的參數(shù)進(jìn)行辨識。選用粒子群算法，利用指數(shù)衰減正弦信號通式形式輔助建立算法模型，并合理設(shè)置了PSO參數(shù)以完成求參工作。

依照上述思路，動態(tài)變形模型參數(shù)的在線辨識流程設(shè)置如下，如圖3所示。

圖3 基于PSO算法的動態(tài)變形模型參數(shù)在線辨識流程Fig.3 Online estimation procedure of dynamic flexure model parameters based on PSO

1）求2套INS在每個采樣周期內(nèi)的角增量輸出差值ΔΘk，并補(bǔ)償 -δΘ1k×Φ項(xiàng)。其中，第1次靜態(tài)變形角Φ可由粗對準(zhǔn)大致給出，之后由變形算法估計(jì)結(jié)果給出；

2）求解 ΔΘk自相關(guān)函數(shù)值 RΔΘ（n）；

3）將 RΔΘ（n）代入粒子群算法，利用式（28）進(jìn)行參數(shù)辨識；

4）將辨識得到的參數(shù)應(yīng)用于角速度匹配算法中進(jìn)行變形角估計(jì)，將所估靜態(tài)變形角Φ代回至步驟1，進(jìn)行下一步解算。

對于PSO算法的參數(shù)設(shè)置，通過嘗試和經(jīng)驗(yàn)，將種群規(guī)模s設(shè)置為120；學(xué)習(xí)因子c1設(shè)為1.3，c2設(shè)為1.7；為保證收斂速度，同時避免陷入局部最優(yōu)，將最大進(jìn)化代數(shù)設(shè)為500；使用時變權(quán)重w（t），以確保PSO算法在全程中都擁有良好的搜索性能。

3 仿真實(shí)驗(yàn)

3.1 仿真條件

設(shè)置載體以5m/s的速度做定速直航行駛，其在3個方向的角運(yùn)行模型如下

式（29）中，幅值 A＝5°、B＝4°、C＝4°，角頻率ω1＝2π/7、ω2＝2π/6、ω3＝2π/8，初始相位 ψ0＝θ0＝γ0＝0。

表1 動態(tài)撓曲變形模型參數(shù)設(shè)置Table 1 Parameter setting of dynamic flexure deformation model

在仿真實(shí)驗(yàn)中，其他參數(shù)設(shè)置如下。

圖4 動態(tài)撓曲變形角曲線Fig.4 Diagram of dynamic deformation angle curve

采樣頻率設(shè)為100Hz，采樣時間為10min。

3.2 仿真結(jié)果

1）考慮無噪聲情況，即動態(tài)撓曲變形信號的SNR→∞。選用20s的仿真數(shù)據(jù)進(jìn)行角增量差值自相關(guān)函數(shù)RΔΘ（n）的計(jì)算，結(jié)果如圖5所示。顯然，2套INS在俯仰、橫滾、航向3個方向的RΔΘ（n）均為指數(shù)衰減正弦信號形式。截取前15s的結(jié)果，代入PSO算法進(jìn)行參數(shù)辨識。對50次獨(dú)立實(shí)驗(yàn)的結(jié)果求取均值，結(jié)果如表2所示。

圖5 無噪聲條件下2套INS在3個軸向上角增量差值的自相關(guān)函數(shù)曲線Fig.5 Autocorrelation function curve of angular increment output difference in three axes without noise

表2 不同信噪比情況下50次獨(dú)立求參的結(jié)果均值Table 2 Mean value of the result of 50 independent parameter estimation experiments under different SNR conditions

2）在實(shí)際情況下，受艦船發(fā)動機(jī)振動等因素的影響，動態(tài)變形信號會帶有高頻噪聲。首先考慮動態(tài)撓曲變形信號的SNR＝20dB的情況，人為地給動態(tài)撓曲信號加入噪聲，并設(shè)置信噪比為20dB，考察在受噪聲干擾情況下PSO算法的參數(shù)辨識能力。將角增量差值自相關(guān)函數(shù)曲線截取前15s結(jié)果，同樣利用PSO算法進(jìn)行50次獨(dú)立求參實(shí)驗(yàn)并取均值，得到的結(jié)果如表2所示。

3）考慮動態(tài)撓曲變形信號的SNR＝10dB的情況。設(shè)置信噪比為10dB，進(jìn)一步考察在高噪聲干擾下的求參效果，最終結(jié)果由表2給出。

將辨識結(jié)果與真實(shí)參數(shù)進(jìn)行比對，不難看出在不同程度的噪聲干擾的情況下，利用PSO算法所辨識出的動態(tài)變形模型參數(shù)均具有較高的精度。將3組參數(shù)均應(yīng)用于角速度匹配Kalman濾波器中進(jìn)行變形角估計(jì)，得到的估計(jì)精度相當(dāng)，這里給出了在SNR＝10dB條件下得到的參數(shù)測量結(jié)果，如圖6所示。

結(jié)合PSO辨識算法之后的船體變形角估計(jì)結(jié)果與真實(shí)值十分接近，且由圖6（c）可以清楚地看到，動態(tài)變形角估計(jì)曲線能夠很好地跟隨真實(shí)值的變化情況。圖7給出了總變形角的估計(jì)誤差曲線，對3個軸向的變形角估計(jì)誤差求均方根（RMS），得到的結(jié)果為：10.11″、9.74″、10.25″。

圖6 船體變形角估計(jì)結(jié)果Fig.6 Result of ship hull deformation angle estimation

圖7 總變形角估計(jì)誤差曲線Fig.7 Estimation error curves of total deformation angle

3.3 結(jié)果分析

仿真結(jié)果證實(shí)，本文提出的基于PSO算法的辨識方法能夠有效辨識出船體動態(tài)變形模型參數(shù)，并且在不同程度的噪聲干擾下均具有較高的辨識精度。

從結(jié)合該參數(shù)辨識算法的船體變形角估計(jì)結(jié)果來看，靜態(tài)變形與動態(tài)變形曲線在3s以內(nèi)即可收斂，且動態(tài)變形角曲線能夠很好地跟隨實(shí)際值，總變形角估計(jì)誤差能夠保證在10″左右，可滿足工程應(yīng)用的需要。

4 結(jié)論

本文通過將2套INS的輸出角增量進(jìn)行匹配，提出可以利用角增量差值進(jìn)行船體動態(tài)變形模型參數(shù)的在線辨識。結(jié)合前人工作，推導(dǎo)了角增量差值自相關(guān)函數(shù)的具體形式，并借助指數(shù)衰減正弦信號通式建立了辨識算法模型。選擇PSO算法作為辨識工具，合理設(shè)計(jì)了算法流程，設(shè)置了相應(yīng)參數(shù)，最終細(xì)化為1套完整的基于角速度匹配并結(jié)合動態(tài)變形模型參數(shù)在線辨識算法的船體變形測量流程。通過仿真結(jié)果，得到如下結(jié)論：

1）本文提出的基于PSO算法的動態(tài)變形模型參數(shù)辨識方案，在不同信噪比條件下均具有良好的辨識性能；

2）結(jié)合該在線辨識方法的船體變形測量算法，所估靜態(tài)變形角與動態(tài)變形角均能夠很好地跟隨真實(shí)值的變化，證明了參數(shù)辨識的可靠性；

3）最終的總變形角估計(jì)誤差在10″左右，且算法簡單有效，具備較高的工程應(yīng)用價值。