趙 雙，徐艷麗，戴樺宇，趙程亮

(1.中國人民解放軍裝備學院 研究生院，北京 101416；2.中國人民解放軍裝備學院 航天裝備系，北京 101416)

基于抽樣理論的改進型星座覆蓋分析方法

趙 雙1，徐艷麗2，戴樺宇1，趙程亮1

(1.中國人民解放軍裝備學院 研究生院，北京 101416；2.中國人民解放軍裝備學院 航天裝備系，北京 101416)

針對傳統(tǒng)網格點方法在分析星座對目標區(qū)域的覆蓋性能時，存在著計算效率不高的問題，利用抽樣理論中借助樣本統(tǒng)計量對總體參數進行估計的思想，對傳統(tǒng)網格點方法加以了改進;通過分析星座對目標區(qū)域內所有網格點的覆蓋情況與網格點經緯度之間的關系，確立了對所有網格點構成的總體采用先分層后隨機的抽樣方式，并分析了層的劃分和樣本量在各層的分配；根據不同指標的統(tǒng)計特征，給出了不同星座覆蓋性能指標的估計方法和不同精度條件下總樣本量的計算過程，并建立了改進網格點法的實施流程；最后利用改進方法計算了仿真時段內北斗衛(wèi)星導航系統(tǒng)對某一地面區(qū)域的平均PDOP值，通過與傳統(tǒng)網格點方法的仿真結果進行比較并對仿真數據進行分析，證明了該改進方法在保證計算精度的同時，有效地提升了計算效率。

網格點法；分層隨機抽樣；總體估計；誤差精度

0 引言

在衛(wèi)星星座設計與研究中，有一個很重要的問題，就是如何分析星座的覆蓋性能。現有的分析方法可以分為解析法和數值法兩種，其中解析法主要是根據衛(wèi)星星座的空間結構，利用幾何方法將覆蓋條件轉化為一系列的約束方程，通過對方程組的求解得到星座的覆蓋性能。典型的解析法包括外接圓方法[1]和覆蓋帶方法[2]，雖然解析法形象直觀、簡潔明了，但是對于由多顆衛(wèi)星組成，星座中包含多種軌道的衛(wèi)星星座而言，由于其空間幾何結構的復雜性，同時又要考慮到攝動、章動對星座構型的影響，使得解析法在解決此類問題時往往力不從心。

因此，在計算復雜構型的星座覆蓋性能時，一般采用數值法，最常用的就是網格點方法[3]。因為其簡單易于實現，便于統(tǒng)計各類覆蓋性能指標，且適用于任意類型的軌道和任意構型的星座而被廣泛應用于航天設計與應用中，但是由于網格點方法本身的特性，網格點方法存在著計算精度與計算效率之間的矛盾：當網格分辨率小時，計算精度高但計算效率低，當網格分辨率大時，計算精度差但計算效率高。

除上述幾種傳統(tǒng)方法之外，在進行星座覆蓋性能計算時，還存在著一些其他的方法。文獻[4]提出了地心法和視函數法。文獻[5]提出了經度條帶法。文獻[6]提出了計算覆蓋率的改進型網格點法，但上述方法都存在著適用范圍不廣，計算精度不夠的缺陷。本文針對目前航天設計與應用中使用最為廣泛的網格點法計算效率低的問題，利用抽樣理論的知識，提出了一種改進型的網格點方法，在保證計算精度的同時，提升了計算效率。

1 網格點方法

利用網格點方法分析星座對某一地面區(qū)域的覆蓋性能時，一般將該區(qū)域按照某種分辨率劃分為一定數量的網格，通過計算星座對每個網格中心點的覆蓋情況來表征整個網格的覆蓋情況，進而表征整個區(qū)域的覆蓋情況[7]。網格點方法的計算步驟如下:

1)劃分網格。網格點的劃分方法按照劃分規(guī)則的不同一般可分為按經緯度劃分、按距離劃分和按面積劃分。

2)計算覆蓋。對于星座中的任意一顆衛(wèi)星，計算并記錄它在整個仿真時段內每個時刻對每一個網格中心點的覆蓋情況。

3)統(tǒng)計結果。統(tǒng)計所有網格中心點在仿真時段內每一時刻的覆蓋情況，根據任務需要，用不同性能指標的統(tǒng)計值(平均值、最大值、最小值)表征星座對該區(qū)域的覆蓋情況。

2 基于抽樣理論的改進型網格點法

2.1 抽樣方法的選擇

抽樣方法是指從N個單元組成的總體中，抽取n個單元組成的樣本的方式。對于網格點方法而言，其中總體是指將目標區(qū)域按照某種分辨率劃分后得到的所有網格的集合，樣本是指從所有網格的集合中通過某種抽樣方法抽取出的若干個網格的集合?；镜某闃臃椒ㄓ泻唵坞S機抽樣、系統(tǒng)抽樣、分層抽樣、多階段抽樣等[8]。

針對星座覆蓋分析時不同網格點會因為地理經緯度的差別，而存在覆蓋情況上的差異，本文利用分層隨機抽樣的方法，對傳統(tǒng)的網格點法加以改進。該方法既保留了簡單隨機抽樣簡單直觀的特性，又保證了樣本中包含有各種特征的抽樣單元，使得樣本的結構與總體的結構相近，從而有效地提高了估計的精度。

2.2 不同星座覆蓋性能指標的估計方法

星座覆蓋性能指標是衡量星座性能最重要的指標，在分析星座對目標區(qū)域的覆蓋性能時，我們采用的指標一般有：覆蓋百分比、N重覆蓋率、平均覆蓋間隙、平均響應時間、平均覆蓋重數、覆蓋總時間、總覆蓋間隙等[9]。除此之外，針對不同功能和類型的星座，往往還采用其它的指標，如導航星座一般采用目標區(qū)域的PDOP值(position dilution of precision)三維位置精度因子，PDOP值不僅反映了覆蓋重數、覆蓋間隙等其他星座覆蓋性能指標的信息，還是衡量衛(wèi)星導航定位精度的關鍵參數[10]。

上述這些性能指標就是由所有網格所構成的總體的某項總體參數，根據這些指標的特點可以將它們劃分到總體均值、總體總值、總體比例中。例如平均覆蓋重數、平均覆蓋間隙、平均響應時間屬于總體均值；覆蓋總時間、總覆蓋間隔屬于總體總值；覆蓋百分比、N重覆蓋率屬于總體比例。與總體參數相對應，在分層隨機抽樣下，利用樣本統(tǒng)計量對總體參數進行估計的方法一般有：均值估計、總值估計、比例估計[11]。

(1)

總值估計。用樣本均值和總體單元數得到總值估計，對于分層隨機樣本，總體總量Yst的估計為：

(2)

2)比例估計。用樣本比例作為總體比例的估計，對于分層隨機樣本，總體比例Pst的估計為：

(3)

式中，ph表示第h層的樣本比例。

2.3 層的劃分和樣本量在各層的分配

2.3.1 分層原則與層數的確定

在分層隨機抽樣的過程中，首先就是要對總體中的所有單元進行分層。對于網格點法而言，最好的分層方法就是根據所有網格點的星座覆蓋性能指標的分布情況，對網格點進行分層，但是在抽樣調查之前，我們并不知道各網格點的覆蓋情況，因此只能根據網格點的其他特征，對網格點進行分層。

對于一個確定的衛(wèi)星星座而言，對地面不同網格點的覆蓋差異，由網格點的地理位置直接決定。圖1顯示了北斗衛(wèi)星導航系統(tǒng)某一時刻對某一地面區(qū)域內不同經緯度的地面點的PDOP值，可以看出隨著經緯度的變化，地面點的PDOP值也在變化，所以本文在設計分層隨機抽樣時，依據網格點經緯度的差異，對網格點進行分層。

圖1 在北斗衛(wèi)星導航系統(tǒng)下某一區(qū)域內不同地面點的PDOP值隨經緯度的變化

在實際抽樣調查中，要保證每個層中都存在樣本單元，而且為了給出估計量方差的無偏估計，一般層數不超過樣本量的一半。但是為了提高估計結果的準確性，一般對總體盡可能多的進行劃分，使層內差異降低。對于分層隨機抽樣，當層數增加到一定程度時，估計精度的收益會非常小。實際分層隨機抽樣中，一般層數不超過6[12]。

2.3.2 樣本量在各層的分配

對于分層隨機抽樣，估計量的精度不僅取決于劃分的層數，還取決于總樣本量在各層的分配方式。常見的樣本量在各層的分配方式有：常數分配、比例分配。常數分配是將總樣本量均勻的分配到各層中，比例分配則是使各層的樣本量與各層的層權成正比，這種分配方式可以使樣本均勻的散布在各層內。因為網格點方法是將地面區(qū)域按照某一分辨率均勻的劃分為一定個數的網格，所以在對地面區(qū)域按照經緯度分層之后，不同層內所包含的網格點與層的面積成正比，這時若采用比例分配的方式會保證抽樣所得樣本的結構與總體的結構更加接近，估計結果的精度更高。

2.4 給定精度條件下總樣本量的確定

在利用分層抽樣和簡單隨機抽樣的方法對傳統(tǒng)的網格點方法加以改進時，目的是為了提高網格點方法的計算效率，但同時也要保證計算結果的精度，而不同精度條件下所對應的總樣本量是不同的。

(4)

(5)

其中：u為標準正太分布的雙側α分位數，絕對誤差d和相對誤差r分別為：

(6)

令wh=nh/n，當估計精度的誤差給定時，分層抽樣中確定樣本量的一般公式為：

(7)

(8)

2.5 改進網格點法的實施流程

1)首先，根據任務要求的分辨率將目標區(qū)域劃分為N個網格點，作為抽樣的總體，分析所關心的覆蓋性能指標屬于總體參數(總體均值、總體比例、總體總量)中的哪一項。

2)其次，根據經緯度將總體劃分為若干層，確定樣本量在各層的分配方式(均勻分配、比例分配等)以及給定精度條件下所對應的總樣本量n。

3)然后，采用簡單隨機抽樣的方式從各層中抽取相應個數的單元作為此次抽樣的樣本，對樣本進行調查獲得樣本的統(tǒng)計量。

4)最后，根據覆蓋性能指標選擇相應的估計方法(均值估計、比例估計、總量估計)，利用已獲得的樣本統(tǒng)計量對總體參數進行估計，得到星座對目標區(qū)域的覆蓋性能指標。其大致流程如圖2所示。

圖2 基于抽樣理論的網格點改進方法的流程示意圖

3 仿真算例

仿真試驗中，星座為北斗衛(wèi)星導航系統(tǒng)，目標區(qū)域為北京地區(qū)(北緯39.4°～41.6°,東經115.7°～117.4°)，假設地面任務要求網格劃分的分辨率為0.005°，仿真時段為2020-01-01 00:00:00～2020-01-02:00:00:00，時間步長為300 s，計算目標區(qū)域的平均PDOP值。

3.1 仿真流程

首先，將目標區(qū)域按照經緯度均勻劃分為6層，如圖3。然后根據給定的誤差精度計算總樣本量和每層中需要抽取的網格數，目標區(qū)域的平均PDOP值屬于總體參數中的總體均值，最后利用總體均值的估計方法對目標區(qū)域的平均PDOP值進行估計。

圖3 將目標區(qū)域按照經緯度分層

3.2 仿真結果

在同樣計算資源和仿真條件下，傳統(tǒng)網格點方法和基于抽樣理論的改進網格點方法的仿真結果見表1、表2。

表1 傳統(tǒng)網格點方法的仿真結果

表2 改進方法的仿真結果

對上述離散的試驗數據進行擬合，繪制出平均PDOP值、理論誤差、計算耗時和效率提升隨抽樣網格數的變化圖像，如圖4。

圖4 平均PDOP值、理論誤差、計算耗時和效率提升隨抽樣網格數的變化圖像

通過上述的仿真結果可以看出，隨著抽樣網格數的增加，樣本的平均PDOP值越來越接近總體的PDOP 值；理論誤差最開始增加而陡然下降而后趨于平緩；計算耗時與效率提升基本與抽樣網格數成正比。如在計算效率提升94.18%的情況下，理論誤差僅為0.001%；隨著樣本量的變化，樣本統(tǒng)計量與總體參數之間的理論誤差始終在0.0003%～0.0040%之間變化，相比傳統(tǒng)方法而言，利用改進的方法在保證計算結果精度的條件下，顯著提高了計算效率。

3.3 誤差分析

根據式(7)可以得到不同樣本量下的絕對誤差和相對誤差：

(9)

本文在對總體分層時，按照經緯度將目標區(qū)域均勻的劃分成6層，樣本量在各層的分配方式為比例分配，可知wh=Wh=1/6，置信水平取95%，查表可得對應的u=1.96，θ為總體的平均PDOP值，利用各層的總體方差，繪制出絕對誤差和相對誤差隨抽樣網格數的變化圖像，如圖5。

圖5 絕對誤差和相對誤差隨抽樣網格數的變化圖像

從圖中曲線可知，絕對誤差和相對誤差的下降速率，隨著抽樣網格數的增加而逐漸減小。起初，絕對始誤差和相對誤差隨著抽樣網格數的增加而急劇下降；然后，當抽樣網格數達到一定數量時，誤差的變化趨于平緩；最后，隨著抽樣網格數的進一步增大，誤差的變化趨近于0。

可以看出當抽樣網格數增加到一定數量時，計算精度的收益會非常小，此時以提升計算精度為目的而犧牲計算效率，不僅會耗費大量的計算時間而且計算精度的提升也微乎其微。所以在實際工程應用中，就存在著一個計算精度和計算效率的平衡問題，本文提出的基于抽樣理論的改進型網格點法，就是尋找一種在保證計算精度的條件下最大程度的提升計算效率的方法。

4 結束語

在利用網格點方法長時間、頻繁地計算星座對不同區(qū)域的覆蓋性能時，為了保證計算結果的精度就需要耗費大量的計算時間，而這在實際工程應用中不僅會影響任務進展還會耗費計算資源。本文針對這一問題提出了基于抽樣理論的網格點改進方法，并結合抽樣理論的相關知識，給出了不同星座覆蓋性能-指標的估計方法和不同精度條件下總樣本量的計算方法,制定了改進型網格點方法的實施流程，最后通過仿真計算證明了該方法的準確性和有效性。

[1] Walker J G. Some circular orbit patterns providing continuous whole earth coverage[J]. Journal of the British Interplanetary Society, 1971, 24(1): 369-384.

[2]Rider L. Optimized polar orbit constellation for redundant earth coverage[J]. Journal of the Astronautical Sciences, 1985, 33: 147-161.

[3]Morrison J J. A system of sixteen synchronous satellites for world-wide navigation and surveillance, TSC-FAA-72-31[R]. [S.I.]: Transportation System Center Cambridge Mass, 1973.

[4]徐 敏，程鳳舟，陳士櫓. 星座覆蓋性能數值仿真[J]. 宇航學報，2000，11, 35: 11-17.

[5]宋志明，戴光明，王茂才，等. 衛(wèi)星星座區(qū)域覆蓋問題的快速仿真算法[A]. 航天控制，2014.10, 32(5): 65-76.

[6]秦睿杰，戴光明，王茂才，等. 一種計算星座區(qū)域覆蓋的高效抽樣網格點法[A]. 計算機應用研究，2015.4, 32(4): 1065-1073.

[7]翁慧慧. 遙感衛(wèi)星對地覆蓋分析與仿真[D]. 鄭州：解放軍信息工程大學, 2006.

[8]金勇進. 抽樣理論與應用[M]. 北京: 高等教育出版社, 2016.

[9]李金昌. 應用抽樣技術[M]. 北京: 科學出版社, 2016.

[10] 倪育德，陳 君，蔚保國，等. 基于STK的BDS星座仿真和性能分析[A]. 計算機測量與控制，2016，24(4): 281-288.

[11]張育林等. 衛(wèi)星星座理論與設計[M]. 北京：科學出版社，2008.

[12]金勇進. 抽樣技術[M]. 北京: 中國人民大學出版社, 2015.

Improved Constellation Coverage Analysis Method Based on Sampling Theory

Zhao Shuang1，Xu Yanli2，Dai Huayu1，Zhao Chengliang1

(1.Department of Graduate Management, Equipment Academy of PLA, Beijing 101416, China 2.Department of Space Equipment, Equipment Academy of PLA, Beijing 101416, China)

In view of the problem of low computational efficiency of traditional grid point method when analyze the constellation coverage performance of target area,the traditional grid point method is improved by using the idea of estimating population parameters through sample statistics in sampling theory.Through analyzing the relationship between the constellation coverage performance of all grid points in the target area and grid points’ latitude and longitude, the stratified random sampling method was established, the layer division method and the distribution way of sample size in each layer are researched.According to the statistical characteristics of different indexes,the improved grid point method gives the way to estimate different coverage performance indexes of constellation and calculate the number of grids with different precision, the process of the improved grid point method is also established. Finally, the average PDOP value of the target region of the Beidou satellite navigation system is calculated through the improved method, by comparing the simulation results with the traditional grid point method, it is proved that the improved method can significantly improve the computational efficiency while ensuring the accuracy of calculation.

grid point method; stratified random sampling; population estimation; precision error

2017-05-07；

2017-05-23。

趙 雙(1991-)，男，湖北洪湖人，碩士研究生，主要從事航天任務分析與設計方向的研究。

徐艷麗(1964-)，女，山西五臺人，副教授，碩士研究生導師，主要從事航天任務分析與設計方向的研究。

1671-4598(2017)12-0162-04

10.16526/j.cnki.11-4762/tp.2017.12.042

TP391.9

A