Leo constellation gnss / inter-satellite link joint autonomous orbit determination method and system based on hierarchical estimation

Abstract

This invention provides a hierarchical estimation-based joint autonomous orbit determination method and system for LEO constellations using GNSS / inter-satellite links. First, GNSS data and inter-satellite link data collected by each satellite are acquired. Then, these data are input into local filters to obtain the state estimates and covariance matrices of each satellite. Next, the state estimates and covariance matrices are input into a main filter to obtain the overall reference error estimate of the constellation. Finally, the overall reference error estimate is fed back to each local filter to correct the state estimates of each satellite, resulting in the final high-precision orbit parameters and clock error parameters for each satellite. This invention decomposes the global orbit determination problem into local estimation and low-dimensional reference error estimation through a hierarchical estimation architecture. Combined with an observation model based on an inter-satellite link closed loop and an adaptive weight adjustment mechanism, it effectively suppresses overall constellation reference drift, significantly reduces computational complexity, and enhances the system's robustness in complex space environments.

Description

Technical Field

[0001] This invention relates to the field of satellite navigation and autonomous orbit determination technology, and in particular to a method and system for joint autonomous orbit determination of LEO constellation GNSS / inter-satellite links based on hierarchical estimation. Background Technology

[0002] With the rise of large low-Earth orbit constellations, their orbit determination missions face enormous challenges. Traditional orbit determination methods relying on ground stations are insufficient to meet the needs of real-time management of tens of thousands of satellites. Using onboard GNSS receivers for orbit determination is the mainstream approach, but when GNSS signals are interfered with or interrupted, orbit determination accuracy drops sharply or even fails.

[0003] Existing technologies attempt to use inter-satellite links for collaborative orbit determination, but relying solely on inter-satellite links (relative measurements) cannot perceive the absolute position and attitude of the constellation in space, causing the constellation to rotate and drift, making it impossible to maintain usable accuracy in the long term. Furthermore, centrally processing GNSS and inter-satellite link observation data from thousands of satellites results in a massive parameter scale, exhibiting the "curse of dimensionality," making real-time or near-real-time calculations difficult. Additionally, the fixed orbit determination mode cannot adaptively respond to dynamic space environments such as changes in GNSS signal quality. Summary of the Invention

[0004] This invention provides a hierarchical estimation-based LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method and system to solve the problem that existing technologies cannot effectively suppress overall reference drift.

[0005] A first aspect of this invention provides a joint autonomous orbit determination method for LEO constellations using GNSS / inter-satellite links based on hierarchical estimation, comprising: Acquire GNSS data and inter-satellite link data collected by each satellite; GNSS data and inter-satellite link data are input into a local filter to obtain the state estimates and covariance matrices of each satellite. The state estimates and covariance matrices of each satellite are input into the main filter to obtain the overall baseline error estimate of the constellation. The overall baseline error estimate is fed back to each local filter to correct the state estimate of each satellite, thus obtaining the final high-precision orbital parameters and clock error parameters of each satellite.

[0006] In one possible implementation, GNSS data and inter-satellite link data are input into a local filter to obtain the state estimates and covariance matrices of each satellite, including: Preprocessing of GNSS observation data and inter-satellite link observation data; The preprocessed inter-satellite link observation data is input into a local filter to obtain the state prediction value and its covariance prediction matrix. The preprocessed GNSS observation data is input into a local filter to correct the state prediction values ​​and their covariance prediction matrix, thereby obtaining the state estimates and covariance matrices of each satellite.

[0007] In one possible implementation, the state estimates and covariance matrices of each satellite are input into the main filter to obtain the overall baseline error estimate of the constellation, including: Based on the state estimates of each satellite, the closure error of the geometric closed loop formed by the inter-satellite links is calculated to construct a baseline error observation. Establish the mathematical relationship between the closure error of the geometric closed loop and the overall rotation / translation parameters of the constellation; The state estimates and covariance matrices of each satellite are used as inputs to the main filter. Based on mathematical relationships and baseline error observations, the overall baseline error of the constellation is estimated. The overall baseline error includes three orbital translation parameters, three rotation parameters, and one clock error baseline parameter.

[0008] In one possible implementation, the overall baseline error estimate is fed back to each local filter to correct the state estimates of each satellite, resulting in the final high-precision orbital parameters and clock error parameters for each satellite, including: The main filter broadcasts the overall baseline error estimate to the local filters of each satellite; Based on the overall baseline error estimate, the state estimates of the local filters of each satellite are corrected. The corrected state estimates are output as the final high-precision orbital parameters and clock error parameters for each satellite.

[0009] In one possible implementation, the method further includes, before acquiring inter-satellite link data: Based on the orbit prediction model, the reachability of links with potential neighboring satellites within a future time window is determined, and an effective neighbor window period is defined. When the effective neighbor window period meets the preset update triggering conditions, a new round of neighbor filtering process is adaptively started; Perform neighbor optimization screening under multi-objective constraints, construct a weighted function with topology ring integrity, communication energy consumption and link quality as optimization objectives, and dynamically adjust the weights of each objective according to the real-time operating status of the constellation; Based on the calculation results of the weighted function, the optimal subset of neighbors is selected from the potential neighbors within the effective neighbor window period; Each satellite establishes a consensus on links by engaging in bidirectional information exchange with candidate neighbor satellites based on the optimal subset of neighbors. At the same time, the overall topology connectivity of the constellation is periodically monitored through a global coordination mechanism. When a node that does not meet the preset connectivity requirements is detected in a local topology, a topology adjustment command is proactively issued to repair it.

[0010] In one possible implementation, when a node is detected in the local topology that does not meet the preset connectivity requirements, a topology adjustment command is proactively issued to repair it, including: When a satellite with connectivity below a preset threshold is detected, a topology adjustment command is sent to force it to add neighbors.

[0011] In one possible implementation, the closure error of the geometric closed loop formed by the inter-satellite links is calculated based on the state estimates of each satellite to construct a baseline error observation, including: A joint observation model coupling the baseline error and the time-varying clock bias is established. In this model, the satellite clock bias is extended from a static parameter to a state vector that includes its time-varying characteristics. The time-varying characteristics are then incorporated into the closed loop error observation equation to quantify the impact of clock bias drift on the baseline error estimation. A dynamic weight is assigned to each geometric closed loop; the dynamic weight is adaptively determined based on the geometric configuration strength and ranging accuracy of the closed loop. Each geometric closed loop difference observation residual is examined, and large residual observations are weighted down. Based on the joint observation model, the dynamic weights of each geometric closed loop, and the state estimates of each satellite, the closure error of the geometric closed loop formed by the inter-satellite links is calculated.

[0012] In one possible implementation, the state estimate includes the satellite position vector, velocity vector, and time-varying clock error parameters; based on the joint observation model, the dynamic weights of each geometric loop, and the state estimates of each satellite, the closure error of the geometric loop formed by the inter-satellite links is calculated, including: Based on the time-varying clock error dynamic equation in the joint observation model, the time-varying clock error parameters of each satellite are substituted into the equation to calculate the static clock error contribution and the time-varying clock error drift integral contribution between each satellite in each closed loop, thus obtaining the total impact of clock error on link ranging. Substitute the dynamic weight of each closed loop into the observation noise variance adjustment formula to perform weighted correction on the original ranging data of each link within the closed loop; The net link ranging value is determined based on the weighted adjusted link ranging data and the total clock bias effect. Based on the link net ranging value, determine the closure error of the geometric closed loop formed by the inter-satellite links.

[0013] A second aspect of this invention provides a LEO constellation GNSS / inter-satellite link joint autonomous orbit determination system based on hierarchical estimation, characterized in that it includes: The data acquisition module is used to acquire GNSS data and inter-satellite link data collected by each satellite; The first calculation module is used to input GNSS data and inter-satellite link data into a local filter to obtain the state estimate and covariance matrix of each satellite. The second calculation module is used to input the state estimates and covariance matrices of each satellite into the main filter to obtain the overall baseline error estimate of the constellation. The feedback correction module is used to feed back the overall reference error estimate to each local filter to correct the state estimate of each satellite, and finally obtain the high-precision orbit parameters and clock error parameters of each satellite.

[0014] A third aspect of the present invention provides an autonomous orbit determination system, characterized in that it includes an onboard computer, a GNSS receiver antenna, an inter-satellite link transceiver antenna, a processing unit, an ion thruster, an onboard clock, and a memory; the processing unit is used to execute the LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation as described in the first aspect above.

[0015] Compared to traditional technologies, this invention provides a hierarchical estimation-based joint autonomous orbit determination method and system for LEO constellations using GNSS / inter-satellite links. First, GNSS and inter-satellite link data collected by each satellite are acquired. Then, the GNSS and inter-satellite link data are input into local filters to obtain the state estimates and covariance matrices of each satellite. Next, the state estimates and covariance matrices of each satellite are input into a main filter to obtain the overall reference error estimate of the constellation. Finally, the overall reference error estimate is fed back to each local filter to correct the state estimates of each satellite, resulting in high-precision orbit parameters and clock error parameters for each satellite. This invention decomposes the global orbit determination problem into local estimation and low-dimensional reference error estimation through a hierarchical estimation architecture. Combined with an observation model based on an inter-satellite link closed loop and an adaptive weight adjustment mechanism, it effectively suppresses overall constellation reference drift, significantly reduces computational complexity, and enhances the system's robustness in complex space environments. Attached Figure Description

[0016] Figure 1 This is a flowchart illustrating the implementation of the LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation provided in this embodiment of the invention. Detailed Implementation

[0017] The embodiments of the present invention will now be described in detail with reference to the accompanying drawings.

[0018] Figure 1 This is a flowchart illustrating the implementation of the LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation, provided in this embodiment of the invention. Figure 1 As shown, the method includes: S110 acquires GNSS data and inter-satellite link data collected by each satellite; S120 inputs GNSS data and inter-satellite link data into a local filter to obtain the state estimates and covariance matrices of each satellite; S130: Input the state estimates and covariance matrices of each satellite into the main filter to obtain the overall baseline error estimate of the constellation; S140 feeds back the overall baseline error estimate to each local filter to correct the state estimate of each satellite, thus obtaining the final high-precision orbital parameters and clock error parameters of each satellite.

[0019] In this embodiment of the invention, each LEO satellite synchronously acquires two types of core observation data through its onboard GNSS receiver and inter-satellite link terminal. For GNSS data acquisition, the satellite receives pseudorange and pseudorange rate observation signals from multiple GNSS navigation satellites in real time, records quality indicators such as the number of visible stars, signal-to-noise ratio (SNR), and geometrical precision factor (GDOP), and simultaneously performs ionospheric and tropospheric delay compensation preprocessing on the raw data to remove obvious outliers. Before acquiring inter-satellite link data, neighbor satellite screening and link establishment must be completed: Based on the simplified orbital model of SGP4 / SDP4, the positional changes of potential neighbors within a 10km radius are predicted within the next 5-10 minutes, and an "effective neighbor window" is determined with an elevation angle ≥5°, relative speed ≤1km / s, and link distance ≤5000km; when the remaining time of the window is ≤2 minutes, a new round of screening is initiated, selecting 3-5 optimal neighbors through a weighted function of "ring integrity - communication energy consumption - link quality" (dynamically adjusting the weight coefficients); each satellite interacts bidirectionally with the candidate neighbors to confirm link consensus, and the main filtering layer monitors the global topology connectivity every 30 seconds, issuing adjustment commands to satellites with insufficient connectivity to forcibly supplement neighbors, ensuring link stability. After the link is established, the satellite collects relative distance data with neighboring satellites through one-way pseudocode ranging technology, records quality parameters such as link delay and packet loss rate, and forms a complete inter-satellite link observation dataset.

[0020] Each satellite's local filtering unit (using UKF or EKF algorithms) processes the preprocessed GNSS and inter-satellite link data in stages. First, data preprocessing is deepened: initial clock errors in the GNSS data are further corrected, and the inter-satellite link data is initially weighted according to link quality indicators. Then, the local filter's time update process is initiated. Based on an orbital dynamics model incorporating Earth's gravity, lunar and solar gravitational perturbations, and solar radiation pressure, combined with time-varying clock error dynamic equations (including clock error, clock error drift rate, and clock error drift acceleration), and inputting the state estimation results from the previous cycle, the system calculates the current state prediction value (including satellite position, velocity, and time-varying clock error parameters) and its covariance prediction matrix. Next, measurement updates are performed: the system's intelligent decision-making module determines the orbit determination environment based on GNSS signal quality (visible satellite count, GDOP value, SNR) and dynamically adjusts the weighting factors of GNSS and inter-satellite link observations. When the GNSS signal is good, GNSS data is used as the main source (weight 0.7-0.8). When the signal weakens, the weight of the inter-satellite link is gradually increased (0.3-0.4). When the signal is completely interrupted, only the inter-satellite link data is used. The weighted observation data is substituted into the filter, and the state prediction value and its covariance prediction matrix are corrected by Kalman gain calculation. Finally, the state estimate value (position vector, velocity vector, time-varying clock error parameter) and its covariance matrix of each satellite are output, and the accurate estimation at the local level is completed.

[0021] The main filter receives the state estimates and their covariance matrices from the local filtering outputs of all satellites, focusing on the estimation and modeling of the constellation's common reference error. First, a reference error observation is constructed: based on the geometric closed loops such as triangles and quadrilaterals formed by inter-satellite links, the position vectors and time-varying clock error parameters are extracted from the state estimates of each satellite. The algebraic sum of all weighted inter-satellite ranging values ​​on the closed loop (i.e., the closure error) is calculated. The closure error is zero when there is no error, but significantly non-zero when there is a reference error or clock drift. Subsequently, a joint observation model was established: the closure error was deeply coupled with the overall constellation baseline error (three orbital translation parameters and three rotation parameters) and the time-varying clock error contribution term. The impact of time-varying clock error on the baseline estimation was quantified by integrally correlating the clock error drift rate with the closure error. Simultaneously, dynamic weights were assigned to each closed loop. These weights were obtained by weighting geometric strength (calculated by side length and interior angles for triangles and by side length and diagonals for quadrilaterals) and ranging accuracy (based on the reciprocal of the link ranging error variance) (geometric strength weight 0.6, ranging accuracy weight 0.4). Abnormal closed loop data were eliminated using the residual test method, and observations with large residuals were downweighted. Finally, the main filter, using the state estimates of each satellite and their covariance matrix as input, combined with the aforementioned joint observation model and dynamic weights, iteratively estimated the overall baseline error of the constellation using the EKF algorithm. The iteration continued until the state estimate change was less than 1e-3m and converged. The final output was a complete baseline error result containing three orbital translation parameters, three rotation parameters, and one clock error baseline parameter.

[0022] The main filter broadcasts the overall reference error estimate to the local filter units of all satellites in the constellation via inter-satellite links. Upon receiving the reference error correction command, each local filter performs state correction according to a pre-defined correction model: substituting orbital translation and rotation parameters into the position and velocity correction formulas to eliminate systematic errors caused by the overall constellation translation and rotation; and fusing clock bias reference parameters with local time-varying clock bias estimates to correct clock bias reference deviations and form a unified clock bias reference system. During the correction process, the local filters synchronously update the covariance matrix of the state estimate and, combined with the error accuracy information fed back by the main filter, optimize the weight allocation strategy for subsequent filtering. After correction, the local filters output the final high-precision state estimation result, where the orbital parameters include the satellite's position vector (x, y, z) and velocity vector (v) in a geocentric rectangular inertial coordinate system. x v The clock error parameters include time-varying clock error values ​​and clock error drift rate after reference correction. The accuracy of the relevant parameters meets the requirements for autonomous operation of the low-Earth orbit mega-constellation. At the same time, the system continuously monitors the GNSS signal quality and dynamically switches the orbit determination mode to ensure the continuity and stability of the subsequent orbit determination process.

[0023] In some embodiments, GNSS data and inter-satellite link data are input into a local filter to obtain the state estimate and covariance matrix of each satellite, including: preprocessing the GNSS observation data and inter-satellite link observation data; inputting the preprocessed inter-satellite link observation data into a local filter to obtain the state prediction value and its covariance prediction matrix; and inputting the preprocessed GNSS observation data into a local filter to correct the state prediction value and its covariance prediction matrix to obtain the state estimate and its covariance matrix of each satellite.

[0024] In this embodiment of the invention, targeted processing is carried out for the two types of data based on their different characteristics to eliminate noise interference and system errors, providing high-quality input for subsequent filtering.

[0025] The focus is on correcting three types of errors that affect orbit determination accuracy. Ionospheric delay is eliminated by differentially analyzing dual-frequency observations (such as L1 / L2 bands) to remove first-order effects; tropospheric delay is corrected by combining the Saastamoinen model with real-time temperature and humidity data collected by onboard sensors; and the initial clock bias between the onboard receiver and the GNSS satellite is initially corrected by referring to the satellite clock bias parameters in the GNSS navigation message.

[0026] A dual screening method of "3σ criterion + trend test" is adopted. First, the mean and standard deviation of all GNSS pseudorange observations are calculated, and extreme data that deviate from the mean by more than 3 times the standard deviation are removed. Then, the trend of change is analyzed on the remaining data. If the pseudorange change rate of a certain GNSS satellite exceeds 0.5 m / s for 5 consecutive sampling periods (1 second per period), it is determined that the signal is interfered with, and all observation data of that satellite are discarded.

[0027] The key indicators of currently available GNSS satellites are statistically analyzed to classify signal status. When the number of visible satellites is ≥6, the geometrical precision factor (GDOP) is ≤3, and the signal-to-noise ratio (SNR) is ≥15dB, it is marked as "Good GNSS signal"; when the number of visible satellites is 4-5, the GDOP is 3-5, and the SNR is 10-15dB, it is marked as "GNSS signal attenuation"; when the number of visible satellites is <4, the GDOP is >5, or the SNR is <10dB, it is marked as "GNSS signal interruption." This marking will be used for subsequent filter weight adjustments.

[0028] Thresholds are set based on parameters (latency, packet loss rate, SNR) recorded by the inter-satellite link terminals to eliminate low-quality data. When the latency is >100 milliseconds, the packet loss rate is >0.1%, or the SNR is <12dB, the link is considered unstable, and the relative ranging data of that link is discarded; only high-quality link data that meets the criteria of "latency ≤100 milliseconds, packet loss rate ≤0.1%, SNR ≥12dB" are retained.

[0029] Further eliminate system biases from high-quality link data. First, correct propagation delay error: calculate the theoretical signal propagation time based on the relative positions between satellites, compare it with the measured delay, and compensate for the hardware delay of the transmitting and receiving antennas; second, correct Doppler frequency shift error: calculate the frequency shift based on the relative velocity of the satellites, compensate for the pseudo-code ranging value, and obtain a more accurate inter-satellite relative distance.

[0030] The corrected ranging data is matched one-to-one with neighboring satellites within the "effective neighbor window" (determined in advance through orbit prediction, meeting the requirements of elevation angle ≥ 5°, relative speed ≤ 1 km / s, and link distance ≤ 5000 km). Each data point is labeled with its neighboring satellite ID, link establishment time, and remaining window duration. If a neighbor's remaining window duration is ≤ 2 minutes, its latest ranging data is retained first to reserve a buffer for subsequent topology updates.

[0031] The local filter (using the UKF or EKF algorithm) uses preprocessed inter-satellite link data as its core, combined with orbital dynamics and clock error models, to complete "time update" (i.e., state prediction): Define the state vector: The state vector contains the satellite's core operational parameters. These include position (x, y, z) and velocity (v) in a geocentric Cartesian inertial coordinate system. x v The system comprehensively covers the dynamic characteristics of orbit and clock difference, including time-varying clock difference parameters (clock difference value, clock difference drift rate, and clock difference drift acceleration).

[0032] Considering the gravitational force of the Earth's center of mass (including the J2 perturbation, i.e. the effect of the Earth's non-spherical symmetry), the gravitational perturbation of the Sun and Moon, and the solar radiation pressure perturbation, calculate the acceleration of the satellite under the three forces to describe the changes in position and velocity (velocity is the rate of change of position, and acceleration is the rate of change of velocity).

[0033] By taking into account the characteristics of spaceborne clocks (such as rubidium clocks), a clock drift time constant (usually 200 seconds) is set to describe the change of clock drift value over time. The clock drift rate decays over time, while a small amount of random noise is added to ensure that the model closely matches the actual clock drift variation.

[0034] The noise level is allocated based on the actual measurement accuracy. The noise variance for the position parameter is set to 0.01 square meters, the velocity parameter to 0.06 square meters per second, the clock error parameter to 0.001 square meters, and the clock error drift rate to 0.01 square meters per second, to avoid excessive or insufficient noise affecting the filtering accuracy.

[0035] Using the state estimation result of the previous filtering cycle as the initial value, and combining it with orbital dynamics and clock error models, the predicted state value at the current moment is calculated through numerical integration (such as the Runge-Kutta method). Simply put, it predicts the current position, velocity, and time-varying clock error parameters based on the satellite's motion patterns and clock error variation trends.

[0036] The covariance matrix reflects the accuracy level of the state estimation. Its calculation requires consideration of the state transition matrix and system noise. The state transition matrix describes the relationship between the previous cycle state and the current predicted state. Adding the influence of system noise yields the covariance prediction matrix for the current moment, thus quantifying the uncertainty of the prediction result.

[0037] The local filter combines GNSS data and dynamically adjusts the observation weights according to the GNSS signal quality to complete the "measurement update" (i.e., correct the prediction results).

[0038] Based on the GNSS signal quality labeling results in the first step, the system assigns weights to GNSS observations and inter-satellite link observations (the sum of the two is 1): When the GNSS signal is "good": GNSS weight accounts for 70%~80%, and inter-satellite link weight accounts for 20%~30%, with GNSS data dominating; When the GNSS signal is "attenuated": GNSS weight drops to 30%~50%, and inter-satellite link weight rises to 50%~70%, balancing the contributions of the two types of data; When the GNSS signal is "interrupted": GNSS weight is 0, and the system relies entirely on inter-satellite link data.

[0039] Meanwhile, GNSS observation data is further weighted according to SNR. Satellites with SNR ≥ 15dB have a weight of 1.0, while the weight of satellites with SNR 10~15dB decreases linearly proportionally (to 0.6 at 10dB), ensuring that high-quality GNSS data can play a greater role.

[0040] The core logic of GNSS pseudorange observations is "inter-satellite distance + clock bias + observation noise". That is, the pseudorange measured by the onboard receiver is equal to the actual distance between the LEO satellite and the GNSS satellite, plus the onboard clock bias (which causes ranging deviation), and then superimposed with a small amount of observation noise.

[0041] Kalman gain is a core parameter for correction strength. Its calculation requires consideration of the state prediction value, covariance prediction matrix, observation sensitivity (i.e., the relationship between the observed value and the state parameter), and observation noise. Simply put, the larger the gain, the stronger the influence of the observed data on the correction prediction result; conversely, the weaker the influence, thus balancing the reliability of prediction and observation.

[0042] First, the "observation residual" (the difference between the measured pseudorange of GNSS and the theoretical pseudorange calculated based on the predicted state) is calculated. Then, the residual is converted into a correction value through Kalman gain to adjust the state prediction value and obtain a more accurate final state estimate (including position, velocity, and time-varying clock error parameters).

[0043] Based on the Kalman gain and observation noise, the covariance matrix is ​​updated to quantify the accuracy level of the corrected state estimate. The smaller the corrected covariance matrix, the more reliable the state estimate result.

[0044] The local filter ultimately outputs the "state estimate" and its "covariance matrix" for each satellite: the state estimate contains high-precision position, velocity, and time-varying clock error parameters; the covariance matrix reflects the estimation accuracy of these parameters. This result will be uploaded to the main filter to estimate the overall constellation reference error; and it will also serve as the initial value for the next filtering cycle, supporting continuous orbit determination calculations.

[0045] In some embodiments, the state estimates and covariance matrices of each satellite are input into the main filter to obtain the overall reference error estimate of the constellation. This includes: calculating the closure error of the geometric closed loop formed by the inter-satellite links based on the state estimates of each satellite to construct a reference error observation; establishing a mathematical relationship between the closure error of the geometric closed loop and the overall rotation / translation parameters of the constellation; using the state estimates and covariance matrices of each satellite as input to the main filter, and estimating the overall reference error of the constellation based on the mathematical relationship and the reference error observation; wherein the overall reference error includes three orbital translation parameters, three rotation parameters, and one clock error reference parameter.

[0046] In this embodiment of the invention, the main filter first utilizes the inter-satellite link topology to construct a "geometric closed loop" and calculate the closure error, transforming the relative measurement data into an "observational basis" that reflects the baseline error. After receiving the state estimates of all satellites (including satellite positions and time-varying clock error parameters), the main filter first sorts out the connection relationships of the inter-satellite links. A "closed topology" consisting of 3 (triangular ring) or 4 (quadrilateral ring) satellites is selected: for example, satellite A-satellite B-satellite C-satellite A forms a triangular ring, and satellite A-satellite B-satellite D-satellite C-satellite A forms a quadrilateral ring. During selection, priority is given to closed loops with "high link quality and uniform satellite distribution" (e.g., triangular ring interior angle ≥ 30°, quadrilateral ring diagonal length difference ≤ 500 km) to ensure the geometric stability of the closed loop and reduce subsequent calculation errors.

[0047] Based on the state estimate (position parameters) of each satellite, the "theoretical relative distance" of each inter-satellite link segment within the closed loop is calculated. For example, in the triangular loop ABCA, the straight-line distance from A to B is calculated based on the positions of satellites A and B; similarly, the straight-line distances from B to C and from C to A are calculated. These distances only reflect the relative positional relationship between satellites and do not include reference errors. The "inter-satellite link measured ranging data" (pre-processed and corrected for link noise) from the local filtering output of each satellite is retrieved and corrected in conjunction with the "time-varying clock error parameters" in the state estimate. Clock errors can cause deviations in measured ranging (for example, if the clock error of satellite A is too large, the measured ranging from A to B will be longer than the actual value). Therefore, the influence of clock error on measured ranging needs to be deducted based on the difference in clock errors between satellites A and B to obtain the "actual ranging after removing clock error interference".

[0048] The "actual ranging after eliminating clock interference" within the same closed loop is algebraically summed along the "link direction": for example, in the triangular loop A→B→C→A, the ranging from A to B is taken as positive, from B to C as positive, and from C to A as negative (or uniformly set as positive / negative in clockwise / counterclockwise direction). The summation result is the "closure error". Ideally, if the constellation has no overall reference error (no rotation, no translation), the closure error should be 0. When the constellation rotates or translates as a whole, the actual ranging of each link within the closed loop will shift as a whole, causing the closure error ≠ 0. For example, if the constellation translates as a whole along the x-axis by 10 meters, the actual ranging of all links within the loop will increase (or decrease) the deviation related to the translation direction, and the final closure error will be significantly different from 0. Therefore, the closure error can be directly used as a "reference error observation" to reflect the degree of reference drift.

[0049] The "rotation parameters" and "translation parameters" in the overall constellation reference error have different effects on the ranging of each link within the closed loop, and the main filter needs to be decomposed separately. When the constellation as a whole translates along a certain direction of x, y, or z (for example, translating Δx along the x-axis), all satellites within the closed loop will move synchronously by the same distance, resulting in the ranging changes of each link within the loop being "equal in magnitude and consistent in direction". For example, in the triangular loop ABCA, A, B, and C synchronously translate Δx along the x-axis. The ranging from A to B, from B to C, and from C to A will all increase (or decrease) by a fixed value related to Δx. Ultimately, these changes will be superimposed in the closure error as "3 × Δx-related deviation" (the triangular loop has 3 links), that is, the closure error is "linearly positively correlated" with the translation parameters.

[0050] When the constellation as a whole rotates around one of the x, y, or z axes (e.g., around the z-axis by an angle Δγ), the changes in link ranging in different directions vary. For example, if satellite A shifts to the right after rotation and satellite B shifts to the left, the ranging from A to B will increase; while the ranging from B to C may decrease. The magnitude of the change is related to the rotation angle Δγ and the initial distance between the satellites (the longer the initial distance, the greater the ranging change for the same rotation angle). The main filter analyzes the geometry of the closed loop (e.g., the angles between each link and the rotation axis) to determine the relationship between the rotation angle Δγ and the ranging change of each link segment, ultimately summarizing this into the closure error, forming a "nonlinear correlation between closure error and rotation parameters" (the larger the rotation angle, the more pronounced the nonlinear change in closure error).

[0051] Based on the aforementioned influence patterns, the main filter establishes a correlation logic of "closure error = translational influence + rotational influence + minor observation noise." For example, for a triangular loop, the magnitude of the closure error equals "the weighted sum of 3 times the translation parameters (such as Δx, Δy, Δz)" plus "the sum of the products of the rotation parameters (such as Δα, Δβ, Δγ) with the link length and the included angle," and then superimposed a small amount of observation noise (residual errors from inter-satellite links). This logic does not require complex formulas; its core is to quantify "how much deviation in closure error is caused by every 1-meter translation and every 0.1° rotation" through geometric analysis, allowing the specific translational and rotational parameters to be deduced from the closure error later.

[0052] The main filter combines local state estimation data, closed-loop difference observations, and correlation logic to calculate the final overall reference error (including 3 translation, 3 rotation, and 1 clock error reference parameters) through filtering iteration.

[0053] The main filter receives the "state estimates" and "covariance matrices" from each satellite. The covariance matrix reflects the accuracy of the state estimates (the smaller the matrix value, the more accurate the position and clock error estimates). Therefore, the main filter assigns weights to the data from different satellites based on the covariance matrix. Satellites with smaller covariance matrices (higher state estimation accuracy) have higher weights in the closed-loop calculation of their position data (e.g., 1.0); satellites with larger covariance matrices (lower state estimation accuracy) have lower weights (e.g., reduced to 0.6-0.8) to prevent low-precision data from interfering with the baseline error estimation.

[0054] The main filter employs algorithms such as EKF (Extended Kalman Filter) to iteratively calculate the baseline error in two steps: prediction and correction. The baseline error estimate from the previous cycle is used as the initial value (if it's the first calculation, the initial value is set to 0, meaning no baseline error is assumed). Considering the stability of the constellation's operation (e.g., the baseline drift of the LEO constellation won't change abruptly in a short time), a preliminary value for the baseline error at the current moment is predicted (e.g., the predicted translation parameter Δx is approximately 1.02 times that of the previous cycle, while the rotation parameter Δγ remains essentially unchanged). The closure error observations are substituted into the mathematical relationship established in the second step to calculate the "theoretical closure error corresponding to the predicted baseline error," which is then compared with the "actually calculated closure error" to obtain the "error residual" (the difference between the theoretical and actual values). Based on the magnitude of the error residual, the predicted baseline error value is adjusted. For example, if the actual closure error is 0.5 meters larger than the theoretical value, it indicates that the predicted translation parameter is too small, and the estimated value of Δx needs to be increased proportionally. Simultaneously, the weighting of the input data ensures that the correction process prioritizes reference to high-precision satellite data. Repeat the "prediction-correction" step until the change in the baseline error estimate between two consecutive steps is ≤0.01 meters (translation parameter) or ≤0.005° (rotation parameter), then the result is considered converged and the iteration stops.

[0055] The final estimated overall baseline error includes three types of parameters: The three orbital translation parameters correspond to the offsets (in meters) of the entire constellation along the x, y, and z axes of the geocentric rectangular inertial coordinate system. For example, Δx = 0.3 meters, Δy = -0.2 meters, and Δz = 0.1 meters represent that the entire constellation is offset 0.3 meters to the right along the x-axis, 0.2 meters to the left along the y-axis, and 0.1 meters up along the z-axis.

[0056] The three rotation parameters correspond to the rotation angles (in degrees or radians) of the constellation around the x, y, and z axes, respectively. For example, Δα=0.05°, Δβ=-0.03°, and Δγ=0.02° represent that the constellation rotates 0.05° clockwise around the x-axis, 0.03° counterclockwise around the y-axis, and 0.02° clockwise around the z-axis.

[0057] A clock bias reference parameter: reflects the clock bias reference deviation of the entire constellation (unit: meters, corresponding to pseudorange deviation). For example, Δl_base=0.04 meters means that the clock bias of all satellites in the constellation has an overall deviation of 0.04 meters relative to absolute time. The clock bias of each satellite needs to be uniformly corrected based on this parameter.

[0058] These parameters will be broadcast from the main filter to the local filters of all satellites to correct the state estimates of each satellite, ultimately achieving high-precision autonomous orbit determination for the entire constellation.

[0059] In some embodiments, the overall reference error estimate is fed back to each local filter to correct the state estimate of each satellite, thereby obtaining the final high-precision orbit parameters and clock error parameters of each satellite. This includes: the main filter broadcasting the overall reference error estimate to the local filters of each satellite; correcting the state estimate of each satellite's local filter based on the overall reference error estimate; and outputting the corrected state estimate as the final high-precision orbit parameters and clock error parameters of each satellite.

[0060] In this embodiment of the invention, after the main filter completes the estimation of the overall constellation reference error (including three orbital translation parameters, three rotation parameters, and one clock bias reference parameter), it needs to synchronize the error information to the local filters of all satellites in the constellation through inter-satellite links to ensure the "uniformity of reference correction".

[0061] The main filter first performs a "validity check" on the overall baseline error estimate. It checks whether each parameter is within a reasonable physical range (e.g., translation parameters are usually ≤1 meter, rotation parameters ≤0.1°, and clock error baseline parameters ≤0.1 meter). If a parameter is out of range, it is judged as an estimation anomaly, and the previous round of baseline error estimation is re-executed (to avoid incorrect correction). After the verification is passed, the seven parameters (Δx, Δy, Δz translation, Δα, Δβ, Δγ rotation angle, Δl_base clock error baseline) are encapsulated in the format of "satellite ID-parameter value-estimation accuracy" (e.g., labeling the covariance value of each parameter to reflect the estimation reliability) to form a standardized "baseline correction instruction".

[0062] The main filter broadcasts a "reference correction command" through the constellation's inter-satellite link topology (preferably using high-bandwidth, low-latency master-slave links): If the main filter is located on the ground, the command is first sent to the "master satellite" (a pre-designated satellite with strong communication capabilities) in the low-Earth orbit constellation, and then the master satellite forwards the command to the surrounding satellites through inter-satellite links, forming a hierarchical broadcast of "ground → master satellite → slave satellite". If the main filter is located on the primary star, it will broadcast directly from the primary star to the entire constellation.

[0063] To ensure successful reception by all satellites, a "retransmission confirmation mechanism" is employed. After receiving the command, each satellite sends a "reception confirmation signal" back to the main filter. The main filter compiles the confirmation information and, for satellites that do not respond (possibly due to temporary link interruptions), retransmits the command through a backup neighbor link until all satellites have received the command (ensuring broadcast coverage ≥ 99.9%).

[0064] After receiving the "reference correction command," the local filters of each satellite first parse the parameter types and values ​​in the command. They distinguish between three types of parameters: translation, rotation, and clock bias references, extracting the corresponding values ​​(e.g., Δx = 0.3 meters, Δα = 0.05°, Δl_base = 0.04 meters). At the same time, they read the estimation accuracy information of each parameter (e.g., the covariance value of Δx = 0.01 square meters, indicating that the estimation of the translation parameter is relatively reliable), providing a basis for the weight allocation during subsequent correction.

[0065] Each satellite's local filter, based on the resolved overall reference error, corrects its previously output state estimate in two steps: "orbit parameter correction" and "clock error parameter correction" (eliminating the systematic error caused by the overall constellation drift).

[0066] The local filter first targets the "position parameters (x, y, z)" and "velocity parameters (v)" in the state estimate. x v "v_z")" respectively eliminates the effects of translation and rotation parameters: Translation parameter correction: Overall constellation translation will cause a deviation in the same direction and magnitude in the position estimates of all satellites. Therefore, the correction is directly achieved by "local position estimate - corresponding translation parameter". For example, if the local filter estimates the satellite position as (x_local, y_local, z_local), and the constellation is translated Δx along the x-axis, Δy along the y-axis, and Δz along the z-axis, then the corrected positions are: x_corr = x_local - Δx, y_corr = y_local - Δy, z_corr = z_local - Δz. (Example: If the local estimate x_local = 10000 km and Δx = 0.3 m, then the corrected x_corr = 10000 km - 0.0003 km = 9999.9997 km, eliminating the 0.3 m deviation caused by the overall translation).

[0067] Rotation parameter correction: The overall rotation of the constellation will cause an angular shift in the satellite positions relative to the geocentric inertial frame. The rotation deviation needs to be calculated and corrected based on the satellite's initial position and rotation angle. The specific logic is as follows: If the constellation rotates around the x-axis by an angle Δα, the y and z coordinates of the satellites will be affected (the x-coordinate remains unchanged). It is necessary to calculate the y and z deviations after rotation based on trigonometric function relationships (e.g., y-direction deviation = z_local × sinΔα, z-direction deviation = y_local × sinΔα), and then use "corrected translation position - rotation deviation" to obtain the final rotation correction position. Similarly, when rotating around the y-axis by an angle Δβ, the x and z coordinates are corrected; when rotating around the z-axis by an angle Δγ, the x and y coordinates are corrected.

[0068] The correction logic for the velocity parameter is consistent with that for the position. The velocity deviation caused by rotation and translation is proportional to the rate of change of the position deviation. The local filter obtains the velocity deviation by "differentiating the position deviation with respect to time", and then completes the velocity correction by "local velocity estimate - velocity deviation".

[0069] Dynamic adjustment of correction weights: The local filter adjusts the correction intensity based on the "parameter estimation accuracy" provided by the main filter. If the covariance of a parameter is small (the estimation is reliable, such as Δx covariance = 0.01 square meters), it is corrected by 100%; if the covariance is large (the estimation reliability is low, such as Δγ covariance = 0.001°²), it is corrected by 80% to 90% (to avoid introducing new errors by unreliable parameters).

[0070] The local filter targets the "time-varying clock error parameters (clock error value, clock error drift rate)" in the state estimate. It utilizes the clock error reference parameter Δl_base to eliminate the overall clock error bias of the constellation. The local clock error estimates of all satellites within the constellation may exhibit an "overall offset" due to initial deviations in their respective clock error models (e.g., all satellite clock errors are 0.04 meters larger than the absolute time). Δl_base precisely quantifies this overall offset. Therefore, the correction logic is: "Local clock error estimate - Clock error reference parameter Δl_base". For example, if the locally estimated clock error value is 0.08 meters and Δl_base = 0.04 meters, then the corrected clock error value = 0.08 meters - 0.04 meters = 0.04 meters, unifying the clock errors of the entire constellation to the same reference (consistent with the absolute time bias). Simultaneously, the correction of the clock error drift rate needs to consider the changing trend of Δl_base. If the previous cycle Δl_base = 0.03 meters and the current Δl_base = 0.04 meters, it means that the overall clock drift rate is 0.01 meters / cycle. The local filter will add this drift rate to its own clock drift rate estimate to ensure the continuity of clock correction.

[0071] After the local filter completes the correction of the orbit and clock error parameters, the correction effect needs to be confirmed through "accuracy verification" before the final high-precision parameters are output to provide a basis for the subsequent operation of the satellite (such as orbit control and communication scheduling).

[0072] The local filter first performs a "consistency check" on the corrected state estimate: comparing the parameter changes before and after correction (e.g., the position change should be consistent with the calculated values ​​of translation + rotation parameters, and the clock error change should be equal to Δl_base). If the deviation exceeds 0.01 meters (position) or 0.001 meters (clock error), it is judged as a correction anomaly, and the correction steps are repeated. The corrected position parameters of neighboring satellites are retrieved, and the inter-satellite distance to its own corrected position is calculated and compared with the measured relative distance of the inter-satellite link (pre-processed). If the difference is ≤0.05 meters, it indicates that the corrected parameters are consistent at the constellation level (no local deviation); if the difference is too large, a "correction anomaly signal" is fed back to the main filter, awaiting a new round of reference error updates.

[0073] After successful verification, the local filter outputs the "corrected state estimate," which is the final high-precision parameter of each satellite, specifically including two types of core parameters: High-precision orbital parameters: position (x_corr, y_corr, z_corr, accuracy ≤ 0.2 meters) and velocity (v) in a geocentric rectangular inertial coordinate system. x _corr、v _corr, v_z_corr (accuracy ≤ 0.001 m / s), these parameters are directly used for satellite orbit control (such as attitude adjustment of ion thrusters) and orbit prediction for ground telemetry and control.

[0074] High-precision clock bias parameters: corrected time-varying clock bias values ​​(accuracy ≤ 0.03 meters) and clock bias drift rate (accuracy ≤ 5e-9 meters / second). These parameters are used for clock bias compensation of onboard GNSS receivers and time synchronization of inter-satellite link communication (ensuring the time accuracy of signal transmission between satellites).

[0075] Meanwhile, the local filter stores the corrected parameters and the corresponding covariance matrix (reflecting the accuracy level) in the on-board memory. On the one hand, this serves as the initial state for the next filtering cycle's "time update." On the other hand, it is uploaded to the ground tracking and control station or constellation management system as needed to support global orbit monitoring and maintenance.

[0076] After outputting the parameters, the local filter will adjust its orbit determination mode synchronously based on the reference error broadcast by the main filter. If the overall reference error is small (e.g., translation ≤ 0.2 meters, rotation ≤ 0.05°), it indicates that the current GNSS signal or inter-satellite link data quality is good, and the original orbit determination mode (e.g., GNSS-dominant mode) is maintained. If the reference error is large (e.g., translation > 0.5 meters), an "environmental anomaly signal" is fed back to the system's intelligent decision-making module, triggering a mode switch (e.g., switching from GNSS-dominant mode to inter-satellite link enhancement mode) to ensure the continued stability of subsequent orbit determination.

[0077] In some embodiments, before acquiring inter-satellite link data, the method further includes: determining the link reachability with potential neighboring satellites within a future time window based on an orbit prediction model, and defining an effective neighbor window period; when the effective neighbor window period meets a preset update triggering condition, adaptively initiating a new round of neighbor screening process; In this embodiment of the invention, SGP4 / SDP4 (Simplified General Perturbations 4) is preferred, as this model can balance prediction accuracy and computational efficiency. For the orbital characteristics of LEO satellites, the model has built-in simplified calculation logic for Earth's non-spherical gravity (J2 term perturbation), lunar and solar gravitational perturbation, and solar radiation pressure perturbation. With limited computing power on the satellite, it can achieve a position prediction accuracy of ≤10 meters within the next 5-10 minutes (meeting the requirements for link reachability determination). If the satellite is equipped with a higher-performance onboard computer, some parameters (such as atmospheric drag correction terms) of the HPOP (High Precision Orbit Propagator) can be superimposed, further improving the prediction accuracy to ≤1 meter and reducing reachability misjudgments caused by orbit prediction errors.

[0078] The local filter first defines the "potential neighbor search range" based on its current orbital position: Spatial range: A spherical region with a radius of 10-20 kilometers centered on itself (the effective communication distance of LEO satellite inter-satellite links is usually 5000 kilometers, but the relative motion of short-range neighbors is smoother, the link stability is higher, and the computational load of subsequent topology optimization can be reduced). Time range: the next 15-30 minutes (covering at least 2 filter cycles to ensure the window period is long enough to support data acquisition); The system accesses the "constellation satellite orbit database" stored on the satellite (which updates the latest orbital parameters of all satellites in real time), filters out satellites that will enter the above-mentioned space range within the next 15-30 minutes, and forms a "potential neighbor candidate list", marking the ID, current orbital position and velocity of each candidate satellite.

[0079] For each satellite in the "potential neighbor candidate list", we determine whether the link is reachable within the future time window based on three dimensions: geometric conditions, communication quality, and obstruction.

[0080] Elevation angle ≥ 5°: The relative positions of satellites are calculated every 10 seconds within the future time window using orbital prediction to ensure that the elevation angle (the angle between the satellite line and the local horizontal plane) between the LEO satellite and potential neighbors is always ≥ 5°. This avoids the link being blocked by the curvature of the Earth (when the elevation angle is too low, the signal needs to penetrate a thicker layer of atmosphere, resulting in severe attenuation). Relative speed ≤ 1 km / s: Calculate the relative speed between the two. If the relative speed exceeds 1 km / s within a certain period of time, the link is prone to disconnection due to excessive relative motion after it is established (signal tracking becomes more difficult and packet loss rate soars). It is determined that the link is unreachable within this period of time. Link distance ≤ 5000 km: The straight-line distance between satellites is calculated based on the predicted position. If it exceeds the maximum communication distance of the inter-satellite link terminal (usually 5000 km), the signal strength is insufficient and it is determined to be unreachable.

[0081] Combining the relative distance between satellites and atmospheric attenuation models (such as ionospheric electron concentration forecast data), the link signal-to-noise ratio (SNR) and latency are predicted within the future time window: if SNR ≥ 12dB (meets the minimum signal requirement for pseudocode ranging) and latency ≤ 100 milliseconds (ensures data real-time performance), the link quality is deemed to be up to standard; if SNR < 12dB or latency > 100 milliseconds within a certain period, even if the geometric conditions are met, the link is deemed unreachable.

[0082] Using onboard solar and lunar position forecast data, we check whether there are any situations where the link is "interfered by strong solar light" (such as the sun shining directly on the link terminal, causing increased reception noise) or "obstructed by the moon / other satellites" within the future time window: if the obstruction time accounts for more than 30% of the total window time, the link stability of that neighbor is judged to be poor and it is excluded from the candidate list.

[0083] For potential neighbors that pass all reachability tests, a "valid neighbor window" is defined. This is a continuous period of time during which the link simultaneously satisfies "geometric reachability, quality compliance, and no obstruction." Specify the start and end times of the window (accurate to the second), and label the key parameters within this time period: neighbor satellite ID, predicted relative position / velocity every 10 seconds, and predicted signal-to-noise ratio / delay value; If a neighbor's reachability time period is split by a brief obstruction (e.g., a 5-minute obstruction in the middle), it is split into two "sub-window periods", and the time range and parameters of each sub-window are marked. The window period information of all valid neighbors is stored in the local "neighbor link management module" and synchronized to the main filter (for global topology monitoring). Subsequent inter-satellite link data collection is only carried out within this window period.

[0084] To avoid data acquisition interruptions due to relative satellite motion and link quality degradation, it is necessary to monitor the status of the effective neighbor window in real time. When the preset trigger conditions are met, a new round of neighbor screening should be initiated to ensure the continuous availability of the link.

[0085] The local filter uses a combination of timer monitoring and real-time parameter detection to determine if the following trigger conditions are met; filtering is initiated if any one of these conditions is satisfied: The remaining duration of the current effective neighbor window is calculated in real time. When the remaining time is ≤2 minutes, filtering is triggered. Because LEO satellites move relatively fast, the 2-minute buffer time ensures that the old link can still provide data before the new link is established, avoiding the gap period of "link disconnection - data loss".

[0086] The system collects the actual operating parameters (signal-to-noise ratio, packet loss rate, and latency) of the current neighbor links in real time. If, for three consecutive sampling periods (1 second per period), the following occurs: SNR decreases from ≥12dB to <10dB, packet loss rate increases from ≤0.1% to >1%, and latency increases from ≤100 milliseconds to >200 milliseconds, the link quality is determined to be unable to meet the orbit determination requirements, and filtering is triggered.

[0087] When a satellite performs orbital adjustments (such as ion thruster activation) or receives a notification from the main filter that a neighboring satellite's orbit is abnormal (such as a neighboring satellite experiencing orbital drift), the original forecast basis for the effective window period becomes invalid, and screening must be triggered immediately to reassess neighbor reachability.

[0088] The main filter sends a "global topology connectivity report" to each satellite every 30 seconds. If the report indicates that the current satellite has fewer than 3 neighbors (which cannot form a closed loop and affect subsequent baseline error estimation) or there is a risk of "topology islands" (insufficient connectivity with other satellites), then a screening process is triggered to add new neighbors to meet the topology requirements.

[0089] When any update trigger condition is detected, the local filter initiates adaptive filtering according to the following steps to ensure that the process does not affect the current data acquisition: First, mark the currently used neighbor link as a "transition link" and continue to collect its data (to avoid interruption), while suspending new data requests for low-quality links (such as links that have triggered quality degradation conditions); Set the filtering priority: prioritize filtering potential neighbors that can form a "triangle / quadrilateral closed loop" with the current satellite (to ensure that the subsequent main filter can calculate the closure error), then select neighbors with low communication energy consumption (short link distance, high signal-to-noise ratio), and finally supplement neighbors that meet basic reachability, so as to ensure that the filtering results take into account both "loop formation requirements" and "energy consumption control".

[0090] The latest orbit prediction model is invoked (updating the current orbit parameters of itself and all satellites), the "potential neighbor search range" is redefined (still within the next 15-30 minutes, with a radius of 10-20 kilometers), and reachability is reassessed in the order of "geometric conditions - communication quality - obstruction conditions", generating a new list of "effective neighbor windows".

[0091] For newly selected valid neighbors, a "link establishment request" is sent to the candidate neighbors through the current "transition link". The request includes the time range of the new window period and its own orbit parameters. Receive the "confirmation response" from the candidate neighbor. If the other party also lists the current satellite as a valid neighbor (two-way confirmation), mark it as a "confirmed neighbor" and negotiate the data transmission protocol (such as sampling frequency and data format) within the window period. Once the new neighbor link is established (usually taking ≤10 seconds), gradually switch data collection to the new link while disconnecting the original low-quality / expiring "transition link" to avoid wasting resources.

[0092] The new list of valid neighbor windows, the IDs of confirmed neighbors, and link parameters are updated to the local "neighbor link management module" and synchronized to the main filter for updating global topology information; Record the triggering reason, screening time (usually ≤30 seconds), and final number of neighbors (3-5) for this screening, and store the "screening log" in the on-board storage for subsequent troubleshooting and algorithm optimization.

[0093] Through the above process, it can be ensured that LEO satellites always have stable and usable neighbor links before inter-satellite link data acquisition, and can adaptively adjust neighbor selection according to dynamically changing orbits and environments, thus providing a guarantee for high-quality acquisition of subsequent inter-satellite link observation data.

[0094] The neighbor optimization screening under multi-objective constraints is performed. A weighted function is constructed with topology ring integrity, communication energy consumption and link quality as optimization objectives, and the weights of each objective are dynamically adjusted according to the real-time operation status of the constellation. Based on the calculation results of the weighted function, the optimal neighbor subset is selected from the potential neighbors within the effective neighbor window period. In this embodiment of the invention, before performing neighbor optimization screening, the three abstract objectives of "topology loop integrity," "communication energy consumption," and "link quality" must be transformed into calculable quantitative indicators to ensure that each objective can be measured by specific parameters. For topology loop integrity, the core is to determine whether potential neighbors can form a triangular or quadrilateral closed loop with the current satellite and existing neighbors. Neighbors that can participate in forming two or more closed loops are assigned an indicator value of 1.0; those that can form one closed loop are assigned 0.8; those that can only establish a link but cannot form a loop are assigned 0.3; and those that cannot form a loop at all are assigned 0.1. This prioritizes retaining key neighbors that support the baseline error estimation. The communication energy consumption indicator is calculated by combining link distance and hardware characteristics. The longer the link distance and the lower the antenna gain, the higher the energy consumption. The theoretical energy consumption value of all potential neighbors is first calculated, and then normalized according to "maximum energy consumption is 1, minimum energy consumption is 0". The closer the indicator value is to 0, the lower the energy consumption, which conforms to the on-board endurance constraints. Link quality indicators need to comprehensively consider signal-to-noise ratio (SNR), latency, and packet loss rate. First, normalize each of the three into sub-indices of 0-1 (e.g., SNR≥15dB is 1.0, ≤10dB is 0.2). Then, weight the values ​​by 50% for SNR, 30% for latency, and 20% for packet loss rate and sum them up. The closer the index value is to 1, the more reliable the link is, and the more accurate the observation data can be guaranteed.

[0095] Based on the three quantified indicators, a linear weighted function is constructed with the objective of "optimal overall performance (minimum function value)", transforming multi-objective optimization into a single numerical comparison. The function is defined as J = ω1 × (1 C)+ω2×E+ω3×(1 Q), where C is the topology ring integrity index, E is the communication energy consumption index, and Q is the link quality index. During the design, C and Q were intentionally converted into (1... C) and (1) The weighted function (Q) unifies the optimization direction of the three indicators (all are "smaller values ​​are better"), avoiding conflicts between objectives. For example, if a neighbor has good ring integrity (C=0.9), low energy consumption (E=0.2), and high quality (Q=0.95), the calculated J value after substituting it into the function will be smaller, representing better overall performance. The core function of the weighted function is to integrate the influence of the three objectives, avoiding the deterioration of other objectives due to the optimization of a single objective (such as selecting only low-energy neighbors but failing to form a ring), and reflecting the priority under different scenarios through weight allocation.

[0096] The allocation of weights ω1 (loop formation), ω2 (energy consumption), and ω3 (quality) is not fixed and needs to be flexibly adjusted according to the real-time operating status of the constellation to ensure that the selection results are suitable for the current needs. When the GNSS signal is good (visible satellites ≥ 6, GDOP ≤ 3), the demand for loop formation is low, so ω1=0.3, ω2=0.4, and ω3=0.3 are set, and energy consumption is prioritized. If GNSS is completely interrupted, the inter-satellite link closed loop needs to be relied upon to maintain the baseline, and the priority of loop formation is the highest, so ω1=0.6, ω2=0.2, and ω3=0.2 are adjusted. From the perspective of energy consumption, when the onboard battery power is ≥ 50%, the normal weight is maintained. When the power drops to 30%-50% (energy consumption is tight), ω2 is increased by 0.1 (e.g., ω2=0.5). When the power is < 30% (energy consumption is urgent), ω2 is increased to 0.6, and short-distance low-energy neighbors are prioritized. If spatial electromagnetic interference (such as solar flares) is detected, the link quality stability will decrease, and ω3 will be increased by 0.1-0.2 (e.g., ω3=0.5) to ensure that the selected link can resist interference and reduce data packet loss.

[0097] After completing indicator quantification, function construction, and weight adjustment, the optimal neighbor subset is selected through a three-step process: "combined calculation - sorting and filtering - loop construction verification." First, all candidate combinations of "3-5 neighbors" are generated (the number satisfies the loop construction requirement while avoiding redundant energy consumption). For example, if there are 5 potential neighbors, all combinations of 3, 4, and 5 neighbors are generated. For each combination, the J-value of each neighbor is calculated and averaged (representing the overall performance of the combination). Combinations are sorted from smallest to largest by average J-value, and the top 3 optimal combinations are initially retained. Next, the loop construction capability is verified. It is checked whether the neighbors within the combination can form at least one closed loop with the current satellite (e.g., 3 neighbors + the current satellite forming a triangle). If a combination cannot form a loop (e.g., there is no overlapping window period between neighbor links), it is eliminated, and the next combination is verified. Finally, the combination with the "smallest average J-value and capable of forming a loop" is determined as the optimal neighbor subset. The IDs and effective window periods of the neighbors within the subset are stored in the onboard management module and synchronized to the main filter, preparing for subsequent bidirectional link consensus establishment and data collection.

[0098] Each satellite establishes a consensus on links by engaging in bidirectional information exchange with candidate neighbor satellites based on the optimal subset of neighbors. At the same time, the overall topology connectivity of the constellation is periodically monitored through a global coordination mechanism. When a node that does not meet the preset connectivity requirements is detected in a local topology, a topology adjustment command is proactively issued to repair it.

[0099] In this embodiment of the invention, after each satellite determines its optimal neighbor subset, it will proactively initiate bidirectional information exchange with candidate neighbors within the subset. The core purpose is to ensure that both parties agree on the feasibility of link establishment through parameter alignment and willingness confirmation, thereby avoiding data transmission failure caused by one-way links. At the start of the interaction, the initiating satellite will send a "link establishment request" to the candidate neighbor through the inter-satellite link. The request includes its own key parameters: current orbital position and forecast data for the next 5 minutes, start and end times of the effective neighbor window, the comprehensive performance value of the candidate neighbor in the weighting function (such as the J value), and the preset protocol for subsequent data transmission (such as the sampling frequency of pseudocode ranging and the data frame format). Upon receiving a request, a candidate neighbor first checks if its optimal neighbor subset includes the initiating satellite. Then, it compares the effective window periods (ensuring an overlap of ≥10 minutes to meet data acquisition requirements) and orbital parameter deviations (relative position change rate ≤1km / s). If all conditions match, it sends an "agreement response" to the initiating satellite, along with its own orbital forecast, window period details, and data transmission protocol confirmation. If there are mismatches (e.g., no window overlap, excessive orbital deviation), it sends a "rejection response" explaining the reason (e.g., "effective window overlap is only 3 minutes, not meeting requirements"). Upon receiving the "agreement response," the initiating satellite reconfirms the consistency of parameters (e.g., window overlap duration, transmission protocol). Once confirmed, it adds the candidate neighbor to its "effective neighbor list" and sends a "consensus confirmation notification." Upon receiving the notification, the candidate neighbor also adds the initiating satellite to its own effective neighbor list. At this point, link consensus is formally established, and both parties simultaneously begin preparations for inter-satellite link data acquisition (e.g., adjusting transmit and receive antenna angles, warming up the ranging module).

[0100] To prevent isolated satellites from becoming "islands" in the constellation topology due to local link disconnections (where a satellite lacks sufficient connectivity with other satellites to participate in global orbit determination), the system establishes a global coordination mechanism through a master filter (usually deployed on the constellation's main satellite or at the ground control center) to periodically monitor topology connectivity. The monitoring cycle is set to 30 seconds (matching the dynamic change rate of LEO satellite orbits, avoiding both excessively frequent monitoring that increases computational burden and long intervals that may lead to omissions). During each monitoring cycle, the master filter sends a "topology information reporting command" to all satellites, requiring them to report their current list of valid neighbors (including neighbor IDs, link status, and remaining window duration) and real-time quality parameters of established links (signal-to-noise ratio, packet loss rate). Upon receiving the command, each satellite packages and uploads this information to the master filter within 5 seconds. After receiving all satellite information, the master filter constructs a global topology graph, presenting the connectivity status of the entire constellation in the form of "nodes (satellites) - edges (valid links)," and performs verification according to preset connectivity requirements. Preset requirements typically include: each satellite has at least 3 effective neighbors (ensuring at least one closed loop can be formed), no isolated nodes (the number of connections between satellites and other satellites is at least 2 to avoid being unable to participate in inter-satellite coordination), and no local small clusters (any two satellites can be connected via ≤3 hops links to ensure that the baseline error correction command can be globally synchronized). The main filter traverses the topology graph, checking the number of neighbors and connection paths of each node one by one, and marking abnormal nodes that do not meet the requirements (such as satellite A having only 1 effective neighbor, satellite B requiring 5 hops connection with other satellites).

[0101] When the main filter detects an abnormal node, it immediately initiates a repair process, restoring topology connectivity through a closed-loop mechanism of "instruction issuance - satellite adjustment - result feedback." First, the main filter generates personalized "topology adjustment instructions" for the abnormal node: If the node has insufficient neighbors (e.g., satellite A has only 1 neighbor), the instruction will include a "recommended candidate neighbor list" (based on the global topology graph, filtering satellites with orbits matching satellite A and overlapping window periods, prioritizing candidates with previously better weighted function values) and adjustment priorities (e.g., "prioritize establishing links with satellites C and D, then consider satellite E"); if the node experiences connectivity failure due to poor link quality (e.g., the packet loss rate of the two neighbor links of satellite B is >1%), the instruction will prompt "prioritize repairing high-quality links (e.g., the link of satellite F, SNR=14dB), and eliminate low-quality links." Upon receiving the instruction, the abnormal node will pause its current low-priority data transmission and prioritize adjustment operations: for recommended candidate neighbors, it will quickly initiate link establishment requests (simplifying some screening steps and directly reusing orbital parameters provided by the main filter). If the recommended neighbor is unreachable (e.g., window period mismatch), it will autonomously expand the search range (increasing the search radius of potential neighbors from 10km to 15km), re-screen, and interact. For low-quality links that need to be removed, it will send a "link termination notification" to the corresponding neighbor, releasing resources for new link establishment. After the abnormal node completes the adjustment, it will report the "repair result" to the main filter (e.g., "Satellites C and D have been added as valid neighbors, current number of neighbors is 4"). Upon receiving the feedback, the main filter will re-check the node's connectivity. If it meets the preset requirements, it will mark it as "repair completed." If it still does not meet the requirements (e.g., the window period for adding a link is too short), it will issue a second round of adjustment instructions until all abnormal nodes have returned to normal connectivity, ensuring global topology stability to support subsequent orbit determination processes.

[0102] In some embodiments, when a node is detected in the local topology that does not meet the preset connectivity requirements, a topology adjustment command is actively issued to repair it, including: when a satellite with connectivity below a preset threshold is detected, a topology adjustment command is sent to force it to add neighbors.

[0103] In this embodiment of the invention, when the main filter discovers through a global coordination mechanism (such as topology information reporting and verification every 30 seconds) that the connectivity (current number of effective neighbors) of a satellite is lower than a preset threshold (usually set to 3 to ensure that the satellite can participate in at least one inter-satellite link closed loop construction and data interaction requirements), it will immediately generate and issue a customized topology adjustment instruction for that satellite. The instruction not only clearly indicates the "number of new neighbors to be added" (e.g., if the current connectivity is 1, at least 2 new neighbors are required), but also includes a "recommended candidate neighbor list" (based on global orbit data and weighted function results, prioritizing satellites with an effective window period overlapping with the satellite's by ≥10 minutes, a relative position change rate ≤1km / s, and a previously superior comprehensive performance J value, while also indicating the priority of each candidate neighbor), and an "operation time limit requirement" (e.g., initiating neighbor interaction within 1 minute to avoid prolonged low connectivity). After receiving the instruction, the satellite will pause its current low-priority inter-satellite data transmission task and prioritize initiating link establishment requests to candidate neighbors according to the recommended list in the instruction (simplifying some preliminary screening steps and directly reusing the orbit forecast and window period data provided by the main filter). To shorten interaction time, if there are unreachable neighbors among the recommended neighbors (e.g., no overlapping window), the search range of potential neighbors is automatically expanded (from the original 10km search radius to 15km). Satellites that meet the basic reachability requirements are re-selected and requests are initiated. At the same time, neighbors with poor link quality (e.g., packet loss rate > 1%, SNR < 10dB) in the original effective neighbor list are removed to free up resources. After the satellite successfully adds a neighbor and the connectivity reaches the preset threshold, it will send a "repair completion notification" to the main filter. The main filter will then perform a second check on its connectivity status. If the requirements are met, the node is marked as having a normal topology. If the requirements are still not met, a second round of adjustment instructions is issued until its connectivity meets the preset standard, ensuring that the satellite can participate normally in the constellation's global orbit determination and data collaboration.

[0104] In some embodiments, the closure error of the geometric closed loop formed by the inter-satellite links is calculated based on the state estimates of each satellite to construct a reference error observation. This includes: establishing a joint observation model coupling the reference error with time-varying clock errors, wherein the joint observation model extends the satellite clock errors from static parameters to a state vector that includes their time-varying characteristics, and incorporates the time-varying characteristics into the closed loop difference observation equation to quantify the impact of clock drift on the reference error estimation; assigning a dynamic weight to each geometric closed loop; the dynamic weight is adaptively determined based on the geometric configuration strength and ranging accuracy of the closed loop; verifying the observation residuals of each geometric closed loop difference and reducing the weight of large residual observations; and calculating the closure error of the geometric closed loop formed by the inter-satellite links based on the joint observation model, the dynamic weights of each geometric closed loop, and the state estimates of each satellite.

[0105] In this embodiment of the invention, to address the interference of satellite clock drift on the estimation of the reference error, a joint observation model coupling the two needs to be established first. Since the traditional static clock bias assumption cannot reflect the characteristics of LEO satellite clock bias changes with temperature and time (e.g., the small drift rate of rubidium clocks), the model expands the clock bias from a single static parameter to a state vector containing time-varying characteristics, specifically encompassing the clock bias value, clock bias drift rate, and drift acceleration, thus fully describing the dynamic changes in the clock bias. Subsequently, this time-varying characteristic is incorporated into the closed-loop error observation equation: when calculating the ranging values ​​of each link within the closed loop, not only is the influence of the reference error caused by the overall rotation and translation of the constellation considered, but also the ranging deviation caused by clock bias drift (e.g., the ranging increment per unit time due to the clock bias drift rate) is quantified. This clarifies the correlation between the closed-loop error and "reference error + time-varying clock bias," ensuring that the final closed-loop error observation reflects both the reference drift and the impact of dynamic clock bias changes, laying the foundation for subsequent accurate estimation.

[0106] Considering the varying reliability of different geometric closed loops, dynamic weights need to be assigned to each loop to ensure that high-quality closed loops play a dominant role in baseline error estimation. The determination of these weights depends on two core factors: First, the geometric strength of the closed loop. For example, a triangular closed loop has significantly higher geometric strength than a quadrilateral one due to its structural stability and strong anti-interference capability. Within the same type of loop, loops with uniform interior angle distribution (e.g., all interior angles of a triangle ≥ 30°) and small differences in side length have higher geometric strength. Second, the ranging accuracy of the closed loop, determined by the measured parameters of the inter-satellite links within the loop. Loops with high link signal-to-noise ratio (≥ 15dB), low latency (≤ 100ms), and low packet loss rate (≤ 0.1%) have higher ranging accuracy and smaller errors. The system adaptively calculates weights based on these two factors, giving higher weights to closed loops with superior geometric configurations and higher ranging accuracy, and vice versa, to prevent errors from inferior loops from interfering with baseline error estimation.

[0107] To further improve the reliability of observations, the observation residuals of each geometric closed loop need to be examined, and outlier data needs to be addressed accordingly. First, the preliminary observation residuals of each closed loop are calculated (i.e., the difference between the measured closed loop difference and the theoretical error-free closed loop difference). A preset threshold (e.g., residuals exceeding three times the standard deviation) is used to determine if a residual is large. These residuals are often caused by anomalies such as electromagnetic interference from inter-satellite links and packet loss in data transmission; directly incorporating them into the calculation would severely affect the accuracy of the baseline error estimation. For observations identified as having large residuals, their impact is reduced by multiplying them by a weighting factor less than 1 (e.g., 0.5-0.8). Finally, combining the previously constructed joint observation model, the dynamic weights of each closed loop, and the residual weighting results, the state estimates of each satellite (position, velocity, and time-varying clock error parameters) are substituted into the values. The ranging values ​​of each link within the closed loop are then algebraically summed in direction to obtain the final closure error. This closure error is the observation that accurately reflects the overall baseline error of the constellation.

[0108] In some embodiments, the state estimate includes satellite position vector, velocity vector, and time-varying clock error parameters. Based on the joint observation model, the dynamic weights of each geometric loop, and the state estimates of each satellite, the closure error of the geometric loop formed by the inter-satellite links is calculated, including: substituting the time-varying clock error parameters of each satellite into the time-varying clock error dynamic equation in the joint observation model to calculate the static clock error contribution term and the time-varying clock error drift integral contribution term between satellites within each loop, obtaining the total influence of clock error on link ranging; substituting the dynamic weights of each loop into the observation noise variance adjustment formula to perform weighted correction on the original ranging data of each link within the loop; determining the net ranging value of the link based on the weighted corrected link ranging data and the total influence of clock error; and determining the closure error of the geometric loop formed by the inter-satellite links based on the net ranging value of the link.

[0109] This invention also provides a hierarchical estimation-based LEO constellation GNSS / inter-satellite link joint autonomous orbit determination system, comprising: The data acquisition module is used to acquire GNSS data and inter-satellite link data collected by each satellite; The first calculation module is used to input GNSS data and inter-satellite link data into a local filter to obtain the state estimate and covariance matrix of each satellite. The second calculation module is used to input the state estimates and covariance matrices of each satellite into the main filter to obtain the overall baseline error estimate of the constellation. The feedback correction module is used to feed back the overall reference error estimate to each local filter to correct the state estimate of each satellite, and finally obtain the high-precision orbit parameters and clock error parameters of each satellite.

[0110] This invention also provides an autonomous orbit determination system, including an onboard computer, a GNSS receiver antenna, an inter-satellite link transceiver antenna, a processing unit, an ion thruster, an onboard clock, and a memory; the processing unit is used to execute the LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation as shown in the above embodiments.

[0111] In this embodiment of the invention, the onboard computer serves as the "control center" of the autonomous orbit determination system, responsible for coordinating the operational sequence of all components, receiving status information from each component, and issuing control commands. The onboard clock provides a unified, high-precision time reference for the system, ensuring time synchronization for GNSS observations, inter-satellite link communication, orbit calculations, and other stages. The memory stores various types of data, including raw GNSS observation data, inter-satellite link ranging data, state estimates and reference error results output by the processing unit, and orbit dynamics model parameters required for the orbit determination method. The three components are closely interconnected: the onboard computer establishes bidirectional communication with the onboard clock and memory via the onboard bus, reads the time signal from the onboard clock in real time to calibrate the timing of each stage, and writes the collected raw data, processing data, and final orbit determination results into the memory. It can also retrieve historical data from the memory for orbit prediction or algorithm optimization.

[0112] The core function of the GNSS receiver antenna is to receive observation signals such as pseudorange and pseudorange rate from GNSS navigation satellites, convert them into electrical signals, and transmit them to the processing unit. The inter-satellite link transceiver antenna has bidirectional communication capabilities, enabling it to receive relative distance data between neighboring satellites and send its own status information (such as position prediction and link requests) to neighboring satellites, also transmitting the received and transmitted electrical signals to the processing unit. The processing unit is the "execution core" of the orbit determination method, and it needs to execute the joint autonomous orbit determination method based on hierarchical estimation (including local filtering, main filtering, and reference error correction processes). Its connection relationship is as follows: the GNSS receiver antenna and the inter-satellite link transceiver antenna are unidirectionally connected to the processing unit through signal interfaces, inputting observation data to the processing unit; the processing unit is bidirectionally connected to the onboard computer, receiving start / adjustment commands (such as data acquisition frequency and filtering period) from the onboard computer on one hand, and feeding back information such as the status estimate of the local filter and the reference error calculation requirements of the main filter to the onboard computer on the other hand, so that the onboard computer can coordinate data upload or command broadcast.

[0113] The ion thruster, acting as an "orbit actuator," adjusts the satellite's orbit based on the orbit determination results. When the final orbit parameters output by the processing unit show that the satellite deviates from the preset orbit, or when the reference error fed back by the main filter needs to be compensated through orbit fine-tuning, the onboard computer generates attitude and orbit control commands based on the orbit determination data and sends them to the ion thruster to drive the satellite to correct its position or velocity. The connection is as follows: the ion thruster is unidirectionally connected to the onboard computer via a control interface and only receives control commands from the onboard computer. Simultaneously, the entire system forms a closed loop: GNSS receiver antenna and inter-satellite link transceiver antenna collect data → the processing unit executes the orbit determination method to generate orbit parameters and reference errors → the onboard computer sends adjustment commands to the ion thruster based on the results → the adjusted orbit state is fed back to the processing unit through the next round of data acquisition, realizing an autonomous orbit determination closed loop of "observation-processing-execution-feedback." The memory stores data from each stage throughout the process, providing support for closed-loop optimization and fault tracing.

[0114] In the above embodiments, the descriptions of each embodiment have their own emphasis. Parts not detailed or described in a particular embodiment can be referred to in the relevant descriptions of other embodiments. Unless otherwise specified or in conflict with logic, the terminology and / or descriptions between different embodiments are consistent and can be referenced interchangeably. Technical features in different embodiments can be combined to form new embodiments based on their inherent logical relationships.

[0115] The above-described embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit it. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention, and should all be included within the protection scope of the present invention.

Claims

1. A joint autonomous orbit determination method for LEO constellation GNSS / inter-satellite links based on hierarchical estimation, characterized in that, include: Acquire GNSS data and inter-satellite link data collected by each satellite; The GNSS data and the inter-satellite link data are input into a local filter to obtain the state estimate and covariance matrix of each satellite. The state estimates and covariance matrices of each satellite are input into the main filter to obtain the overall baseline error estimate of the constellation. The overall baseline error estimate is fed back to each local filter to correct the state estimate of each satellite, thus obtaining the final high-precision orbital parameters and clock error parameters of each satellite.

2. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 1, characterized in that, The GNSS data and the inter-satellite link data are input into a local filter to obtain the state estimates and covariance matrices of each satellite, including: Preprocessing of GNSS observation data and inter-satellite link observation data; The preprocessed inter-satellite link observation data is input into a local filter to obtain the state prediction value and its covariance prediction matrix. The preprocessed GNSS observation data is input into a local filter to correct the state prediction values ​​and their covariance prediction matrix, thereby obtaining the state estimates and covariance matrices of each satellite.

3. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 1, characterized in that, The state estimates and covariance matrices of each satellite are input into the main filter to obtain the overall baseline error estimate of the constellation, including: Based on the state estimates of each satellite, the closure error of the geometric closed loop formed by the inter-satellite links is calculated to construct a baseline error observation. Establish the mathematical relationship between the closure error of the geometric closed loop and the overall rotation / translation parameters of the constellation; The state estimates and covariance matrices of each satellite are used as inputs to the main filter. Based on the mathematical relationships and the baseline error observations, the overall baseline error of the constellation is estimated. The overall baseline error includes three orbital translation parameters, three rotation parameters, and one clock error baseline parameter.

4. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 1, characterized in that, The overall baseline error estimate is fed back to each local filter to correct the state estimate of each satellite, resulting in the final high-precision orbital parameters and clock error parameters for each satellite, including: The main filter broadcasts the overall baseline error estimate to the local filters of each satellite; Based on the overall baseline error estimate, the state estimates of the local filters of each satellite are corrected. The corrected state estimates are output as the final high-precision orbital parameters and clock error parameters for each satellite.

5. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 1, characterized in that, Before acquiring inter-satellite link data, the method further includes: Based on the orbit prediction model, the reachability of links with potential neighboring satellites within a future time window is determined, and an effective neighbor window period is defined. When the effective neighbor window period meets the preset update triggering conditions, a new round of neighbor filtering process is adaptively started; Perform neighbor optimization screening under multi-objective constraints, construct a weighted function with topology ring integrity, communication energy consumption and link quality as optimization objectives, and dynamically adjust the weights of each objective according to the real-time operating status of the constellation; Based on the calculation results of the weighting function, the optimal subset of neighbors is selected from the potential neighbors within the effective neighbor window period; Each satellite establishes a consensus on links by engaging in bidirectional information exchange with candidate neighbor satellites based on the optimal subset of neighbors. At the same time, the overall topology connectivity of the constellation is periodically monitored through a global coordination mechanism. When a node that does not meet the preset connectivity requirements is detected in a local topology, a topology adjustment command is proactively issued to repair it.

6. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 5, characterized in that, When a node is detected in the local topology that does not meet the preset connectivity requirements, a topology adjustment command is proactively issued to repair it, including: When a satellite with connectivity below a preset threshold is detected, a topology adjustment command is sent to force it to add neighbors.

7. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 3, characterized in that, Based on the state estimates of each satellite, the closure error of the geometric closed loop formed by the inter-satellite links is calculated to construct a baseline error observation, including: A joint observation model coupling the baseline error and the time-varying clock bias is established. The joint observation model extends the satellite clock bias from static parameters to a state vector that includes its time-varying characteristics, and incorporates the time-varying characteristics into the closed loop error observation equation to quantify the impact of clock bias drift on the baseline error estimation. A dynamic weight is assigned to each geometric closed loop; the dynamic weight is adaptively determined based on the geometric configuration strength and ranging accuracy of the closed loop. Each geometric closed loop difference observation residual is examined, and large residual observations are weighted down. Based on the joint observation model, the dynamic weights of each geometric closed loop, and the state estimates of each satellite, the closure error of the geometric closed loop formed by the inter-satellite links is calculated.

8. The LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation according to claim 7, characterized in that, The state estimates include satellite position vectors, velocity vectors, and time-varying clock error parameters; based on the joint observation model, the dynamic weights of each geometric loop, and the state estimates of each satellite, the closure error of the geometric loop formed by the inter-satellite links is calculated, including: Based on the time-varying clock error dynamic equation in the joint observation model, the time-varying clock error parameters of each satellite are substituted into the equation to calculate the static clock error contribution and the time-varying clock error drift integral contribution between each satellite in each closed loop, so as to obtain the total influence of clock error on link ranging. Substitute the dynamic weight of each closed loop into the observation noise variance adjustment formula to perform weighted correction on the original ranging data of each link within the closed loop; The net link ranging value is determined based on the weighted adjusted link ranging data and the total clock bias effect. Based on the net ranging value of the link, determine the closure error of the geometric closed loop formed by the inter-satellite links.

9. A joint autonomous orbit determination system for LEO constellation GNSS / inter-satellite links based on hierarchical estimation, characterized in that, include: The data acquisition module is used to acquire GNSS data and inter-satellite link data collected by each satellite; The first calculation module is used to input the GNSS data and the inter-satellite link data into a local filter to obtain the state estimate and covariance matrix of each satellite; The second calculation module is used to input the state estimates and covariance matrices of each satellite into the main filter to obtain the overall baseline error estimate of the constellation. The feedback correction module is used to feed back the overall reference error estimate to each local filter to correct the state estimate of each satellite, so as to obtain the final high-precision orbit parameters and clock error parameters of each satellite.

10. An autonomous orbit determination system, characterized in that, It includes an onboard computer, a GNSS receiver antenna, an inter-satellite link transceiver antenna, a processing unit, an ion thruster, an onboard clock, and a memory; the processing unit is used to execute the LEO constellation GNSS / inter-satellite link joint autonomous orbit determination method based on hierarchical estimation as described in any one of claims 1-8 above.

Images