CN113189619B - Low-orbit constellation phase maintaining parameter estimation method - Google Patents

Low-orbit constellation phase maintaining parameter estimation method Download PDF

Info

Publication number
CN113189619B
CN113189619B CN202110362735.9A CN202110362735A CN113189619B CN 113189619 B CN113189619 B CN 113189619B CN 202110362735 A CN202110362735 A CN 202110362735A CN 113189619 B CN113189619 B CN 113189619B
Authority
CN
China
Prior art keywords
orbit
semi
ping
gen
relative
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110362735.9A
Other languages
Chinese (zh)
Other versions
CN113189619A (en
Inventor
李新刚
裴胜伟
齐彧
黄华
佟金成
侯凤龙
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Academy of Space Technology CAST
Original Assignee
China Academy of Space Technology CAST
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Academy of Space Technology CAST filed Critical China Academy of Space Technology CAST
Priority to CN202110362735.9A priority Critical patent/CN113189619B/en
Publication of CN113189619A publication Critical patent/CN113189619A/en
Application granted granted Critical
Publication of CN113189619B publication Critical patent/CN113189619B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13—Receivers
    • G01S19/14—Receivers specially adapted for specific applications
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42—Determining position
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/15—Correlation function computation including computation of convolution operations
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
    • Y02D30/00—Reducing energy consumption in communication networks
    • Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Theoretical Computer Science (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Computational Mathematics (AREA)
  • Algebra (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Computing Systems (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

A low orbit constellation phase maintaining parameter estimation method belongs to the technical field of satellite orbit control. According to the invention, the orbit semi-major axis Ping Gen is indirectly estimated by using the latitude amplitude angle measured value with higher precision without directly using the orbit semi-major axis measured value output by the satellite-borne GNSS orbit determination, so that the orbit determination precision of the existing GNSS receiver is not improved, the high precision requirement of satellite constellation maintenance on the orbit semi-major axis can be met, and the configuration of a double-frequency GNSS receiver with higher price on satellites is avoided; according to the method, polynomial fitting is carried out after the relative latitude amplitude angle measured value is converted into the flat root, then the change relation between the relative latitude amplitude angle flat root and the track semi-long axis flat root is utilized, polynomial fitting is carried out on the change trend of the track semi-long axis flat root, and the fitting coefficient is calculated, so that the track semi-long axis Ping Gen at any moment is estimated indirectly, the semi-long axis flat root precision is superior to 1m, and constellation phase maintenance can be completed based on a unilateral drift ring control strategy under the condition that any semi-long axis attenuation rate is met.

Description

Low-orbit constellation phase maintaining parameter estimation method
Technical Field
The invention relates to a low-orbit constellation phase maintaining parameter estimation method, and belongs to the technical field of satellite orbit control.
Background
The low-rail constellation needs to remain relatively configuration stable during operation, preventing changes in coverage characteristics to the ground. Limited by orbit determination accuracy and orbit control accuracy, the initial orbit parameters of satellites in orbit deviate from the design parameters, so that the relative orbit parameters among satellites are changed, and the relative configuration of the constellation is gradually destroyed; in addition, satellites are affected by various environmental perturbation forces during long-term operation, such as earth gravitational field perturbation, solar-lunar gravitational perturbation, weather resistance perturbation, light pressure perturbation and the like, and the perturbation forces cause the orbit parameters of the satellites to change and damage the constellation configuration. Therefore, a series of constellation phase maintaining tasks are required to be performed during the in-orbit operation to ensure that the constellation configuration meets the design requirement, namely, the relative latitude and amplitude angles between satellites and the relative ascending intersection point right-angle between different track surfaces are required to be maintained within a certain design range.
The phase maintaining strategy between co-orbit satellites generally takes a constellation reference orbit as a phase maintaining target, changes the relative orbit angular velocity of the satellite relative to the reference orbit by adjusting the semi-long axis Ping Gen of the satellite relative to the reference orbit, indirectly controls the relative latitude amplitude Ping Gen of the satellite relative to the reference orbit, and maintains the relative latitude amplitude flat root in a preset control box, namely, at any moment, the relative latitude amplitude flat root always does not exceed the boundary (relative latitude amplitude flat root maximum threshold) of the control box.
The more accurate the semi-long axis flat root control of the satellite relative to the reference orbit, the longer the period of constellation phase maintenance, the smaller the semi-long axis flat root that needs to be adjusted, and the less propellant is consumed.
In addition, in order to reduce the carrying cost, the low-orbit constellation generally adopts a one-arrow multi-satellite batch launching deployment scheme, so that the satellite outboard layout is very compact, the thrusters can be often only installed on one deck of the satellite, and because the main task of the thrusters during long-term on-orbit operation of the satellite is constellation phase maintenance, under the condition of limited layout, the thrusters are generally installed on the-X plane of the satellite (the +X plane in the flight speed direction in the case of three-axis on-orbit orientation of the satellite is regulated), and the thrusters are output along the +X direction to generate tangential thrusts along the flight direction. In order not to affect inter-satellite links and ground communication services during constellation phase maintenance, the same orbit phase maintenance strategy generally adopts single-side drift ring control, a satellite adjusts a semi-long axis only when the relative latitude amplitude root of a relative reference orbit reaches the right boundary of a control box, and the relative latitude amplitude root of the satellite gradually drifts to the right boundary of the control box after drifting to the left boundary of the control box to the minimum in the process that the semi-long axis gradually decreases under the action of atmospheric resistance, so that the single-side drift ring is formed. Because the solar activity affects the atmospheric density, the attenuation speed of the semi-long axis also changes correspondingly, when the attenuation speed of the semi-long axis is smaller, the excessive adjustment quantity of the semi-long axis flat root of the satellite relative to the reference orbit can cause the latitude amplitude angle Ping Gen of the satellite relative to the reference orbit to drift out of the left boundary of the control box, and the phase maintaining precision among the co-orbit satellites is poor.
At present, a low orbit constellation satellite generally realizes autonomous orbit determination on a satellite through a satellite-borne GNSS navigation receiver, the root mean square of a measurement error of a semi-long axis of an orbit is generally more than 10m, and if the semi-long axis flat root adjustment quantity is directly calculated according to the measurement value, the condition that the amplitude angle of relative latitude exceeds the boundary of a control box is extremely easy to occur, so that the control strategy of a unilateral drift ring is invalid.
Disclosure of Invention
The invention solves the technical problems that: the method has the advantages that the defect of the prior art is overcome, the low orbit constellation phase maintaining parameter estimation method is provided, the orbit semi-long axis is estimated indirectly through the relative latitude amplitude angle measurement value, the estimation error of the orbit semi-long axis flat root can be reduced to be within 1m, and the requirement of the satellite for constellation phase maintaining based on a unilateral drift ring control strategy is met.
The technical scheme of the invention is as follows: a low-orbit constellation phase-preserving parameter estimation method comprises the following steps:
determining a constellation reference orbit and a constellation phase maintaining parameter;
performing relative latitude amplitude angle flat root polynomial fitting to obtain a relative reference orbit latitude amplitude angle Ping Gen;
calculating the coefficients of the opposite semi-major axis flat root polynomials;
the relative semi-major axis flat root deviation is estimated.
Further, the determining constellation reference track specifically includes: the low-orbit constellation nominal orbit semi-long axis flat root is taken as a reference orbit semi-long axis Ping Gen a R The reference orbit extrapolation only considers the gravitational perturbation, extrapolates any time t, and the orbit latitude amplitude root is u R (t) track semi-major axis Ping Gen mean
Figure BDA0003006243470000031
Is kept unchanged, i.e.)>
Figure BDA0003006243470000032
Further, determining constellation phase maintaining parameters, specifically: for the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as: Δu (t) =u (t) -u R (t); the semi-long axis flat root of the track is a (t), and the semi-long axis flat root of the track relative to the reference track is defined as delta a (t) =a (t) -a R 。
Further, the method for performing the relative latitude amplitude angle flat root polynomial fitting comprises the following steps: under the condition of no orbit control, obtaining instantaneous orbit parameters of a satellite by GNSS orbit determination, calculating latitude argument Ping Gen u (t) of the satellite by using conversion of the instantaneous orbit parameters and the average orbit parameters, and using Deltau (t) =u (t) -u R (t) calculating a relative reference orbit latitude argument Ping Gen u (t); approximating the relative reference orbit latitude argument Ping Gen u (t) to Deltau (t) ≡k with a cubic polynomial u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 ；t 0 Start time of trackless control, k u0 、k u1 、k u2 、k u3 Is a coefficient.
Further, the calculating the coefficients of the relative semi-major axis flat root polynomial comprises the following steps:
approximating the relative reference orbit semi-major axis Ping Gen a (t) to a quadratic polynomial
Figure BDA0003006243470000033
Wherein t is 0 Is the starting moment of trackless control;
according to the reference orbit semi-major axis Ping Gen a R And
Figure BDA0003006243470000034
calculating coefficient K a The method comprises the steps of carrying out a first treatment on the surface of the Wherein μ is a gravitational coefficient;
by means of
Figure BDA0003006243470000035
Calculating the coefficient Δa (t 0 )、
Figure BDA0003006243470000036
And->
Figure BDA0003006243470000037
Further, the estimation is biased with respect to the minor axis Ping Gen, specifically: from the determination of
Figure BDA0003006243470000038
Calculating the starting time t of the trackless control 0 To the end time t of the trackless control f An estimate of the relative reference orbit semi-major axis Ping Gen a (t) at any time in between.
A low-rail constellation phase-preserving parameter estimation system, comprising:
a first module for determining a constellation reference orbit and a constellation phase maintaining parameter;
the second module is used for carrying out relative latitude amplitude angle flat root polynomial fitting to obtain a relative reference orbit latitude amplitude angle Ping Gen;
the third module is used for calculating the polynomial coefficient of the flat root relative to the semi-major axis;
and a fourth module estimating a flat root deviation from the semi-major axis.
Further, the determining constellation reference track specifically includes: the low-orbit constellation nominal orbit semi-long axis flat root is taken as a reference orbit semi-long axis Ping Gen a R The reference orbit extrapolation only considers the gravitational perturbation, extrapolates any time t, and the orbit latitude amplitude root is u R (t) track semi-major axis Ping Gen mean
Figure BDA0003006243470000046
Is kept unchanged, i.e.)>
Figure BDA0003006243470000047
Further, determining constellation phase maintaining parameters, specifically: for the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as: Δu (t) =u (t) -u R (t); the semi-long axis flat root of the track is a (t), and the semi-long axis flat root of the track relative to the reference track is defined as delta a (t) =a (t) -a R ；
Further, the method for performing the relative latitude amplitude angle flat root polynomial fitting comprises the following steps: under the condition of no orbit control, obtaining instantaneous orbit parameters of a satellite by GNSS orbit determination, calculating latitude argument Ping Gen u (t) of the satellite by using conversion of the instantaneous orbit parameters and the average orbit parameters, and using Deltau (t) =u (t) -u R (t) calculating a relative reference orbit latitude argument Ping Gen u (t); approximating the relative reference orbit latitude argument Ping Gen u (t) to Deltau (t) ≡k with a cubic polynomial u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 ；t 0 Start time of trackless control, k u0 、k u1 、k u2 、k u3 Is a coefficient;
further, the calculating the coefficients of the relative semi-major axis flat root polynomial comprises the following steps:
approximating the relative reference orbit semi-major axis Ping Gen a (t) to a quadratic polynomial
Figure BDA0003006243470000041
Wherein t is 0 Is the starting moment of trackless control;
according to the reference orbit semi-major axis Ping Gen a R And
Figure BDA0003006243470000042
calculating coefficient K a The method comprises the steps of carrying out a first treatment on the surface of the Wherein μ is a gravitational coefficient;
by means of
Figure BDA0003006243470000043
Calculating the coefficient Δa (t 0 )、
Figure BDA0003006243470000044
And->
Figure BDA0003006243470000045
The estimation is biased with respect to the semi-major axis Ping Gen, specifically: from the determination of
Figure BDA0003006243470000051
Calculating the starting time t of the trackless control 0 To the end time t of the trackless control f An estimate of the relative reference orbit semi-major axis Ping Gen a (t) at any time in between.
A computer readable storage medium storing a computer program which when executed by a processor performs the steps of the low-rail constellation phase-preserving parameter estimation method.
A low-rail constellation phase-preserving parameter estimation device comprising a memory, a processor and a computer program stored in said memory and executable on said processor, said processor implementing the steps of said one low-rail constellation phase-preserving parameter estimation method when executing said computer program.
Compared with the prior art, the invention has the advantages that:
(1) According to the invention, the orbit semi-major axis Ping Gen is indirectly estimated by using the latitude amplitude angle measured value with higher precision without directly using the orbit semi-major axis measured value output by the satellite-borne GNSS orbit determination, so that the orbit determination precision (typical value is 10 m) of the existing GNSS receiver is not improved, the high precision requirement of satellite constellation maintenance on the orbit semi-major axis can be met, and the double-frequency GNSS receiver with higher configuration price on satellites is avoided;
(2) According to the method, polynomial fitting is carried out after the relative latitude amplitude angle measured value is converted into the flat root, then the change relation between the relative latitude amplitude angle flat root and the track semi-long axis flat root is utilized, polynomial fitting is carried out on the change trend of the track semi-long axis flat root, and the fitting coefficient is calculated, so that the track semi-long axis Ping Gen at any moment is estimated indirectly, the semi-long axis flat root precision is better than 1m, and the constellation phase retention can be completed based on a unilateral drift ring control strategy under the condition of meeting any semi-long axis attenuation rate;
(3) According to the invention, polynomial fitting is carried out on the relative latitude amplitude angle measured value and the semi-long axis flat root, so that the semi-long axis Ping Gen estimation can be carried out on orbit autonomously, the satellite full-autonomous constellation phase maintenance is realized, and the ground operation and control operation is reduced.
Drawings
FIG. 1 is a flow chart of the method of the present invention;
FIG. 2 is a graph of the relative reference orbit semi-major axis flat root variation using a single-sided drift ring control strategy in accordance with the present invention;
FIG. 3 is a graph of the relative reference orbit latitude amplitude root variation using a single-sided drift ring control strategy in accordance with the present invention;
FIG. 4 is a graph of the variation of the root-mean-square of the latitude amplitude of the trackless control segment relative to a reference orbit (GNSS output and theoretical value) according to the present invention;
FIG. 5 is a graph of the semi-long axis flat root variation curve (GNSS output and theoretical value) of the trackless control segment relative to a reference track according to the present invention;
FIG. 6 is a graph of the semilong axis flat root variation curve (theoretical and estimated) of a trackless control segment of the present invention relative to a reference track;
FIG. 7 is a plot of error of the semi-major axis flat root estimate of the trackless control segment relative to a reference track according to the present invention relative to a theoretical value.
Detailed Description
In order to better understand the technical solutions described above, the following detailed description of the technical solutions of the present application is provided through the accompanying drawings and specific embodiments, and it should be understood that the specific features of the embodiments and embodiments of the present application are detailed descriptions of the technical solutions of the present application, and not limit the technical solutions of the present application, and the technical features of the embodiments and embodiments of the present application may be combined with each other without conflict.
The following describes in further detail a low-rail constellation phase-preserving parameter estimation method provided in the embodiments of the present application with reference to the accompanying drawings of the specification, and the specific implementation manner may include (as shown in fig. 1 to 7):
in the scheme provided by the embodiment of the application, the relative semi-major axis Ping Gen is expressed as a quadratic polynomial according to the perturbation characteristic of the semi-major axis Ping Gen of the relative reference orbit, the relative latitude root is fitted by using a high-precision relative latitude root obtained by GNSS, the quadratic polynomial coefficient of the relative semi-major axis root is determined by using the equivalent relation between the relative latitude root of the relative reference orbit and the relative semi-major axis Ping Gen, and finally the relative semi-major axis root at different moments is estimated by using the quadratic polynomial. The method comprises the following specific steps:
1) Defining constellation reference trajectories
The nominal orbit semi-long axis flat root of the low orbit constellation is taken as a reference orbit semi-long axis Ping Gen and is defined as a R When the reference orbit extrapolation only considers the gravitational perturbation, the orbit latitude amplitude root is u when the arbitrary time t is extrapolated R (t) track semi-major axis Ping Gen mean
Figure BDA0003006243470000061
Is kept unchanged, i.e.)>
Figure BDA0003006243470000062
2) Defining constellation phase retention parameters
For the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as:
Δu(t)＝u(t)-u R (t) (1)
the semi-long axis flat root of the track is a (t), and the definition of the semi-long axis flat root of the track relative to the reference track is as follows:
Δa(t)＝a(t)-a R (2)
3) Relative latitude argument flat root polynomial fitting
Under the condition of no orbit control, the instantaneous orbit parameters of the satellite can be obtained through GNSS orbit determination, the latitude breadth angle Ping Gen u (t) of the satellite can be calculated by using a conversion formula of the instantaneous orbit parameters and the average orbit parameters, and the relative reference orbit latitude breadth angle Ping Gen u (t) can be calculated by using the formula (1).
The relative reference orbit latitude argument Ping Gen u (t) can be approximated by a cubic polynomial:
Δu(t)≈k u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 (3)
4) Polynomial coefficient calculation relative to semi-major axis flat root
For low-orbit constellations, the semi-long axis root-flattening relative to the reference orbit is mainly affected by atmospheric resistance, and in the case of trackless control, it can be approximated by a taylor's expanded quadratic polynomial:
Figure BDA0003006243470000071
wherein: Δa (t) 0 ) At t is Δa (t) 0 The value of the time of day,
Figure BDA0003006243470000072
at t is Δa (t) 0 First derivative of time of day>
Figure BDA0003006243470000073
At t is Δa (t) 0 The second derivative of time.
The rate of change of the relative reference orbit latitude argument Ping Gen u (t) can be expressed approximately as:
Figure BDA0003006243470000074
due to a (t) and a R The phase difference is small, and the expression is developed by using the Taylor series:
Figure BDA0003006243470000075
wherein the method comprises the steps of
Figure BDA0003006243470000081
Substituting the formula (4) into the formula (6) can obtain:
Figure BDA0003006243470000082
the derivation of (3) is as follows:
Figure BDA0003006243470000083
comparison of formulas (8) and (9) can be made:
Figure BDA0003006243470000084
from the above, the calculation steps of the relative semi-major axis flat root polynomial coefficients are as follows:
a) Let t be 0 From time to t f At the moment, the satellite is in a constant orbit state, the instantaneous orbit parameter of the satellite can be obtained through GNSS orbit determination, the orbit latitude breadth angle Ping Gen u (t) of the satellite can be calculated by using a conversion formula of the instantaneous orbit parameter and the average orbit parameter, and the relative reference orbit latitude breadth angle Ping Gen u (t) can be calculated by using (1).
b) By (3) and t 0 From time to t f The polynomial coefficient k can be obtained by performing polynomial fitting for three times on the value of Deltau (t) at the moment u0 、k u1 、k u2 And k u3 。
c) According to the reference orbit semi-major axis Ping Gen a R And (7) calculating the coefficient K a 。
d) Calculating the coefficient delta a (t) of the formula (4) by using the formula (10) 0 )、
Figure BDA0003006243470000085
And->
Figure BDA0003006243470000086
5) Deviation DA estimation relative to semi-major axis Ping Gen
Calculating t by using the formula (4) 0 From time to t f An approximation of the relative reference orbit semi-major axis Ping Gen a (t) at any time between times.
The specific implementation process of the invention is as follows:
1) Defining constellation reference trajectories
The nominal orbit semi-long axis flat root of the low orbit constellation is taken as a reference orbit semi-long axis Ping Gen and is defined as a R When the reference orbit extrapolation only considers the gravitational perturbation, the orbit latitude amplitude root is u when the arbitrary time t is extrapolated R (t) track semi-major axis Ping Gen mean a R Is maintained unchanged, i.e. a R ＝a R 。
2) Defining constellation phase retention parameters
For the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as:
Δu(t)＝u(t)-u R (t) (1)
the semi-long axis flat root of the track is a (t), and the definition of the semi-long axis flat root of the track relative to the reference track is as follows:
Δa(t)＝a(t)-a R (2)
3) Relative latitude argument flat root polynomial fitting
Under the condition of no orbit control, the instantaneous orbit parameters of the satellite can be obtained through GNSS orbit determination, the latitude breadth angle Ping Gen u (t) of the satellite can be calculated by using a conversion formula of the instantaneous orbit parameters and the average orbit parameters, and the relative reference orbit latitude breadth angle Ping Gen u (t) can be calculated by using the formula (1).
The relative reference orbit latitude argument Ping Gen u (t) can be approximated by a cubic polynomial:
Δu(t)≈k u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 (3)
4) Polynomial coefficient calculation relative to semi-major axis flat root
For low-orbit constellations, the semi-long axis root-flattening relative to the reference orbit is mainly affected by atmospheric resistance, and in the case of trackless control, it can be approximated by a taylor's expanded quadratic polynomial:
Figure BDA0003006243470000091
wherein: Δa (t) 0 ) At t is Δa (t) 0 The value of the time of day,
Figure BDA0003006243470000092
at t is Δa (t) 0 First derivative of time of day>
Figure BDA0003006243470000093
At t is Δa (t) 0 The second derivative of time.
The rate of change of the relative reference orbit latitude argument Ping Gen u (t) can be expressed approximately as:
Figure BDA0003006243470000094
due to a (t) and a R The phase difference is small, and the expression is developed by using the Taylor series:
Figure BDA0003006243470000101
wherein the method comprises the steps of
Figure BDA0003006243470000102
Substituting the formula (4) into the formula (6) can obtain:
Figure BDA0003006243470000103
the derivation of (3) is as follows:
Figure BDA0003006243470000104
comparison of formulas (8) and (9) can be made:
Figure BDA0003006243470000105
from the above, the calculation steps of the relative semi-major axis flat root polynomial coefficients are as follows:
a) By (3) and t 0 From time to t f The polynomial coefficient k can be obtained by performing polynomial fitting for three times on the value of Deltau (t) at the moment u0 、k u1 、k u2 And k u3 。
b) According to the reference orbit semi-major axis Ping Gen a R And (7) calculating the coefficient K a 。
c) Calculating the coefficient delta a (t) of the formula (4) by using the formula (10) 0 )、
Figure BDA0003006243470000106
And->
Figure BDA0003006243470000107
5) Deviation DA estimation relative to semi-major axis Ping Gen
Calculating t by using the formula (4) 0 From time to t f An approximation of the relative reference orbit semi-major axis Ping Gen a (t) at any time between times.
Examples:
assume that the nominal orbit semi-long axis flat root of the low orbit constellation is 7578.137km, the eccentricity is Ping Gen 0.001562, the inclination angle is flat root 60 degrees, the near-place amplitude angle is 90 degrees, and the true near-place angle is 270 degrees. The phase maintaining precision of the same orbit satellite requires +/-0.1 degrees, the boundary of a horizontal root control box of the latitude amplitude angle relative to a reference orbit is +/-0.1 degrees, the root mean square of the measurement error of the GNSS orbit determination semi-long axis is 10m, and the root mean square of the measurement error of the corresponding latitude amplitude angle is 0.0001 degrees.
And a single-side drift ring control strategy is adopted, a semi-long axis flat root change curve of the relative reference track is shown in figure 2, and a flat root change curve of the latitude amplitude angle of the relative reference track is shown in figure 3 within 280 days. When the relative reference orbit horizontal root is close to the right boundary of the control box, the relative reference orbit semi-long axis horizontal root is required to be controlled to be close to a target value through orbit maneuver, the target value is in the range of 10 m-14 m, and because the orbit semi-long axis error root mean square output by the satellite-borne GNSS orbit determination is 10m, the error root mean square of the relative reference orbit semi-long axis horizontal root is also 10m, the actual relative reference orbit semi-long axis horizontal root control result is between 0 and 24m, and when the relative semi-long axis horizontal root control quantity exceeds 14m, the relative latitude horizontal root floats out of the left boundary of the control box, and finally the unilateral drift ring control strategy is invalid.
By using the method herein, taking the section 1 uncontrolled orbit as an example, the time from day 0 to day 50, the measured value and the theoretical value change curve of the latitude angle Ping Gen u (t) of the relative reference orbit obtained according to the GNSS orbit determination are shown in fig. 4, the measured value and the theoretical value change curve of the semi-long axis Ping Gen a (t) of the relative reference orbit obtained according to the GNSS orbit determination are shown in fig. 5, it can be seen that the difference between the GNSS output of Δu (t) and the theoretical value is very small, and the difference between the GNSS output of Δa (t) and the theoretical value is very large.
The coefficient K can be calculated according to equation (7) a = 0.0009455 °/day/m.
The third polynomial fitting is performed by using the formula (4) to obtain k u0 ＝0.0622°，k u1 = -0.0094 °/day, k u2 ＝1.3793×10 -4 Degree/day 2 And k u3 ＝3.1113×10 -7 Degree/day 3 The method comprises the following steps:
Δu(t)≈0.0622 0 -0.0094(t-t 0 )+1.3793×10 -4 (t-t 0 ) 2 +3.1113×10 -7 (t-t 0 ) 3 (11)
from the expression (10), the coefficient Δa (t) of expression (3) can be calculated 0 )＝9.9002m，
Figure BDA0003006243470000111
Tianhe->
Figure BDA0003006243470000112
Namely:
Δa(t)≈9.9002m-0.2918(t-t 0 )-0.0010(t-t 0 ) 2 (12)
from day 0 to day 50, the change curves of the approximate value and the theoretical value with respect to the half major axis Ping Gen a (t) of the reference track are calculated according to the formula (12) and are shown in fig. 6, and the error change curve of the approximate value and the theoretical value with respect to the half major axis Ping Gen a (t) of the reference track is shown in fig. 7. As can be seen from the figure, with the method herein, the approximation error with respect to the reference orbit semi-major axis Ping Gen a (t) is better than ±0.7m, satisfying the on-orbit implementation of the single-sided drift ring control strategy.
It will be apparent to those skilled in the art that various modifications and variations can be made in the present application without departing from the spirit or scope of the application. Thus, if such modifications and variations of the present application fall within the scope of the claims and the equivalents thereof, the present application is intended to cover such modifications and variations.
What is not described in detail in the present specification is a well known technology to those skilled in the art.

Claims (7)

1. The low-orbit constellation phase maintaining parameter estimation method is characterized by comprising the following steps of:
determining a constellation reference orbit and a constellation phase maintaining parameter;
performing relative latitude amplitude angle flat root polynomial fitting to obtain a relative reference orbit latitude amplitude angle Ping Gen;
calculating the coefficients of the opposite semi-major axis flat root polynomials;
estimating a relative semi-major axis Ping Gen deviation;
the determining constellation reference track specifically comprises the following steps: the low-orbit constellation nominal orbit semi-long axis flat root is taken as a reference orbit semi-long axis Ping Gen a R The reference orbit extrapolation only considers the gravitational perturbation, extrapolates any time t, and the orbit latitude amplitude root is u R (t) track semi-major axis Ping Gen mean
Figure FDA0004088852800000012
Is kept unchanged, i.e.)>
Figure FDA0004088852800000013
The constellation phase maintaining parameter is determined, specifically: for the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as: Δu (t) =u (t) -u R (t); the semi-long axis flat root of the track is a (t), and the semi-long axis flat root of the track relative to the reference track is defined as delta a (t) =a (t) -a R ；
The method for carrying out the relative latitude amplitude angle flat root polynomial fitting comprises the following steps: under the condition of no orbit control, obtaining instantaneous orbit parameters of a satellite by GNSS orbit determination, calculating latitude argument Ping Gen u (t) of the satellite by using conversion of the instantaneous orbit parameters and the average orbit parameters, and using Deltau (t) =u (t) -u R (t) calculating a relative reference orbit latitude argument Ping Gen u (t); approximating the relative reference orbit latitude argument Ping Gen u (t) to Deltau (t) ≡k with a cubic polynomial u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 ；t 0 Start time of trackless control, k u0 、k u1 、k u2 、k u3 Is a coefficient.
2. A method of estimating phase-preserving parameters of a low-rail constellation according to claim 1, wherein said calculating coefficients of a parallel root polynomial with respect to a major axis comprises the steps of:
approximating the relative reference orbit semi-major axis Ping Gen a (t) to a quadratic polynomial
Figure FDA0004088852800000011
Wherein t is 0 Is the starting moment of trackless control;
according to the reference orbit semi-major axis Ping Gen a R And
Figure FDA0004088852800000021
calculating coefficient K a The method comprises the steps of carrying out a first treatment on the surface of the Wherein,,μ is the gravitational coefficient;
by means of
Figure FDA0004088852800000022
Calculating the coefficient Δa (t 0 )、
Figure FDA0004088852800000023
And->
Figure FDA0004088852800000024
3. The method according to claim 1, wherein the estimation is biased with respect to a semi-major axis Ping Gen, specifically: from the determination of
Figure FDA0004088852800000025
Calculating the starting time t of the trackless control 0 To the end time t of the trackless control f An estimate of the relative reference orbit semi-major axis Ping Gen a (t) at any time in between.
4. A low-rail constellation phase-preserving parameter estimation system, comprising:
a first module for determining a constellation reference orbit and a constellation phase maintaining parameter;
the second module is used for carrying out relative latitude amplitude angle flat root polynomial fitting to obtain a relative reference orbit latitude amplitude angle Ping Gen;
the third module is used for calculating the polynomial coefficient of the flat root relative to the semi-major axis;
a fourth module that estimates a deviation from the median axis Ping Gen;
the determining constellation reference track specifically comprises the following steps: the low-orbit constellation nominal orbit semi-long axis flat root is taken as a reference orbit semi-long axis Ping Gen a R The reference orbit extrapolation only considers the gravitational perturbation, extrapolates any time t, and the orbit latitude amplitude root is u R (t) track semi-major axis Ping Gen mean
Figure FDA0004088852800000027
Is kept unchanged, i.e.)>
Figure FDA0004088852800000026
The constellation phase maintaining parameter is determined, specifically: for the actual running orbit of satellite affected by earth attraction perturbation, solar-lunar three-body attraction perturbation, atmospheric resistance perturbation and light pressure perturbation, the orbit latitude amplitude angle flat root is u (t), and the relative reference orbit latitude amplitude angle flat root is defined as: Δu (t) =u (t) -u R (t); the semi-long axis flat root of the track is a (t), and the semi-long axis flat root of the track relative to the reference track is defined as delta a (t) =a (t) -a R ；
The method for carrying out the relative latitude amplitude angle flat root polynomial fitting comprises the following steps: under the condition of no orbit control, obtaining instantaneous orbit parameters of a satellite by GNSS orbit determination, calculating latitude argument Ping Gen u (t) of the satellite by using conversion of the instantaneous orbit parameters and the average orbit parameters, and using Deltau (t) =u (t) -u R (t) calculating a relative reference orbit latitude argument Ping Gen u (t); approximating the relative reference orbit latitude argument Ping Gen u (t) to Deltau (t) ≡k with a cubic polynomial u0 +k u1 (t-t 0 )+k u2 (t-t 0 ) 2 +k u3 (t-t 0 ) 3 ；t 0 Start time of trackless control, k u0 、k u1 、k u2 、k u3 Is a coefficient.
5. The low-rail constellation phase preserving parameter estimation system of claim 4, wherein said calculating relative semi-major axis flat root polynomial coefficients comprises the steps of:
approximating the relative reference orbit semi-major axis Ping Gen a (t) to a quadratic polynomial
Figure FDA0004088852800000031
Wherein t is 0 Is the starting moment of trackless control;
according to the reference orbit semi-major axis Ping Gen a R And
Figure FDA0004088852800000032
calculating coefficient K a The method comprises the steps of carrying out a first treatment on the surface of the Wherein μ is a gravitational coefficient;
by means of
Figure FDA0004088852800000033
Calculating the coefficient Δa (t 0 )、
Figure FDA0004088852800000034
And->
Figure FDA0004088852800000035
The estimation is biased with respect to the semi-major axis Ping Gen, specifically: from the determination of
Figure FDA0004088852800000036
Calculating the starting time t of the trackless control 0 To the end time t of the trackless control f An estimate of the relative reference orbit semi-major axis Ping Gen a (t) at any time in between.
6. A computer readable storage medium storing a computer program, which when executed by a processor performs the steps of the method according to any one of claims 1 to 3.
7. A low-rail constellation phase-preserving parameter estimation device comprising a memory, a processor and a computer program stored in said memory and executable on said processor, characterized by: the processor, when executing the computer program, performs the steps of the method of any one of claims 1 to 3.
CN202110362735.9A 2021-04-02 2021-04-02 Low-orbit constellation phase maintaining parameter estimation method Active CN113189619B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110362735.9A CN113189619B (en) 2021-04-02 2021-04-02 Low-orbit constellation phase maintaining parameter estimation method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110362735.9A CN113189619B (en) 2021-04-02 2021-04-02 Low-orbit constellation phase maintaining parameter estimation method

Publications (2)

Publication Number Publication Date
CN113189619A CN113189619A (en) 2021-07-30
CN113189619B true CN113189619B (en) 2023-05-09

Family

ID=76975004

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110362735.9A Active CN113189619B (en) 2021-04-02 2021-04-02 Low-orbit constellation phase maintaining parameter estimation method

Country Status (1)

Country Link
CN (1) CN113189619B (en)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114006646B (en) * 2021-09-27 2023-09-29 中国人民解放军战略支援部队航天工程大学 Track control frequency analysis method and device for maintaining Walker constellation configuration
CN114254262B (en) * 2021-11-22 2023-04-18 浙江大学 Method and device for maintaining autonomous configuration of heterogeneous quality ratio satellite constellation and electronic equipment
CN114852375B (en) * 2022-03-24 2024-08-30 北京控制工程研究所 Formation satellite relative orbit change estimation method and estimation device
CN115396002B (en) * 2022-05-11 2024-09-17 航天行云科技有限公司 Control method, device and processing system for constellation phase of low-orbit satellite
CN115806061B (en) * 2022-11-10 2023-05-09 北京航天驭星科技有限公司 Modeling method, model and acquisition method of satellite relative phase maintaining strategy model
CN115636111B (en) * 2022-12-21 2023-03-28 北京航天驭星科技有限公司 Phase difference maintaining method, system, device, and medium
CN117421532B (en) * 2023-12-18 2024-02-27 北京航天驭星科技有限公司 Method, system, equipment and medium for calculating attenuation rate of medium-dip angle low orbit satellite

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108055069A (en) * 2017-12-11 2018-05-18 中国人民解放军战略支援部队航天工程大学 Low rail communication and navigation enhancing hybrid constellation maintain control feature modeling and control method
CN111541477A (en) * 2019-11-25 2020-08-14 航天科工空间工程发展有限公司 Method and device for suppressing internal frequency interference of low-orbit constellation system
CN111591469A (en) * 2020-03-03 2020-08-28 航天科工空间工程发展有限公司 Low-orbit constellation system phase keeping method, system, equipment and storage medium

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20130002484A1 (en) * 2011-07-03 2013-01-03 Daniel A. Katz Indoor navigation with gnss receivers

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108055069A (en) * 2017-12-11 2018-05-18 中国人民解放军战略支援部队航天工程大学 Low rail communication and navigation enhancing hybrid constellation maintain control feature modeling and control method
CN111541477A (en) * 2019-11-25 2020-08-14 航天科工空间工程发展有限公司 Method and device for suppressing internal frequency interference of low-orbit constellation system
CN111591469A (en) * 2020-03-03 2020-08-28 航天科工空间工程发展有限公司 Low-orbit constellation system phase keeping method, system, equipment and storage medium

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
基于参考轨道的Walker星座相对相位保持策略;胡松杰等;《空间控制技术与应用》;20101031;第36卷(第05期);第45-49页 *
基于约化相对轨道拟平根数的长期稳定高精度卫星编队导航技术;杨盛庆等;《空间控制技术与应用》;20170228;第43卷(第01期);第30-35页 *

Also Published As

Publication number Publication date
CN113189619A (en) 2021-07-30

Similar Documents

Publication Publication Date Title
CN113189619B (en) Low-orbit constellation phase maintaining parameter estimation method
US5528502A (en) Satellite orbit maintenance system
US10532829B2 (en) Orbit transfer method for a spacecraft using a continuous or quasi-continuous thrust and embedded driving system for implementing such a method
US8099186B2 (en) Satellite navigation using long-term navigation information and autonomous orbit control
CN110429974B (en) Fast alignment method and device based on regression orbit constellation
US8620496B2 (en) Systems and method of controlling a spacecraft using attitude sensors
US20080105788A1 (en) Methods and apparatus for node-synchronous eccentricity control
CN109542112B (en) Fixed time convergence anti-interference control method for return flight of vertical take-off and landing reusable rocket
US6439507B1 (en) Closed-loop spacecraft orbit control
CN111045457B (en) Optical axis pointing adjustment method based on satellite-borne remote sensing instrument
CN112486196B (en) Aircraft rapid trajectory optimization method meeting strict time and position constraints
Blandino et al. Feasibility for orbital life extension of a CubeSat in the lower thermosphere
CN109760854B (en) Volume-controllable inflatable antenna and unfolding volume control method thereof
JP2000001200A (en) Improved method and system for determining posture of space craft
US5411227A (en) Satellite thruster uncertainty estimation in transition mode
CN117674965A (en) Constellation configuration maintaining method, device, terminal equipment and storage medium
EP0932549B1 (en) Autonomous solar torque management
Dilssner A note on the yaw attitude modeling of BeiDou IGSO-6
RU2721813C1 (en) Autonomous collocation method in geostationary orbit
US6283415B1 (en) Simplified yaw steering method for satellite antenna beam control
US7185858B2 (en) Spacecraft gyro calibration system
Leomanni et al. An adaptive groundtrack maintenance scheme for spacecraft with electric propulsion
US6457679B1 (en) Method for maintaining the position of geostationary satellites
Chavez et al. Relative-orbit element estimation for satellite navigation and guidance
Chao Semiautonomous stationkeeping of geosynchronous satellites

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant