阅读说明：本技术 一种导航卫星定轨跟踪站优选算法 (Optimal selection algorithm for orbit determination tracking station of navigation satellite ) 是由 张睿 涂锐 卢晓春 韩军强 范丽红 张鹏飞 洪菊 刘金海 王星星 于 2020-06-16 设计创作，主要内容包括：本发明涉及卫星轨道计算技术领域,具体是涉及一种导航卫星定轨跟踪站优选算法,利用指标对不同跟踪站分布优劣以及定轨跟踪站数量饱和度进行定量分析,为各导航系统卫星定轨跟踪站选取以及核心基准站选址提供技术支持；本发明可在无地面跟踪站观测数据的前提下对导航卫星定轨地面跟踪站的位置以及数量进行优选,相较于已有的跟踪站分布优选算法,提高了计算效率,加入了跟踪站数量的优选算法,同时可以兼容多类观测值的融合分析；而且,本发明设计了的优选算法普遍适用于不同导航系统,同时适用于多类观测值融合处理时的定轨跟踪站优选。(The invention relates to the technical field of satellite orbit calculation, in particular to a navigation satellite orbit determination tracking station optimization algorithm, which utilizes indexes to carry out quantitative analysis on the distribution advantages and disadvantages of different tracking stations and the number saturation of orbit determination tracking stations, and provides technical support for the selection of the satellite orbit determination tracking stations and the site selection of a core reference station of each navigation system; the method can optimize the positions and the number of the orbit determination ground tracking stations of the navigation satellite on the premise of no observation data of the ground tracking stations, improves the calculation efficiency compared with the existing tracking station distribution optimization algorithm, adds the optimization algorithm of the number of the tracking stations, and can be compatible with the fusion analysis of various types of observation values; moreover, the optimization algorithm designed by the invention is generally suitable for different navigation systems and is also suitable for the optimization of the orbit determination tracking station during the fusion processing of multiple types of observed values.)

1. A navigation satellite orbit determination tracking station optimization algorithm is characterized by comprising the following steps:

s1: determination of distributed positions of tracking stations

S11: selecting a plurality of basic tracking station networks;

s12: verifying whether each satellite and each observation station in the basic tracking station network meet the observation requirements one by one, and establishing an observation equation for the observation meeting the requirements;

s13: when the optimal distribution of the ground tracking stations in the inter-satellite link state is analyzed, an inter-satellite range observation equation is established epoch by epoch according to an inter-satellite link routing algorithm;

s14: dividing grid points of 10 degrees multiplied by 10 degrees according to longitude and latitude in a global reasonable site selectable area to obtain m grid points;

s15: taking the three-dimensional coordinates of each grid point as the three-dimensional coordinates of the tracking station, verifying whether each satellite and each observation station in the basic tracking station network meet the observation requirements or not epoch by epoch, and establishing an observation equation for the observation meeting the requirements;

s16: combining the observation equations established in the steps S12, S13 and S15 aiming at each grid point to obtain a combined observation equation, and then calculating the DOP value of the kinetic parameter corresponding to each grid point; at the moment, m tracking station distribution schemes are obtained globally, wherein the tracking station distribution scheme with the minimum dynamic parameter DOP value is the optimal distribution of the tracking stations when 1 station is added;

s17: adding the selected 1 tracking station in the step S16 into a basic tracking station network, and screening out the optimal distribution when adding 2 tracking stations by using the processes from the step S2 to the step S6;

s18: obtaining the optimal distribution when the tracking stations are added into the n tracking stations after the steps from S12 to S17 are circulated for n times;

s2: tracking the number of stations for optimization

S21: obtaining DPDP value sequence (d) corresponding to optimal distribution scheme when n tracking stations corresponding to different arc segments are obtained₁，d₂，…，d_n)；

S22: averaging the DPDPDPOP values corresponding to the nth tracking station of all the arc segments to obtain an average value d of the convergence of the DPDPDP values of different arc segments_ave；

S23: and subtracting the adjacent DPDPDPDP values to obtain a change sequence of the slope of the DPDP value sequence, namely:

se_i＝d_i+1-d_i(i＝1，2，…，n-1) (14)

s24: calculating each DPDPDP value (i, d)_i) And DPDP value (n, d) when joining n tracking stations_n) The slope change between, i.e.:

s25: sn is to_iAnd se_iAnd correspondingly making differences between the two sequences, comparing the difference values with a threshold value, and taking the number of the newly added tracking stations as the optimal result of the number of the orbit determination tracking stations when the corresponding difference value of the number of the newly added tracking stations is smaller than the threshold value.

2. The method as claimed in claim 1, wherein in step S12, for the GNSS observation of the ground tracking station, in the case of studying only the distribution structure of the stations ignoring various errors, the observation model is represented as:

in the formula: rho represents the geometric distance between the satellite and the survey station, and i represents the satellite number; (x, y, z) represents the three-dimensional coordinates of the tracking station, t represents a time parameter, F_X、F_Y、F_ZOrbital integral functions, M, representing the three coordinate directions of the satellite, respectivelyⁱSatellite dynamics parameters representing an initial epoch of a satellite;

linearizing equation (3) to obtain:

in the formula:

。

3. the navigation satellite orbital determination tracking station optimization algorithm according to claim 2, wherein in the step S15, for each grid point, observation equations are established for each satellite satisfying the observation requirements and each station in the basic tracking station network based on the formula (4).

4. The method as claimed in claim 1, wherein in step S13, for the GNSS observation of the ground tracking station in the inter-satellite link state, the observation model of the inter-satellite ranging can be expressed as follows, under the condition that only the station distribution structure is studied by ignoring various errors:

in the formula: g_X、G_Y、G_ZOrbital integral function, M, representing three coordinate directions of satellite j, respectively^jRespectively representing the satellite dynamic parameters of the initial epoch of the satellite j;

linearize equation (8) to give:

in the formula:represents the initial value of the kinetic parameter at the initial moment of satellite j,

5. the method as claimed in claim 1, wherein in step S16, the calculation method of the DOP value of the dynamic parameter corresponding to each grid point is as follows:

based on the formula (4) and the formula (9), in the observation equations established by the two types of observation values, the parameters to be estimated are the kinetic parameters of each satellite, so that the observation equations established by the two types of observation values can be combined to obtain the observation equation:

in the formula: v representing the correction of the observed value and l representing the observation equationA constant term, a denotes the coefficient part of the observation equation,representing orbital dynamics parameters

From this, a co-factorial array of kinetic parameters can be obtained:

Q＝(A^TPA)^-1(12)

the DOP value of the satellite dynamics parameter is the arithmetic square root of the sum of diagonal elements corresponding to the dynamics parameter in the covariance matrix, and then the following parameters are obtained:

。

Technical Field

The invention relates to the technical field of satellite orbit calculation, in particular to an optimal selection algorithm for a fixed orbit tracking station of a navigation satellite.

Background

For any spacecraft, precise orbit determination is one of the primary tasks of stable operation of the spacecraft. For a Global Navigation Satellite System (GNSS), satellite orbit errors are one of main error sources of observation, and the navigation positioning accuracy of a user is directly influenced by the orbit accuracy, so the orbit accuracy of each navigation system is always a hot spot problem concerned by each navigation system operation department and vast users. Due to the high-precision orbit determination requirement of the GNSS navigation satellite, in practical application, an observation equation is mainly established by utilizing GNSS measurement data collected by a ground tracking station, and satellite dynamic parameters are solved by a parameter estimation method. In the resolving process, the geometric strength of the tracking station directly influences the resolving strength of the observation equation and then influences the final resolving precision of the orbit parameters, so that the distribution of the ground tracking station plays an important role in the satellite orbit determination of the navigation system.

The traditional orbit determination tracking station selection method is to select a certain number of tracking stations which are uniformly distributed on the whole world according to experience, but related methods for describing the advantages and disadvantages of the distribution of a certain tracking station through a quantitative method are few, and particularly for the current Beidou satellite navigation system (BDS), the satellite constellation comprises three types of satellites, the satellite constellation is not uniformly distributed on the whole world, so that the traditional method of uniformly selecting the orbit determination tracking stations on the whole world by depending on the experience is not suitable for the satellite orbit determination of the BDS.

Aiming at the problems, the invention carries out research on the optimal selection algorithm of the orbit determination tracking station of the navigation satellite, carries out quantitative analysis on the distribution advantages and disadvantages of different tracking stations and the number saturation of the orbit determination tracking station by using indexes, and provides technical support for the selection of the orbit determination tracking station of the satellite and the site selection of a core reference station of each navigation system.

Disclosure of Invention

In order to achieve the above purpose, the invention provides a navigation satellite orbit determination tracking station optimization algorithm, which utilizes indexes to carry out quantitative analysis on the advantages and disadvantages of the distribution of different tracking stations and the number saturation of the orbit determination tracking stations, and provides technical support for the selection of the satellite orbit determination tracking stations and the site selection of a core reference station of each navigation system, and the specific technical scheme is as follows:

firstly, obtaining the DOP value of the kinetic parameter

In the parameter estimation process, the DOP value of a parameter is often used to evaluate the accuracy of the parameter solution. In the navigation satellite orbit determination calculation based on the observation data of the ground tracking station, the DOP value of the kinetic parameter can be obtained by the following method.

For GNSS observations of ground tracking stations, the observation model can be simply expressed as:

in the formula: rho represents the geometric distance between the satellite and the survey station, and i represents the satellite number; (X)ⁱ，Yⁱ，Zⁱ) Representing the three-dimensional coordinates of the satellite, (x, y, z) representing the three-dimensional coordinates of the tracking station, and t representing a time parameter.

Based on the dynamic orbit determination method, the position of the satellite in each epoch can be represented by the dynamic parameters of the satellite in the initial epoch. Thus, the satellite position of a satellite at any time can be expressed as:

in the formula: f_x、F_Y、F_ZOrbit respectively representing three coordinate directions of satelliteIntegral function, MⁱA satellite dynamics parameter representing an initial epoch of a satellite.

The formula (2) is taken into the formula (1) to obtain:

linearizing equation (3) to obtain:

in the formula:initial values, p, representing kinetic parameters of the satellite at the initial time₀Can be expressed as:

aiming at other observation information, an observation equation established by the observation information and an observation equation established by GNSS observation information of a ground tracking station can be fused, and a dynamic parameter DOP value is analyzed together. In the case of studying only the geometry distribution structure ignoring various errors, the observation model for inter-satellite ranging can be simply expressed as:

in the formula: (X)^j，Y^j，Z^j) Respectively represent the three-dimensional coordinates of the satellite j, and the satellite position of the satellite j at any time can be respectively represented as

In the formula: g_X、G_Y、G_ZRespectively representing three coordinate squares of satellite jOrbital integral function of direction, M^jRespectively representing the satellite dynamics parameters of the initial epoch of the satellite j. By substituting formula (6) with formula (2) and formula (7), the compound is obtained

Linearize equation (8) to give:

in the formula:represents the initial value of the kinetic parameter at the initial moment of satellite j,

can be expressed as:

based on the formula (4) and the formula (9), in the observation equations established by the two types of observation values, the parameters to be estimated are the kinetic parameters of each satellite, so that the observation equations established by the two types of observation values can be combined. An observation equation is established based on the formula (4) and the formula (9), and the observation equation can be obtained

In the formula: v denotes the number of corrections of the observed value, l denotes a constant term of the observation equation, a denotes a coefficient part of the observation equation,

representing the orbital dynamics parameters. The co-factor matrix of the kinetic parameters can be expressed as

Q＝(A^TPA)^-1(12)

The DOP value of the satellite dynamics parameter is the arithmetic square root of the sum of diagonal elements corresponding to the dynamics parameter in the covariance matrix, i.e.

For different distribution schemes of the tracking stations, the distribution scheme of the tracking stations with smaller dynamic parameter DOP value is better.

Second, optimal selection algorithm of orbit determination tracking station of navigation satellite

1. Determination of distributed positions of tracking stations

11. And selecting a plurality of basic tracking station networks.

12. And verifying whether each satellite and each observation station in the basic tracking station network meet the observation requirements (if the observation altitude angle is more than 5 degrees) one by one, and establishing an observation equation for the observation meeting the requirements based on the formula (4).

13. When the optimal distribution of the ground tracking stations in the inter-satellite link state is analyzed, an inter-satellite ranging observation equation is established epoch by epoch based on the formula (9) according to an inter-satellite link routing algorithm (the algorithm is applied in engineering and is not deeply discussed in the patent).

14. In a global reasonable site selectable area, grid points of 10 degrees multiplied by 10 degrees are divided according to longitude and latitude to obtain m grid points.

15. And (3) taking the three-dimensional coordinate of each grid point as the three-dimensional coordinate of the tracking station, verifying whether each satellite and each observation station in the basic tracking station network meet the observation requirements (if the observation altitude angle is more than 5 degrees) or not by epoch, and establishing an observation equation for the observation meeting the requirements based on the formula (4).

16. And (3) combining the observation equations established in the steps (12), (13) and (15) aiming at each grid point to obtain a combined observation equation, and then calculating the DOP value of the kinetic parameter corresponding to each grid point according to the formulas (11), (12) and (13). At this time, m kinds of tracking station distribution schemes are obtained globally, wherein the tracking station distribution scheme with the minimum dynamic parameter DOP value is the optimal distribution of the tracking stations when 1 station is added.

17. And (4) adding the 1 selected tracking station in the step (16) into a basic tracking station network, and screening out the optimal distribution when the 2 tracking stations are added by utilizing the processes from the step (2) to the step (6).

18. And (5) obtaining the optimal distribution when the n tracking stations are added after n times of circulation of the steps 12 to 17.

2. Tracking the number of stations for optimization

21. Obtaining DPDP value sequence (d) corresponding to optimal distribution scheme when n tracking stations corresponding to different arc segments are obtained₁，d₂，…，d_n)。

22. Averaging the DPDPDPOP values corresponding to the nth tracking station of all the arc segments to obtain an average value d of the convergence of the DPDPDP values of different arc segments_ave。

23. And (3) carrying out subtraction on adjacent DPDPDPOP values, namely subtracting the previous DPOP value from the next DPOP value to obtain a change sequence of the slope of the DPOP value sequence, namely:

se_i＝d_i+1-d_i(i＝1，2，…，n-1) (14)

24. calculating each DPDPDP value (i, d)_i) And DPDP value (n, d) when joining n tracking stations_n) The slope change between, i.e.:

25. sn is to_iAnd se_iThe two sequences are correspondingly differenced, and the difference value is compared with a threshold (such as d)_aveAnd/n), when the corresponding difference value of the number of the newly increased tracking stations is smaller than the threshold value, taking the number of the newly increased tracking stations as the optimal result of the number of the orbit determination tracking stations.

Compared with the existing selection method of the distribution of the tracking stations, the method has the beneficial effects that:

(1) the method can optimize the positions and the number of the orbit determination ground tracking stations of the navigation satellite on the premise of no observation data of the ground tracking stations, improves the calculation efficiency compared with the existing tracking station distribution optimization algorithm, adds the optimization algorithm of the number of the tracking stations, and can be compatible with the fusion analysis of various types of observation values.

(2) The invention designs an optimal selection algorithm of the orbit determination tracking station of the navigation satellite, which is generally suitable for different navigation systems and is also suitable for optimal selection of the orbit determination tracking station during fusion processing of a plurality of types of observed values.

Drawings

FIG. 1 is a flow chart of a preferred algorithm of the tracking station of the present invention;

FIG. 2 is a distribution diagram of DPOP values when 1 ground station is added on the premise of only utilizing ground station tracking in the present invention;

FIG. 3 is a distribution diagram of DPOP values when 1 ground station is added on the premise of jointly fixing orbit by using the ground station and the inter-satellite link data in the invention;

FIG. 4 is a prophetic ground station optimization map for ground station tracking only according to the present invention;

FIG. 5 is a proposed lower ground station optimal distribution map for joint orbit determination using ground station and inter-satellite link data in accordance with the present invention;

fig. 6 shows DPDOP value sequences corresponding to the number of each added station when the present invention is observed jointly using only a ground station and using a ground station and an inter-satellite link.

Detailed Description

To further illustrate the manner in which the present invention is made and the effects achieved, the following description of the present invention will be made in detail and completely with reference to the accompanying drawings.