王歆 李鑫冉

(1中國科學(xué)院紫金山天文臺 南京 210008) (2中國科學(xué)院空間目標(biāo)與碎片重點(diǎn)實(shí)驗(yàn)室 南京 210008) (3中國科學(xué)技術(shù)大學(xué) 合肥 230026)

基于核密度估計(jì)的極短弧定軌的分布估計(jì)方法?

王歆1,2?李鑫冉1,2,3

(1中國科學(xué)院紫金山天文臺 南京 210008) (2中國科學(xué)院空間目標(biāo)與碎片重點(diǎn)實(shí)驗(yàn)室 南京 210008) (3中國科學(xué)技術(shù)大學(xué) 合肥 230026)

采用分布估計(jì)算法,通過建立解空間的概率模型,實(shí)現(xiàn)了一種測角資料的極短弧初軌計(jì)算方法.在概率模型建立中,采用非參數(shù)核密度估計(jì),無需對分布進(jìn)行任何假設(shè).不同于遺傳算法、粒子群算法等進(jìn)化算法,方法不僅考慮解的適值優(yōu)劣,同時考慮了解的整體特性.根據(jù)基于實(shí)測數(shù)據(jù)的數(shù)值計(jì)算表明:方法在沒有任何約束條件情況下對于一般觀測精度仍可獲得有效解.

航天器,天體力學(xué):定軌,方法:數(shù)值,方法:統(tǒng)計(jì)

1 引言

空間目標(biāo)編目的主要目的是了解在軌的工作衛(wèi)星和空間碎片的準(zhǔn)確軌道,從而可推算它們的位置,有效地支持各類航天活動.為了能夠獲取準(zhǔn)確軌道必須獲得足夠弧長和數(shù)量的觀測資料,然而隨著空間目標(biāo)數(shù)量的不斷增加,采用面向目標(biāo)的跟蹤觀測方式已不能滿足需求.近年來空間目標(biāo)觀測逐步轉(zhuǎn)變?yōu)槊嫦蚩沼虻挠^測模式,目前我國低軌道空間目標(biāo)光學(xué)測角資料超過80%來自于這種觀測模式,是我國自主編目的主要數(shù)據(jù)來源.這種方式在提高觀測量的同時帶來了新的問題,對于每個目標(biāo)每次采集的弧長是非常短的,往往只有數(shù)十秒,甚至十幾秒,這么短的資料是不足以獲取足夠準(zhǔn)確的軌道的,需要將許多短弧段關(guān)聯(lián)在一起共同實(shí)現(xiàn)準(zhǔn)確的軌道確定.通過極短弧定軌獲取弧段有效信息、提高關(guān)聯(lián)效率就成為一項(xiàng)重要的研究內(nèi)容.即使在跟蹤模式下,如能以較短弧長獲得有效軌道顯然也是有助于提高觀測效率的[1].

短弧定軌問題是一個經(jīng)典問題,主要困難在于問題本身的病態(tài)性,弧段越短這種病態(tài)性越突出[2].文獻(xiàn)[3]將遺傳算法(Genetic A lgorithm,GA)用于解決天基極短弧定軌問題,文獻(xiàn)[4]采用了相同方法研究了一般的地基短弧定軌問題,但在GA算子選擇上和文獻(xiàn)[3]不同,這些工作克服了經(jīng)典方法不能獲得合理結(jié)果的困難,但解的偏差和彌散度仍相當(dāng)大.文獻(xiàn)[5-6]采用完全不同的優(yōu)化變量、適值函數(shù)和算子給出了一種計(jì)算方案,在一定約束條件下縮小了解的偏差和彌散度,初步滿足了后續(xù)應(yīng)用的需求.文獻(xiàn)[7]初步探索了采用粒子群算法(Particle Swarm Optim ization,PSO)求解極短弧定軌問題的方法,得到了和文獻(xiàn)[5-6]相似的結(jié)果.上述工作雖然各自不同,但都采用了進(jìn)化計(jì)算(Evolutionary A lgorithm),文獻(xiàn)[8]則針對進(jìn)化計(jì)算不同于經(jīng)典方法的計(jì)算流程,提出了具有一般意義的資料處理策略,通過穩(wěn)健估計(jì)實(shí)現(xiàn)了極高崩潰點(diǎn)的定軌.這些工作初步建立起了采用進(jìn)化計(jì)算解決極短弧定軌的框架,但都需要針對具體問題進(jìn)行方法、算子以及參數(shù)和約束條件的選擇,而對于編目而言,面對的是大批量目標(biāo),因此過多的約束和參數(shù)并不是最佳選擇.

本文從另一個角度做了探索,初軌計(jì)算的實(shí)踐告訴我們,對于初軌計(jì)算而言不能以殘差均方根(RMS)最小作為唯一判別.在精度范圍內(nèi),短弧定軌問題是多解的,這些解的RMS十分接近,而且往往RMS最小的解并不是實(shí)際最優(yōu)解[9].之前所有工作都著眼于尋找RMS最小的單個解,而忽視了解的整體情況.本文引入了分布估計(jì)算法(Estimation of Distribution A lgorithm,EDA)[10],EDA是目前國際上進(jìn)化計(jì)算研究的熱點(diǎn)之一,與GA、PSO的不同在于,EDA立足于對優(yōu)勢解群體的整體刻畫,而并不局限于個別優(yōu)勢個體.

EDA僅是一個算法框架,本文基于核密度估計(jì)(Kernel Density Estim ation)實(shí)現(xiàn)了一種具體的算法,采用非參數(shù)方法,無需對分布進(jìn)行任何假設(shè),構(gòu)造了具體的極短弧定軌方法,并結(jié)合中國科學(xué)院光學(xué)觀測網(wǎng)的實(shí)際數(shù)據(jù)情況做了計(jì)算驗(yàn)證.

2 分布估計(jì)算法

EDA是進(jìn)化計(jì)算領(lǐng)域新興起的一類隨機(jī)優(yōu)化算法,是當(dāng)前國際進(jìn)化計(jì)算領(lǐng)域的研究熱點(diǎn),其概念最早在1996年提出并得到了迅速發(fā)展.EDA是GA和統(tǒng)計(jì)學(xué)習(xí)的結(jié)合,通過統(tǒng)計(jì)學(xué)習(xí)的手段建立解空間內(nèi)個體分布的概率模型,然后對概率模型隨機(jī)采樣產(chǎn)生新的群體,如此反復(fù)進(jìn)行,實(shí)現(xiàn)群體的進(jìn)化.EDA中沒有傳統(tǒng)的交叉、變異等遺傳操作,是一種全新的進(jìn)化模式[11].和其他方法一樣,EDA用于解決如下優(yōu)化問題:

其中x=(x1,x2,···,xn)為實(shí)數(shù)域R上的n維解向量.和分別為第i維的上下界. EDA的流程如圖1.EDA首先隨機(jī)生成N個解作為初始種群XG={x1,x2,···,xN}, G=0為進(jìn)化代數(shù);然后依據(jù)適值優(yōu)劣選擇P個優(yōu)勢群體;接著根據(jù)優(yōu)勢群體構(gòu)造概率模型;由概率模型通過采樣,生成G+1代群體XG+1,如此迭代進(jìn)化.與GA相比,EDA方法在選擇優(yōu)勢群體后,對于優(yōu)勢群體內(nèi)個體在建立概率模型時是等權(quán)的,因此概率模型中考慮了優(yōu)勢群體的數(shù)密度,是解空間整體特性和適值優(yōu)劣的綜合反映,這個特點(diǎn)非常適合極短弧定軌問題.不同的EDA實(shí)現(xiàn)主要區(qū)別在于概率模型建立的方法不同.

與GA一樣,由于采用二進(jìn)制編碼比較利于概率模型的建立,EDA起初也是針對二進(jìn)制編碼在離散域上的優(yōu)化問題,后來在此基礎(chǔ)上發(fā)展到采用實(shí)數(shù)編碼處理連續(xù)域問題.由于連續(xù)域概率模型的復(fù)雜性給設(shè)計(jì)有效的分布估計(jì)算法增加了難度,因此連續(xù)域EDA的發(fā)展相對緩慢.目前在連續(xù)域概率模型建立時基本都采用了Gauss分布假設(shè),由于這種處理對于許多工程問題并不合適,在具體問題中直方圖分布也常被使用,但直方圖分布不連續(xù)光滑,也有較大局限性[12].

圖1 分布估計(jì)算法的流程圖Fig.1 T he flow chart of EDA

3 初軌計(jì)算的分布估計(jì)方法

3.1 變量和編碼的選擇

文獻(xiàn)[3]選擇觀測首末時刻的斜距作為優(yōu)選變量,這種選擇使優(yōu)選變量需要和首末觀測資料共同得到軌道狀態(tài)量,特別是當(dāng)首末資料出現(xiàn)較大誤差時,無法得到合理結(jié)果.文獻(xiàn)[5-6]針對這點(diǎn)選擇了和觀測量不相關(guān)的歷元時刻t0的3個Kepler根數(shù)(a,e,M0)作為優(yōu)選變量,在保持參數(shù)空間維數(shù)較低的同時,使得結(jié)果不依賴于具體的某個觀測量.這種選擇對于近圓軌道是合適的,但當(dāng)偏心率稍大時,每個優(yōu)化變量都符合上下界約束條件的情況下,仍會出現(xiàn)a(1?e)

可見其中x1為近地點(diǎn)高度,確定了x1的下界,每個變量都符合上下界條件時,解都是合理的軌道.x和(a,e,M0)的逆變換是簡單的.

對于定軌問題采用二進(jìn)制編碼使得優(yōu)化空間過大,因此仍然采用了實(shí)數(shù)編碼.

3.2 適值函數(shù)

采用和文獻(xiàn)[5-6]相同的適值函數(shù),令觀測量為{ti,αi,δi,i=1,2,···,N},其中αi和δi分別為赤經(jīng)和赤緯.由歷元時刻t0的x經(jīng)過簡單轉(zhuǎn)換即可得到相應(yīng)的a、e、M0,任意時刻ti的平近點(diǎn)角M為:

此時目標(biāo)斜距ρi為:

其中Ri=|Ri|,Ri為測站的地心位置矢量,Zi為目標(biāo)天頂距.由測量方程可得:

對任意兩個觀測時刻tk、tj,tk>tj,分別得到(rk,rj),定義

適值函數(shù)為:

3.3 初值的生成

初始種群生成和一般進(jìn)化計(jì)算采用的方法完全相同,對每個維在上下界取值范圍內(nèi)根據(jù)均勻分布隨機(jī)選取.盡管取值范圍并不影響計(jì)算的開展,但根據(jù)已知信息減小取值范圍縮小求解空間顯然是有益的.考慮到極短弧軌道計(jì)算中涉及觀測時間很短,目標(biāo)的地心距變化不會太大,根據(jù)圓軌道假設(shè),任選兩個觀測時刻tk、tj,令rj=rk=r,由(4)–(7)式可得?fkj,將

作為目標(biāo)函數(shù)可優(yōu)選得到r,這是個一維優(yōu)化問題[13].顯然x1不能大于r,同時考慮空間目標(biāo)最低軌道高度情況,可確定x1的取值范圍為[1.03ae,r].事實(shí)上遠(yuǎn)地點(diǎn)地心距必須大于r,但考慮到計(jì)算誤差,實(shí)際取值時需要對近地點(diǎn)略微放大一些,而對于遠(yuǎn)地點(diǎn)略微縮小一些,因此沒有直接采用a(1?e)和a(1+e)作為優(yōu)選變量,避免了遠(yuǎn)地點(diǎn)近于近地點(diǎn)的不可行情況.這個步驟僅是參考,所以計(jì)算要求不高,一般選取首尾觀測資料計(jì)算,如考慮穩(wěn)健性因素可對每對觀測資料做同樣的優(yōu)選,根據(jù)得到的全部r選取參考值.

x2的取值范圍對近地空間目標(biāo)取[0,4ae]已涵蓋了全部范圍,x3∈[0,2π).

3.4 概率模型的建立方法

EDA中概率模型是根據(jù)優(yōu)勢群體建立的,確定優(yōu)勢群體規(guī)模后可采用GA的各種選擇方法,由于EDA對于優(yōu)勢群體內(nèi)個體是同等對待的,多采用截?cái)噙x擇或錦標(biāo)賽選擇.這兩種方法都只涉及適值大小比較,無需對適值進(jìn)行標(biāo)定,算法大為簡化.采用截?cái)噙x擇,選擇適值最好的P個個體作為優(yōu)勢群體.利用第G代的優(yōu)勢群體XGE進(jìn)行概率建模.

建模中將各維作為獨(dú)立變量考慮,因此問題轉(zhuǎn)換為3個1維建模問題.令XGE中某維數(shù)據(jù)為{y1,y2,···,yP}.不作任何分布假設(shè),采用核密度估計(jì)方法可給出任意y的概率密度D(y)[14]:

其中h為帶寬,K(·)為核函數(shù),采用常用的Epanechnikov核函數(shù)[15]:

(9)式即概率模型,令從第G代得到的優(yōu)勢群體建立的概率模型記為DEG,在概率模型更新過程中采用阻尼方式,第G+1代的概率模型為:

其中α為學(xué)習(xí)效率,反映了新信息的吸收率,也是舊信息的遺忘率.

對于核密度估計(jì)而言,帶寬h的選擇至關(guān)重要,在EDA中概率密度函數(shù)是不斷收斂的,因此不宜采用固定帶寬,根據(jù)Silverman大拇指法則(rule of thumb)動態(tài)選擇帶寬:

其中σ為建立核密度估計(jì)樣本的標(biāo)準(zhǔn)差.

3.5 概率模型的采樣方法

從概率模型DG中采樣已無實(shí)質(zhì)困難,對于核密度估計(jì),Silverman給出了一種無需計(jì)算密度本身的采樣方法[16],令數(shù)據(jù)樣本為{y1,y2,···,yP},采樣步驟如下:(1)從[1,P]中均勻選擇一個隨機(jī)整數(shù)v;(2)從K(·)中隨機(jī)采樣得到?;(3)采樣z=yv+?h.重復(fù)上述步驟即可得到任意樣本容量的采樣,新樣本和原樣本的概率密度函數(shù)一致.上述方法將一個非參數(shù)分布的采樣問題,轉(zhuǎn)換為了一個解析的參數(shù)分布采樣問題,僅需要從核函數(shù)采樣.

對于Epanechnikov核,有一個極其快速的采樣方法[16]:生成3個[0,1)均勻分布的隨機(jī)數(shù)v1、v2和v3,則采樣?為:

將(13)式代入上述采樣過程,可快速實(shí)現(xiàn)從任意樣本的核密度估計(jì)中采樣.

3.6 終止條件

進(jìn)化計(jì)算方法一般采用最大進(jìn)化代數(shù)Gmax和連續(xù)C代適值沒有進(jìn)化作為終止條件, EDA也可采用同樣的終止條件.考慮到EDA是對解的整體刻畫,因此考慮概率密度收斂作為終止條件,由(12)式可知,當(dāng)樣本σ很小時,h趨于0,從而解不再更新.σ反映了概率密度的收斂,以樣本數(shù)據(jù)在3個維度上的標(biāo)準(zhǔn)差之和σ1+σ2+σ3

EDA在進(jìn)化過程中和GA一樣對父代是破壞性的,GA中采用精英保留策略使得適值最好的解始終進(jìn)入子代.在EDA的計(jì)算過程中記錄下適值最好的結(jié)果x best作為最優(yōu)解;由于EDA考慮了解的整體特征,當(dāng)計(jì)算終止時,σ趨于0,說明解都集中在一點(diǎn),記這個解為x prob.對于每次計(jì)算同時記錄下這兩個解.

3.7 (i,?,ω)的求解

由優(yōu)選得到歷元t0時刻的(a,e,M0)可求得ti時刻的地心位置矢量ri,在此基礎(chǔ)上求解完整的6個軌道根數(shù)已無任何困難.考慮到極短弧定軌精度要求不高,由每對(ri,rj)可得到一組(i,?,ω),取其平均(,,)作為結(jié)果,考慮穩(wěn)健性用中值替代均值作為結(jié)果.

4 數(shù)值驗(yàn)證

采用MATLAB1h ttp://www.m athw orks.com編寫了程序,并選取中國科學(xué)院空間目標(biāo)與碎片光學(xué)觀測網(wǎng)中的一圈實(shí)測資料進(jìn)行了數(shù)值驗(yàn)證,設(shè)備精度優(yōu)于5′′,事后驗(yàn)證該圈實(shí)際精度1′、采樣頻率1 Hz,采用了前10 s弧段做計(jì)算驗(yàn)證.根據(jù)圓軌道假設(shè)計(jì)算得到該弧段的地心距約為1.13ae,因此優(yōu)化變量的取值范圍為x1∈[1.03ae,1.15ae]、x2∈[0,4.0ae]以及x3∈[0,2π],對于a和e并未利用任何先驗(yàn)信息.由于EDA需要對概率密度進(jìn)行估計(jì),過少數(shù)量不利于模型的建立,因此種群大小大于其他進(jìn)化算法,參數(shù)取N=100、P=30、α=0.1、Tσ=10?6.

圖2給出了一次計(jì)算中適值的收斂過程,可見算法全局搜索能力強(qiáng),進(jìn)化開始收斂很快,30代后適值已非常接近最終結(jié)果.

初軌計(jì)算中(a,e)是最為關(guān)心的,將(x1,x2)轉(zhuǎn)換為(a,e)后,圖3給出了概率密度的變化情況.圖中顏色越紅表示概率密度越大,每行采用了相同的尺度,行與行之間尺度不同.圖中可看出進(jìn)化開始,解集范圍很大,峰值也不明顯;隨著進(jìn)化過程,解集范圍不斷縮小,峰值也越來越顯著.

與其他進(jìn)化計(jì)算一樣,隨機(jī)數(shù)對結(jié)果會產(chǎn)生影響,通過MATLAB中rng(shuffle)命令采用不同隨機(jī)數(shù)種子計(jì)算10次,結(jié)果列于表1,對于每次計(jì)算依次給出了x best和x prob,分別用A和B表示.表2給出了這些結(jié)果的統(tǒng)計(jì)信息,包括全部結(jié)果的均值(M ean)、中值(M ed)和標(biāo)準(zhǔn)差(σ),以及A組和B組的分組結(jié)果,POD表示由多天多站精密定軌得到的結(jié)果作為參考標(biāo)準(zhǔn).

表1 不同隨機(jī)數(shù)種子下的定軌結(jié)果Tab le 1 The resu lts of o rb it determ ination with d ifferen t random seed s

從表1中結(jié)果可以看出:每次計(jì)算結(jié)果雖不同,同一次計(jì)算結(jié)果中A組和B組的結(jié)果也并不完全一致,但RMS都很小并且非常接近,從尋優(yōu)角度可視作等價,隨機(jī)數(shù)生成對計(jì)算并無實(shí)質(zhì)影響,本文方法是有效的.獲得的軌道結(jié)果不同則反映了極短弧定軌本身的困難.從表2的統(tǒng)計(jì)結(jié)果看,各組解的均值或中值非常接近,也都能滿足需求,和文獻(xiàn)[5-6]相比,在減少約束條件情況下,解的偏差更小,標(biāo)準(zhǔn)差相當(dāng),說明關(guān)于解的整體刻畫對于極短弧定軌有很好的效果,在實(shí)踐中可采用多次計(jì)算的均值或中值作為參考解.

對于極短弧定軌而言,測量誤差影響更為顯著,文獻(xiàn)[3]強(qiáng)調(diào)了極高精度觀測對于極短弧定軌的必要性,即使是2′′的觀測誤差,對結(jié)果造成的影響也是巨大的.文獻(xiàn)[5-6]采用參數(shù)化的Bootstrap方法[17]考察了測量誤差對求解的影響,采用相同方法,根據(jù)設(shè)備精度對觀測資料增加5′′的隨機(jī)誤差后采用本文方法定軌,以不同隨機(jī)數(shù)種子重復(fù)50次,解的統(tǒng)計(jì)結(jié)果見表3.

表2 定軌的統(tǒng)計(jì)結(jié)果Tab le 2 The statistical resu lts of o rb it determ ination

圖2 適值的收斂過程Fig.2 T he convergence p rogress of the fitness

從表3中結(jié)果可見:和文獻(xiàn)[3]一樣,觀測精度下降后,定軌受到大幅度影響.從中值和均值看,A組受到很大影響,而B組影響較小.說明以概率密度收斂位置作為解比適值最小解對于實(shí)際問題而言更優(yōu),也再次說明考察解的整體情況是十分有利的.圖4給出了初軌計(jì)算最為關(guān)心的a的概率密度圖,包括了全體數(shù)據(jù)以及A組和B組的結(jié)果.圖中可見都是重尾形狀,解釋了表3中中值結(jié)果優(yōu)于均值結(jié)果的原因.兩組解的彌散程度都大幅度增加,相比較而言,B組相對集中,但兩組的概率密度函數(shù)峰值的位置卻很接近,分別為7228.928274 km和7218.417037 km,都在可接受范圍內(nèi),和沒有增加誤差情況下相當(dāng).高精度觀測對極短弧定軌顯然是十分有利的,采用本文方法在一般精度情況下仍可獲取有效解,以概率密度最大作為評估選取標(biāo)準(zhǔn)解則更為穩(wěn)健.

圖3 概率密度收斂過程Fig.3 T he convergence p rogress of the p robab ility density

圖4 半長徑a的概率密度Fig.4 The p robability density of sem i-m a jor ax is a

表3 增加隨機(jī)誤差后的結(jié)果Tab le 3 T he resu lts of orb it d eterm ination with rand om noises

5 結(jié)論與討論

綜上所述,本文方法和其他進(jìn)化計(jì)算方法類似,可以克服經(jīng)典方法對于極短弧定軌的困難,相比于GA、PSO等進(jìn)化計(jì)算方法,在減少約束條件、一般觀測精度的情況下仍可獲取令人滿意的結(jié)果,主要得益于算法不僅考慮解的適值大小,同時考慮了解的整體分布情況.不同于其他進(jìn)化算法,算法本質(zhì)上不再尋找單個解,而是尋找解的區(qū)域,不斷縮小解的區(qū)域最終獲得結(jié)果,直觀上看,就是尋找優(yōu)勢解最為稠密的區(qū)域.種群大小是決定進(jìn)化計(jì)算效率的主要因素之一,EDA采用的種群數(shù)量要大于其他進(jìn)化算法,因此計(jì)算效率較低,然而對于極短弧而言,通常觀測數(shù)據(jù)不多,在當(dāng)前計(jì)算能力下已可滿足實(shí)踐應(yīng)用的需求.本文不僅是引入一種不同的計(jì)算方法,更重要的是從優(yōu)勢解稠密性去考慮極短弧定軌問題,拓展了研究思路.

參考文獻(xiàn)

[1]Bennett J C,Sang J,Sm ith C,et a l.A dSpR,2015,55:617

[2]吳連大,賈沛璋.天文學(xué)報(bào),1997,38:288

[3]Ansalone L,Cu rti F.Ad SpR,2013,52:477

[4]H inagaw a H,Yam aoka H,Hanada T.A d SpR,2014,53:532

[5]李鑫冉,王歆.天文學(xué)報(bào),2016,57:66

[6]Li X R,W ang X.ChA&A,2017,41:undeterm ined

[7]李鑫冉,王歆.飛行器測控學(xué)報(bào),2015,34:540

[8]李鑫冉,王歆.天文學(xué)報(bào),2016,57:181

[9]賈沛璋,吳連大.天文學(xué)報(bào),1998,39:337

[10]Larra?naga P,Lozano J A.Estim ation of D istribu tion A lgorithm s:A New Too l for Evo lu tionary Com pu tation.New York:Sp ringer Science+Business M ed ia,2002:3-26

[11]周樹德,孫增圻.自動化學(xué)報(bào),2007,33:113

[12]王圣堯,王凌,方晨,等.控制與決策,2012,27:961

[13]吳連大.人造衛(wèi)星和空間碎片的軌道和探測.北京:中國科學(xué)技術(shù)出版社,2011:158-191

[14]Parzen E.Anna ls of M athem atica l Statistics,1962,33:1065

[15]Epanechn ikov V A.Theory of P robab ility and Its A pp lications,1969,14:153

[16]Silverm an B W.Density Estim ation for Statistics and Data Analysis.London:Chapm an and Hall, 1986

[17]E fron B,T ibshiran i R.A n In troduction to the Bootstrap.New York:Chapm an&Ha ll,1993:53-35

Estim ation of D istribu tion A lgorithm for Initial O rbit Determ ination of Too-short-arc Based on K ernel Density Estim ation

WANG Xin1,2LIXin-ran1,2,3

(1 Pu rp le M oun tain Observatory,Chinese Academ y of Scien ces,Nan jing 210008) (2 K ey Labo rato ry fo r Space O bjec t and D ebris O bserva tion,Pu rp le M oun tain O bserva to ry,Chinese Academ y of Scien ces,Nan jing 210008) (3 Un iversity of Scien ce and Techno logy of China,Hefei 230026)

A new approach of the initial orbit determ ination for too-short-arc with angular measurements is im plemented by building the probabilistic model in the solution space with the estimation of distribution algorithm.without any assumption of distribution,the non-parametric kernel density estimation is employed in themodel building.Themethod,unlike other evolutionary algorithms,such as genetic algorithm and particle swarm optim ization,considers the fitness as well as the characteristic of solution space.Numerical check with real observations indicates that without any constraints,the proposed technique has a good performance for themeasurements of general accuracy.

space vehicles,celestialmechanics:orbit determ ination,methods:numerical,methods:statistical

P135;

A

10.15940/j.cnki.0001-5245.2016.06.005

2016-01-27收到原稿,2016-02-29收到修改稿

?國家自然科學(xué)基金項(xiàng)目(11373072)資助

?wangxin@pm o.ac.cn