地球臨邊觀測衛(wèi)星的姿態(tài)角仿真分析方法

朱永生 胡海鷹 鄭珍珍 蘇瑞豐 張科科

(上海微小衛(wèi)星工程中心， 上海 201203)

?

地球臨邊觀測衛(wèi)星的姿態(tài)角仿真分析方法

朱永生 胡海鷹 鄭珍珍 蘇瑞豐 張科科

(上海微小衛(wèi)星工程中心， 上海 201203)

針對地球臨邊觀測衛(wèi)星的姿態(tài)角計算問題，給出一種應(yīng)用STK軟件進(jìn)行仿真的分析方法。利用STK軟件的連接接口聯(lián)立STK軟件，創(chuàng)建仿真場景；在二體模型下，求得臨邊觀測衛(wèi)星姿態(tài)角的解析解，并作為下一步精確迭代求解的初值；在WGS84地球模型下，調(diào)用STK軟件中的高精度軌道外推模型(HPOP)，并結(jié)合推導(dǎo)的迭代公式，求得衛(wèi)星姿態(tài)的高精度數(shù)值解。以20 km和100 km的臨邊觀測為例，仿真計算姿態(tài)角，并將計算結(jié)果與二體模型方法的進(jìn)行對比。結(jié)果表明，文章提出的方法可有效解決地球臨邊觀測衛(wèi)星的姿態(tài)角計算問題，且具有收斂迅速和計算精度高的優(yōu)勢。

地球臨邊觀測；衛(wèi)星姿態(tài)角；高精度計算

1 引言

地球臨邊觀測是指利用衛(wèi)星遙感器對準(zhǔn)地球某一特定高度大氣的切線方向進(jìn)行觀測，獲得觀測視場內(nèi)一個大氣層條帶里的散射輻射，并配合掃描對所要觀測大氣的高度進(jìn)行覆蓋探測。臨邊觀測是地球大氣探測的主要方式之一，由于其空間覆蓋范圍大、垂直分辨率高且無需特定目標(biāo)輔助，近來受到越來越多關(guān)注[1-3]。進(jìn)行臨邊觀測時，衛(wèi)星一般通過姿態(tài)機(jī)動使遙感器指向不同高度大氣的切線方向，因此衛(wèi)星的姿態(tài)角計算是臨邊觀測必須要解決的問題。然而，衛(wèi)星運動復(fù)雜[4-6]，且地球形狀不規(guī)則[7]，使姿態(tài)角的精確計算成為一個難點。

目前，地球臨邊觀測衛(wèi)星的姿態(tài)角計算，主要是假設(shè)地球為理想球形的二體模型，即使考慮了地球扁率攝動等對衛(wèi)星運動的影響，在計算姿態(tài)角時仍要認(rèn)為地球半徑一致[8-9]。這種方法雖然計算簡單，但是計算精度差。本文充分利用STK軟件中的高精度軌道外推(HPOP)模型[10]和WGS84地球模型，結(jié)合MATLAB軟件的編程計算能力，提出了地球臨邊觀測衛(wèi)星姿態(tài)角[11]的仿真分析方法。該方法可求出高精度數(shù)值解[12]，從而解決利用臨邊觀測切點高度求解衛(wèi)星姿態(tài)角問題，為臨邊觀測衛(wèi)星姿態(tài)控制提供參考目標(biāo)；還可應(yīng)用于有探測結(jié)果時，對應(yīng)時間和姿態(tài)角測量值計算出探測高度。

2 動力學(xué)建模及方法設(shè)計

2.1 運動方程

假設(shè)近地衛(wèi)星飛行過程中只受到地球引力作用，且地球是半徑Re為6 378.137 km的球體(理想情況)，則在二體模型下描述衛(wèi)星的運動方程為

(1)

式中：r為衛(wèi)星位置矢量；r為衛(wèi)星地心距；μ為地球引力常數(shù)。

求解式(1)，可得

(2)

(3)

式中：a為軌道半長軸；i為軌道傾角；Ω為升交點赤經(jīng)；ω為近地點幅角；M為平近點角；fr，ft，fn分別為徑向、橫向和法向攝動加速度；p為半通徑。

2.2 臨邊觀測模型

地球臨邊觀測三維示意如圖1所示。假設(shè)衛(wèi)星初始姿態(tài)對地定向，遙感器視軸方向沿本體坐標(biāo)系+yB軸，則通過調(diào)整衛(wèi)星姿態(tài)滾動角便可調(diào)整臨邊觀測高度，二維示意如圖2所示。其中：P1為在地球為理想球形假設(shè)下，遙感器視軸與任務(wù)要求觀測高度Δh處大氣的切點，對應(yīng)的姿態(tài)滾動角為φ1，與地心的距離為r1；P2為地球模型選擇WGS84模型時，遙感器視軸與任務(wù)要求觀測高度Δh處大氣的切點，對應(yīng)的姿態(tài)滾動角為φ2，與地心的距離為r2。

圖1 臨邊觀測三維示意Fig.1 3-D view of limb observation

由圖2可知，臨邊觀測衛(wèi)星姿態(tài)滾動角同任務(wù)要求觀測高度存在如下關(guān)系，可用來求得φ1和φ2。

(4)

式中：φ為姿態(tài)滾動角。

圖2 臨邊觀測二維示意Fig.2 2-D view of limb observation

另外，結(jié)合圖1和圖2還可以看出：衛(wèi)星繞zB軸轉(zhuǎn)動，設(shè)置不同的偏航角，可以觀測不同地理位置的大氣，但對觀測大氣的高度影響不大。因此，本文在分析地球臨邊觀測衛(wèi)星的姿態(tài)角時，只考慮滾動角，不考慮偏航角。

在大數(shù)據(jù)背景下，一個強(qiáng)大的智能平臺是支撐網(wǎng)站急需前行的關(guān)鍵。目前，基于大數(shù)據(jù)技術(shù)的智能平臺有許多，比如智能搜索引擎、智能化交互共享平臺等。因此，大數(shù)據(jù)背景下的高校圖書館也需依靠智能平臺為其各項服務(wù)提供支持，利用云計算等眾多智能技術(shù)將高校圖書館原本分散的服務(wù)進(jìn)行統(tǒng)一，使大數(shù)據(jù)時代的高校圖書館服務(wù)真正地與用戶聯(lián)系在一起，為其提供更加優(yōu)質(zhì)的信息資源與服務(wù)。

2.3 方法介紹

由式(4)可以看出，臨邊觀測的滾動角只與衛(wèi)星地心距r、地球半徑Re和觀測高度Δh有關(guān)。r可通過調(diào)用STK軟件中的HPOP模型求得；Δh由任務(wù)需求決定；由于地球形狀不規(guī)則，無法確定具體的Re值，因此無法直接利用式(4)求解滾動角φ。本文根據(jù)衛(wèi)星動力學(xué)及臨邊觀測幾何關(guān)系，給出了滾動角數(shù)值解求法，具體流程見圖3。在地球為理想球形假設(shè)的基礎(chǔ)上(Re=6 378.137 km)，聯(lián)立式(2)和式(4)，求得在二體模型下的滾動角解析解φ0，并以φ0作為初值，通過STK軟件建模，在HPOP模型下迭代求出滿足精度要求的滾動角數(shù)值解。其中：t為預(yù)設(shè)誤差；Δhi為高度計算誤差；i為迭代次數(shù)；迭代公式φi=f(φi-1)由臨邊觀測的幾何關(guān)系推導(dǎo)得出，具體推導(dǎo)過程如下。

圖3 滾動角計算流程Fig.3 Flow for computing roll angle

由圖2的幾何關(guān)系得

(5)

式中：Δh1為滾動角φ1時WSG84地球模型下的觀測高度；Δh2為滾動角φ2時WSG84地球模型下的觀測高度；Re1和Re2分別為切點P1和P2所在地理位置對應(yīng)的地球?qū)嶋H半徑。

令Re1=Re2，故

(6)

對于大氣切點高度，可以通過STK/VectorGeometryTool工具創(chuàng)建切點，再通過STKReport命令直接輸出[10]。

2.4 誤差分析

地球臨邊觀測衛(wèi)星的姿態(tài)角計算誤差主要來自衛(wèi)星姿態(tài)控制精度、遙感器安裝精度、地球形狀模型和軌道外推模型。其中：衛(wèi)星姿態(tài)控制精度和遙感器安裝精度屬于衛(wèi)星自身系統(tǒng)誤差，直接影響姿態(tài)角計算結(jié)果；地球形狀模型和軌道外推模型屬于外部計算誤差，分別影響大氣高度值和衛(wèi)星預(yù)報位置，繼而影響姿態(tài)角的計算結(jié)果。若考慮姿態(tài)控制精度σ1，遙感器安裝精度σ2，地球形狀模型導(dǎo)致的大氣高度誤差σ3，軌道外推模型導(dǎo)致的衛(wèi)星徑向位置預(yù)報誤差σ4，則推導(dǎo)得到這些誤差導(dǎo)致的計算誤差Δσ1，Δσ2，Δσ3，Δσ4分別為

(7)

式(7)中，Δσ4是在地球為理想球形的假設(shè)下得到的，后續(xù)可考慮采用更精確的地球形狀模型來分析軌道外推模型精度引起的計算誤差。

3 仿真分析

3.1 計算方法比較

假設(shè)衛(wèi)星的軌道為700 km高的太陽同步軌道，對100 km、50 km、20 km高度地球臨邊大氣進(jìn)行觀測，迭代誤差設(shè)置為5.00×10-3km，觀測時間為2016-01-01T00:00:00(UTC)，本文方法及二體模型方法的計算結(jié)果如表1所示，其中STK軟件驗證高度指利用STK軟件HPOP模型結(jié)合WSG84地球模型求得的大氣高度。

由表1不難看出，在理想球形假設(shè)的二體模型方法下，計算誤差超過了3.00 km；而本文方法的計算誤差控制在預(yù)定誤差5.00×10-3km以內(nèi)，兩種方法的計算誤差相差約600倍。同時，利用本文方法，每次的迭代次數(shù)均控制在5次以內(nèi)，收斂迅速。由此可見，本文方法計算精度高，收斂迅速，可有效解決地球臨邊觀測衛(wèi)星的姿態(tài)角計算問題。

表1 兩種方法計算結(jié)果比較

3.2 單軌100 km高度大氣臨邊觀測

衛(wèi)星軌道為700 km高的太陽同步軌道，初始姿態(tài)對地定向，考慮衛(wèi)星1個軌道周期內(nèi)(5 922.26 s)對100 km高度的大氣進(jìn)行臨邊觀測，預(yù)設(shè)誤差為5.00×10-3km。以二體模型方法求得的解析解25.32°作為初值，利用本文方法求得衛(wèi)星滾動角的變化如圖4所示，對應(yīng)的切點大氣實際高度如圖5所示，衛(wèi)星遙感器在1個軌道周期內(nèi)切過的大氣軌跡(地球固連坐標(biāo)系下)如圖6所示。

圖4 滾動角變化曲線Fig.4 Graph of roll angle change

圖5 切點大氣實際高度變化曲線Fig.5 Corresponding altitude change of tangent atmosphere

圖6 切過的大氣軌跡示意Fig.6 Tracing points of tangent atmosphere

由圖4可知，衛(wèi)星滾動角隨時間呈正弦曲線變化。在1個軌道周期內(nèi)，滾動角變化曲線出現(xiàn)2個波峰和2個波谷。其中：波峰出現(xiàn)的時刻分別為第979 s和第4065 s，此時衛(wèi)星星下點正好處于南北極點處；波谷出現(xiàn)的時刻分別為第2553 s和第5449 s，此時衛(wèi)星星下點正好處于赤道附近。這主要是因為地球在赤道處的半徑最大，在兩極的半徑較小。因此，觀測同樣高度的臨邊，在兩極時衛(wèi)星姿態(tài)向地心方向偏得最嚴(yán)重，在赤道時向地心方向偏得最小。

由圖5可知，利用本文方法得到的切點大氣高度變化范圍為99.994～100.005 km。對比二體模型方法的計算誤差(如圖7所示)，顯然本文方法計算精度較高。

圖7 兩種方法計算誤差對比Fig.7 Calculation error comparison of two methods

3.3 全天候20 km高度大氣臨邊觀測

衛(wèi)星軌道為700 km高的太陽同步軌道，初始姿態(tài)對地定向，考慮衛(wèi)星1 d對20 km高度大氣臨邊觀測。利用本文方法得到衛(wèi)星滾動角變化曲線和切過的大氣高度變化曲線，分別如圖8和圖9所示。為了更直觀了解20 km高度大氣臨邊觀測情況，圖10和圖11給出了衛(wèi)星1 d內(nèi)切過的20 km高度大氣軌跡的三維圖(地球固連坐標(biāo)系下)和二維圖。

圖8 1 d內(nèi)滾動角變化示意Fig.8 Graph of roll angle change in one day

圖9 1 d內(nèi)切過的大氣高度變化示意Fig.9 Altitude change of tangent atmosphere in one day

圖10 1 d內(nèi)切過的大氣軌跡三維圖Fig.10 3-D view of tracing points of tangent atmosphere in one day

圖11 1 d內(nèi)切過的大氣軌跡二維圖Fig.11 2-D view of tracing points of tangent atmosphere in one day

從圖8可以看出，衛(wèi)星滾動角呈周期性變化，且變化周期與衛(wèi)星軌道周期幾乎一致。這主要是由地球形狀和衛(wèi)星周期性運動造成的。當(dāng)衛(wèi)星星下點緯度接近0°時，觀測到的是赤道上空的大氣，此時地球半徑最大，故衛(wèi)星在赤道上空進(jìn)行臨邊觀測時滾動方向上的姿態(tài)機(jī)動量較小。

地球1 d內(nèi)繞自轉(zhuǎn)軸轉(zhuǎn)動360°，因此衛(wèi)星全天候的臨邊觀測可以保證掃過整個經(jīng)度上的大氣。另外，700 km高的太陽同步軌道衛(wèi)星，單軌覆蓋的緯度范圍為[81.84°S，81.83°N]。因此，利用本文方法，在1個回歸周期內(nèi)可以覆蓋經(jīng)度[180°E，180°W]、緯度[81.84°S，81.83°N]的臨邊大氣。

3.4 誤差分析

假設(shè)衛(wèi)星姿態(tài)控制精度為0.005°，遙感器安裝精度為0.000 1°，地球形狀模型引起的大氣高度誤差為0.05 km，軌道外推模型導(dǎo)致的衛(wèi)星徑向位置預(yù)報誤差為0.10 km，則由式(7)計算得到20 km和100 km地球臨邊觀測時的姿態(tài)角計算誤差，見表2。

表2 誤差分析結(jié)果

由表2可以看出：當(dāng)?shù)厍蛐螤钅Ｐ途葍?yōu)于0.05 km，軌道預(yù)報精度優(yōu)于0.10 km時，與衛(wèi)星姿態(tài)控制精度和遙感器安裝精度相比，它們引起的計算誤差很小。因此，采用本文方法引起的計算誤差主要來源于衛(wèi)星自身的姿態(tài)控制精度和遙感器安裝精度。

4 結(jié)束語

本文針對地球臨邊觀測衛(wèi)星的姿態(tài)角計算問題，提出了一種利用STK軟件的高精度分析方法，并給出了20 km和100 km臨邊觀測的計算結(jié)果。基于仿真分析結(jié)果，可以對實際臨邊觀測任務(wù)中的衛(wèi)星姿態(tài)導(dǎo)引律進(jìn)行驗證評估，也可以為觀測大氣的區(qū)域規(guī)劃提供依據(jù)；同時，根據(jù)測軌數(shù)據(jù)及姿態(tài)實測數(shù)據(jù)，計算出實際的高度，可以對任務(wù)結(jié)果進(jìn)行驗證評估。需要注意的是，在本文的仿真算例中，僅考慮了衛(wèi)星沿一個方向(向地心方向)的姿態(tài)機(jī)動，在實際任務(wù)中，為了提高觀測效率，要求衛(wèi)星姿態(tài)在其他方向(偏航)也要進(jìn)行機(jī)動。這將是后續(xù)研究的重點，通過引進(jìn)優(yōu)化算法使衛(wèi)星在臨邊觀測任務(wù)中實現(xiàn)高精度和高效率的觀測。

References)

[1]Fleig A, Rault D F. Ozone vertical profiles in the upper troposphere and stratosphere from the OMPS limb sensor[C]//Proceedings of SPIE Remote Sensing 2011. Bellingham，WA：SPIE，2011

[2]Rault D F, Spurr R. The OMPS limb profile instrument two-dimensional retrieval algorithm [C]//Proceedings of SPIE Remote Sensing 2010. Bellingham，WA：SPIE，2010

[3]Dittman M G, Leitch J W, Chrisp M, et al. Limb broad-band imaging spectrometer for the NPOESS ozone mapping and profile suite (OMPS)[C]//Proceedings of SPIE 4814， Earth Observing System Ⅶ，120. Bellingham，WA：SPIE， 2002: 120-130

[4]劉林.航天器軌道理論[M].北京：國防工業(yè)出版社，2000

Liu Lin. Spacecraft orbit theory [M]. Beijing: National Defense Industry Press, 2000 (in Chinese)

[5]楊嘉墀.航天器軌道動力學(xué)與控制[M].北京：中國宇航出版社，1995

Yang Jiachi. Spacecraft orbit dynamics and control [M]. Beijing: China Astronautics Press, 1995 (in Chinese)

[6]Kyle T, Alfriend H S, Dong W G. Gravitational perturbations nonlinearity and circular orbit assumption effects on formation flying control strategies [J]. Advances in the Astronautical Sciences, 2000, 104(6):139-158

[7]于鳳軍. 地球自轉(zhuǎn)、形狀和重力加速度隨緯度的變化[J]. 大學(xué)物理，2013，32(6):14-17

Yu Fengjun. The rotation and shape of earth and the gravity acceleration’s variety with latitude [J]. College Physics，2013，32(6): 14-17 (in Chinese)

[8]Janz S J, Hilsenrath E, Flittner D E, et al. Rayleighscattering attitude sensor[C]//Proceedings of SPIE’s 1996 International Symposium on Optical Science, Engineering and Instrumentation. Bellingham，WA：SPIE， 1996: 146-153

[9]Hilsenrath E, Bhartia P K. Global Environmental Monitoring Satellite (GEMS): a pathfinder to measure gas constituents and aerosols in the earth’s lower atmosphere [C]//Proceedings of SPIE’s 1996 International Symposium on Optical Science, Engineering and Instrumentation. Bellingham，WA：SPIE，1996：224-330

[10]楊穎, 王琦. STK在計算機(jī)仿真中的應(yīng)用[M].北京: 國防工業(yè)出版社，2005

Yang Ying, Wang Qi. STK applied in computer simulation [M]. Beijing: National Defense Industry Press, 2005 (in Chinese)

[11]章仁為.衛(wèi)星軌道姿態(tài)動力學(xué)與控制[M].北京: 北京航空航天大學(xué)出版社，1998

Zhang Renwei. Satellite orbit and attitude dynamics and control [M]. Beijing: Beihang University Press, 1998 (in Chinese)

[12]Curtis F Gerald, Patrick O Wheatley.應(yīng)用數(shù)值分析[M].呂淑娟，譯.北京：機(jī)械工業(yè)出版社，2006

Curtis F G, Patrick O W. Applied numerical analysis [M]. Lv Shujuan, translated. Beijing: China Machine Press, 2006 (in Chinese)

(編輯：夏光)

Simulation Analysis Method on Satellite Attitude Angle in Earth Limb Observation

ZHU Yongsheng HU Haiying ZHENG Zhenzhen SU Ruifeng ZHANG Keke

(Shanghai Microsatellite Engineering Center, Shanghai 201203， China)

A simulation analysis method using STK is presented to calculate the attitude of a satellite for limb observation. The STK simulation scenario is created by connecting STK software. The two-body model is utilized to perform an analytical solution, which provides a reasonable initial guess for subsequent accurate iterative calculation. The WGS84 model is used to improve the initial solution by calling the HPOP (high precision orbit propagator) and the iterative formula. Taking the 20km and 100km limb observation for example, both the proposed method and the two-body method are employed to calculate the attitude angle. And the results show that the method proposed in this paper has faster convergence and higher accuracy.

earth limb observation； satellite attitude angle； high precision calculation

2014-05-19；

2015-06-18

上海市科委揚帆計劃(14YF1413600)

朱永生，男，碩士，從事航天器動力學(xué)與控制研究工作。Email：juesin@163.com。

V412.4

A

10.3969/j.issn.1673-8748.2015.04.007