A long-term prediction method for navigation satellite orbits

Abstract

The present invention relates to a method for long-term orbit prediction of a navigation satellite in a medium-orbit orbit, comprising: S100, obtaining a solar pressure parameter time series and a satellite orbit initial state sequence through orbit fitting; S200, analyzing the solar pressure parameter time series, establishing a pressure parameter model, and obtaining pressure parameter model coefficients; S300, performing long-term orbit prediction using the satellite orbit initial state sequence and the pressure parameter model coefficients; and S400, comparing the predicted orbit with a precise orbit to obtain an accuracy assessment of the predicted orbit. The orbit prediction method of the present invention can achieve high-precision long-term orbit predictions applicable to 30 to 90 days.

Description

Technical Field

[0001] The present invention relates to the field of positioning technology, and in particular to a method for long-term prediction of a navigation satellite orbit. Background Art

[0002] Navigation satellites use predicted broadcast ephemeris to provide positioning services to users on the Earth's surface and near-Earth, thereby achieving high-precision delivery of spatial references. Navigation satellite orbit forecasts are a crucial component of navigation service data processing. Currently, the update cycle for navigation satellite broadcast ephemeris is hourly. For example, GPS uses daily injections every two hours, while Beidou uses hourly injections. Therefore, orbit forecasts in navigation service data processing only focus on short-term orbit prediction accuracy, typically on a daily or even hourly basis. The three-dimensional position error of short-term orbit forecasts is on the decimeter level. Currently, demand for long-term, high-precision predictions of navigation satellite orbits is low. Typically, long-term predictions of navigation satellite orbits are primarily used for rapid signal acquisition by ground receivers and constellation orbit maintenance, requiring kilometer-level orbit accuracy.

[0003] To enhance the autonomy of the space segment and user terminals, research on long-term orbit predictions tailored to navigation service needs has been conducted. O. Montenbruck and other researchers at the German Aerospace Center (DLR) studied data processing methods for 14-day orbit predictions. The paper proposed using two-day arc-length rectangular orbit data for orbit fitting and prediction. The estimated dynamic parameters included the initial orbital state (position and velocity) and three parameters of the ECOM-1 solar pressure model. The orbit results obtained using this orbit prediction method were compared with the precise orbits. The radial orbit error (RMS) on the 14th day was better than 2 meters, and the orbit signal-in-space accuracy (SISREorb) was better than 30 meters. This method provides a feasible data processing method for 14-day orbit prediction. Essentially, it uses short-arc (2-day) dynamic characteristics to model satellite motion over a longer period. Evaluating orbit accuracy shows that the radial orbit and signal-in-space accuracy decrease rapidly from one to 14 days. Considering that the on-orbit dynamics of navigation satellites are influenced not only by perturbations such as Earth's gravity and the gravitational pull of the sun and moon, which are relatively easy to accurately model, but also by significant influences from solar radiation pressure, it is difficult to model the long-term variations in solar radiation pressure using the aforementioned method and relatively short arc segments, making it impossible to achieve orbit predictions for longer periods (over 14 days). There is currently no effective solution for long-term orbit predictions with accuracies on the order of hundreds of meters, such as for 30, 60, or 90 days. Summary of the Invention

[0004] In order to solve the problem of low orbit prediction accuracy of more than 14 days in the above-mentioned existing technology, the present invention provides a long-term orbit prediction method for navigation medium orbit satellites (MEO), which can achieve high-precision long-term orbit prediction for 30 to 90 days.

[0005] In order to achieve the above-mentioned object of the invention, the technical solution of the present invention is:

[0006] The present invention provides a method for long-term prediction of a navigation satellite orbit, comprising:

[0007] S100, obtaining a solar pressure parameter time series and a satellite orbit initial state sequence through orbit fitting;

[0008] S200, analyzing the solar light pressure parameter time series, establishing a solar light pressure parameter model, and obtaining solar light pressure parameter model coefficients;

[0009] S300, performing long-term orbit prediction using the satellite orbit initial state sequence and the light pressure parameter model coefficient;

[0010] S400: Compare the predicted orbit with the precise orbit to obtain an accuracy assessment of the predicted orbit.

[0011] According to one aspect of the present invention, the S100 includes:

[0012] S110, using the five-parameter ECOM model to describe the solar pressure perturbation;

[0013] S120, using the post-precision ephemeris as a pseudo-observation value, sliding back one day to generate a set of observations with a certain arc length;

[0014] S130, using a dynamic fitting method to obtain solar pressure parameters and initial state of the satellite orbit;

[0015] S140 , repeating S120 to S130 to obtain the solar pressure parameter time series and satellite orbit initial state sequence for each day of the year.

[0016] According to one aspect of the present invention, the initial state of the satellite orbit includes an initial position and an initial velocity of the satellite orbit.

[0017] According to one aspect of the present invention, the initial state of the satellite orbit is obtained by fitting a 10-day arc length orbit.

[0018] According to one aspect of the present invention, the S200 includes:

[0019] S201, using Fourier series and a method based on the angle between the sun and the satellite orbital plane to establish a solar pressure parameter model for each of the solar pressure parameter time series;

[0020] S202, eliminating the data during the Earth's shadow based on the angle between the Sun and the satellite's orbital plane;

[0021] S203. Fitting the light pressure parameter model coefficients using the least square method based on the eliminated data.

[0022] According to one aspect of the present invention, in S201, the solar pressure parameter model established by the Fourier series method for the solar pressure parameter time series adopts the following second-order Fourier series model:

[0023]

[0024] Among them, a0 is a constant term, a n is the coefficient of the n-order cosine periodic term, b n is the coefficient of the nth-order sinusoidal periodic term, P is the orbital period, and the orbital periods are years and half a year respectively.

[0025] According to one aspect of the present invention, in S201, based on the angle β between the sun and the satellite orbital plane, the calculation is performed using the following formula:

[0026]

[0027] in, is the satellite position, is the satellite speed, is the sun position;

[0028] In S201, the solar pressure parameter model established for the solar pressure parameter time series using a method based on the angle between the sun and the satellite orbital plane is:

[0029] F(β)=A+B sinβ+C / sinβ+D cosβ

[0030] Among them, A is the constant term, B is the coefficient of the solar altitude angle sine period term, C is the coefficient of the solar altitude angle sine inverse period term, and D is the coefficient of the solar altitude angle chord period term.

[0031] According to one aspect of the present invention, the S300 includes:

[0032] S301, reading the solar pressure parameter model coefficients and the solar pressure parameters obtained by orbit fitting, performing orbit integration on the initial state of the satellite orbit using a dynamic model, and obtaining an orbit prediction result for the first day;

[0033] S302: Read the solar pressure parameter model coefficient and the orbit forecast result of the previous day, calculate the solar pressure parameter, determine the initial state of the orbit of the current day based on the orbit forecast result of the previous day, and perform orbit integration on the initial state of the orbit of the current day using a dynamic model;

[0034] S303: Repeat S302 until all orbit predictions are completed.

[0035] According to one aspect of the present invention, the precise orbit in S400 is obtained by precise ephemeris orbit fitting.

[0036] According to one aspect of the present invention, in S400, the accuracy assessment of the predicted orbit includes: Cartesian coordinate accuracy assessment and orbit plane orientation accuracy assessment.

[0037] The three-dimensional rectangular coordinates of the predicted orbit and the precise orbit are compared to obtain the Cartesian coordinate accuracy assessment, as follows:

[0038] dX(t)=X pre (t)-X sp3 (t)

[0039] Among them, X pre (t) is the three-dimensional rectangular coordinate of the predicted orbit, X sp3 (t) is the three-dimensional rectangular coordinate of the precision track, and dX(t) is the difference between the two;

[0040] The three-dimensional rectangular coordinates of the predicted orbit and the precise orbit are converted into the orbital plane coordinate system and then compared to obtain the orbital plane orientation accuracy assessment, as follows:

[0041] dRTN(t)=RNT pre (t)-RTN sp3 (t)

[0042]

[0043] RNT pre (t)=(e x ,e y ,e z )X pre (t)

[0044] RNT sp3 (t)=(e x ,e y ,e z )X sp3 (t)

[0045] in, is the satellite position, is the sun's position, (e x ,e y ,e z ) is the rotation matrix for transforming the three-dimensional rectangular coordinate system into the orbital plane coordinate system, RNT pre (t) is the orbital plane coordinate of the predicted orbit, RNT sp3 (t) is the orbital plane coordinate of the precision orbit, and dRNT(t) is the difference between the two.

[0046] Compared with the prior art, the present invention has the following beneficial effects:

[0047] According to the solution of the present invention, compared with the existing technology, the long-term prediction method of MEO orbit takes into account the time domain characteristics of solar pressure, performs time domain characteristic analysis on the solar pressure parameters of MEO satellites, obtains the annual variation characteristics of the solar pressure parameters of MEO satellites, that is, the annual variation law of the solar pressure parameters, and reduces the error of the solar pressure parameter dynamic model by modeling the solar pressure parameters, thereby effectively improving the accuracy of long-term orbit prediction. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] To more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be derived from these drawings without inventive effort.

[0049] Figure 1 A flowchart schematically illustrates an implementation method for long-term prediction of a navigation satellite orbit provided by an embodiment of the present invention;

[0050] Figure 2 The specific implementation process of a method for long-term prediction of a navigation satellite orbit provided by an embodiment of the present invention is schematically shown;

[0051] Figure 3 Schematically showing the comparison results of the solar pressure parameter time series calculated by the BeiDou C19 satellite using the Fourier series and the β-angle-based solar pressure parameter model provided by the embodiment of the present invention and the solar pressure parameter time series estimated by orbit fitting;

[0052] Figure 4 The following schematically shows the orbit prediction result of the BeiDou C19 satellite provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0053] The description of the embodiments in this specification should be combined with the corresponding drawings, which should be considered a complete part of this specification. In the drawings, the shapes and thicknesses of the embodiments may be exaggerated and indicated for simplicity or convenience. Furthermore, the various structural components in the drawings will be described separately. It is worth noting that components not shown in the drawings or not described in words are known to those of ordinary skill in the art.

[0054] The description of the embodiments herein and any references to directions and orientations are for ease of description only and are not to be construed as limiting the scope of the present invention. The following description of the preferred embodiments may involve combinations of features, which may exist independently or in combination. The present invention is not specifically limited to the preferred embodiments. The scope of the present invention is defined by the claims.

[0055] Considering that the on-orbit dynamic characteristics of navigation satellites are affected by the perturbations of the Earth's gravity, the Sun's gravity, and the Moon's gravity, which are relatively easy to accurately model, they are also greatly affected by solar pressure. Existing methods basically use 2-day short-arc dynamic characteristics for modeling to predict the satellite motion status for a longer period of time. However, it is difficult to model the long-term variation of solar pressure using shorter arc segments, and it is difficult to predict orbits with an arc length of more than 14 days and obtain orbit products with higher accuracy. To this end, the present invention proposes a long-term MEO orbit prediction method that can be applied for 30 to 90 days and takes into account the annual variation characteristics of solar pressure, and effectively solves the problem of low accuracy of orbit prediction for more than 14 days in the existing technology.

[0056] like Figure 1 and Figure 2 As shown, the processing steps of the MEO orbit long-term prediction method disclosed in this embodiment include the following:

[0057] S100. Obtaining a solar pressure parameter time series and a satellite orbit initial state sequence through orbit fitting.

[0058] In step S100, long-term orbit products, such as IGS post-precision ephemeris, are used to perform orbit fitting on a daily basis to estimate the orbital elements and solar pressure parameters of each MEO satellite, and obtain a time series of solar pressure parameters for each day of a certain arc length, such as a whole year.

[0059] According to one embodiment of the present invention, the specific implementation process of obtaining the solar pressure parameter time series and the MEO satellite orbit initial state sequence through orbit fitting in step S100 includes: S110, using the five-parameter ECOM (EmpiricalCODE Orbit Model) model to describe the solar pressure perturbation; S120, using the post-precision ephemeris as pseudo-observation values, generating a set of observations with a certain arc length each day; S130, using a dynamic fitting method to obtain the solar pressure parameters and the MEO satellite orbit initial state; S140, repeating steps S120 to S130 to obtain the solar pressure parameter time series and the MEO satellite orbit initial state sequence for each day throughout the year. Specifically, the MEO satellite orbit initial state includes the initial position and initial velocity of the MEO satellite orbit, which are obtained from the 10-day arc length orbit fitting.

[0060] S200: Analyze the solar pressure parameter time series, establish a solar pressure parameter model, and obtain solar pressure parameter model coefficients. The model coefficients are used for subsequent orbit prediction.

[0061] According to one embodiment of the present invention, the specific implementation process of analyzing the solar light pressure parameter time series, establishing the solar light pressure parameter model, and obtaining the solar light pressure parameter model coefficients in step S200 includes: S201, using Fourier series and a method based on the angle between the sun and the MEO satellite orbital plane to establish the solar light pressure parameter model for the solar light pressure parameter time series; S202, eliminating data during the Earth's shadow based on the angle between the sun and the MEO satellite orbital plane (e.g., setting a elimination threshold of 10°); and S203, fitting the solar light pressure parameter model coefficients using the least squares method based on the eliminated data. This can reduce the error of the solar light pressure parameter dynamic model.

[0062] Specifically, in step S201, the solar pressure parameter model established by the Fourier series method for the solar pressure parameter time series adopts the following second-order Fourier series model:

[0063]

[0064] Among them, a0 is a constant term, a n is the coefficient of the n-order cosine periodic term, b n is the coefficient of the nth-order sinusoidal periodic term, P is the annual period of the orbital node, and the annual periods of the orbital node are years and half a year respectively.

[0065] Based on the angle β (i.e., solar altitude angle) between the sun and the MEO satellite orbital plane, it is calculated using the following formula:

[0066]

[0067] in, is the MEO satellite position, is the MEO satellite speed, is the sun position;

[0068] In step S201, the solar pressure parameter model established for the solar pressure parameter time series is as follows:

[0069] F(β)=A+B sinβ+C / sinβ+D cosβ

[0070] Among them, A is the constant term, B is the coefficient of the solar altitude angle sine period term, C is the coefficient of the solar altitude angle sine inverse period term, and D is the coefficient of the solar altitude angle chord period term.

[0071] According to one embodiment of the present invention, the basis for removing the data during the Earth's shadow is that the angle β between the sun and the MEO satellite orbital plane is 14.5 degrees.

[0072] S300. Use the initial state sequence of the MEO satellite orbit and the light pressure parameter model coefficient to perform long-term orbit prediction.

[0073] According to one embodiment of the present invention, the specific implementation process of using the MEO satellite orbit initial state sequence and the light pressure parameter model coefficients for long-term orbit prediction in step S300 includes: S301, reading the light pressure parameter model coefficients and the solar light pressure parameters obtained by orbit fitting, using the dynamic model to perform orbit integration on the MEO satellite orbit initial state, and obtaining the orbit prediction result for the first day; S302, reading the light pressure parameter model coefficients and the orbit prediction result of the previous day, obtaining the solar light pressure parameters by calculation, determining the initial state of the orbit of the day according to the orbit prediction result of the previous day, and using the dynamic model to perform orbit integration on the initial state of the orbit of the day; S303, repeating step S302 until all orbit predictions are completed.

[0074] The long-term orbit prediction process described above first reads the initial orbit state and the solar pressure parameter model coefficients. The initial orbit state of a MEO satellite is derived from a 10-day arc-length orbit fit, which is consistent with the long-term orbit prediction based on the solar pressure parameters of the 10-day arc-length orbit fit. The solar pressure parameter model coefficients are derived from the above modeling results. Based on the satellite list and dynamic configuration specified by the configuration parameters, the long-arc orbit dynamics fit is performed to complete the long-term orbit prediction.

[0075] S400: Compare the predicted orbit with the precise orbit to obtain an accuracy assessment of the predicted orbit.

[0076] According to one embodiment of the present invention, the specific implementation process of using the MEO satellite orbit initial state sequence and the light pressure parameter model coefficients to perform long-term orbit prediction in step S400 includes: comparing the predicted orbit with the three-dimensional rectangular coordinates of the precise orbit to obtain a Cartesian coordinate accuracy assessment, specifically as follows:

[0077] dX(t)=X pre (t)-X sp3 (t)

[0078] Among them, X pre (t) is the three-dimensional rectangular coordinate of the predicted orbit, X sp3 (t) is the three-dimensional rectangular coordinate of the precision track, and dX(t) is the difference between the two;

[0079] The three-dimensional rectangular coordinates of the predicted orbit and the precise orbit are converted into the orbital plane coordinate system and then compared to obtain the orbital plane orientation accuracy evaluation, as follows:

[0080] dRTN(t)=RNT pre (t)-RTN sp3 (t)

[0081]

[0082] RNT pre (t)=(e x ,e y ,e z )X pre (t)

[0083] RNT sp3 (t)=(e x ,e y ,e z )X sp3 (t)

[0084] in, is the MEO satellite position, is the sun's position, (e x ,e y ,e z ) is the rotation matrix for transforming the three-dimensional rectangular coordinate system into the orbital plane coordinate system, RNT pre (t) is the orbital plane coordinate of the predicted orbit, RNT sp3 (t) is the orbital plane coordinate of the precision orbit, and dRNT(t) is the difference between the two.

[0085] Specifically, the above precise orbit is obtained by precise ephemeris orbit fitting.

[0086] Figure 3 Comparison of the solar pressure parameter time series for the BeiDou C19 satellite using the Fourier series model and the β-angle model, respectively, with those estimated using orbital fitting. The agreement between the two indicates model validity. The black circles represent the solar pressure parameters estimated using a 7-day arc-length orbital fitting, while the black dots represent the values ​​calculated using the solar pressure model. The five sub-figures on the left represent the values ​​calculated using the Fourier series model, and the five sub-figures on the right represent the values ​​calculated using the β-angle model. The horizontal axis on the left is time, and the horizontal axis on the right is solar altitude.

[0087] from Figure 3 As can be seen from the figure, the solar pressure parameters obtained by orbital fitting exhibit distinct periodic variations. Both the Fourier series model and the β-angle model can be used to model the solar pressure parameters. The relationship between the fitted parameters and the β-angle indicates that the fitting of the solar pressure parameters becomes unstable when the solar altitude angle is small. To address this situation, neither model considers solar altitudes below a certain threshold (e.g., 10°).

[0088] Based on the light pressure parameters obtained from the light pressure parameter model, the orbit is predicted for 90 days starting from February 4, 2021, and compared with the model of using 10-day arc length to solve the light pressure parameters to predict the 90-day orbit in the same period.

[0089] Compared with the 10-day arc length model, the present invention improves the accuracy of orbit prediction results for most MEO satellites. Figure 4 This is the C19 satellite orbit prediction result. Figure 4 The horizontal axis is days, and the vertical axis is the difference between the orbit prediction result and the precise ephemeris. The left column shows the orbit prediction results using the β-angle model to calculate solar pressure, while the right column shows the orbit prediction results using Fourier series modeling. 10d represents the 90-day orbit prediction result using a 10-day arc length solution to calculate solar pressure parameters. As can be seen from the figure, both the Fourier series modeling and the β-angle model can produce highly accurate orbit predictions, and both significantly improve orbit accuracy compared to the predictions using the 10-day arc length solution to calculate solar pressure parameters.

[0090] The serial numbers of the above-mentioned steps involved in the method of the present invention do not mean the order of execution of the method. The execution order of each step should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiment of the present invention.

[0091] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A method for long-term prediction of a navigation satellite orbit, comprising: S100, obtaining a solar pressure parameter time series and a satellite orbit initial state sequence through orbit fitting; S200, analyzing the solar light pressure parameter time series, establishing a solar light pressure parameter model, and obtaining solar light pressure parameter model coefficients, wherein S200 includes: S201, using Fourier series and a method based on the angle between the sun and the satellite orbital plane to establish a solar pressure parameter model for each of the solar pressure parameter time series; S202, eliminating the data during the Earth's shadow based on the angle between the Sun and the satellite's orbital plane; S203, fitting the light pressure parameter model coefficients using the least squares method based on the eliminated data; S300, using the satellite orbit initial state sequence and the light pressure parameter model coefficient to perform long-term orbit prediction, S300 includes: S301, reading the solar pressure parameter model coefficients and the solar pressure parameters obtained by orbit fitting, performing orbit integration on the initial state of the satellite orbit using a dynamic model, and obtaining an orbit prediction result for the first day; S302: Read the solar pressure parameter model coefficient and the orbit forecast result of the previous day, calculate the solar pressure parameter, determine the initial state of the orbit of the current day based on the orbit forecast result of the previous day, and perform orbit integration on the initial state of the orbit of the current day using a dynamic model; S303, repeating S302 until all orbit predictions are completed; S400: Compare the predicted orbit with the precise orbit to obtain an accuracy assessment of the predicted orbit. The accuracy assessment of the predicted orbit includes: Cartesian coordinate accuracy assessment and orbit plane orientation accuracy assessment. The three-dimensional rectangular coordinates of the predicted orbit and the precise orbit are compared to obtain the Cartesian coordinate accuracy assessment, as follows: dX(t)=X pre (t)-X sp3 (t) Among them, X pre (t) is the three-dimensional rectangular coordinate of the predicted orbit, X sp3 (t) is the three-dimensional rectangular coordinate of the precision track, and dX(t) is the difference between the two; The three-dimensional rectangular coordinates of the predicted orbit and the precise orbit are converted into the orbital plane coordinate system and then compared to obtain the orbital plane orientation accuracy assessment, as follows: dRTN(t)=RNT pre (t)-RTN sp3 (t) RNT pre (t)=(and x ,and y ,and z )X pre (t) RNT sp3 (t)=(and x ,and y ,and z )X sp3 (t) in, is the satellite position, is the sun's position, (e x ,e y ,e z ) is the rotation matrix for transforming the three-dimensional rectangular coordinate system into the orbital plane coordinate system, RNT pre (t) is the orbital plane coordinate of the predicted orbit, RNT sp3 (t) is the orbital plane coordinate of the precision orbit, and dRNT(t) is the difference between the two.

2. The method according to claim 1, characterized in that The S100 includes: S110, using the five-parameter ECOM model to describe the solar pressure perturbation; S120, using the post-precision ephemeris as pseudo-observation values, and sliding back one day at a time to generate a set of observations with a certain arc length; S130, using a dynamic fitting method to obtain solar pressure parameters and initial state of the satellite orbit; S140 , repeating S120 to S130 to obtain the solar pressure parameter time series and satellite orbit initial state sequence for each day of the year.

3. The method according to claim 2, characterized in that The satellite orbit initial state includes the initial position and initial velocity of the satellite orbit.

4. The method according to claim 2, characterized in that The initial state of the satellite orbit is obtained by fitting a 10-day arc length orbit.

5. The method according to claim 1, wherein In S201, the light pressure parameter model established for the solar light pressure parameter time series using the Fourier series method adopts the following second-order Fourier series model: Among them, a0 is a constant term, a n is the coefficient of the n-order cosine periodic term, b n is the coefficient of the nth-order sinusoidal periodic term, P is the orbital period, and the orbital periods are years and half a year respectively.

6. The method according to claim 1, characterized in that In S201, based on the angle β between the sun and the satellite orbital plane, the following formula is used for calculation: in, is the satellite position, is the satellite speed, is the sun position; In S201, the solar pressure parameter model established for the solar pressure parameter time series using a method based on the angle between the sun and the satellite orbital plane is: F(β)=A+Bsinβ+C / sinβ+Dcosβ Among them, A is the constant term, B is the coefficient of the solar altitude angle sine period term, C is the coefficient of the solar altitude angle sine inverse period term, and D is the coefficient of the solar altitude angle chord period term.

7. The method according to claim 1, characterized in that The precise orbit in S400 is obtained by precise ephemeris orbit fitting.

Images