Two-line element generation method adaptive to analytic propagation model and orbit prediction method

Abstract

The invention provides a two-line element generation method and orbit prediction method adaptive to an analytic propagation model, and the method comprises the steps: obtaining the historical osculating orbit elements at each epoch moment, and carrying out the orbit element optimization step of a plurality of iteration rounds on each historical osculating orbit element through employing the analytic propagation model, obtaining a target osculating orbital element number at each epoch moment, enabling the historical osculating orbital element number to be the same as the target osculating orbital element number, and performing inversion to obtain an average orbital element number; and extracting each average semi-major axis of each average orbital element, fitting each average semi-major axis by adopting a least square method to obtain a semi-major axis attenuation rate, inverting the semi-major axis attenuation rate to obtain an atmospheric resistance parameter, and equivalently converting the atmospheric resistance parameter into a standard resistance term in two rows of elements to obtain the two rows of elements. According to the method, the problems of inconsistent models, insensitive residual errors, unable physical explanation of parameters and the like in the two-line element generation process can be solved, and the ephemeris extrapolation orbit prediction precision is improved.

Description

Technical Field

[0001] This application relates to the field of satellite positioning, and in particular to a two-row root generation method and orbit prediction method adapted to an analytical propagation model. Background Technology

[0002] Two-Line Element Set (TLE) has been the internationally accepted format for orbital elements since the 1960s. Its concise two 80-character ASCII line structure is easy for computers to read and write, and can efficiently represent the orbital state of low Earth orbit (LEO) satellites. TLE includes six mean Keplerian elements and one atmospheric drag term (…). The data epochs are used to determine the orbital decay caused by aerodynamic drag and solar radiation pressure. In practical applications, the drag term is a freely tuned parameter that considers orbital decay caused by phenomena such as aerodynamic drag and solar radiation pressure. The SGP4 (Simplified General Perturbations 4) model is a closed perturbation propagation algorithm based on the average orbital elements input. TLE and SGP4 have become the de facto standards for global satellite cataloging and orbit prediction.

[0003] Currently, common methods for automatically generating TLE parameters mainly revolve around the inversion and optimization of instantaneous orbital elements (i.e., the actual observed state at a certain epoch) to average orbital elements. Then, combined with analytical orbital dynamics models such as SGP4, the average orbital element parameters are continuously adjusted to achieve standardized generation and efficient conversion of orbital parameters. However, existing methods rely on simplified analytical models such as SGP4, which have inherent limitations in terms of mechanical modeling accuracy and long-term propagation capability. Furthermore, many existing methods depend on traditional finite difference or manually derived Jacobian matrices, increasing algorithmic complexity and error risk, and limiting efficient adaptive optimization of high-dimensional parameter spaces and complex dynamic models, making it difficult to guarantee numerical convergence and stability. Simultaneously, traditional TLE cannot seamlessly integrate with J2 / J4, semi-analytical perturbation models, or fully perturbation analytical models. This leads to a series of practical obstacles in ephemeris generation, long-term prediction, and multi-model fusion, such as incompatible parameter formats, uncertain accuracy, and low data conversion efficiency. Moreover, the accuracy of current orbital prediction algorithms is insufficient, with errors exhibiting exponential divergence. Especially in commercial aerospace scenarios involving high-frequency telemetry and control, cataloging updates, and large-scale satellite management, the lack of a unified, high-precision input bridge between models severely limits the engineering application and promotion of these high-order dynamic methods. Summary of the Invention

[0004] In view of this, embodiments of this application provide a two-row root generation method and an orbit prediction method adapted to the analytical propagation model, so as to eliminate or improve one or more defects existing in the prior art.

[0005] The first aspect of this application provides a method for generating two-row roots that adapts to an analytical propagation model, the method comprising: The process involves obtaining the historical tangent orbital elements corresponding to each epoch, and then using a pre-defined analytical propagation model to perform multiple iterations of orbital element optimization steps on each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, ensuring that the historical tangent orbital elements corresponding to each epoch are the same as their respective target tangent orbital elements. The average orbital elements corresponding to each target tangent orbital element are then derived. The orbital element optimization steps include: in the analytical propagation model, constructing a cost function for each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, calculating the residual vector and Jacobian matrix for the current iteration, and updating the target tangent orbital elements for the current iteration using the LM method. The historical tangent orbital elements, the target tangent orbital elements, and the average orbital elements are all represented as Kepler six-element orbitals. Extract each average semi-major axis corresponding to each of the average orbital elements, fit each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, invert the semi-major axis decay rate to obtain the atmospheric drag parameter, and convert the atmospheric drag parameter into the standard drag term in the two rows of elements to obtain the two rows of elements.

[0006] In some embodiments of this application, before the step of performing multiple iterations of orbital element optimization on each of the historical coincidence orbital elements using a preset analytical propagation model, the method further includes: The satellite's velocity and position vectors are converted into Kepler orbital six roots, and the initial drag parameters of the satellite are calculated based on the Knudsen number and drag acceleration. The Kepler orbital six roots include the orbital semi-major axis, orbital eccentricity, orbital inclination, right ascension of the ascending node, argument of perigee, and mean perigee. The initial drag parameters are a combination of atmospheric drag coefficient, the equivalent surface-to-mass ratio of the drag experienced by the satellite, and the average atmospheric density at the satellite's location.

[0007] In some embodiments of this application, before the step of performing multiple iterations of orbital element optimization on each of the historical coincidence orbital elements using a preset analytical propagation model, the method further includes: Based on the Kepler orbital elements and the initial drag parameters, a residual function is constructed and the Jacobian matrix is ​​obtained using the difference quotient method; wherein, the Jacobian matrix is ​​represented as the partial derivative matrix of the residual function with respect to the orbital elements.

[0008] In some embodiments of this application, the calculation of the residual vector for the current iteration includes: The parameters corresponding to each of the said tangent orbital elements are perturbed, and the perturbed parameters are input into the orbital propagation model to obtain the residual vector.

[0009] In some embodiments of this application, updating the target tangent orbital elements in the current iteration using the LM method includes: The damping factor is added to form the parameter update formula and the parameter update increment is calculated. The increment is then added to the number of tangent orbital elements in the current iteration to obtain the input data for the next iteration. If the residual norm corresponding to the residual vector in the current iteration decreases, the damping factor is reduced. If the parameter update increment is less than a preset threshold, or the residual norm is less than a set threshold, or the number of iterations reaches a preset threshold, then the iteration stops.

[0010] In some embodiments of this application, the step of extracting each average semi-major axis corresponding to each of the average orbital elements, fitting each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, and inverting the semi-major axis decay rate to obtain atmospheric drag parameters includes: Extract the average semi-major axis corresponding to each of the average orbital elements to form a dataset, and construct a physical model of how the average semi-major axis changes over time. Based on the dataset, the least squares method is used to fit each of the average semi-major axes to obtain the semi-major axis attenuation rate. The atmospheric drag parameter is calculated by substituting the semi-major axis attenuation rate into the relationship between the semi-major axis attenuation rate and the average orbital elements.

[0011] A second aspect of this application provides an orbit prediction method, the method comprising: Based on the two rows of roots, higher-order band harmonics and fan harmonics are introduced, and a spatiotemporal dynamic atmospheric model is used for orbit prediction. The two rows of roots are obtained by the method described in the first aspect above.

[0012] A third aspect of this application provides an electronic device including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the method described in the first aspect above, and / or to implement the method described in the second aspect above.

[0013] A fourth aspect of this application provides a computer-readable storage medium having a computer program stored thereon that, when executed by a processor, implements the method described in the first aspect above, and / or implements the method described in the second aspect above.

[0014] The fifth aspect of this application provides a computer program product comprising a computer program that, when executed by a processor, implements the method described in the first aspect above, and / or implements the method described in the second aspect above.

[0015] This application provides a two-row root generation method adapted to an analytical propagation model. The method includes: obtaining the historical tangent orbital roots corresponding to each epoch; performing multiple iterations of orbital root optimization steps on each historical tangent orbital root using a preset analytical propagation model to obtain the target tangent orbital roots corresponding to each epoch, ensuring that the historical tangent orbital roots corresponding to each epoch are the same as their respective target tangent orbital roots; and reversing the average orbital roots corresponding to each target tangent orbital root. The orbital root optimization step includes: constructing, for each historical tangent orbital root in the analytical propagation model... The cost function is used to obtain the target tangent orbital elements corresponding to each epoch, and the residual vector and Jacobian matrix of the current iteration are calculated. The target tangent orbital elements of the current iteration are updated using the LM method. The historical tangent orbital elements, the target tangent orbital elements, and the average orbital elements are all represented as Kepler orbital six elements. The average semi-major axis corresponding to each average orbital element is extracted, and the average semi-major axis is fitted using the least squares method to obtain the semi-major axis decay rate. The semi-major axis decay rate is inverted to obtain the atmospheric drag parameter. The atmospheric drag parameter is equivalently converted into the standard drag term in the two rows of elements to obtain the two rows of elements. By designing an orbital parameter inversion and conversion mechanism for multi-model propagation environments, the Jacobi construction method based on the difference quotient method and the LM (Levenberg-Marquardt) nonlinear least squares joint inversion method can flexibly adapt to various orbital dynamics models such as SGP4, J2, and J4 without changing the propagator structure, significantly improving the adaptability, convergence efficiency, and system stability of orbital inversion. Furthermore, a multi-point semi-major axis time series fitting method is designed to accurately extract orbital energy decay characteristics and invert physically consistent atmospheric drag parameters, improving the accuracy of drag modeling. This achieves high-precision automatic conversion and standardized compatibility of parameters between different orbital dynamics models, improving the physical consistency of orbital ephemeris propagation in different models and the accuracy of ephemeris extrapolation orbit prediction. The proposed method possesses high-precision and highly versatile automatic orbital element conversion capabilities and physical parameter inversion mechanisms, applicable to various orbital service scenarios such as satellite cataloging management, long-term constellation evolution control, and terminal frequency shift pre-compensation, providing highly reliable orbital prediction support for low-Earth orbit communication, remote sensing, and navigation missions.

[0016] Additional advantages, objectives, and features of this application will be set forth in part in the description which follows, and will in part become apparent to those skilled in the art upon review of the following description, or may be learned by practice of the application. The objectives and other advantages of this application can be realized and obtained by means of the structures specifically pointed out in the specification and drawings.

[0017] Those skilled in the art will understand that the purposes and advantages that can be achieved with this application are not limited to those specifically described above, and that the above and other purposes that this application can achieve will be more clearly understood from the following detailed description. Attached Figure Description

[0018] The accompanying drawings, which are included to provide a further understanding of this application and form part of this application, do not constitute a limitation thereof. The components in the drawings are not drawn to scale but are merely for illustrating the principles of this application. For ease of illustration and description of certain parts of this application, corresponding portions in the drawings may be enlarged, i.e., may appear larger relative to other components in an exemplary device actually manufactured according to this application. In the drawings: Figure 1 This is a flowchart illustrating a two-row root generation method adapted to the analytic propagation model in one embodiment of this application.

[0019] Figure 2 This is a flowchart illustrating the LM method in one embodiment of this application.

[0020] Figure 3 This is a flowchart illustrating the method for generating atmospheric drag parameters in one embodiment of this application.

[0021] Figure 4 This is a flowchart illustrating the trajectory prediction method in an application example of this application. Detailed Implementation

[0022] To make the objectives, technical solutions, and advantages of this application clearer, the application will be further described in detail below with reference to the embodiments and accompanying drawings. Here, the illustrative embodiments and their descriptions are used to explain this application, but are not intended to limit it.

[0023] It should also be noted that, in order to avoid obscuring this application with unnecessary details, only the structures and / or processing steps closely related to the solution according to this application are shown in the accompanying drawings, while other details that are not closely related to this application are omitted.

[0024] It should be emphasized that the term "including / comprises" as used herein refers to the presence of a feature, element, step, or component, but does not exclude the presence or addition of one or more other features, elements, steps, or components.

[0025] It should also be noted that, unless otherwise specified, the term "connection" in this article can refer not only to a direct connection, but also to an indirect connection involving an intermediary.

[0026] In the following description, embodiments of the present application will be illustrated with reference to the accompanying drawings. In the drawings, the same reference numerals represent the same or similar parts, or the same or similar steps.

[0027] As satellite networks evolve towards lower latency and higher bandwidth, satellite communication, with its ability to provide wireless access anytime and anywhere, is poised to become an effective supplement to terrestrial systems, thus attracting widespread research and application. However, due to significant differences in channel characteristics and system configurations between low Earth orbit (LEO) satellite systems and terrestrial 5G New Radio (5G NR) systems, which utilize a new OFDM-based global air interface design, such as larger Doppler frequency offsets, propagation delays, and dynamic variations, terrestrial 5G NR systems cannot be directly applied to LEO satellite communication. Considering the periodic nature of LEO satellite orbits, if the satellite's position in its operational orbit can be predicted in advance, combined with the geographical location of the user equipment (UE), channel parameters such as Doppler frequency offset and propagation delay can be estimated in advance. This will greatly simplify the integration design of 5G and LEO satellite systems.

[0028] The average elements describe the orbit's inclination, right ascension of the ascending node, eccentricity, perigee angle (w), mean perigee angle, and mean velocity. In practice, the drag term is a freely tuned parameter that considers orbital decay caused by aerodynamic drag and solar radiation pressure. Current mainstream methods for automatically generating TLE parameters mainly revolve around the inversion and optimization of instantaneous orbital elements (i.e., the actual observed state at a given epoch) to average orbital elements. These methods use the Keplerian elements of the space target or the observed position and velocity as initial values, combined with analytical orbital dynamics models such as SGP4 / SDP4, and continuously adjust the average orbital element parameters through iterative optimization algorithms (such as Newton's method, finite difference method, least squares fitting, etc.) to ensure a high degree of consistency between the model's output state and the target state at a specified epoch. When the deviation between the model output and the observed state converges to a predetermined threshold, the resulting average orbital elements serve as the core input for TLE, thereby achieving standardized generation and efficient conversion of orbital parameters.

[0029] To address existing techniques for automatically generating mean Kepler elements based on orbital state inversion, a new research proposes an automated calculation method and system for converting orbital Kepler elements into two-row elements, applicable to the standardized generation of orbital parameters for space targets. This method first utilizes the SGP4 / SDP4 analytical model, iteratively optimizing it to ensure a high degree of consistency between the Kepler elements at a specified epoch and the orbital parameters in the two-row elements. Subsequently, based on multi-time-period extrapolated ephemeris, the method employs least-squares fitting to further correct the semi-major axis and drag parameters in the TLE (Tencent-Like-Leader) model, improving the physical consistency of the parameters and the accuracy of long-term forecasts. The system architecture includes an iterative calculation module, a fitting calculation module, and a correction calculation module, achieving efficient automatic conversion of orbital elements and generation of parameter sets.

[0030] A study proposes an automatic conversion method and system for satellite instantaneous elements to mean orbital elements, suitable for standardized expression and rapid prediction of satellite orbital states. This method first obtains the satellite's instantaneous Keplerian elements and calculates its position and velocity in the J2000 coordinate system based on these elements. Then, it converts to the True Equator Mean Equinox (TEME) coordinate system as the initial state for orbital conversion. Subsequently, iterative optimization is performed using an SGP4 or SDP4 analytical model, repeatedly converting, fitting, and correcting the satellite's position and velocity. Error parameters are assessed in real time to adjust the perturbation factor and step size until a convergent mean orbital element is output as the final result.

[0031] A study proposes a two-row element generation method based on model error compensation, applicable to the automatic generation and optimization of orbital parameters for space targets. This method first transforms the observed position and velocity vectors from the geocentric equatorial inertial coordinate system to the true equatorial mean year coordinate system, and further converts them into instantaneous orbital elements as initial values ​​for two-row element optimization. Subsequently, using the initial values ​​in conjunction with orbital propagation models such as SGP4, multi-period ephemeris simulations are performed on the two-row elements. By comparing the deviations between the predicted position and velocity and the actual observations, the two-row elements are iteratively corrected using the least squares method and regularization factors. The algorithm effectively improves the physical consistency of the final TLE parameters and the long-term orbit prediction accuracy by introducing a smoothing factor and an error compensation matrix.

[0032] A study has proposed a system and method for the automatic generation of extended satellite ephemeris data, applicable to long-term ephemeris acquisition by satellite terminal equipment. This method collects historical orbital parameter data, such as TLE (Time Leakage Time), constructs a future orbital parameter prediction model, and stores this model on the terminal device. Users can use this model to generate orbital parameters (such as predicted TLE) at any desired epoch, and then input this model into a standard orbital predictor to achieve efficient estimation of satellite position or overpass time.

[0033] A study has proposed a method to extend the two-planetary ephemeris prediction period for geostationary Earth Orbit (GEO) communication satellites, applicable to high-precision prediction of the orbital position of GEO satellites. This method first acquires the TLE data of the target satellite, extracting the average number of orbits and various Keplerian orbital parameters. Second, it corrects the average number of orbits to more closely approximate the actual geosynchronous motion, and accordingly adjusts key orbital parameters such as the semi-major axis and mean anomaly angle. Finally, it uses the corrected orbital parameters to calculate the satellite's three-dimensional spatial coordinates at any future time.

[0034] The SGP4 model is open-source and easy to implement, enabling rapid orbit prediction and meeting the situational awareness, collision warning, and service coverage needs of most commercial satellites. Furthermore, TLE+SGP4 has been widely integrated into mainstream aerospace simulation software and global space target databases, serving as the foundation for commercial satellites to conduct international cooperation, compliance registration, and external information dissemination. The SGP4 model integrates data on Earth's non-spherical gravity (…). Item), index atmospheric drag (through (Approximate fitting), supporting the direct calculation of orbital position and velocity from two rows of root data. In typical applications, SGP4 maintains positioning errors at epochs within the meter range and can sustain an error growth rate of 1–3 km per day, thus becoming a core component of systems such as aerospace monitoring, situation assessment, and space debris tracking. However, in addition to SGP4, researchers are constantly proposing new analytical or semi-analytical extrapolation methods.

[0035] To address existing analytical models and extrapolation methods for orbital dynamics, a new research proposes a high-precision semi-analytical method and system for extrapolating orbital elements based on a second-order high-precision J2 perturbation approach. This method is suitable for long-term, high-precision prediction of near-Earth satellite orbits. The method employs singularity-free semi-versus-factor elements, using the average orbital elements as input. Through Fourier series expansion and second-order perturbation theory, it analytically models the influence of the Earth's gravitational field J2 term and its higher-order perturbations on orbital evolution. The system achieves high-precision extrapolation of orbital elements and three-dimensional positions for satellites with arbitrary eccentricities and inclinations by jointly compensating for short-period terms, long-term perturbation terms, and higher-order average motions.

[0036] A study proposed the J2Propagator model, a closed-loop analytical propagation method based on the first-order Earth oblateness perturbation (J2 term). This model uses the average orbital elements as input to directly calculate the long-term drift effects of the right ascension of the ascending node (RAAN) and the argument of perigee (ω). Its processing excludes atmospheric drag, solar and lunar gravitational pull, or other higher-order perturbations, resulting in extremely fast propagation speeds and numerical stability, making it highly suitable for preliminary orbit design and dynamic analysis of sun-synchronous orbits.

[0037] Building upon this foundation, research has proposed the J4Propagator model. This model introduces a second-order contribution from J2 and a first-order perturbation term from J4, based on the first-order J2 model, and continuously uses closed-form expressions to rapidly predict long-term perturbations. By enhancing the modeling of higher-order non-spherical features of the Earth's gravitational field, J4Propagator can more accurately describe orbital elements—especially small variations in the semi-major axis and minor drifts in the nodal / perigee angles—compared to J2, exhibiting lower cumulative errors during medium- to long-term propagation (e.g., tens of days), making it suitable for orbital propulsion analysis and orbital maintenance planning.

[0038] A study proposes a fully perturbation-based analytical orbit propagation method based on the average orbital elements, applicable to satellite orbit design and long-term extrapolation with arbitrary eccentricity and inclination. This method takes the singularity-free average orbital elements as input and systematically constructs analytical series solutions encompassing nine types of perturbations, including Earth's non-spherical gravity, the gravitational pull of the Sun and Moon, Earth deformation, post-Newtonian effects, solar radiation pressure, and atmospheric drag, forming a complete perturbation model for orbital evolution. Through series expansion and averaging techniques, as well as joint compensation for short-period and long-period perturbations, the system accurately outputs orbital elements and three-dimensional position changes, achieving high-precision propagation over days or even tens of days, while eliminating singularity problems under conditions of small Kepler element eccentricity and high inclination.

[0039] Existing automatic generation techniques for mean Kepler elements rely on simplified analytical models such as SGP4, which have inherent limitations in terms of mechanical modeling accuracy and long-term propagation capability. SGP4's mechanical perturbation terms only simplify lower-order terms of the Earth's gravitational field (such as J2, J3, and J4) and atmospheric drag; complex perturbations such as higher-order perturbations, real atmospheric changes, and solar radiation pressure are often handled using empirical parameters. The simultaneous fitting of parameters leads to insufficient physical interpretation of the drag term. Furthermore, these methods typically use single-epoch or short-term windows for parameter fitting, failing to fully utilize multi-period, high-frequency observation data. This results in limited error control for medium- and long-term ephemeris predictions, especially during long-duration satellite missions or when orbital environments change drastically, easily leading to cumulative deviations in the track direction and orbital altitude. Simultaneously, existing methods largely rely on traditional finite difference methods or manually derived Jacobian matrices, increasing algorithmic complexity and error risk, and limiting efficient adaptive optimization of high-dimensional parameter spaces and complex dynamic models, making it difficult to guarantee numerical convergence and stability.

[0040] The fundamental limitation of existing analytical or semi-analytical orbital dynamics models lies in their strong dependence on initial input conditions such as mean orbital elements and atmospheric drag parameters. However, different models have significant differences in the definition and extraction methods of "mean orbital elements," and the TLE scheme itself is not a strictly physical mean orbital element; its drag term... Often adjusted based on experience, traditional TLE lacks physical interpretability that can be directly applied to higher-order models. Therefore, it struggles to seamlessly integrate with J2 / J4, semi-analytical perturbation models, or fully perturbation analytical models. This leads to a series of practical obstacles in ephemeris generation, long-term prediction, and multi-model fusion, including incompatible parameter formats, uncertain accuracy, and low data conversion efficiency. Especially in commercial aerospace high-frequency telemetry, tracking, and command (TT&C), cataloging updates, and large-scale satellite management scenarios, the lack of a unified, high-precision input bridge between models severely limits the engineering application and promotion of these higher-order dynamic methods. Furthermore, current orbit prediction algorithms such as SGP4 suffer from insufficient accuracy and exponential error divergence, necessitating a new orbit prediction algorithm that can compensate for SGP4's short-term accuracy and long-term error divergence while maintaining computational efficiency.

[0041] Therefore, to address the issues of model inconsistency, residual insensitivity, and uninterpretable parameters in the two-row root generation process, this application provides a two-row root generation method, orbit prediction method, electronic device, computer-readable storage medium, and computer program product adapted to the analytical propagation model, respectively, to improve the accuracy of ephemeris extrapolation orbit prediction and enhance prediction performance and generalization applicability in various scenarios such as constellation extrapolation, orbit maintenance, and frequency pre-compensation.

[0042] The following examples will provide a detailed description.

[0043] Based on this, embodiments of this application provide a two-row root generation method for an adapted analytic propagation model that can be executed by a two-row root generation device adapted to the analytic propagation model. See [link to relevant documentation]. Figure 3 The method specifically includes the following: Step 100: Obtain the historical tangent orbital elements corresponding to each epoch. Utilize a preset analytical propagation model to perform multiple iterations of orbital element optimization on each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, ensuring that the historical tangent orbital elements corresponding to each epoch are the same as their respective target tangent orbital elements. Then, derive the average orbital elements corresponding to each target tangent orbital element. The orbital element optimization step includes: in the analytical propagation model, constructing a cost function for each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, calculating the residual vector and Jacobian matrix for the current iteration, and updating the target tangent orbital elements for the current iteration using the LM method. The historical tangent orbital elements, the target tangent orbital elements, and the average orbital elements are all represented as Kepler six-element orbitals.

[0044] It should also be noted that the preset analytical propagation model can be SGP4, J2, or J4, etc., and the propagation model used in step 100 must be consistent with the orbit extrapolation model to ensure physical consistency and engineering coupling between the parameter inversion process and the extrapolation mechanism, and to avoid error propagation caused by model mismatch. The residual vector and Jacobian matrix will be described in subsequent embodiments. The LM method is a parameter estimation method for solving nonlinear least squares problems.

[0045] Specifically, the acquisition of spatial targets at multiple epochs The set of historical trajectories with the given tangent, where N represents the number of epochs. Represents the set of all moments, from Until Represented in set form as: Each of them For the first The historical tangent orbit element vectors at each moment can be obtained from ground-based radar measurements, on-board orbit estimations, simulation results, or other high-precision data sources.

[0046] For each of the above historical tangent orbital elements The difference quotient-driven LM method is called to update the target tangent orbital elements in the current iteration, and the corresponding average orbital element solutions are derived by inversion, as shown in the following formula: in, This represents the orbital state inversion operator. This inversion process ensures that: Among them, Model( ) represents any analytical or semi-analytic orbit extrapolation algorithm.

[0047] The cost function uses the average number of orbital elements as the optimization variable, minimizing the difference between the corresponding propagation trajectory and the target state. The formula for the cost function is as follows: in, The average orbital element vector, The cost function is the residual vector under undisturbed conditions. It is a nonlinear least squares optimization problem, and the objective is to minimize the cost function.

[0048] Step 200: Extract each average semi-major axis corresponding to each of the average orbital elements, fit each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, invert the semi-major axis decay rate to obtain the atmospheric drag parameter, and convert the atmospheric drag parameter into the standard drag term in the two rows of elements to obtain the two rows of elements.

[0049] Specifically, based on the convergent average orbital elements, further analysis is performed on the extrapolated ephemeris, and the tangent orbit is fitted at multiple points with the average orbital elements and key parameters such as the semi-major axis obtained from the extrapolation of various analytical models. Atmospheric drag parameters are obtained through joint inversion using the least squares method. and convert it to its equivalent form The standard resistance term is given by the following formula: in, This represents the average atmospheric density at the perigee of the satellite's orbit.

[0050] As described above, the two-row root generation method adapted to analytical propagation models provided in this application, by designing an orbital parameter inversion and conversion mechanism for multi-model propagation environments, and combining the Jacobi construction method based on the difference quotient method with the LM nonlinear least squares joint inversion method, can flexibly adapt to various orbital dynamics models such as SGP4, J2, and J4 without changing the propagator structure, significantly improving the adaptability, convergence efficiency, and system stability of orbital inversion. Furthermore, a multi-point semi-major axis time series fitting method is designed to accurately extract orbital energy decay characteristics and invert atmospheric drag parameters with physical consistency, improving the accuracy of drag modeling. This achieves high-precision automatic conversion and standardized compatibility of parameters between different orbital dynamics models, thereby improving the physical consistency of orbital ephemeris propagation in different models and the accuracy of ephemeris extrapolation orbit prediction. The method proposed in this application has high-precision and highly versatile automatic orbital root conversion capability and physical parameter inversion mechanism, applicable to various orbital service scenarios such as satellite cataloging management, long-term constellation evolution control, and terminal frequency shift pre-compensation, and can provide highly reliable orbital prediction support for low-Earth orbit communication, remote sensing, navigation, and other tasks.

[0051] To further improve the adaptability, convergence efficiency, and system stability of orbit inversion, in the two-row root generation method adapted to the analytical propagation model provided in this application embodiment, before step 100, the following content is specifically included: Step 010: Convert the satellite velocity vector and position vector into Kepler orbital six roots, and calculate the initial drag parameters of the satellite based on the Knudsen number and drag acceleration; wherein, the Kepler orbital six roots include the orbital semi-major axis, orbital eccentricity, orbital inclination, right ascension of the ascending node, argument of perigee, and mean perigee angle, and the initial drag parameters represent a combination of atmospheric drag coefficient, equivalent surface-to-mass ratio of drag borne by the satellite, and average atmospheric density at the satellite's location.

[0052] Specifically, based on the satellite's velocity and position vectors at each epoch, the Keplerian orbit six-root number of the satellite at that epoch can be determined. The formula for the Keplerian orbit six-root number is as follows: in, Indicates the semi-major axis of the track. Indicates the orbital eccentricity. Indicates the orbital inclination angle. Indicates the right ascension of the ascending node of the orbit. Indicates the argument of perigee. It represents the angle of approach.

[0053] remember and These are the satellite's position and velocity vectors in the J2000 geocentric inertial coordinate system, respectively. The size and shape of the elliptical orbit, as well as the satellite's position relative to its periapsis within the orbital plane, can be determined using the three roots a, e, and M, as shown in the following formula: in, Indicates the angle closer to the point. Represents the satellite's position vector. Represents the satellite's velocity vector. Indicates the satellite's instantaneous velocity. Let x represent the Earth's gravitational constant, and let x, y, and z represent the three components of the position vector on the three-dimensional coordinate axes. , , These represent the three components of the velocity vector along the three-dimensional coordinate axes. The three orientation roots for the orbital plane are also represented. , , It can be determined by three unit vectors in the perigee direction, the semi-path direction, and the angular momentum direction, denoted as […]. The formula is as follows: in, , Let P and Q represent the components of vectors P and Q along the z-axis, respectively. , , These represent the components of vector R along the three-dimensional coordinate axes x, y, and z, respectively. Therefore, for any given satellite measurement velocity and position, it can be losslessly converted into the Kepler orbital six-root number.

[0054] It should also be noted that the Knudsen number K is defined as Where λ represents the mean free path of gas molecules, and L is the characteristic scale of the spacecraft, the corresponding division is as follows: The formula for drag acceleration is as follows: Where CD represents the drag coefficient, and CD = 2.2 ± 0.2. S represents the equivalent cross-sectional area of ​​the satellite, and S / m represents the surface-to-mass ratio of the satellite. This indicates the atmospheric density at the location of the satellite. This represents the satellite's velocity vector relative to the atmosphere. and These are the velocity vectors of the satellite and the atmosphere relative to the Earth's center of mass, respectively. This represents the satellite's area. The drag coefficient is then... The formula is as follows: in, This represents the average atmospheric density at the perigee of the satellite's orbit. This is a combined parameter that includes the atmospheric drag coefficient, the equivalent surface-to-mass ratio of the satellite's drag resistance, and the average atmospheric density at the selected point. These three parameters are generally difficult to define precisely and cannot be separated during orbit determination; therefore, they appear as a combined parameter. In the absence of prior satellite information, this can be derived from the general case of CD=2.2. , Take initial value .

[0055] To further improve the adaptability, convergence efficiency, and system stability of orbit inversion, the two-row root generation method adapted to the analytical propagation model provided in this application embodiment includes the following specific content before step 100 and after step 010: Step 020: Construct a residual function based on the Kepler orbital six elements and the initial drag parameters, and obtain the Jacobian matrix using the difference quotient method; wherein, the Jacobian matrix is ​​represented as the partial derivative matrix of the residual function with respect to the orbital elements.

[0056] It should also be noted that the residual function is composed of the number of tangent orbital elements and the number of average orbital elements. The tangent orbital is reduced to a linear superposition of short-period terms of various orders of the average orbital, as shown in the following formula: in, This represents the average orbital elements. The following terms are short-period expansions of the power series of the perturbation force's minor parameters. The minor parameters are the ratio of the perturbation force's influence to the gravitational acceleration of the Earth's particles. The short-period first-order term representing the small parameters of the perturbation force often includes the mean anterior angle. Since it is a periodic function, it has a short-period term. This represents the short-period second-order term of the small parameter of the perturbation force.

[0057] orbital function It is a complex six-variable nonlinear function, with each dimension corresponding to one of the six Kepler orbital roots, as shown in the following formula: in, The average orbital elements of the initial epoch. This represents the linear change in the angle of approach due to the motion of two bodies without perturbation. , Let represent the long-term linear changes of the average orbit at each iteration, and t represent the target epoch time. Indicates the initial epoch time. This represents the various long-period variations of the average orbit. This indicates the number of orbital elements that tangent at the initial epoch. , Equations represent the short-period variations of the tangent orbit. This represents the orbital element of a long-period term in the target epoch. This represents the orbital element of the first long-period term at the initial epoch. A single data point is used to calculate the Time-Likelihood Expiration (TLE), as shown in the following formula: Where SGP4(·) represents the SGP4 algorithm, which uses the orbital six-root number of the epoch time contained in the common TLE and The output includes the target epoch, and the satellite's position and velocity vectors at that epoch. TLE indicates a common TLE. This indicates the epoch time (i.e., the same epoch as in TLE). Let be the position and velocity vectors at each epoch, respectively. Then, for any analytical or semi-analytical orbit extrapolation algorithm, the residual function formula is as follows: in, Indicates the target time The elements of the kissing track, The output is based on the Kepler orbital six-roots method, which converts the velocity and position vectors into Kepler orbital six-roots values ​​using step 010. Let the average orbital element vector be the vector to be optimized, and the formula is as follows: Specifically, the Jacobian matrix is ​​defined as the matrix of partial derivatives of the residual function with respect to the optimization parameters, as shown in the following formula: Each element Represents the residual function of the first The partial derivatives of each component with respect to the j-th optimization parameter of the orbital elements.

[0058] To further improve the adaptability, convergence efficiency, and system stability of orbit inversion, in a two-row root generation method adapted to an analytical propagation model provided in this application embodiment, step 100 of the two-row root generation method for the adapted analytical propagation model, specifically includes the following: The parameters corresponding to each of the said tangent orbital elements are perturbed, and the perturbed parameters are input into the orbital propagation model to obtain the residual vector.

[0059] Specifically, for each orbital element parameter (such as semi-major axis) eccentricity (etc.), apply a very small perturbation. That is, it can be adopted Modify the parameter value by a factor of 1 or other multiples to form the perturbed parameter vector, as shown in the following formula: The perturbated parameters Substituting the values ​​into the orbital propagation model for state propagation, we obtain the perturbed residual vector: The approximate value of the partial derivative of the j-th column can be obtained using the difference quotient formula: This yields the residual vector under the undisturbed state. .

[0060] To improve the scalability and engineering deployment efficiency of the algorithm, and to ensure stable operation of the system under various modeling accuracy requirements and task backgrounds, this application provides a method for generating two-row roots that adapts to the analytical propagation model. In step 100 of the method, see... Figure 2 The method of updating the target tangent orbital elements in the current iteration using the LM method further includes the following: Step 110: Add the damping factor to form the parameter update formula and calculate the parameter update increment. Add the increment to the number of tangent orbital elements in the current iteration round to obtain the input data for the next iteration round. If the residual norm corresponding to the residual vector in the current iteration round decreases, then reduce the damping factor. Step 120: If the parameter update increment is less than a preset threshold, or the residual norm is less than a set threshold, or the number of iterations reaches a preset threshold, then stop the iteration.

[0061] Specifically, the parameter update formula of the LM method consists of the Jacobian matrix, residual vector, and damping factor of the current iteration, as shown in the following formula: in, This represents the Jacobian matrix for the current iteration round. Let represent the residual vector, and λ represent the damping factor, where λ > 0.

[0062] The increment obtained from the above calculation This value is superimposed on the current orbital element estimate and used as the input for the next iteration. The specific formula is as follows: If the residual norm decreases significantly in the new round, the damping factor λ can be reduced to accelerate convergence; if the residual increases, λ should be increased.

[0063] After each iteration, determine whether one of the following criteria is met: Rule 1: Small parameter update magnitude: ; Criterion 2: The residual norm is sufficiently small: ; Criterion 3: Reach the maximum number of iterations: .

[0064] If the condition is met, the iteration ends and the current average orbital element is output. As the optimal solution.

[0065] To improve the accuracy of drag modeling and achieve high-precision automatic conversion and standardized compatibility of parameters between different orbital dynamics models, thereby enhancing the physical consistency of orbital ephemeris propagation across different models, this application provides a method for generating two-row roots for an adaptive analytical propagation model. In step 200 of this method, see... Figure 3 The step of extracting the average semi-major axis corresponding to each of the average orbital elements, fitting each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, and inverting the semi-major axis decay rate to obtain the atmospheric drag parameter, also specifically includes the following: Step 210: Extract the average semi-major axis corresponding to each of the average orbital elements to form a dataset, and construct a physical model of the average semi-major axis changing over time. Step 220: Based on the dataset, fit each of the average semi-major axes using the least squares method to obtain the semi-major axis attenuation rate; Step 230: Substitute the semi-major axis attenuation rate into the relationship between the semi-major axis attenuation rate and the average orbital elements, and calculate the atmospheric drag parameter.

[0066] It should also be noted that the physical model consists of the average semi-major axis and the semi-major axis decay rate.

[0067] Specifically, after reversing the average orbital elements, the average semi-major axis component at each time step is extracted from them. The dataset is formed as follows: in, , ... These represent different times. Due to the orbital energy decay trend caused by atmospheric drag, the average semi-major axis is constructed. A physical model that varies with time t, with a linear decay model for the mean semi-major axis: Based on dataset The initial average semi-major axis of the parameters is solved using the least squares fitting method. With semi-major axis attenuation rate This minimizes the sum of squared residuals between the predicted and observed orbits: The solution is given directly from the closed form. Let: Fitting parameter vector The solution is: Substituting the known relationship between the semi-major axis attenuation rate γ and the average orbital elements from the perturbation theory into the following equation: in, , denoted as: in, This represents a combination function of various perturbation coefficient terms. , It is the flattening of the Earth's geometry. The distance from the geocentric point of perigee. This indicates the atmospheric density elevation at the perigee of the satellite's orbit. , , These are Bessel functions of the first kind, representing 0th, 1st, and 2nd order dummy variables, respectively, where n is the average angular velocity of the satellite. The diurnal variation factor of atmospheric density. The gravitational constant of Earth, The average orbital elements output from step 100 are determined. Atmospheric drag parameters are then derived. .

[0068] To further illustrate the above embodiments, this application also provides a specific application example of a two-row root generation method adapted to an analytical propagation model. This method relates to technical fields such as space orbital dynamics, aerospace telemetry and control, and satellite orbital data processing. It possesses high-precision, highly versatile automatic orbital root conversion capabilities and a physical parameter inversion mechanism, making it suitable for various orbital service scenarios such as satellite cataloging management, long-term constellation evolution control, and terminal frequency shift pre-compensation. It can provide highly reliable orbit prediction support for low-Earth orbit communication, remote sensing, navigation, and other missions.

[0069] To achieve high-precision automatic conversion and standardized compatibility of parameters between different orbital dynamics models, and to improve the physical consistency and accuracy of orbital ephemeris propagation across different models, the following key issues need to be addressed. Currently, TLE+SGP4, as the mainstream standard for global commercial spaceflight and space monitoring, has fundamentally different definitions of mean orbital elements and drag parameters compared to models such as J2 / J4 and fully perturbation analytical models. This difference hinders the interoperability of ephemeris data, parameter formats, and physical meanings in multi-model collaboration and large-scale engineering automation, thus restricting the development of applications such as multi-source orbital data fusion, long-term mission forecasting, and high-precision situational awareness. (1) Design of a multi-point inversion mechanism based on difference quotient and nonlinear least squares optimization: In view of the limitations of existing inversion methods that rely on manual derivation or fixed finite differences and are difficult to extend to multidimensional parameter space and complex dynamic models, it is urgent to develop a multi-point / multi-model joint inversion mechanism based on automatic approximation of difference quotient and nonlinear least squares optimization, so as to achieve efficient convergence and numerical stability of the orbital parameter inversion process, and support cross-model, high-precision, and engineering-oriented automated ephemeris generation and sharing.

[0070] (2) Based on physical consistency and model convertibility Parameter inversion mechanism design: How to directly invert drag parameters with clear physical meaning based on the principle of physical consistency. Effectively avoids the problems in the traditional TLE / SGP4 system The uncertainty brought about by parameter empirical fitting provides a reliable basis for the automatic conversion of parameters between different models, ensuring the accuracy of ephemeris extrapolation and the standardized sharing of parameters.

[0071] The two-row root generation method for the adapted analytic propagation model provided in this application example specifically includes the following: Step 1: Based on the satellite's instantaneous orbital elements, mass, cross-sectional area, and drag coefficient, the satellite's spatial state at a given epoch is initially calculated using a standard Keplerian dynamics model and a selected analytical propagation algorithm (such as SGP4, J2 / J4, etc.). This yields the target state and reference orbit information for subsequent parameter optimization.

[0072] Step 2: For the spatial state of the satellite at a given epoch, a nonlinear least squares residual function is constructed. The average orbital elements after model propagation of the tangent orbital elements are compared with the target state. The difference quotient method is used to automatically numerically approximate the Jacobian matrix for the orbital analytical propagation algorithm, efficiently calculating the sensitivity of all Kepler orbital six elements without requiring analytical differential derivation. The propagation model used in this step is consistent with the final orbital extrapolation model, ensuring physical consistency and engineering coupling between the parameter inversion process and the extrapolation mechanism, and avoiding error propagation caused by model mismatch.

[0073] Step 3: Based on nonlinear optimization algorithms such as Levenberg-Marquardt, the average orbital elements are adjusted iteratively through multiple rounds. In each iteration, the Jacobian matrix is ​​updated using the difference quotient method to ensure flexible adaptation across multiple models (SGP4 / J2 / J4, etc.). The error between the model output and the target state is continuously reduced until it converges to a set threshold, outputting an average orbital element with strong physical consistency.

[0074] Step 4: Based on the convergent average orbital elements, further analysis is performed on the extrapolated ephemeris. The intersecting orbits are then fitted at multiple points with the intersecting orbital elements, semi-major axis, phase, and other key parameters obtained from the extrapolation of various analytical models. Atmospheric drag parameters are then retrieved using the least squares method. And convert it to an equivalent form as needed. Standard resistance terms are used to ensure that the results are compatible with the TLE and SGP4 standard systems.

[0075] The specific implementation process for each step is given below: Step 1: (1-1) Obtain the six roots of the satellite’s Kepler orbit based on the satellite’s velocity vector and position vector.

[0076] Based on the satellite's velocity and position vectors at each epoch, the six roots of the satellite's Keplerian orbit at that epoch can be determined. (Note: The original text contains a typo and can be left as is.) and Let be the position vector and velocity vector of the satellite in the J2000 geocentric inertial coordinate system, respectively. in, It is the semi-major axis of the track. It is the orbital eccentricity. It is the track inclination angle. It is the right ascension of the ascending node of the orbit. It is the perigee argument of the orbit. It is the angle of approach. Consider the following system of nonlinear equations. The three roots a, e, and M determine the size and shape of the elliptical orbit and the satellite's position relative to its periapsis within the orbital plane. These are quantities within the orbital plane and can be calculated using the following formula: in: in, It's closer to the point angle. It is the satellite's position vector. It is the satellite's velocity vector. It is the satellite's instantaneous velocity. It is the Earth's gravitational constant, and x, y, and z represent the three components of the position vector on the three-dimensional coordinate axes, respectively. , , These represent the three components of the velocity vector along the three-dimensional coordinate axes. The three orientation roots for the orbital plane are also represented. , , It can be determined by three unit vectors in the perigee direction, the semi-path direction, and the angular momentum direction, denoted as […]. It can be calculated using the following formula: in, , Let P and Q represent the components of vectors P and Q along the z-axis, respectively. , , These represent the components of vector R along the three-dimensional coordinate axes x, y, and z, respectively. Thus, for any given satellite measurement velocity and position, the application examples of this application can losslessly convert them into Kepler orbital six-root numbers and save them as initialization parameters.

[0077] (1-2) Calculate the initial atmospheric drag parameters of the satellite.

[0078] The force conditions experienced by artificial Earth satellites during atmospheric flight are extremely complex, making it a crucial research topic in hypersonic gas dynamics. The satellite force problem involves continuous medium flow in the lower atmosphere, free molecular flow in the upper atmosphere, and transitional flows in between, as well as complex factors such as neutral gases, ionized gases, and multi-component mixtures. Regarding high-speed, high-altitude gas dynamics, the basic assumption for a neutral atmospheric state in engineering applications generally adopts the classification based on the Knudsen number K, proposed by Qian Xuesen in 1946. K is defined as... Here, λ is the mean free path of gas molecules, and L is the characteristic scale of the spacecraft. The corresponding division is as follows: Artificial Earth satellites typically operate at altitudes above 200 km, and their characteristic scales are not very large. Therefore, the atmospheric state in this airspace is characterized by free molecular flow, and the damping effect manifests as atmospheric drag. The corresponding drag acceleration is expressed as follows: Where CD is the drag coefficient, CD = 2.2 ± 0.2, and S / m is the satellite's surface-to-mass ratio (S is the equivalent cross-sectional area in terms of drag). It is the atmospheric density at the location of the satellite. . is the satellite's velocity vector relative to the atmosphere. and These are the velocity vectors of the satellite and the atmosphere relative to the Earth's center of mass, respectively. It is the equivalent cross-sectional area of ​​the satellite. This refers to the satellite's area. Therefore, the drag coefficient... It can be represented as: in, This represents the average atmospheric density at the perigee of the satellite's orbit. In reality, it is a combination of parameters, including the atmospheric drag coefficient, the equivalent surface-to-mass ratio of the satellite's drag resistance, and the average atmospheric density at the selected point. These three parameters are generally difficult to define precisely and cannot be separated during orbit determination; therefore, they appear as a combination of parameters. In the absence of prior satellite information, the application example of this application can be derived from the general case of CD=2.2. , Take initial value .

[0079] Step 2: (2-1) Construction of residual function based on average orbital elements and tangent orbital elements.

[0080] In the field of existing orbital dynamics analysis and space target orbit prediction, mainstream analytical propagation models such as SGP4, DSST, and Vinti all use mean orbital elements as input parameters. Mean orbital elements refer to a set of elements characterizing the average orbital state obtained by averaging various periodic perturbation terms (such as the Earth's non-spherical gravitational field, atmospheric drag, and three-body perturbations) over a certain period under perturbation conditions. Compared with osculating orbital elements (instantaneous orbital parameters at real time), mean orbital elements eliminate the influence of high-frequency, short-period perturbations, making them more suitable for long-term propagation and prediction of analytical models.

[0081] in, The first term represents the average orbital elements. The subsequent terms are short-period expansions of the power series of the perturbation force's minor parameters, where the minor parameters represent the ratio of the perturbation force's influence to the gravitational acceleration of the Earth's particles. The short-period first-order term representing the small parameters of the perturbation force often includes the mean anterior angle. Since it is a periodic function, it has a short-period term. This represents the short-period second-order term of the perturbation force's minor parameter. As described in the above equation, the tangent orbit is reduced to a linear superposition of short-period terms of various orders of the average orbit. The calculation of each short-period term also depends on the average orbit as prior input. Orbit function It is a complex six-variable nonlinear function, with each dimension corresponding to one of the six Kepler orbital roots, which can be reduced to: As shown in the above formula, where, The average orbital elements of the initial epoch. This represents the linear change in the angle of approach due to the motion of two bodies without perturbation. , Let represent the long-term linear changes of the average orbit at each iteration, and t represent the target epoch time. Indicates the initial epoch time. This represents the various long-period variations of the average orbit. This indicates the number of orbital elements that tangent at the initial epoch. , Equations represent the short-period variations of the tangent orbit. This represents the orbital element of a long-period term in the target epoch. This represents the number of long-period orbital elements at the initial epoch. Calculate the initial average orbital... The main terms in TLE can be obtained by subtracting the initial short-period terms of each order from the initial snip orbit. However, the calculation of the short-period terms depends on the average orbit input as prior information, and the average orbit cannot be obtained directly.

[0082] Furthermore, the procedures used by the organizations that release TLEs are not publicly available, making the creation of SGP4-compatible TLEs a challenge. Researchers have made numerous attempts to bypass standard procedures for obtaining TLEs to utilize independently collected measurement data. Fitting methods aim to directly reverse engineer the actual fitting process using known SGP4 algorithms. While these fitting methods can provide accurate results, research indicates that they require collecting multiple days of measurement data to achieve a good fit. In some cases, such as low-power CubeSats, data collection over such long timeframes may be impossible due to power budget constraints. Therefore, an alternative approach that can take a single data point and use it to compute TLEs is valuable. The single-point method is designed to convert a single state measurement into an SGP4-compatible TLE. As shown in the equation above, this is accomplished by separating the perturbations added to the orbit by SGP4 and removing them during the TLE generation process. The underlying process is abstracted into the following equation.

[0083] Here, SGP4(·) refers to the SGP4 algorithm, which uses the orbital six-root number of the epoch time contained in the common TLE and The output includes the target epoch, and the satellite's position and velocity vectors at that epoch. TLE refers to the common TLE. It refers to the epoch time (that is, the same epoch as in TLE). Let be the position and velocity vectors at each epoch. Then, for any analytical or semi-analytical orbit extrapolation algorithm Model( The residual function can be expressed as: in, Refers to the target time The kissing track element, Model( The output is also in Kepler orbital six-root form. In this application example, the method in step 1 can also be used to convert the velocity and position vectors into Kepler orbital six-root forms. The mean orbital element vector to be optimized includes: (2-2) Calculation of Jacobian Matrix Based on Difference Quotient Method To achieve nonlinear least squares optimization, it is necessary to calculate the Jacobian matrix of the residual function with respect to the optimization variable (orbital elements) to characterize the sensitivity of the objective function to each parameter. This application innovatively employs the difference quotient method to automate and efficiently calculate the Jacobian matrix, avoiding complex analytical derivations.

[0084] The Jacobian matrix is ​​defined as the matrix of partial derivatives of the residual function with respect to the optimization parameters: Each element Represents the residual function of the first The partial derivatives of each component with respect to the j-th optimization parameter of the orbital elements.

[0085] For each orbital element parameter (such as semi-major axis) eccentricity (etc.), apply a very small perturbation. (For example (Multiply the parameter value by 1), forming the perturbed parameter vector: The perturbated parameters Substituting the values ​​into the orbital propagation model for state propagation, we obtain the perturbed residual vector: The approximate value of the partial derivative of the j-th column can be obtained using the difference quotient formula: in: This represents the residual vector under undisturbed conditions. The above process is performed one by one for each orbital parameter to be optimized, thereby quickly and efficiently obtaining all columns of the Jacobian matrix.

[0086] Step 3: (3-1) Construct the cost function.

[0087] This problem belongs to nonlinear least squares optimization, and the objective is to minimize the following objective function: This function uses the average number of orbital elements. To optimize variables, the cost function minimizes the difference between their corresponding propagation trajectories and the target state. This cost function has a clear geometric meaning: it seeks a set of trajectory elements in the state space that generate a trajectory curve that is "closest" to the target trajectory at all sampling times. The function is general in form, adaptable to SGP4, and also supports higher-order J2 / J4 analytical propagators or semi-analytical models, providing a unified and clear solution foundation for subsequent nonlinear optimization algorithms such as Levenberg-Marquardt. Through continuous iterative descent of the cost function, the optimal average trajectory elements that are physically consistent and satisfy trajectory dynamics constraints can be gradually approximated.

[0088] (3-2) Iteration and convergence mechanism for average orbital element optimization based on Levenberg-Marquardt After constructing the cost function and calculating the current residual and Jacobian matrix, the Levenberg-Marquardt (LM) algorithm, based on nonlinear least squares, iteratively approximates the optimal average orbital root solution. This method completes parameter updates and convergence checks for each iteration through the following three sub-steps: LM parameter update formula: The update direction of the LM algorithm is given by the following analytical expression: in, Let Jacobian matrix be the value at the current parameter point. Let λ be the residual vector, and λ>0 be the damping factor. This parameter adjustment and update strategy flexibly switches between the Gauss-Newton method (λ→0) and the gradient descent method (λ→∞), taking into account both convergence speed and stability.

[0089] In the application examples of this application, It is not obtained through analytical differentiation, but rather through numerical approximation using the difference quotient method in step 2. That is, a small perturbation is introduced into each orbital element dimension, and the changes in the state propagation results are observed to estimate the approximate partial derivatives of the target state with respect to the orbital parameters, thus constructing a complete Jacobian matrix.

[0090] This method of numerically approximating the Jacobian matrix based on the difference quotient method has significant advantages. Firstly, it avoids the tedious process of analytical differentiation for complex orbit propagation models, greatly simplifying the model adaptation process. Secondly, this method is highly versatile and flexible across propagation models; whether it's the SGP4 or J2 analytical propagator, or more complex J4 or multi-perturbation models, the corresponding Jacobian matrix can be directly obtained through numerical perturbation. This decoupling from the propagation model not only improves the algorithm's scalability and engineering deployment efficiency but also ensures stable operation of the system under various modeling accuracy requirements and task backgrounds, exhibiting good modularity and cross-platform compatibility.

[0091] Parameter update and damping adjustment mechanism: The increment obtained from the above calculation This value is superimposed on the current orbital element estimate and used as the input for the next iteration. If the residual norm decreases significantly in the new round, it indicates that the current search direction is reasonable, and the damping factor λ can be appropriately reduced to accelerate convergence; conversely, if the residual increases, λ should be increased to approach the gradient descent path and improve stability.

[0092] Convergence condition determination and iteration termination: After each iteration, the system checks whether it meets any of the following convergence criteria: Rule 1: Small parameter update magnitude: Criterion 2: The residual norm is sufficiently small: Criterion 3: Reach the maximum number of iterations: Once the conditions are met, the optimization process terminates and the current optimal solution is output. This method serves as the input of the average orbital elements for the analytical orbital extrapolation model. Leveraging the difference quotient method for rapid Jacobian calculation and the powerful convergence performance of the LM algorithm, it can flexibly adapt to multiple propagation models (such as SGP4, J2, J4, etc.) to achieve efficient and physically consistent TLE inversion generation.

[0093] Step 4: In step 3, this application example uses any given model (such as SGP4, J2, or J4) as the dynamics propagator. By constructing a residual function and employing the Levenberg-Marquardt nonlinear least squares algorithm, a set of average orbital elements that ensure the matching orbit is completely consistent with the observed state at the target epoch is generated iteratively. This method, combined with the difference quotient method to calculate the Jacobian matrix, has advantages such as no analytical differentiation required, adaptability to various propagation models, and robustness to singularities, and can effectively generate orbital parameters with high physical consistency.

[0094] However, it should be noted that this step can only ensure that the tangent trajectory is precisely aligned with the target state at a single epoch, and cannot directly generate atmospheric drag parameters (such as...). or This is because the main effect of atmospheric drag is manifested in the long-term dissipation of orbital energy, with its influence on the orbit primarily manifested in the slow decay of the semi-major axis, while its instantaneous impact on short-period terms (such as the mean anterior angle and the right ascension of the ascending node) is relatively small. Therefore, within the local fitting framework of step 3, the drag term is not sensitive to the contribution of the residuals, making it difficult to achieve effective identification.

[0095] Therefore, to address this deficiency, step 4 requires further analysis of the orbital energy decay trend and inversion of drag-related parameters through multi-point extrapolation and least-squares fitting, ensuring compatibility with the SGP4 model. This two-stage strategy enables the system to achieve both local fit and long-term physical parameter estimation, providing a solid foundation for high-precision TLE generation.

[0096] (4-1) Input the number of tangent orbital elements in multiple epochs and invert the average orbital elements point by point.

[0097] This step first obtains information about the space target at multiple epochs. The set of tangent orbital elements, where N represents the number of epochs. Represents the set of all moments, from Until Represented in set form as: Each of them For the first The six-dimensional tangent orbital element vectors at each moment, these input data can come from ground-based radar measurements, on-board orbit estimation, simulation results or other high-precision data sources.

[0098] For each of the above-mentioned number of kissing track elements The Levenberg-Marquardt least squares iterative algorithm driven by the difference quotient proposed in the application example of this application is invoked to derive the corresponding average orbital root solutions, denoted as: in This represents the orbital state inversion operator. This inversion process ensures that: After inverting the average orbital elements, the average semi-major axis component at each time step is extracted. For subsequent modeling and analysis: in, , ... These represent different times.

[0099] (4-2) Least square fitting inversion of decay rate parameters.

[0100] Considering the tendency of orbital energy decay due to atmospheric drag, construct the average semi-major axis. Physical models that change over time generally assume that the mean semi-major axis follows a linear decay model: Based on observation dataset The initial average semi-major axis of the parameters is solved using the least squares fitting method. With semi-major axis attenuation rate This minimizes the sum of squared residuals between the predicted and observed orbits: This problem can be directly solved using a closed-form solution. Let: Fitting parameter vector The solution is: (4-2) Atmospheric drag parameters are inverted based on theoretical models.

[0101] Substituting the known relationship between the semi-major axis decay rate γ and the Kepler-mean orbital elements at the current point from the perturbation theory into the following equation: in: Recorded as: in, This represents a combination function of various perturbation coefficient terms. , It is the flattening of the Earth's geometry. The distance from the geocentric point of perigee. This indicates the atmospheric density elevation at the perigee of the satellite's orbit. , , These are Bessel functions of the first kind, representing 0th, 1st, and 2nd order dummy variables, respectively, where n is the average angular velocity of the satellite. The diurnal variation factor of atmospheric density. The gravitational constant of Earth, The average orbital elements output from step 3 are used to determine the atmospheric drag parameters. These parameters are then estimated by substituting the values ​​back into the equation. Atmospheric drag parameters It can also be converted to TLE. Parameters, drag parameters in the SGP4 model It is a pseudo-drag parameter, essentially the result of adjusting the satellite's ballistic coefficient for atmospheric density. This application example demonstrates... The estimate yields: In traditional Time-Like Exception (TLE) fitting, B is usually used as an empirical adjustment parameter, obtained through repeated fitting of residuals. Its physical meaning is unclear, and it often results in negative or unreasonable values. In contrast, the application example in this application fits the semi-major axis attenuation trend at multiple points, yielding a result with good physical consistency. Then, by converting it to B using the above formula, not only is its compatibility with the SGP4 model preserved, but the physical interpretability of the parameters and the accuracy of long-term forecasts are also significantly improved. This method achieves the unification of drag terms between the SGP4 / TLE system and higher-order dynamic models, laying the foundation for the interoperability and sharing of orbital data among multiple models.

[0102] Based on the embodiments and / or application examples of the two-row root generation method for the above-mentioned adaptive analytical propagation model, in order to improve its prediction performance and generalization applicability in multiple scenarios such as constellation extrapolation, orbit maintenance, and frequency pre-compensation, and to improve the accuracy of ephemeris extrapolation orbit prediction, this application also provides an embodiment of an orbit prediction method, which specifically includes the following: Based on the two rows of roots, higher-order band harmonics and fan harmonics are introduced, and a spatiotemporal dynamic atmospheric model is used for orbit prediction. The two rows of roots are obtained by the two-row root generation method of the above-mentioned adapted analytical propagation model.

[0103] It should also be noted that the higher-order zonal harmonics are composed of the satellite's distance from the Earth's center, gravity coefficient, latitude, orbital inclination, latitudinal argument, and the Earth's mean equatorial radius, while the fan harmonics are composed of gravity coefficient, latitude, and relative geographical longitude. Specifically, the formula for the harmonic terms of Earth's non-uniform gravitational field is as follows: in, It is the distance of the satellite from the center of the Earth. (n=2, 3, ...) represents the zonal gravity coefficients. For latitude, sin =sin sinµ and Let µ be the orbital inclination angle. +f represents the latitude argument. It is the perigee argument of the orbit. It is the true anterior angle, where re is the Earth's average equatorial radius. (·) denotes a Legendre polynomial of the first kind with an order of n.

[0104] The formula for the J22 sector harmonic of Earth's non-uniform gravitational field is as follows: in, Indicates the main harmonic term and ,and and For harmonic coefficients, This represents a second-order quadratic associated Legendre polynomial. It is the geographical latitude, and , Indicates relative geographical longitude and = , Indicates geographical longitude. For harmonic coefficients and The geographical longitude of the determined direction of the Earth's equatorial "axis of symmetry" , and It is the harmonic coefficient.

[0105] To accurately model the decay of orbital energy caused by non-conservative forces, the influence of atmospheric drag must be considered in low Earth orbit satellite orbits. A spatiotemporal dynamic atmospheric model is adopted, and the formula is as follows: Where S is the equivalent cross-sectional area, m is the satellite mass, and CD is the drag coefficient. Where is atmospheric density, and v is the satellite's velocity vector. It is the atmospheric rotational velocity vector. It is the daily average density of the reference ellipsoid passing through the orbital perigee. The diurnal variation factor of atmospheric density. Let r be the angle between the satellite's geocentric radius and the direction of the diurnal density peak (corresponding to the maximum atmospheric density), and let r be the distance from the satellite to the Earth's center. Let be the distance from the point where the satellite's radial vector intersects the reference ellipsoid to the Earth's center. This is the atmospheric density elevation.

[0106] To further illustrate the above embodiments, this application also provides a specific application example of the trajectory prediction method. Specifically, it includes the following: After generating a high-precision Time Limit Exceeded (TLE) suitable for the analytical propagation model using the two-row root generation method described above, this application example proposes a high-precision orbit prediction algorithm based on this TLE as the initial condition. This algorithm introduces higher-order band harmonic terms (such as J2, J3, J4, J22, etc.) and a spatiotemporal dynamic atmospheric model during the analytical propagation process, effectively weakening the exponential divergence trend of errors in long-term extrapolation, making the error grow slowly and approximately linearly. Simultaneously, by introducing the short-period principal term J22, the rapid oscillation term of the satellite orbit is compensated more accurately, significantly improving the orbital position accuracy in the short term. Compared to the traditional SGP4 model, this application example can significantly extend the high-precision prediction time while maintaining computational efficiency, and can further improve the overall orbit extrapolation stability and reliability through rolling updates of the TLE. The specific implementation process is as follows: This application proposes a method to improve the short- and long-term accuracy of orbit prediction while minimizing perturbation. The two-row root generation method adapted to the analytical propagation model described above can be used to transform the TLE from the common TLE into a TLE suitable for the prediction method of this application, and then orbit prediction can be performed.

[0107] First, considering the perturbation effect of band harmonics, the band harmonic term of the Earth's non-uniform gravitational field can be expressed as: in, It is the distance of the satellite from the center of the Earth. (n=2, 3, ...) represents the zonal gravity coefficients. For latitude, sin =sin sinµ, Let u be the orbital inclination angle. +f represents the latitude argument. It is the perigee argument of the orbit. It is the true anterior angle, where re is the Earth's average equatorial radius. (·) denotes a Legendre polynomial of the first kind with an order of n.

[0108] The principal term with harmonics is J2, which reflects the flattening effect of a non-spherical gravitational potential. In low Earth orbit satellites, its magnitude can reach that of point mass gravity. This force is on the order of magnitude larger than any other perturbing force, the largest of which is the gravitational force of a point mass. The magnitude is several times that of the J2 term. Therefore, not only should the first-order perturbation of J2 be considered, but also the second-order perturbation of J2. In band harmonics, to ensure accuracy, it is also necessary to consider the J3 and J4 terms, which mainly affect long-term evolution, so as to reduce the divergence of prediction errors.

[0109] Besides the band harmonic J2, the largest magnitude harmonic in low Earth orbit satellite orbit is the fan harmonic J22. To ensure a certain level of accuracy, the application examples in this application need to consider its perturbation effects. The fan harmonic J22 under the non-uniform gravitational field of Earth can be expressed as: in, Indicates the main harmonic term and ,and and For harmonic coefficients, This represents a second-order quadratic associated Legendre polynomial. It is the geographical latitude, and , Indicates relative geographical longitude and = , Indicates geographical longitude. For harmonic coefficients and The geographical longitude of the direction of the Earth's equatorial "axis of symmetry" has been determined. , and It is the harmonic coefficient.

[0110] To accurately model the decay of orbital energy caused by nonconservative forces, the effect of atmospheric drag must be considered in low Earth orbit satellite orbits, which can be expressed as: Where S is the equivalent cross-sectional area, m is the satellite mass, and CD is the drag coefficient. Where is atmospheric density, and v is the satellite's velocity vector. It is the atmospheric rotational velocity vector. It is the daily average density of the reference ellipsoid passing through the orbital perigee. The diurnal variation factor of atmospheric density. Let r be the angle between the satellite's geocentric radius and the direction of the diurnal density peak (corresponding to the maximum atmospheric density), and let r be the distance from the satellite to the Earth's center. Let be the distance from the point where the satellite's radial vector intersects the reference ellipsoid to the Earth's center. This is the atmospheric density elevation.

[0111] This is a spatiotemporal dynamic atmospheric model that considers the diurnal variation of atmospheric density driven by sunlight and latitudinal gradient changes. It constructs a spatiotemporal dynamic atmospheric drag model, integrating variations in atmospheric pressure elevation, wind disturbances caused by Earth's rotation, and aspherical atmospheric distribution, thus achieving spatiotemporal dynamic evolution modeling of orbital drag acceleration. The β value can be directly substituted into this spatiotemporal dynamic atmospheric drag model for high-precision orbit prediction. The overall process of the orbit prediction method is as follows: Figure 4shown.

[0112] This application also provides an electronic device, which may include a processor, a memory, a receiver, and a transmitter. The processor is used to execute the two-row root generation method for the adapted analytical propagation model mentioned in the above embodiments, and / or the orbit prediction method. The processor and the memory can be connected via a bus or other means, taking a bus connection as an example. The receiver can be connected to the processor and the memory via wired or wireless means.

[0113] The processor can be a central processing unit (CPU). The processor can also be other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, or combinations of the above types of chips.

[0114] Memory, as a non-transitory computer-readable storage medium, can be used to store non-transitory software programs, non-transitory computer-executable programs, and modules, such as the program instructions / modules corresponding to the two-row root generation method of the adapted analytical propagation model in the embodiments of this application. The processor executes various functional applications and data processing by running the non-transitory software programs, instructions, and modules stored in the memory, thereby implementing the two-row root generation method of the adapted analytical propagation model in the above method embodiments.

[0115] The memory may include a program storage area and a data storage area. The program storage area may store the operating system and applications required for at least one function; the data storage area may store data created by the processor, etc. Furthermore, the memory may include high-speed random access memory and non-transitory memory, such as at least one disk storage device, flash memory device, or other non-transitory solid-state storage device. In some embodiments, the memory may optionally include memory remotely located relative to the processor, which can be connected to the processor via a network. Examples of such networks include, but are not limited to, the Internet, corporate intranets, local area networks, mobile communication networks, and combinations thereof.

[0116] The one or more modules are stored in the memory, and when executed by the processor, the two-row root generation method adapted to the parsing propagation model in the implementation embodiment is executed.

[0117] In some embodiments of this application, the user equipment may include a processor, a memory, and a transceiver unit. The transceiver unit may include a receiver and a transmitter. The processor, memory, receiver, and transmitter may be connected via a bus system. The memory is used to store computer instructions, and the processor is used to execute the computer instructions stored in the memory to control the transceiver unit to send and receive signals.

[0118] As one implementation method, the functions of the receiver and transmitter in this application can be implemented by transceiver circuits or dedicated transceiver chips, and the processor can be implemented by dedicated processing chips, processing circuits or general-purpose chips.

[0119] As another implementation approach, the server provided in this application embodiment can be implemented using a general-purpose computer. That is, the program code implementing the processor, receiver, and transmitter functions is stored in memory, and the general-purpose processor implements the processor, receiver, and transmitter functions by executing the code in memory.

[0120] This application also provides a computer-readable storage medium storing a computer program thereon. When executed by a processor, the computer program implements the steps of the aforementioned two-row root generation method for the adapted analytical propagation model, and / or the steps of the orbit prediction method. The computer-readable storage medium can be a tangible storage medium, such as random access memory (RAM), main memory, read-only memory (ROM), electrically programmable ROM, electrically erasable programmable ROM, registers, floppy disks, hard disks, removable storage disks, CD-ROMs, or any other form of storage medium known in the art.

[0121] This application also provides a computer program product, specifically comprising a computer program that, when executed by a processor, implements the steps of the two-row root generation method for the adapted analytical propagation model mentioned in the foregoing embodiments, and / or the steps of the orbit prediction method.

[0122] Those skilled in the art will understand that the exemplary components, systems, and methods described in conjunction with the embodiments disclosed herein can be implemented in hardware, software, or a combination of both. Whether implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application. When implemented in hardware, it can be, for example, electronic circuits, application-specific integrated circuits (ASICs), appropriate firmware, plug-ins, function cards, etc. When implemented in software, the elements of this application are programs or code segments used to perform the required tasks. The programs or code segments can be stored on a machine-readable medium or transmitted over a transmission medium or communication link via data signals carried on a carrier wave.

[0123] It should be clarified that this application is not limited to the specific configurations and processes described above and shown in the figures. For the sake of brevity, detailed descriptions of known methods are omitted here. In the above embodiments, several specific steps are described and shown as examples. However, the method process of this application is not limited to the specific steps described and shown. Those skilled in the art can make various changes, modifications, and additions, or change the order of steps, after understanding the spirit of this application.

[0124] In this application, features described and / or illustrated for one embodiment may be used in the same or similar manner in one or more other embodiments, and / or combined with or in place of features of other embodiments.

[0125] The above description is merely a preferred embodiment of this application and is not intended to limit this application. Various modifications and variations can be made to the embodiments of this application by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A method for generating two-row roots adapted to an analytical propagation model, characterized in that, The method includes: The process involves obtaining the historical tangent orbital elements corresponding to each epoch, and then using a pre-defined analytical propagation model to perform multiple iterations of orbital element optimization steps on each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, ensuring that the historical tangent orbital elements corresponding to each epoch are the same as their respective target tangent orbital elements. The average orbital elements corresponding to each target tangent orbital element are then derived. The orbital element optimization steps include: in the analytical propagation model, constructing a cost function for each historical tangent orbital element to obtain the target tangent orbital elements corresponding to each epoch, calculating the residual vector and Jacobian matrix for the current iteration, and updating the target tangent orbital elements for the current iteration using the LM method. The historical tangent orbital elements, the target tangent orbital elements, and the average orbital elements are all represented as Kepler six-element orbitals. Extract each average semi-major axis corresponding to each of the average orbital elements, fit each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, invert the semi-major axis decay rate to obtain the atmospheric drag parameter, and convert the atmospheric drag parameter into the standard drag term in the two rows of elements to obtain the two rows of elements.

2. The two-row root generation method for adapting the analytical propagation model according to claim 1, characterized in that, Before performing multiple iterations of orbital element optimization on each of the historical matching orbital elements using a preset analytical propagation model, the method further includes: The satellite's velocity and position vectors are converted into Kepler orbital six roots, and the initial drag parameters of the satellite are calculated based on the Knudsen number and drag acceleration. The Kepler orbital six roots include the orbital semi-major axis, orbital eccentricity, orbital inclination, right ascension of the ascending node, argument of perigee, and mean perigee. The initial drag parameters are a combination of atmospheric drag coefficient, the equivalent surface-to-mass ratio of the drag experienced by the satellite, and the average atmospheric density at the satellite's location.

3. The two-row root generation method for adapting the analytical propagation model according to claim 2, characterized in that, Before performing multiple iterations of orbital element optimization on each of the historical matching orbital elements using a preset analytical propagation model, the method further includes: Based on the Kepler orbital elements and the initial drag parameters, a residual function is constructed and the Jacobian matrix is ​​obtained using the difference quotient method; wherein, the Jacobian matrix is ​​represented as the partial derivative matrix of the residual function with respect to the orbital elements.

4. The two-row root generation method for adapting the analytical propagation model according to claim 1, characterized in that, The calculation of the residual vector for the current iteration includes: The parameters corresponding to each of the said tangent orbital elements are perturbed, and the perturbed parameters are input into the orbital propagation model to obtain the residual vector.

5. The two-row root generation method for adapting the analytical propagation model according to claim 1, characterized in that, The step of updating the target snipe orbital elements in the current iteration using the LM method includes: The damping factor is added to form the parameter update formula and the parameter update increment is calculated. The increment is then added to the number of tangent orbital elements in the current iteration to obtain the input data for the next iteration. If the residual norm corresponding to the residual vector in the current iteration decreases, the damping factor is reduced. If the parameter update increment is less than a preset threshold, or the residual norm is less than a set threshold, or the number of iterations reaches a preset threshold, then the iteration stops.

6. The two-row root generation method for adapting the analytical propagation model according to claim 1, characterized in that, The process of extracting the average semi-major axis corresponding to each of the average orbital elements, fitting each of the average semi-major axes using the least squares method to obtain the semi-major axis decay rate, and inverting the semi-major axis decay rate to obtain atmospheric drag parameters includes: Extract the average semi-major axis corresponding to each of the average orbital elements to form a dataset, and construct a physical model of how the average semi-major axis changes over time. Based on the dataset, the least squares method is used to fit each of the average semi-major axes to obtain the semi-major axis attenuation rate. The atmospheric drag parameter is calculated by substituting the semi-major axis attenuation rate into the relationship between the semi-major axis attenuation rate and the average orbital elements.

7. A trajectory prediction method, characterized in that, include: Based on the two rows of roots, higher-order band harmonics and fan harmonics are introduced, and a spatiotemporal dynamic atmospheric model is used for orbit prediction. The two rows of roots are obtained by the method described in any one of claims 1-6.

8. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the two-row root generation method for the adapted analytical propagation model as described in any one of claims 1 to 6, and / or implements the orbit prediction method as described in claim 7.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that, When executed by a processor, the computer program implements the two-row root generation method for the adapted analytical propagation model as described in any one of claims 1 to 6, and / or implements the orbit prediction method as described in claim 7.

10. A computer program product, characterized in that, It includes a computer program, which, when executed by a processor, is used to implement the two-row root generation method for the adaptive analytical propagation model as described in any one of claims 1 to 6, and / or to implement the orbit prediction method as described in claim 7.

Images