基于小波分析和密度估計的紅移測速導(dǎo)航研究

摘 要：基于光譜紅移天文自主導(dǎo)航思路以及天文學(xué)領(lǐng)域求解紅移值的方法，提出一種用小 波分析和密度估計相結(jié)合的方法測飛行器紅移值，進(jìn)而求解速度的測速導(dǎo)航方法。首先采用小波 變換濾除噪聲并進(jìn)行光譜歸一化，然后采用整體閾值和局部閾值相結(jié)合的方法，提取特征譜線， 用密度估計法求出紅移值，最后結(jié)合天文測角信息，計算出飛行器相對于慣性坐標(biāo)系的速度。采 用實測的太陽光譜數(shù)據(jù)，進(jìn)行了紅移模擬與解算。實驗證明該方法在高信噪比的情況下，計算速 度快，求解結(jié)果精度高達(dá)96%，為基于光譜紅移的天文導(dǎo)航提供一種新的速度計算方法。

關(guān)鍵詞：天文導(dǎo)航；測速；紅移；小波分析；密度估計

中圖分類號：V249.32+9 文獻(xiàn)標(biāo)識碼：A 文章編號：1673-5048（2014）06-0003-05

ResearchonMeasurementofSpectrumRedshiftinVelocityNavigation UsingWaveletAnalysisandDensityEstimationMethod

WANGYong，ZHAOYan，YANGKui

（SchoolofInstrumentationScienceandOptoelctronicsEngineering，BeihangUniversity，Beijing100191，China）

Abstract：Basedonideaofspectrumredshiftastronomicalautonomousnavigationandmethodforastron omysolvingaredshift，presentsavelocitycelestialnavigationmethod，whichmeasuresavehicleredshiftsu singthemethodofwaveletanalysiscombinedwithdensityestimation，andsolvesthevehicle’sspeed.First， wavelettransformisusedtofilterthenoiseandnormalizespectrum.Thenthecharacteristicspectrallinesare extractedbycombiningtheglobalthresholdandlocalthreshold，andtheredshiftvalueiscalculatedusingden sityestimationmethod.Finallycombiningwiththeastronomicalangleinformation，thespeedofthevehiclerel ativetoinertialcoordinatesystemiscalculated.Usingtheobserveddataofthesolarspectrum，simulationand calculationoftheredshiftareconducted.TheexperimentsprovethatinthehighSNR，themethodhasfeatures offastcalculation，highprecisionofsolvingresults，theaccuracyratereached96%，andprovidesanewspeed calculationmethodforcelestialnavigationbasedonspectralredshift.

Keywords：celestialnavigation；velocitymeasurement；redshift；waveletanalysis；densityestimation

0 引 言

深空探測是指發(fā)射航天器在深空對地外天體或空間進(jìn)行探測的活動，是當(dāng)今世界航天活動的 重要領(lǐng)域。常用的自主導(dǎo)航方法是天文光學(xué)自主 導(dǎo)航，其導(dǎo)航性能受制于目標(biāo)天體的數(shù)量和距離、 空間環(huán)境、航天器測量精度等多種因素，且不能直 接測速，為此提出了基于光譜紅移的航天器自主 導(dǎo)航方法[1]。

天然光源的光譜包含了天體相對于航天器的 速度信息，可通過光譜紅移量結(jié)合航天器姿態(tài)信息與天文測角參數(shù)，快速解算出航天器在慣性空 間中的速度與位置信息，從而實現(xiàn)深空探測高精 度天文光學(xué)測角測速組合自主導(dǎo)航。這是一種幾 何導(dǎo)航定位方法，不需要復(fù)雜的軌道動力學(xué)，原理 簡潔清晰，計算量小。該方法的核心就是求解出飛 行器相對于參考天體的紅移值。

目前尚未有直接研究基于光譜紅移的天文自 主導(dǎo)航的論文報告，但是天文學(xué)家采用光譜紅移 測量天體運行速度的研究可以追溯到1859年。 1868年英國天文學(xué)家威廉·哈金斯利用多普勒效 應(yīng)，首次從譜線的微小位移測出了天狼星的視向 速度[2]。而最早出現(xiàn)的紅移自動測量技術(shù)則是 Tonry和Davis的交叉相關(guān)法[3]。此方法目前被認(rèn) 為是最成功的紅移自動測量方法，后來澳大利亞 天文臺的Glazebrook在此基礎(chǔ)上提出了基于主分 量（PCA）的PCAZ方法[4]。但兩種方法運算量都偏 大。國內(nèi)的相關(guān)工作者有周虹的基于Hough變換 和神經(jīng)網(wǎng)絡(luò)求紅移的方法[5]；邱波的偽三角法[6]。 但是這些方法只能處理發(fā)射線天體，不適合用于 對吸收譜天體度的測量。段福慶提出的密度估計 法[7]又比較依賴譜線的提取結(jié)果。

本文結(jié)合文獻(xiàn)[8]的天文導(dǎo)航思路以及天文學(xué) 領(lǐng)域求解紅移值的方法，提出一種用小波分析和 密度估計相結(jié)合的方法測量飛行器紅移值，進(jìn)而 求解速度的導(dǎo)航方法。

1 測速原理

1.1 夫瑯和費線

天體光譜是由連續(xù)譜、譜線和各種噪聲組成 的，譜線是由天體中的各種原子、分子等在連續(xù)譜 基礎(chǔ)上吸收或輻射能量所體現(xiàn)出的特征[9]。以太 陽光譜為例，德國物理學(xué)家約瑟夫·夫瑯和費發(fā) 現(xiàn)太陽光譜中有一系列暗特征線，后來被稱之為 夫瑯和費線[10]。這些主要的特征線用字母A到K 標(biāo)注，較弱的則用其他字母。后來基爾霍夫和羅伯 特·本生確認(rèn)了每一條譜線所對應(yīng)的化學(xué)元素， 并推論在太陽光譜中的暗線是由太陽上層的那些 元素吸收造成的，有些被觀察到的特征線則是在 地球大氣層中的氧分子造成的[11]。主要的夫瑯和 費線和對應(yīng)的元素參見表1。

在這些譜線中，C，F(xiàn)，G，h，H，K，L，N等吸 收線比較明顯，可以作為測量紅移值的參考特征。

1.2 光波多普勒效應(yīng)（紅移與藍(lán)移）

除此之外，天體光譜的一個重要參數(shù)紅移值， 蘊藏著天體的運行速度信息。當(dāng)天體以很高的速度背離地球使得觀測到的譜線比靜止的譜線波長， 表現(xiàn)為譜線朝紅端移動一段距離，叫做紅移現(xiàn)象。 反之，則會出現(xiàn)藍(lán)移現(xiàn)象。

紅移的大小由“紅移值”衡量，紅移值用z表 示，定義為

式中：v1，v2，v3為各天體在慣性系中的速度矢量，可以通過星歷確定；u1，u2，u3為慣性系中各天體 指向航天器位置矢量的單位矢量，可由太陽敏感 器或星敏感器測得。

至此可以建立關(guān)于速度矢量、位置矢量的狀 態(tài)估計方程組：

給定初值后，求解方程組（6）可得航天器在慣 性坐標(biāo)系中的速度矢量vp，知道時間后進(jìn)行積分 即可以求得位置矢量rp。

2 基于小波分析和密度估計的紅移值求解

測出參考天體的光譜紅移值是實現(xiàn)光譜測速 的重要環(huán)節(jié)，本文側(cè)重于研究光譜紅移值的測量。 一般包括光譜預(yù)處理、特征譜線提取、特征匹配和 紅移值解算三個步驟。

2.1 光譜預(yù)處理

天體光譜的預(yù)處理包括噪聲去除和連續(xù)譜歸 一化[12]。

其中：b代表天光背景，與成像視野里的天體背景 相關(guān)；i代表設(shè)備噪聲，它在成像階段可以得到很 好的抑制，它們的殘余和其他噪聲源引起的噪聲 混合在一起，其表現(xiàn)形式是隨機(jī)白噪聲以及波長 不固定的干擾噪聲，可以采用小波軟閾值法去 除[13]；x代表未知源。一般情況下噪聲主要來自i。

常用的去噪方法包括高斯濾波、中值濾波、小 波閾值法[14]。鑒于收集到的光譜信號中主要的噪 聲是高斯白噪聲，而小波閾值法的優(yōu)點就是可以 很好地去除高斯白噪聲，所以采用小波分層閾值 降噪的方法，具體闡述如下：

其中：dm，dt等是前j階中系數(shù)最大的k個小波系 數(shù)。閾值確定采用的是Birge-Massart法則：給定 一個指定的分解層數(shù)j，對j+1以及更高層的所有 系數(shù)保留；對于第i（i=1～j）層保留絕對值最大的 k個小波系數(shù)，其余的置0。

連續(xù)譜的存在使得譜線的強(qiáng)度不能真實地反映 出來，所以必須使光譜中緩慢變化的部分歸為零。 通常采用的方法是原始光譜減去連續(xù)光譜，該過程 叫做連續(xù)譜歸一化。常采用多項式逼近法、中值濾 波、小波變換等方法來進(jìn)行連續(xù)光譜的擬合[15]。多 項式逼近實際是多項式擬合，涉及到基向量和最小 二乘，計算量大。中期濾波是一個典型的非線性濾 波技術(shù)，雖簡單易行，但對于吸收譜，會由于吸收帶 的影響造成中值位置與真實連續(xù)譜不一致。

一般采用的小波變換方法其缺點在于：生成 的譜線光譜中仍然有噪聲和部分較弱的吸收譜線 混雜在一起，不利于后期譜線特征的提取，因此本 文在譜線光譜的基礎(chǔ)上再次進(jìn)行光譜歸一化，得 到新的譜線光譜（由于譜線光譜中連續(xù)譜線已經(jīng)基 本去除，嚴(yán)格地說所進(jìn)行的操作并不叫做“光譜歸 一化”）。假定原始的譜線光譜成為一次光譜，新 得到的光譜稱為二次光譜。具體步驟如下：

2.2 特征譜線提取

已知太陽光譜中有很多夫瑯和費線，理論上 只要準(zhǔn)確測定了同一條吸收譜線在靜止和運動狀 態(tài)所對應(yīng)的頻率值f0和f以及夾角θ，參考公式 （2），就可以精確解算某一時刻的徑向速度vcosθ， 在實際操作中提取多條譜線以減少計算誤差。本文決定提取表1中的吸收強(qiáng)度比較大的譜線，也即 是C，F(xiàn)，G，h，H，K，L，N。

吸收譜線多為局部極小值，若只采用局部閾 值會導(dǎo)致某些譜線稀疏的地方產(chǎn)生偽譜線。而只 采用整體閾值又會導(dǎo)致吸收譜線強(qiáng)度弱的地方譜 線丟失。因此本文采用局部閾值和整體閾值相結(jié) 合的方法[16]。具體步驟如下：

Step1：設(shè)置局部閾值t，提取局部極小值點。

Step2：設(shè)置整體閾值T，提取低于整體閾值的 所有點（本文提取吸收譜線，故取低于整體閾值的 譜線點）。

2.3 特征匹配和紅移值解算

對于觀測光譜來說，每一時刻確定的速度v唯 一地和一個z值對應(yīng)。靜止的吸收譜線被提取時就 變成一個離散的數(shù)列{λi，i=1，2，…，M}。假設(shè)與 z值對應(yīng)的被提取的吸收譜線的數(shù)列為{λ′i，i=1， 2，…，M}。根據(jù)公式（1）可以知道λ=（1+z）λ0。 由此可知，理想狀態(tài)下兩個離散數(shù)列中的對應(yīng)吸 收譜線必然滿足此種關(guān)系。按照ki=λ′i/λi-1，計 算可以得出一個數(shù)據(jù)集k={ki，i=1，2，…，M· N}。很顯然，在k=ki處數(shù)列的點比較密集。因此 只要找到數(shù)據(jù)集中最密集的點域，而后進(jìn)行均值 化，就可以得到估計的紅移值k。

3 驗證實驗

為了測試算法的可行性，采用太陽光譜數(shù)據(jù)： ASTMG173-03ReferenceSpectraDerivedfromSM ARTSv.2.9.2[17]。該數(shù)據(jù)來自于美國材料與試驗 協(xié)會，采集于大氣表層外，很好地避開了天光線的 干擾，很類似于已在深空飛行的航天器獲得的光 譜數(shù)據(jù)。把此數(shù)據(jù)作為標(biāo)準(zhǔn)數(shù)據(jù)，借助MATLAB， 為數(shù)據(jù)添加信噪比不同的高斯白噪聲。分別對其 進(jìn)行不同的紅移值模擬，而后采用以上提出的方 法進(jìn)行紅移值的測定。

閾值降噪時，采用層數(shù)為三的多尺度分解，具 有較好的效果，處理結(jié)果如圖1。從圖中可以看出 該方法可以成功地濾除噪聲，并對吸收譜線造成 較小的影響。

擬合譜線光譜時，多尺度分解的層數(shù)選擇4， 處理結(jié)果如圖2所示。

為了突出譜線而對譜線光譜再進(jìn)行一次光譜 歸一化，處理前后的對比結(jié)果如圖3所示。從圖中 看出低頻和噪聲信號又一次降低，高強(qiáng)度的吸收 譜線更加突出。

本文采用局部閾值和整體閾值相結(jié)合的方法 來提取譜線。經(jīng)實驗調(diào)試整體閾值T采用1.2rms， 而局部閾值t采用局部均方根（rms）。由于譜線的 寬度多為2～3nm，本文采用的數(shù)據(jù)分辨率為0.5 ～1nm，所以局部窗口選擇2.5～5nm。實驗證明 該方法有較好的效果（如圖4所示）。

采用密度估計法測100組信噪比不低于20， 模擬紅移值在1.0001～1.01之間的數(shù)據(jù)，可知正 確率達(dá)到96%。其中一組的測量結(jié)果如圖5所示 （圖5（a）橫軸為k值從小到大的排序，縱軸代表所 測得k值，越密集的地方越有可能是紅移值。圖5（b）為統(tǒng)計k值最密集的點的放大圖，也就是測 得的紅移值）。舉例時光譜紅移的模擬值為k=1. 01，均值化時的點域半徑σ=0.5（也是光譜數(shù)據(jù)的 分辨率）密度估計法計算結(jié)果也是k=1.01。可見 本方法在實踐中是可行的。

4 結(jié) 論

從實驗結(jié)果可以得知，該方法可以僅僅依靠 一些高效的算法實現(xiàn)快速、準(zhǔn)確地解求紅移值。如 果再結(jié)合角度信息，利用公式（4）～（6）便可以解 求飛行器速度。本文在提取吸收譜線時，只針對8 條特別強(qiáng)的吸收譜線做處理。從理論和算法上證 明該方法的可行性，并且首次采用兩次連續(xù)譜歸 一化方法，使得歸一化后的譜線特征更為明顯。這 對于吸收譜線比較明顯的天體紅移測量具有新的 借鑒意義。為適應(yīng)將來更高精度的天文導(dǎo)航，還可 以先進(jìn)行譜線擬合再進(jìn)行匹配以獲得更高精度的 測量值。

參考文獻(xiàn)：

[1]張偉，陳曉，尤偉，等.光譜紅移自主導(dǎo)航新方法[J]. 上海航天，2013，30（2）：32-38.

[2]HugginsW.FurtherObservationsontheSpectraofSomeoftheStarsandNebulae，withanAttempttoDetermine TherefromWhetherTheseBodiesareMovingtowardsor fromtheEarth，AlsoObservationsontheSpectraofthe SunandofCometII，1868[J].PhilosophicalTransactions oftheRoyalSocietyofLondon，1868，158：529-564.

[3]TonryJ，DavisM.ASurveyofGalaxyRedshifts.I.Data ReductionTechniques[J].TheAstronomicalJournal， 1979，84（10）：1511-1525.

[4]GlazebrookK，OfferAR，DeeleyK.AutomaticRedshift DeterminationbyUseofPrinicipalComponentAnalysis. Ⅰ.Fundamentals[J].TheAstrophysicalJournal，1998， 492（1）：98-109.

[5]周虹，黃凌云，羅曼麗.一種基于Hough變換和神經(jīng) 網(wǎng)絡(luò)的分層類星體識別方法[J].電子科學(xué)學(xué)刊， 2000，22（4）：529-535.

[6]邱波，胡占義，趙永恒.一種快速求紅移和證認(rèn)譜線 的新方法———偽三角法[J].光譜學(xué)與光譜分析， 2002，22（4）：695-698.

[7]段福慶，吳福朝，羅阿理，等.用于紅移測量的基于密 度估計的模板匹配法[J].光譜學(xué)與光譜分析，2005， 25（11）：1895-1898.

[8]張偉，方寶東，陳曉，等.航天器光譜紅移自主導(dǎo)航方 法[P].中國專利：CN103206955A，2013-07-17.

[9]段福慶，吳朝福，羅阿理，等.基于均值漂移實現(xiàn)非射 線天體光譜的譜線自動提取[J].光譜學(xué)與光譜分析， 2005，25（11）：1884-1888.

[10]趙梅芳，吳潮，羅阿理，等.基于K近鄰方法的窄線 與寬線活動星系核的自動光譜分類[J].天文學(xué)報， 2007，48（1）：1-10.

[11]JenkinsFA，WhiteHE.FundamentalsofOptics[M]. McGrawHillSciense/Engineering/Math，1981：18- 20.

[12]劉蓉，段福慶，劉三陽，等.基于知識的紅移測量和 譜線證認(rèn)方法[J].電子與信息學(xué)報，2006，28（1）： 76-79.

[13]趙瑞珍，羅阿理.天體光譜信號的連續(xù)譜歸一化新方 法[J].光譜學(xué)與光譜分析，2006，26（3）：587-590.

[14]趙瑞珍，胡占義，胡紹海.天體光譜信號去噪的小波 域復(fù)合閾值新算法[J].光譜學(xué)與光譜分析，2007，27 （8）：1644-1647.

[15]劉蓉，段福慶，劉三陽.基于特征約束的譜線自動提 取[J].光譜學(xué)與光譜分析，2007，27（8）：1648- 1652.

[16]劉蓉，劉三陽，趙瑞珍.正常星系光譜的一種譜線自 動提取方法[J].光譜學(xué)與光譜分析，2006，26（3）： 583-586.

[17]ASTMG173-03ReferenceSpectraDerivedfromSMAR TSv.2.9.2[OL/E].http：//rredc.nrel.gov/solar/ spectra/am1.5/ASTMG173/.