魏宗康，高榮榮，趙 友

（北京航天控制儀器研究所，北京 100854）

慣性制導(dǎo)工具誤差分離精度對(duì)于提升彈道導(dǎo)彈命中精度有著舉足輕重的意義。在利用彈道導(dǎo)彈遙外測(cè)數(shù)據(jù)分離慣性測(cè)量系統(tǒng)誤差系數(shù)時(shí)，可采用遙外測(cè)速度誤差作為觀測(cè)量，其優(yōu)點(diǎn)是速度誤差反映了加速度計(jì)組合和陀螺儀組合的測(cè)量誤差，另外一個(gè)優(yōu)點(diǎn)是建立速度環(huán)境函數(shù)矩陣后可直接通過解方程求解誤差系數(shù)，過程中沒有微分計(jì)算。

基于速度環(huán)境函數(shù)的遙外測(cè)誤差分離時(shí)，首先要確定誤差模型的結(jié)構(gòu)。彈道導(dǎo)彈飛行軌跡的特點(diǎn)是慣性系統(tǒng)只工作于主動(dòng)段，而不能全姿態(tài)任意方位都存在大過載或大機(jī)動(dòng)，這就決定了選擇的誤差模型結(jié)構(gòu)矩陣中部分系數(shù)之間具有相關(guān)性，而最小二乘法對(duì)于強(qiáng)相關(guān)結(jié)構(gòu)矩陣有較低的適應(yīng)能力，會(huì)導(dǎo)致分離的誤差系數(shù)偏離真值較大。為此，如何在強(qiáng)相關(guān)條件的約束下，實(shí)現(xiàn)基于彈道導(dǎo)彈遙外測(cè)數(shù)據(jù)的慣性測(cè)量系統(tǒng)誤差系數(shù)的精確分離是一個(gè)難題。為解決該問題，文獻(xiàn)[1]采用靈敏度理論對(duì)制導(dǎo)工具誤差進(jìn)行分離，通過比較兩個(gè)誤差系數(shù)的靈敏度曲線變化趨勢(shì)的一致性來判斷兩個(gè)誤差系數(shù)的相關(guān)性，但是此方法的缺點(diǎn)是沒有具體的量化指標(biāo)。文獻(xiàn)[2]和文獻(xiàn)[3]分別提出了主成分估計(jì)、偏最小二乘估計(jì)、嶺估計(jì)法、粒子群算法等不同的估計(jì)方法，這些方法的核心思想是在模型不降階的前提下估計(jì)出系數(shù)的近似解，但這些方法帶來的問題就是估計(jì)有偏差，沒有從根本上解決慣性制導(dǎo)工具誤差系數(shù)的精確分離問題。文獻(xiàn)[4]采用顯著性分析方法對(duì)石英加速度計(jì)的二次項(xiàng)誤差系數(shù)進(jìn)行分離，但是并沒有考慮誤差系數(shù)之間的相關(guān)性問題，分離方法缺乏嚴(yán)謹(jǐn)性。

針對(duì)該問題本文提出了一種把相關(guān)性檢驗(yàn)和顯著性檢驗(yàn)分階段結(jié)合起來的慣性制導(dǎo)工具誤差分離方法，重點(diǎn)不是在于估計(jì)結(jié)構(gòu)矩陣相關(guān)時(shí)的每個(gè)系數(shù)，而是估計(jì)不同相關(guān)系數(shù)之間組合后的值，采用該方法后分離的結(jié)果更準(zhǔn)確，也較符合工程實(shí)際。

1 慣性測(cè)量系統(tǒng)高階誤差模型

加速度計(jì)組合的高階誤差模型分別為[5-6]：

式(1)中，k0x、k0y、k0z為x、y、z加速度計(jì)的零偏，單位為g；δkx、δky、δkz為x、y、z加速度計(jì)的線性度誤差；δKax、δKay、δKaz為x、y、z加速度計(jì)的線性度不對(duì)稱性誤差；k yx、k zx、k xy、k zy、k xz、kyz為x、y、z加速度計(jì)的安裝誤差角，單位為rad；K2x、K2y、K2z為x、y、z加速度計(jì)的二次項(xiàng)誤差系數(shù)，單位為g/g2；δK2x、δK2y、δK2z為x、y、z加速度計(jì)的奇二次項(xiàng)誤差系數(shù)，單位為g/g2；K xxy、K xxz、K xyz、Kyxy、Kyxz、Kyyz、Kzxy、Kzxz、Kzyz為x、y、z加速度計(jì)的交叉耦合項(xiàng)誤差系數(shù)，單位為g/g2；K3x、K3y、K3z分別為x、y、z加速度計(jì)的三次項(xiàng)誤差系數(shù)，單位為g/g3。

陀螺儀組合的高階誤差模型為[7-8]

式(2)中，DFx、DFy、DFz為x、y、z陀螺儀的常值漂移，單位為°/h；D1x、D1yD1z、D2x、D2y、D2z、D3x、D3y、D3z為x、y、z陀螺儀的一次項(xiàng)漂移，單位為°/h/g；D4x、D4y、D4z、D5x、D5y、D5z、D6x、D6y、D6z為x、y、z陀螺儀的二次項(xiàng)漂移，單位為°/h/g2；D7x、D7y、D7z、D8x、D8y、D8z、D9x、D9y、D9z為x、y、z陀螺儀的交叉耦合二次項(xiàng)漂移，單位為°/h/g2。

2 遙外測(cè)速度誤差

對(duì)于在彈道導(dǎo)彈主動(dòng)段關(guān)機(jī)點(diǎn)就結(jié)束工作的慣性系統(tǒng)來說，由于時(shí)間較短，可在誤差分析時(shí)忽略各反饋回路的影響，簡(jiǎn)化的導(dǎo)航誤差傳播流程見圖1。

圖1 基于地球坐標(biāo)系的導(dǎo)航誤差傳播開環(huán)流程圖Fig.1 Open-loop flow chart of navigation error transmission based on earth coordinate system

采用上述簡(jiǎn)化誤差模型的優(yōu)點(diǎn)是可直接對(duì)姿態(tài)微分方程、速度微分方程和位置微分方程進(jìn)行積分，對(duì)應(yīng)列寫出各項(xiàng)慣性儀表誤差系數(shù)的環(huán)境函數(shù)矩陣，采用最小二乘法直接求解誤差系數(shù)，該求解方法相對(duì)簡(jiǎn)單[9]。圖1 中，

Aε——由陀螺漂移引起姿態(tài)角速度誤差時(shí)的變換矩陣；

Aφ——由姿態(tài)角誤差引起速度誤差變化量時(shí)的變化矩陣；

根據(jù)速度誤差傳播方程即可列寫出把陀螺儀、加速度計(jì)各項(xiàng)誤差系數(shù)對(duì)應(yīng)的速度環(huán)境函數(shù)[10]。速度誤差方程的模型結(jié)構(gòu)為：

上式采用矩陣形式可寫為：

其中，Y表示速度誤差矩陣；C表示速度環(huán)境函數(shù)矩陣；X表示慣性器件的誤差系數(shù)。

3 主成分估計(jì)的缺陷

當(dāng)結(jié)構(gòu)矩陣C不為維列滿秩時(shí)，則不能用最小二乘法求解[11-12]。所謂主成分估計(jì)就是把信息矩陣Φ=CTC中較小的特征值近似簡(jiǎn)化為0，然后求解X的過程。

當(dāng)C不為列滿秩時(shí)，Φ的特征值λ1、λ2…λm中至少有一個(gè)為0。為分析方便，設(shè)λ1、λ2、…、λm為從小到大排序后的結(jié)果，則D=diag（λ1,λ2…λm）可寫為：

其中，ΛB為對(duì)角線所有元素非零的對(duì)角陣。

根據(jù)矩陣D對(duì)應(yīng)的正交變換矩陣表示為：

則有：

由于PA為信息矩陣Φ=CTC的特征值0 對(duì)應(yīng)的特征向量，有

對(duì)上式可解出：

由此可知，當(dāng)C不為列滿秩時(shí)，求解方程Y=CX等價(jià)于求解方程：

一種求解X的應(yīng)用較廣但不準(zhǔn)確的方法是，令：

可得到：

因此，采用主成分估計(jì)時(shí)，仍然需要非相關(guān)處理。只有在對(duì)結(jié)構(gòu)矩陣進(jìn)行降維處理后，才可求得部分系數(shù)的解。主成分辨識(shí)法的另外一個(gè)缺點(diǎn)是特征值大并不意味著對(duì)應(yīng)誤差系數(shù)也顯著，因此，主成分估計(jì)方法無(wú)法識(shí)別誤差系數(shù)在實(shí)際使用中的最優(yōu)估計(jì)問題[13]。

4 結(jié)合相關(guān)性和顯著性的辨識(shí)方法

針對(duì)結(jié)構(gòu)矩陣不為列滿秩的情形，本節(jié)的研究思路則是對(duì)模型適當(dāng)降階以求出各系數(shù)不同組合的精確解[14]。

4.1 結(jié)構(gòu)矩陣列向量相關(guān)時(shí)的系數(shù)計(jì)算

在結(jié)構(gòu)矩陣不為列滿秩時(shí)，則必然有相應(yīng)的列向量相關(guān)。在求解系數(shù)時(shí)，可先找到結(jié)構(gòu)矩陣中最相關(guān)的列后，去掉其中的一列來對(duì)該結(jié)構(gòu)矩陣進(jìn)行簡(jiǎn)化[15]。

下面給出相關(guān)性定義，設(shè)有向量Ci（i=1,2…m）為結(jié)構(gòu)矩陣的列矩陣，當(dāng)Ci全不為零時(shí)，若有某些非零標(biāo)量αi（i=1,2…m）滿足：

則稱向量C1、C2···Cm線性相關(guān)。

若僅當(dāng)：

時(shí)式(15)才成立，則稱向量C1、C2···Cm線性獨(dú)立，或不相關(guān)。

為分析方便，設(shè)向量C1、C2···Cj-1、Cj、Cj+1、···Cm線性相關(guān)，且αj非零，則必然存在某些非零標(biāo)量rj,i（i=1,2…j-1，j+1…m）滿足：

把Cj代入線性方程Y=CX中，有：

4.2 辨識(shí)具體步驟

顯著性檢驗(yàn)只是解決了因子與輸出的關(guān)系，但不能解決因子之間的相關(guān)性問題；相關(guān)性檢驗(yàn)只是解決了因子之間的相關(guān)性問題，而沒有解決因子與輸出的顯著性問題。下面試圖把顯著性檢驗(yàn)和相關(guān)性檢驗(yàn)結(jié)合起來以辨識(shí)系統(tǒng)的系數(shù)。

在結(jié)構(gòu)矩陣中存在相關(guān)的兩個(gè)列時(shí)，采用離線最小二乘法時(shí)信息矩陣CTC的行列式將近似為0，此時(shí)分離的系數(shù)值將嚴(yán)重偏離真實(shí)值，甚至無(wú)法求解。雖然采用顯著性檢驗(yàn)可以簡(jiǎn)化模型的結(jié)構(gòu)，但這種簡(jiǎn)化也存在一個(gè)缺點(diǎn)，就是在結(jié)構(gòu)矩陣奇異或近似奇異時(shí)有可能把兩個(gè)相關(guān)的變量都作為顯著項(xiàng)而保留到模型結(jié)構(gòu)[16]。

為解決奇異線性系統(tǒng)的離線系數(shù)辨識(shí)問題，在本節(jié)提出一種綜合相關(guān)性和顯著性的辨識(shí)方法。就是對(duì)一個(gè)考慮到各種因子的線性系統(tǒng)模型結(jié)構(gòu)，先對(duì)該模型結(jié)構(gòu)進(jìn)行相關(guān)性檢驗(yàn)并簡(jiǎn)化模型結(jié)構(gòu)，直至結(jié)構(gòu)矩陣的所有列之間互不相關(guān)；然后，對(duì)簡(jiǎn)化模型進(jìn)行顯著性檢驗(yàn)并在此基礎(chǔ)上進(jìn)一步簡(jiǎn)化模型結(jié)構(gòu)，直至最終的模型結(jié)構(gòu)為列滿秩，且對(duì)應(yīng)的系數(shù)都顯著。

具體步驟為：

(1) 求出結(jié)構(gòu)矩陣的基向量，并以所有的基向量構(gòu)成基矩陣；

(2) 對(duì)結(jié)構(gòu)矩陣中除了基矩陣之外的列向量，求出用基矩陣表示的系數(shù)，這些系數(shù)用于構(gòu)建基矩陣對(duì)應(yīng)的新的待估計(jì)狀態(tài)變量；在求解過程中采用顯著性檢驗(yàn)，把不顯著的系數(shù)直接置為零；

(3) 根據(jù)基矩陣和觀測(cè)向量，求解出新的狀態(tài)變量的值；在求解過程中采用顯著性檢驗(yàn)，把不顯著的狀態(tài)變量直接置為零。

5 采用高階誤差模型的工具誤差系數(shù)分離

遙外測(cè)速度誤差如圖2，可以看出，速度誤差隨著時(shí)間積累。

圖2 遙外測(cè)速度誤差的擬合結(jié)果Fig.2 Fitting result of remote measurement velocity error

先對(duì)高階誤差模型進(jìn)行相關(guān)性檢驗(yàn)[9-10]，在環(huán)境函數(shù)矩陣中選取46 個(gè)列作為基，這46 個(gè)列對(duì)應(yīng)的誤差系數(shù)分別為ΔT、k0x、δKax、kyx、Kxxy、Kxyz、K3x、k0y、δKay、kxy、kzy、δK2y、Kyxy、Kyyz、K3y、k0z、δkz、δKaz、kyz、Kzxy、Kzyz、DFx、DFy、DFz、D1x、D1y、D1z、D2x、D2y、D2z、D3y、D3z、D4x、D4z、D5x、D5y、D5z、D6x、D6z、D7x、D7y、D7z、D8x、D8y、D8z、D9z。由這46 個(gè)基構(gòu)成基矩陣C′，另外18 個(gè)系數(shù)δkx、kzx、K2x、δK2x、Kxxz、δky、K2y、Kyxz、kxz、K2z、δK2z、Kzxz、K3z、D3x、D4y、D6y、D9x、D9y對(duì)應(yīng)的列向量都可表示為基的線性組合，這些列向量構(gòu)成矩陣C″。

C″中的每一列都可用C′中基的線性組合表示，設(shè)m=64、d=46，用矩陣方程描述為

式中，p1=3、p2=6、p3=7、p4=8、p5=10、p6=14、p7=18、p8=21、p9=27、p10=29、p11=30、p12=32、p13=34、p14=44、p15=48、p16=54、p17=62、p18=63；q1=1、q2=2、q3=4、q4=5、q5=9、q6=11、q7=12、q8=13、q9=15、q10=16、q11=17、q12=19、q13=20、q14=22、q15=23、q16=24、q17=25、q18=26、q19=28、q20=31、q21=33、q22=35、q23=36、q24=37、q25=38、q26=39、q27=40、q28=41、q29=42、q30=43、q31=45、q32=46、q33=47、q34=49、q35=50、q36=51、q37=52、q38=53、q39=55、q40=56、q41=57、q42=58、q43=59、q44=60、q45=61、q46=64。

矩陣R中各元素的值可通過采用最小二乘法求得，但需強(qiáng)調(diào)的是，經(jīng)過顯著性檢驗(yàn)后有些元素值可置為零，從而簡(jiǎn)化方程。比如，

圖3 δkx 對(duì)應(yīng)的列向量元素及其擬合殘差Fig.3 Contrast of δkxvalue and its’ fitting residual error in different conditions

圖4 δky 對(duì)應(yīng)的列向量元素及其擬合殘差Fig.4 Contrast of δkyvalue and its’ fitting residual error in different conditions

在求解出矩陣R中的各元素后，考慮到相關(guān)性后的基于遙外測(cè)速度差的因果關(guān)系如圖5 所示。

圖5 結(jié)構(gòu)矩陣列相關(guān)時(shí)的遙外測(cè)速度差因果關(guān)系圖Fig.5 The causality diagram of remote measurement velocity error when the structure matrix columns are correlated

重新定義一組新的誤差系數(shù)：

基于新的誤差系數(shù)的遙外測(cè)速度差因果關(guān)系如圖6 所示，其結(jié)構(gòu)矩陣為列滿秩，用線性函數(shù)進(jìn)行描述見式(28)。

圖6 結(jié)構(gòu)矩陣列滿秩時(shí)的遙外測(cè)速度差因果關(guān)系圖Fig.6 The causality diagram of remote measurement velocity error based on new error coefficient

再采用最小二乘法并經(jīng)顯著性檢驗(yàn)后，保留顯著項(xiàng)后的遙外測(cè)速度差因果關(guān)系如圖7 所示，表達(dá)式為：

圖7 保留顯著項(xiàng)后的遙外測(cè)速度差因果關(guān)系圖Fig.7 The causality diagram of remote velocity error retaining significance

最終得到顯著的慣性器件誤差模型為：

除上述系數(shù)顯著外，其余各項(xiàng)誤差系數(shù)都不顯著。誤差系數(shù)不顯著的原因分為三種，第一種原因是在地面上零次項(xiàng)和一次項(xiàng)誤差系數(shù)都已進(jìn)行了補(bǔ)償，第二種原因是有些高階誤差項(xiàng)在飛行條件下不能充分激勵(lì)，第三種原因部分誤差只是數(shù)學(xué)描述而并不存在明確的物理意義。

圖8 遙外測(cè)速度誤差補(bǔ)償殘差曲線Fig.8 Residual curve of remote velocity error compensation

6 結(jié)論

本文給出了一種結(jié)合相關(guān)性檢驗(yàn)和顯著性檢驗(yàn)的提高慣性制導(dǎo)精度的方法，通過引進(jìn)環(huán)境函數(shù)矩陣到相關(guān)檢驗(yàn)中，把相互之間關(guān)聯(lián)的制導(dǎo)工具誤差系數(shù)進(jìn)行整合，重構(gòu)出新的遙外測(cè)速度誤差模型。整合后的新系數(shù)對(duì)應(yīng)的環(huán)境函數(shù)矩陣為列滿秩，從而可精確求解誤差系數(shù)，克服了主成分估計(jì)、嶺估計(jì)等方法不能精確求解的缺點(diǎn)。把顯著性檢驗(yàn)引入制導(dǎo)工具誤差系數(shù)分離過程中，有利于簡(jiǎn)化模型，也有利于分析各制導(dǎo)工具誤差系數(shù)與測(cè)量值之間的本質(zhì)特性，并且能夠大幅降低模型的維數(shù)，克服了結(jié)構(gòu)矩陣維數(shù)過多的問題，具有簡(jiǎn)單快捷、容易實(shí)現(xiàn)的優(yōu)點(diǎn)。將采用本文提出的方法分離得到的誤差系數(shù)代入遙外測(cè)速度誤差模型，計(jì)算得到的擬合殘差大部分在±0.2 m/s 之內(nèi)，計(jì)算結(jié)果表明本文所提出的方法可以有效、準(zhǔn)確地分離慣性制導(dǎo)工具誤差系數(shù)，在慣性制導(dǎo)精度提升領(lǐng)域具有很好的應(yīng)用價(jià)值。