快速最大似然估計用于脈沖星時間延遲測量

江若玫,佘 顏

(1.西安交通大學城市學院,陜西西安 710018; 2.西安電子科技大學微電子學院,陜西西安 710071)

快速最大似然估計用于脈沖星時間延遲測量

江若玫1,佘 顏2

(1.西安交通大學城市學院,陜西西安 710018; 2.西安電子科技大學微電子學院,陜西西安 710071)

針對X射線脈沖星時間延遲測量方法存在的測量精度不高、計算量大等問題,提出了一種基于快速最大似然估計時間延遲測量方法.首先通過歷元折疊提取脈沖星脈沖輪廓,在此基礎上建立似然函數(shù),進而通過求解最大似然函數(shù)得到脈沖輪廓時間延遲量;通過縮小標準脈沖輪廓相位間隔對快速最大似然估計進行改進,采用粗略估計與精細估計相結(jié)合的方法精確計算時間延遲量.以脈沖星實測數(shù)據(jù)為基礎,從觀測時間和運算量對快速最大似然方法與其他方法進行仿真對比實驗,結(jié)果表明,文中方法在保證運算效率基本不變的前提下,明顯提高了時間延遲測量的精確度.

X射線脈沖星;脈沖輪廓;歷元折疊;時間延遲量;快速最大似然估計

由于脈沖星信號具有高度的周期穩(wěn)定性和不易被干擾等特性,使得脈沖星導航成為空間導航領域的研究熱點.在基于X射線脈沖星導航系統(tǒng)(X-ray Pulsar-based NAVigation,XPNAV)中,脈沖光子從航天器到太陽系質(zhì)心的到達時間差(Time Difference Of Arrival,TDOA)是XPNAV系統(tǒng)的基本觀測量,其測量精度決定了整個XPNAV系統(tǒng)的精度,從而決定了航天器測速、定位和定時的精度[1-2].TDOA的測量精度主要取決于脈沖星測量脈沖輪廓與標準脈沖輪廓之間的時間延遲量的計算精度[3-4].

最直接的時間延遲量計算方法,可通過在時域中比較測量脈沖輪廓和標準脈沖輪廓的相似度來進行計算,這就是傳統(tǒng)的基于時域的互相關算法的基本思路.但是,該方法的測量精度與脈沖輪廓的采樣率有關.為克服采樣率的限制,研究者提出了經(jīng)典的基于頻域的Taylor快速傅里葉變換(Fast Fourier Transform,FFT)算法[4],其對累積脈沖輪廓進行離散傅里葉變換,將時域問題轉(zhuǎn)換到頻域再進行計算,從而使其測量精度可不受累積脈沖輪廓采樣率的影響,但脈沖輪廓信噪比(Signal to Noise Ratio,SNR)對該方法影響較大,并且迭代搜索測量時間延遲量會使得計算量增大.文獻[5]利用歷元折疊來重新獲得脈沖星光子強度函數(shù),并用該函數(shù)估計時間延遲量,提出了基于非線性最小均方差的時間延遲測量算法.隨后,文獻[6]利用脈沖星脈沖光子到達時間進行數(shù)學建模,利用光子到達時間概率函數(shù)構(gòu)建最大似然準則對X射線脈沖星時間延遲量進行測量.文獻[5-6]的兩種方法均是非凸優(yōu)化問題,為避免局部最優(yōu)解,必須通過直接搜索得到時間延遲量,而直接搜索的過程需要在計算復雜度和測量精度之間進行權(quán)衡:若提高測量精度,則計算復雜度也會提高;若計算復雜度降低,則測量精度又會降低.

為解決上述問題,筆者提出了一種基于快速最大似然(Fast Maximum Likelihood,FML)估計的時間延遲量測量方法,并通過縮小標準脈沖輪廓的相位間隔對其進行了改進.該方法首先利用歷元折疊提取的脈沖輪廓建立似然函數(shù),通過求最大似然函數(shù)得到X射線脈沖星脈沖輪廓時間延遲量;再對快速最大似然估計算法進行改進,縮小標準脈沖輪廓的相位間隔,采用粗略估計與精細估計相結(jié)合的策略求得時間延遲量.對真實數(shù)據(jù)的仿真結(jié)果表明,在運算量允許的范圍內(nèi),文中方法能夠提高時間延遲的測量精度.

1　快速最大似然估計的時間延遲測量方法

1.1脈沖星時間延遲測量的基本原理

對于某一觀測歷元t,假設脈沖星標準脈沖輪廓模型函數(shù)為s(t),滿足E[s(t)]=0,則測量脈沖輪廓p(t)可用時間延遲量為τ的標準脈沖輪廓s(t-τ)來表示[1],即

其中,a為測量脈沖輪廓與標準脈沖輪廓的常數(shù)偏差;b為脈沖輪廓尺度變換因子;τ為時間延遲量,或時間原點的平移量;g(t)為X射線背景和探測器本底噪聲,并假設為加性隨機噪聲.

在實際數(shù)據(jù)處理過程中,離散的測量脈沖輪廓可用離散的標準脈沖輪廓表示為

由式(2)可知,測量脈沖輪廓和標準脈沖輪廓歸一化后,僅相差i個采樣間隔,因此,測量脈沖輪廓相對于標準脈沖輪廓的時間延遲量τ為

為提取脈沖輪廓,所有光子都必須按照單個脈沖的光子到達時間對齊并排列,這種所有測量光子排列成一個脈沖輪廓的數(shù)據(jù)處理過程稱為歷元折疊,即歷元折疊是在期望的脈沖周期[0,1]或[0,2π]內(nèi)同步平均所有光子能量.對于某顆脈沖星,假設觀測時間為Tobs(Tobs?P),包含NP個脈沖周期,即Tobs≈NPP,這里P代表脈沖星周期.在觀測時間內(nèi),星載探測器探測并記錄了n個脈沖光子到達時間,光子到達探測器的原子時為τ0至τn-1.光子到達探測器的原子時轉(zhuǎn)換為太陽系質(zhì)心坐標時為t0至tn-1.將觀測時間分為NP個周期,將第2,3,…,N個周期的脈沖光子按照到達時間的順序疊加到第1個周期上.再將第1個脈沖周期分為m個相等的時間間隔Δt,每一個間隔代表一個脈沖輪廓曲線的“相位盒”,統(tǒng)計相位盒內(nèi)的脈沖光子數(shù),這樣一些相位盒就會擁有比較多的光子,表明具有很強的光子輻射流量;而一些相位盒擁有的光子數(shù)量就比較少,從而產(chǎn)生一條可識別的脈沖輪廓曲線.

1.2基于快速最大似然估計的時間延遲測量方法

現(xiàn)有的基于最大似然估計的時間延遲計算算法,從概率統(tǒng)計學的角度直接測量得到光子到達時間集,對脈沖光子到達時間建模,建立似然函數(shù),再利用最大似然準則估計求得時間延遲量.但隨著觀測時間增長,探測器探測到的光子數(shù)越來越多,直接利用光子到達時間計算時間延遲量將要付出高昂的計算代價.基于上述原因,文中提出一種基于FML時間延遲測量方法,先將探測到的光子按周期進行歷元折疊,得到脈沖輪廓信息;然后將測量脈沖輪廓與標準脈沖輪廓進行比較,即構(gòu)造似然函數(shù);最后對標準脈沖輪廓進行不同時間的移動,計算對應的似然函數(shù)值,記錄得到最大似然函數(shù)值所對應的標準脈沖輪廓移動時間,即可作為時間延遲量.

在觀測時間間隔(t0,tf)內(nèi),觀測時間為Tobs=tf-t0,其中,Tobs≈NPP,P代表脈沖星周期.將探測到的光子按周期進行歷元折疊,再將每個相位間隔內(nèi)的光子數(shù)量進行標準化處理,即得到經(jīng)驗脈沖輪廓.第i個相位間隔內(nèi)時間標記的光子密度函數(shù)可表示為

其中,NP表示觀測時間包含的脈沖星周期個數(shù);Nb表示將脈沖周期時間劃分為Nb個等長相位間隔,每個相位間隔的長度為Tb;cj(i)表示在觀測時間內(nèi),第j個脈沖周期中第i個相位間隔內(nèi)探測到的光子數(shù).

測量脈沖輪廓第i個相位間隔內(nèi)的光子數(shù)可表示為

同理,通過插值運算可將標準脈沖輪廓表示為第1個周期內(nèi)各相位間隔所對應的一系列脈沖幅值,即

考慮到脈沖輪廓具有周期性,采用快速最大似然估計時,對標準脈沖輪廓進行循環(huán)移位操作,將測量脈沖輪廓與標準脈沖輪廓進行如下運算,可得到時間延遲量的代價函數(shù)[7],即

去掉和τ無關的項,則有

采用一維的網(wǎng)格搜索法,對第1個周期內(nèi)的每個相位間隔進行計算,記錄得到最大似然函數(shù)的相位值,即為時間延遲量:

可以注意到,由式(10)求得的時間延遲量比直接使用脈沖光子到達時間求得的時間延遲量運算量小,但測量精度依賴于脈沖輪廓的時間分辨率,這是由m的取值所決定的,為了緩解這一問題,文中進一步提出了改進的方法.

1.3改進的快速最大似然時間延遲估計

如果直接采用擴大相位間隔個數(shù)的標準脈沖輪廓與測量脈沖輪廓進行比較,對標準脈沖輪廓進行不同相位移動的次數(shù)就要增加,這樣大幅度地增大了計算量.為此,文中采用先粗略估計、后精確估計的時間延遲測量方法對其進行改進.

對標準脈沖輪廓進行插值采樣:

其中,m0為縮小相位間隔后的標準脈沖輪廓的相位間隔個數(shù);Nb為測量脈沖輪廓的相位間隔個數(shù).

通過求測量脈沖輪廓與插值采樣后標準脈沖輪廓的最大似然函數(shù)值得到粗估計的相位延遲量,即

其中,m=0,1,2,…,Nb-1,此時得到的相位延遲量只是一個粗略的估計,真實的相位延遲量應在mpos附近,或者就是mpos,即在[mpos-m0/m,mpos+m0/m]范圍內(nèi).為得到更精確的時間延遲量,在粗略估計得到的相位延遲量附近進行精確搜索,這樣計算量不會大幅度增大,但測量精度又有所提高.

在[mpos-m0/m,mpos+m0/m]范圍內(nèi),對標準脈沖輪廓進行不同相位移動,然后記錄似然函數(shù)值,找到最大似然函數(shù)值對應的信號移動相位,即可求出時間延遲量.對標準脈沖輪廓移相后再進行插值采樣,則有

求標準脈沖輪廓與測量脈沖輪廓的似然函數(shù).當似然函數(shù)最大時,對應的信號移動時間即為時間延遲量:

2　實驗結(jié)果與分析

2.1實測數(shù)據(jù)樣本來源及選擇

文中采用真實數(shù)據(jù)來驗證所提算法的有效性,并將其與經(jīng)典的Taylor FFT算法進行比較.以X射線脈沖星B0531+21為例[4],提取脈沖輪廓的RXTE觀測數(shù)據(jù)號為40805-01-05-000,該觀測號的有效觀測時間T=15 730.508 077 s.同時提取相位間隔個數(shù)m=100,200,400,800,1 000時的累積脈沖輪廓作為標準脈沖輪廓.文中的Taylor FFT算法和FML算法標準脈沖輪廓與測量脈沖輪廓的相位間隔個數(shù)相同,改進后的FML算法中標準脈沖輪廓的相位間隔個數(shù)是測量脈沖輪廓的10倍.這里,仿真平臺所用計算機的配置為Inter雙核2.2 GHz CPU,2 GBDDR2內(nèi)存.

在脈沖輪廓周期折疊過程中,以FTOOLS軟件的t1=50 849.000 731 64(MJD)作為標準脈沖輪廓周期折疊的起始歷元;為便于進行時間延遲測量精度計算,選擇與t1相近的另一時刻t2=50 849.000 731 74 (MJD)作為起始歷元來獲取測量脈沖輪廓.為比較測量誤差,需計算測量脈沖輪廓相對于標準脈沖輪廓的時間延遲量真實值,即

其中,%為求余運算符,p為脈沖星周期,其值取自Princeton脈沖星星歷數(shù)據(jù)庫.

2.2觀測時間及采樣間隔對測量精度的影響

從觀測時間及采樣間隔的角度出發(fā),測試觀測時間和采樣間隔對時間延遲測量精度的影響.圖1是觀測時間分別在500～2 500 s,10～100 s,相位間隔個數(shù)m=100,200,400,800,1 000時,頻域Taylor FFT算法、FML算法和改進的FML算法得到的累積脈沖輪廓時間延遲量的測量誤差絕對值.改進的FML算法采用的標準脈沖輪廓的相位間隔個數(shù)依次為m0=1 000,2 000,4 000,8 000,10 000.

由圖1(a)和(b)可知,頻域Taylor FFT算法受噪聲影響較為明顯,當觀測時間一定時,時間延遲量的測量誤差絕對值隨著相位間隔個數(shù)m的增大,即信噪比的減小而逐漸減小,這是因為頻域Taylor FFT算法是采用步長迭代方法全局搜索得到時間延遲量,迭代步數(shù)為1 000m,當m增大時,迭代步數(shù)增多,測量誤差就減小.Taylor FFT算法隨著觀測時間的增長,測量誤差也逐漸減小;由圖1(c)和(d)可知,FML算法受相位間隔的影響較小,測量誤差絕對值比Taylor FFT算法的測量誤差絕對值小;由圖1(e)和(f)可知,改進的FML算法的測量誤差絕對值明顯減小,且測量誤差絕對值幾乎不受觀測時間的影響.

2.3運算量的比較

頻域Taylor FFT算法、FML算法和改進的FML算法分別進行一次時間延遲測量所需要的時間如表1所示.由表1可知,FML算法和改進的FML算法,隨著相位間隔個數(shù)的增大,耗時增長.這是因為隨著相位間隔個數(shù)的增大,求似然函數(shù)值的次數(shù)就越多.改進的FML算法比FML算法的耗時要長,但耗時并沒有增大很多,而測量精度卻提高了很多.FFT算法的迭代步長是m/1 000.FFT算法計算量比較大,這是因為FFT算法采用的是迭代搜索的方法在全局范圍內(nèi)搜索時間延遲量,迭代步數(shù)越多,計算量越大.

圖1　不同觀測時間下不同算法的時間延遲量測量誤差絕對值

表1　3種算法耗時對比s

3　結(jié)束語

為提高脈沖星累積脈沖輪廓時間延遲量的測量精度,筆者提出了一種快速最大似然估計的時間延遲測量方法,并通過縮小標準脈沖輪廓的相位間隔對其進行改進.該方法首先通過歷元折疊提取脈沖星脈沖輪廓,然后建立似然函數(shù),最后通過求最大似然函數(shù)得到X射線脈沖星脈沖輪廓時間延遲量;再采用縮小標準脈沖輪廓相位間隔對快速最大似然估計進行改進,能夠使標準脈沖輪廓包含更多的相位信息,通過粗略估計與精細估計相結(jié)合的方法求得時間延遲量.仿真結(jié)果表明,文中方法具有以下優(yōu)點:測量精度高于頻域FFT算法;FML算法的運算耗時明顯比Taylor FFT算法的少;改進的FML算法測量精度比FML算法的測量精度明顯提高,運算耗時略有增大,但相對測量精度的提高幅度,該計算量的增大可忽略.

[1]SHEIKH S I.The Use of Variable Celestial X-ray Sources for Spacecraft Navigation[D].Maryland:University of Maryland,2005.

[2]豐大軍,許錄平,宋詩斌,等.利用X射線脈沖星的星間差分聯(lián)合定位方法[J].西安電子科技大學學報,2014,41 (1):81-86. FENG Dajun,XU Luping,SONG Shibin,et al.Differential Co-localization Method Using the X-ray Pulsar[J].Journal of Xidian University,2014,41(1):81-86.

[3]劉秀平,李小平,孫海峰,等.利用Bayesian多尺度模型的X射線脈沖星信號估計[J].西安電子科技大學學報,2013,40(4):90-95. LIU Xiuping,LI Xiaoping,SUN Haifeng,et al.Estimation of the X-ray Pulsar Signal Based on the Bayesian MultiscaleModel[J].Journal of Xidian University,2013,40(4):90-95.

[4]蘇哲,許錄平,王婷,等.一種新的脈沖星累積脈沖輪廓時間延遲測量算法[J].宇航學報,2011,32(6):1256-1261. SU Zhe,XU Luping,WANG Ting,et al.A New Time Delay Measurement Algorithm for Pulsar Accumulated Pulse Profile[J].Journal of Astronautics,2011,32(6):1256-1261.

[5]EMADZADEH A A,GOLSHAN A R,SPEYER J L.Consistent Estimation of Pulse Delay for X-ray Pulsar Based Relative Navigation[C]//Joint 48th IEEE Conference on Decision and Control,and 28th Chinese Control Conference. Piscataway:IEEE,2009:1488-1493.

[6]EMADZADEH A A,SPEYER J L.On Modeling and Pulse Estimation of X-ray Pulsars[J].IEEE Transactions on Signal Processing,2010,58(9):4484-4495.

[7]RINAURO S,COLONNESE S,SCARANO G.Fast Near-maximum Likelihood Phase Estimation of X-ray Pulsars[J]. Signal Processing,2013,93(1):326-331.

(編輯:齊淑娟)

Method of time delay measurement based on fast maximum likelihood for pulsar pulse profile

JIANG Ruomei1,SHE Yan2
(1.City College,Xi’an Jiaotong Univ.,Xi’an 710018,China; 2.School of Microelectronics,Xidian Univ.,Xi’an 710071,China)

The time delay of the pulsar integrated pulse profile relative to the standard pulse profile is one of the important observations in an X-ray pulsar-based navigation system,with the measuring accuracy directly affecting the the pulsar-based navigation system.In order to improve the measuring accuracy of the time delay of the pulsar integrated pulse profile and reduce the computation complexity,a method called fast maximum likelihood(FML) estimation is proposed in this paper,which first distills the pulse profile of the X-ray pulsar by the epoch folding. Then,the likelihood function is estimated to obtain the time delay of the pulsar integrated pulse profile.This paper also improves the fast maximum likelihood time delay estimation by narrowing the phase interval of the standard pulse profile and combines rough estimation with accurate estimation to effectively compute the time delay.Real pulsar data observed by the Rossi X-ray Timing Explorer(RXTE)are used to validate the effectiveness of the proposed method from the aspects of the observed time and computational complexity.Compared with the Taylor FFT algorithm,the measurement accuracy of the proposed method can be significantly improved while keeping the total computational complexity.

X-ray pulsar;pulse profile;epoch folding;time delay;fast maximum likelihood estimation

V448.2

A

1001-2400(2016)05-0147-06

10.3969/j.issn.1001-2400.2016.05.026

2015-07-16 網(wǎng)絡出版時間:2015-12-10

江若玫(1979-),女,講師,E-mail:13991918588@139.com.

佘 顏(1988-),女,西安電子科技大學碩士研究生,E-mail:sheyan51@163.com.

網(wǎng)絡出版地址:http://www.cnki.net/kcms/detail/61.1076.TN.20151210.1529.052.html