諶炳漢， 左崢嶸,*， 夏魯瑞

(1. 華中科技大學(xué)自動(dòng)化學(xué)院多譜信息處理技術(shù)國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室， 武漢 430074； 2. 航天工程大學(xué)太空安全研究中心， 北京 101416)

目標(biāo)的遙感測(cè)量可有效獲取目標(biāo)運(yùn)動(dòng)特性，有助于對(duì)目標(biāo)飛行試驗(yàn)效果進(jìn)行評(píng)估。目標(biāo)的軌跡重建是遙感測(cè)量系統(tǒng)獲取目標(biāo)信息的關(guān)鍵技術(shù)，目標(biāo)的空間軌跡及運(yùn)動(dòng)參數(shù)通常需要雙平臺(tái)觀察，并采用雙目視覺方法進(jìn)行重建，這提高了觀測(cè)系統(tǒng)的復(fù)雜度，例如，天基紅外預(yù)警系統(tǒng)的高軌部分包括4顆靜止軌道衛(wèi)星和2顆大橢圓軌道衛(wèi)星，低軌部分由21～28顆小型衛(wèi)星組成[1]，其部署和維護(hù)成本極大，因此，在實(shí)際觀測(cè)系統(tǒng)中，常常只進(jìn)行單平臺(tái)的觀測(cè)[2]。這對(duì)單平臺(tái)觀測(cè)條件下的目標(biāo)運(yùn)動(dòng)特性反演提出了迫切的需求。研究人員對(duì)單星條件下的目標(biāo)運(yùn)動(dòng)參數(shù)估計(jì)方法進(jìn)行了較多的研究，參數(shù)估計(jì)方法主要包括基于彈道模板的彈道參數(shù)估計(jì)方法[3-4]和基于目標(biāo)運(yùn)動(dòng)特性建模的參數(shù)估計(jì)方法[5]。為提高估計(jì)精度，研究人員進(jìn)一步將目標(biāo)動(dòng)力特征的先驗(yàn)知識(shí)應(yīng)用于單星彈道導(dǎo)彈參數(shù)估計(jì)中[6]，但目前基于單星的目標(biāo)運(yùn)動(dòng)特性反演方法仍存在反演精度較低的問題，該類方法還不足以支撐目標(biāo)參數(shù)測(cè)量的需求。因此，本文結(jié)合實(shí)際應(yīng)用中所能獲得的目標(biāo)先驗(yàn)知識(shí)，研究具有更高估計(jì)精度的目標(biāo)特性反演方法。

1 單平臺(tái)目標(biāo)運(yùn)動(dòng)特性反演模型

在單平臺(tái)下目標(biāo)軌跡重建思路為：對(duì)獲取的序列圖像進(jìn)行目標(biāo)的提取，確定目標(biāo)在圖像上的位置，獲得不同成像時(shí)刻的目標(biāo)圖像位置。對(duì)每一個(gè)目標(biāo)位置，根據(jù)對(duì)應(yīng)的成像參數(shù)，計(jì)算目標(biāo)視線的方向矢量，從成像器所在位置出發(fā)，考慮大氣折射，對(duì)目標(biāo)視線的光學(xué)路徑進(jìn)行反向追蹤，追蹤路徑與發(fā)射面的交點(diǎn)即為目標(biāo)的空間位置，并通過多幀圖像反演獲得的序列目標(biāo)空間位置進(jìn)行擬合處理，修正目標(biāo)空間位置，反演目標(biāo)運(yùn)動(dòng)參數(shù)。計(jì)算過程中所需的發(fā)射面參數(shù)可根據(jù)經(jīng)驗(yàn)設(shè)定，也可以由其他計(jì)算方式獲得。單平臺(tái)下成像目標(biāo)運(yùn)動(dòng)特性反演的總流程如圖1所示。

通過目標(biāo)的提取，獲得成像時(shí)刻在相機(jī)成像的第k幀圖像上目標(biāo)的位置p=(xk,yk)T，圖像行數(shù)r及列數(shù)c，以及像元視場(chǎng)角φpixel，進(jìn)而可求出該位置對(duì)應(yīng)于目標(biāo)相對(duì)于光軸的偏角(αx,αy)為

目標(biāo)視線方向矢量示意圖如圖2所示,Ose-XseYseZse為傳感器坐標(biāo)系。根據(jù)光軸的偏角，可計(jì)算得到由相機(jī)發(fā)出并最終到達(dá)目標(biāo)的光線在傳感器成像坐標(biāo)系下的初始方向向量[7-8]Tc=(sinαxcosαy,sinαy,cosαxcosαy)。

根據(jù)如圖3所示的成像投影相關(guān)坐標(biāo)系的關(guān)系，可將先驗(yàn)獲得的發(fā)射點(diǎn)位置及發(fā)射面轉(zhuǎn)換到地心坐標(biāo)系下，同時(shí)將目標(biāo)視線矢量和發(fā)射面方程轉(zhuǎn)換到統(tǒng)一的地心坐標(biāo)系，坐標(biāo)系變換涉及傳感器坐標(biāo)系、平臺(tái)坐標(biāo)系、軌道坐標(biāo)系、地心坐標(biāo)系等坐標(biāo)系之間的轉(zhuǎn)換[9-10]。

參見圖3，成像時(shí)刻tk時(shí)的坐標(biāo)系轉(zhuǎn)換矩陣為

圖1 實(shí)現(xiàn)方案總流程圖Fig.1 General flowchart of implementation scheme

圖2 目標(biāo)視線方向矢量示意圖Fig.2 Schematic diagram of target line-of-sight direction vector

圖3 成像投影相關(guān)坐標(biāo)系的關(guān)系示意圖Fig.3 Schematic diagram of relationship between image projection related coordinate systems

(1)

式中：Mearth2orbit為地心坐標(biāo)系Oe-XeYeZe到軌道坐標(biāo)系Oo-XoYoZo的轉(zhuǎn)換矩陣；Morbit2sat為軌道坐標(biāo)系Oo-XoYoZo到平臺(tái)坐標(biāo)系Os-XsYsZs的轉(zhuǎn)換矩陣；Msat2sensor為平臺(tái)坐標(biāo)系Os-XsYsZs到傳感器坐標(biāo)系Ose-XseYseZse的轉(zhuǎn)換矩陣。圖3中的ept和αpt分別為傳感器指向在Os-XsYsZs下的俯仰角和方位角。

根據(jù)目標(biāo)視線方向及平臺(tái)位置可得到過目標(biāo)的射線方程為

(2)

式中：(xs,ys,zs)為傳感器的空間位置；(tx,ty,tz)為射線的法向量；(x,y,z)為目標(biāo)點(diǎn)P的空間位置,可由以下聯(lián)立方程解得

(3)

式中：(qx,qy,qz)為平面法向量。

2 考慮大氣折射影響的單平臺(tái)目標(biāo)空間位置反演

由于受到地心引力的作用，地表大氣密度分布不均。大量的實(shí)驗(yàn)研究表明，大氣垂直方向比水平方向的變化要大1～3個(gè)數(shù)量級(jí)[11]。因此在實(shí)際研究大氣折射對(duì)光線傳播的影響時(shí)，常常忽略大氣水平方向的變化，將大氣視為水平方向均勻分布。將大氣視為球面分層，每層大氣均勻分布并且折射率相同，從而大氣折射率可簡(jiǎn)化成僅隨海拔高度h而變化的量，即N=N(h)[11-14]。

如圖4所示，設(shè)相機(jī)C在地心坐標(biāo)系中坐標(biāo)為(xc,yc,zc)，目標(biāo)P在地心坐標(biāo)系中坐標(biāo)為(xp,yp,zp)，地心Oe坐標(biāo)為(0,0,0)。相機(jī)、目標(biāo)、地心3 點(diǎn)組成的平面坐標(biāo)系Oe-XpYpZp定義如下：平面坐標(biāo)系的原點(diǎn)與地心直角坐標(biāo)系的原點(diǎn)重合，平面OeCP與平面OeXeYe的交線OeXp為平面坐標(biāo)系的Xp軸，Yp軸垂直平面OeCP，則由右手螺旋定則確定Zp軸。平面坐標(biāo)系Oe-XpYpZp到地心坐標(biāo)系Oe-XeYeZe的轉(zhuǎn)換矩陣為

(4)

式中：Mearth2plan為地心坐標(biāo)系Oe-XeYeZe到平面坐標(biāo)系Oe-XpYpZp的轉(zhuǎn)換矩陣。

為了方便計(jì)算光在每層中傳輸路徑，相機(jī)、目標(biāo)、地心平面內(nèi)的點(diǎn)坐標(biāo)和向量可轉(zhuǎn)化到平面坐標(biāo)系下，因?yàn)閅p軸坐標(biāo)始終為0，所以將坐標(biāo)記為二維平面坐標(biāo)(xp,zp)的形式[15]。

當(dāng)光線由第i層射向第i+1層時(shí)，根據(jù)第i+1層入射角，第i層大氣折射率，第i+1層大氣折射率，計(jì)算第i+1層出射角。如圖5所示，當(dāng)光線由第i層射向第i+1層時(shí)，第i+1層光線的入射方向向量vi+1，入射點(diǎn)坐標(biāo)Ti+1=(xplan(i+1),zplan(i+1))，第i+1層入射角λi+1，第i層大氣折射率ni，第i+1層大氣折射率ni+1，每層的地心夾角dθ已知。光線在每層中傳輸路徑的計(jì)算原理如下：

圖4 地心坐標(biāo)系與平面坐標(biāo)系之間的關(guān)系圖Fig.4 Relationship between geocentric coordinate system and plane coordinate system

圖5 相機(jī)、目標(biāo)和地心組成的平面內(nèi)，光在每層大氣中的傳輸路徑圖Fig.5 In plane composed of camera, target and center of the earth, transmission path of light in each layer of atmosphere

根據(jù)大氣球面分層Snell定律，由入射角λi+1、大氣折射率ni、ni+1，求得第i+1層出射角φi+1。

大氣折射率采用Hopfield模型進(jìn)行計(jì)算，利用該模型可以求取任意海拔高度h的大氣折射率n=N(h)，Hopfield模型考慮了大氣含有水汽的情況，將大氣折射率分為干、濕2項(xiàng)[16-17]。

根據(jù)第i+1層入射點(diǎn)法線向量和出射角計(jì)算出射方向的單位向量。由第i+1層入射點(diǎn)法線向量Qi+1=[xplan(i+1),zplan(i+1)]和出射角φi+1，利用向量?jī)?nèi)積求得出射方向的單位向量ue=[ux,uz]，具體求解方法為

(5)

根據(jù)第i+1層的地心夾角、出射角、其距地心的距離以及出射方向的單位向量，計(jì)算該層出射光線的長(zhǎng)度，進(jìn)而得到其出射光線方向向量，即第i+2層入射光線方向向量。在△OTi+1Ti+2中，根據(jù)正弦定理，利用公式：

(6)

求得第i+1層出射光線的長(zhǎng)度Li+1，進(jìn)而求得第i+1層出射光線方向向量Li+1ue，該方向向量即為第i+2層的入射光線方向向量vi+2，其中第i+1層的地心夾角為dθ、第i+1層距地心的距離為Hi+l。

由第i+1層入射點(diǎn)法線向量和第i+2層入射光線方向向量，計(jì)算得到第i+2層光線入射點(diǎn)法線向量。進(jìn)而根據(jù)向量間的關(guān)系Qi+2-Qi+1=Li+1ue，求得第i+2層光線入射點(diǎn)Ti+2法線向量為

Qi+2=Qi+1+Li+1ue=(xplan(i+2),zplan(i+2))

(7)

根據(jù)第i+1層光線的入射點(diǎn)坐標(biāo)和入射點(diǎn)向量計(jì)算第i+1層光線的交點(diǎn)；判斷是否達(dá)到最后一層。在每層中求解射線與發(fā)射面的交點(diǎn)，已知發(fā)射點(diǎn)F(Xf,Yf,Zf)和發(fā)射面法向量fαj=[μαj,ναj,γαj]，可使用點(diǎn)法式求解交點(diǎn)K(xk,yk,zk)，當(dāng)?shù)趇+1層光線的入射點(diǎn)Ti+1和第i+1層中交點(diǎn)Ki+l之間距離小于第i+1層光線的入射點(diǎn)Ti+1和第i+2層光線的入射點(diǎn)Ti+2之間距離時(shí)，交點(diǎn)Ki+l(xk,yk,zk)即為目標(biāo)點(diǎn)。

3 目標(biāo)運(yùn)動(dòng)參數(shù)估計(jì)

(8)

(9)

在用二階多項(xiàng)式擬合時(shí)，為了保證足夠高的擬合精度，采用分段擬合的方式來(lái)減小擬合誤差。圖6為對(duì)實(shí)際軌跡采用上述擬合方式時(shí)，位置平均誤差和徑向速度平均誤差曲線。

考慮到實(shí)際檢測(cè)數(shù)據(jù)點(diǎn)數(shù)，僅僅比對(duì)兩者擬合點(diǎn)數(shù)在10～30范圍內(nèi)的情況?？梢钥闯?，在擬合點(diǎn)數(shù)為10～30范圍內(nèi)，位置平均誤差隨擬合點(diǎn)數(shù)增加而增加，而徑向速度平均誤差隨擬合點(diǎn)數(shù)增加而減小。發(fā)生上述情況的原因是擬合時(shí)3個(gè)分量各自分別擬合，而實(shí)際中軌跡的3個(gè)分量存在一定的耦合，隨著點(diǎn)數(shù)增多位置擬合誤差會(huì)增大，而擬合點(diǎn)數(shù)增多時(shí)擬合出的曲線形狀與實(shí)際更相似，所以速度誤差會(huì)降低。取折中的情況，采用擬合點(diǎn)數(shù)為20的分段擬合。

圖6 位置及徑向速度平均誤差隨擬合點(diǎn)數(shù)變化曲線Fig.6 Curves of location and radial velocity average error changing with number of fitting points

4 估計(jì)參數(shù)的迭代優(yōu)化

先驗(yàn)給出的發(fā)射面參數(shù)經(jīng)常存在誤差，且該誤差對(duì)重建精度有極大影響，因此需要在重建過程中同時(shí)優(yōu)化發(fā)射面參數(shù)，以提高重建精度。目標(biāo)運(yùn)動(dòng)在短時(shí)間內(nèi)可近似為勻加速直線運(yùn)動(dòng)，發(fā)射面準(zhǔn)確時(shí)即理想情況下，所重建的運(yùn)動(dòng)加速度應(yīng)為一常量。若發(fā)射面參數(shù)不準(zhǔn)確即發(fā)射面發(fā)生偏轉(zhuǎn)，目標(biāo)將變?yōu)樽兗铀龠\(yùn)動(dòng)。故在對(duì)重建目標(biāo)點(diǎn)進(jìn)行分段擬合時(shí)，理論上每一段的加速度較之于前一段的加速度差值其絕對(duì)值和應(yīng)接近于零，而對(duì)于變加速運(yùn)動(dòng)，發(fā)射面偏轉(zhuǎn)越大，差值越大。

(10)

對(duì)應(yīng)于發(fā)射面偏角αj，可反演出相應(yīng)的目標(biāo)運(yùn)動(dòng)參數(shù)(x,y,z,v,a,da)aj，分別為位置坐標(biāo)(x,y,z)，速度v，加速度a和加速度差da，則最優(yōu)目標(biāo)運(yùn)動(dòng)參數(shù)為

(11)

圖7為發(fā)射點(diǎn)經(jīng)緯高為(E97.306°,N41.090°,6 369.853 km)，發(fā)射面法向量為(-0.137,0.604,0.784)的軌跡在不同偏角情況下的反演加速度差值和的變化曲線，記錄加速度變化平均值的最小值以及相對(duì)應(yīng)的發(fā)射面偏角、發(fā)射點(diǎn)位置、發(fā)射面法向量。將此條件下得到的軌跡點(diǎn)作為最后的反演軌跡。

圖7 發(fā)射面優(yōu)化分析圖Fig.7 Diagram of launch plane optimization analysis

5 實(shí)驗(yàn)及分析

在發(fā)射點(diǎn)和發(fā)射面等先驗(yàn)信息準(zhǔn)確情形下，不考慮成像和檢測(cè)的影響因素，單平臺(tái)定位的位置平均誤差標(biāo)準(zhǔn)差和徑向速度平均誤差與指向誤差標(biāo)準(zhǔn)差的關(guān)系如圖8、圖9所示。

從圖8可以看出，單平臺(tái)定位時(shí)，分段擬合的方法在指向誤差極小的情形下會(huì)帶來(lái)誤差增大，但是指向誤差標(biāo)準(zhǔn)差在5″內(nèi)時(shí)能夠有效地抑制指向誤差引起的定位誤差，并且誤差在500 m以下，有足夠空間考慮成像和檢測(cè)的誤差。圖9顯示指向誤差標(biāo)準(zhǔn)差變化范圍為(0″,50″)時(shí)，徑向速度平均誤差與指向誤差的標(biāo)準(zhǔn)差成正相關(guān)，徑向速度平均誤差小于60 m/s。

圖8 單平臺(tái)定位的位置平均誤差隨指向誤差標(biāo)準(zhǔn)差變化曲線Fig.8 Curves of single platform positioning’s location average error changing with pointing error standard deviation

圖9 單平臺(tái)定位的徑向速度平均誤差隨指向誤差標(biāo)準(zhǔn)差變化曲線Fig.9 Curve of single platform positioning’s radial velocity average error changing with pointing error standard deviation

由單平臺(tái)下檢測(cè)到的目標(biāo)像面位置重建三維軌跡是典型的病態(tài)問題，先驗(yàn)知識(shí)發(fā)射點(diǎn)位置和發(fā)射面指向顯的尤為重要，其中發(fā)射點(diǎn)位置能夠精確獲得，故發(fā)射面指向是影響單平臺(tái)定位誤差的關(guān)鍵因素。

從圖10可以看出,即使是根據(jù)已知的準(zhǔn)確軌跡擬合出來(lái)的發(fā)射面都不是十分準(zhǔn)確，而在實(shí)際目標(biāo)飛行過程中，受到環(huán)境影響因素太多，發(fā)射面不能保持準(zhǔn)確，并會(huì)產(chǎn)生很大的誤差。從圖中可以看出在發(fā)射面偏轉(zhuǎn)500″，也就是0.139°，在實(shí)際中可看做是較小的誤差，但是目標(biāo)重建位置誤差已經(jīng)在1 km以上。故單平臺(tái)定位在實(shí)際應(yīng)用中存在很大的難題也就是發(fā)射面準(zhǔn)確度的問題，其定位誤差遠(yuǎn)大于雙平臺(tái)定位誤差。

在空間目標(biāo)點(diǎn)位置確定的情況下，影響定位精度的另一個(gè)重要因素就是平臺(tái)位置誤差，在高空中平臺(tái)位置的變化給定位帶來(lái)很大的影響，分別對(duì)平臺(tái)位置誤差添加均值為0，不同標(biāo)準(zhǔn)差的隨機(jī)誤差，統(tǒng)計(jì)分析平臺(tái)位置誤差給定位精度帶來(lái)的影響，圖11為平臺(tái)位置誤差給目標(biāo)重建位置和速度帶來(lái)的誤差。

從圖11可以看出平臺(tái)位置相較于其他影響因素，對(duì)目標(biāo)重建速度誤差帶來(lái)的誤差影響是很大的。在平臺(tái)位置誤差為420 m時(shí)，重建速度誤差已經(jīng)比較大，而重建位置誤差符合定位要求。

圖10 單平臺(tái)定位的位置誤差隨發(fā)射面誤差變化曲線Fig.10 Curve of single platform positioning’s location error changing with launch plane error

圖11 平臺(tái)位置誤差對(duì)目標(biāo)重建位置誤差及速度誤差的影響Fig.11 Effect of platform location error on target reconstruction location error and speed error

6 結(jié) 論

本文針對(duì)單平臺(tái)下的觀測(cè)數(shù)據(jù)進(jìn)行空間軌跡重建是典型的病態(tài)問題，提出單平臺(tái)成像目標(biāo)運(yùn)動(dòng)特性反演方法，得到：

1) 本文提出從平臺(tái)位置出發(fā)，反向追蹤，將大氣劃分多層，計(jì)算在大氣傳輸路徑，最終與發(fā)射面相交得到目標(biāo)位置，并通過迭代優(yōu)化發(fā)射面和分段擬合等方法抑制重建誤差，此方法能夠有效反演單平臺(tái)目標(biāo)運(yùn)動(dòng)特性。

2) 通過本文提出的單平臺(tái)目標(biāo)運(yùn)動(dòng)特性反演方法，目標(biāo)反演位置誤差在200 m以下，速度誤差在60 m/s以下，能夠用作實(shí)際場(chǎng)景。

3) 對(duì)重建出的原始數(shù)據(jù)進(jìn)行分段擬合，能夠有效的抑制重建的位置誤差及速度誤差。

4) 在各影響因素中，發(fā)射面的精確度的影響最大，較小的偏轉(zhuǎn)就會(huì)帶來(lái)很大的反演誤差，而通過發(fā)射面的迭代優(yōu)化，可以較為有效地抑制小范圍發(fā)射面的偏轉(zhuǎn)帶來(lái)誤差。