Satellite ground point calculation method based on SGP4 orbit model

Abstract

This invention discloses a method for calculating the sub-satellite point based on the SGP4 orbital model, relating to the field of satellite navigation and orbit prediction technology. The polar spherical projector of this invention eliminates longitude jumps in latitude >75°, ensuring continuous trajectory output for scenarios such as polar scientific expeditions and glacier monitoring. By directly mapping three-dimensional radial information to generate latitude, it avoids the iterative convergence problem of the ellipsoidal model, ensuring natural symmetry of the North and South Poles coordinates and continuous longitude domain. Second-order differential compensation of the mean anomaly angle, combined with dynamic correction of the Earth's rotational angular velocity, strictly aligns the UTC timestamp with the orbital position. It suppresses the cumulative error generated by traditional linear interpolation at perigee. It enhances adaptability to complex perturbation environments; atmospheric drag vectorization compensation captures real-time aerodynamic direction through relative velocity cross terms, combined with continuous processing of the density index model, significantly improving the orbital attenuation prediction accuracy of low-Earth orbit satellites during geomagnetic storms.

Description

Technical Field

[0001] This invention relates to the field of satellite navigation and orbit prediction technology, and in particular to a method for calculating the satellite nadir point based on the SGP4 orbit model. Background Technology

[0002] Satellite nadir point calculation is a core technology in satellite navigation, Earth observation and other fields, and it is based on TLE data and SGP4 orbit model to achieve position prediction. In scenarios such as polar scientific research and climate monitoring, the positioning needs of high latitude regions (latitude > 75°) are becoming increasingly prominent, such as monitoring the Arctic shipping route or tracking glacier changes. Such applications require that the nadir point output can stably cover the global range, including the area near the poles.

[0003] When a satellite enters a high-latitude region, existing calculation methods use a standard ellipsoidal projection model, such as WGS84, for coordinate transformation. However, due to the abrupt curvature changes in the Earth's geometric model near the poles, the longitude calculation results may exhibit unexpected jumps. These jumps manifest as drastic fluctuations in the output longitude value within a short period, affecting the continuity and reliability of the trajectory. In high-latitude missions, such as when a satellite passes over an Antarctic research station, positional deviations may interfere with real-time data correlation.

[0004] Some solutions mitigate this problem by optimizing the ellipsoidal projection algorithm, such as by introducing an adaptive coordinate transformation strategy. These methods maintain computational stability at normal latitudes, but still rely on the same mathematical model framework near the poles. Some implementations attempt to add latitude threshold judgments to switch between different projection methods, but fail to completely eliminate the singularity effect. Summary of the Invention

[0005] In view of the aforementioned existing problems, the present invention is proposed.

[0006] This invention provides a satellite nadir point calculation method based on the SGP4 orbit model to solve the problems of existing methods causing longitude jumps in high-latitude regions due to ellipsoidal projection singularities, and insufficient time synchronization and perturbation compensation accuracy, which affect the reliability of scenarios such as polar monitoring.

[0007] To solve the above-mentioned technical problems, the present invention provides the following technical solution:

[0008] This invention provides a method for calculating the sub-satellite point based on the SGP4 orbital model, which includes:

[0009] Step S1: Obtain the satellite's TLE orbit data and target UTC time;

[0010] Step S2: Calculate the satellite's position in the Earth-fixed coordinate system based on the SGP4 orbit model;

[0011] Step S3: When the absolute value of the satellite latitude is greater than or equal to the preset threshold, perform the polar spherical projection operation.

[0012] Step S4: When the absolute value of the satellite latitude is less than the preset threshold, latitude and longitude conversion is performed using the WGS84 ellipsoid model.

[0013] Step S5: Output the corrected latitude and longitude coordinates of the nadir point and the corresponding UTC time.

[0014] As a preferred embodiment of the satellite nadir point calculation method based on the SGP4 orbit model described in this invention, wherein: in step S3, the polar spherical projection operation includes:

[0015] Calculate the sign of the satellite's geocentric position vector component along the Z-axis;

[0016] Elevation is determined based on satellite geocentric distance and Earth's polar radius;

[0017] Output longitude based on the arctangent function of the X and Y coordinates in the Earth-fixed coordinate system;

[0018] The latitude value is generated based on the geometric relationship between elevation and polar radius, including dividing the elevation by a preset polar radius scaling factor, and multiplying the scaled elevation ratio by π / 2 and then multiplying by the sign of the Z-axis component.

[0019] During the polar spherical projection operation in step S3, when the satellite's approximate latitude... When entering the polar projection branch, first obtain the geocentric coordinates:

[0020] , , ,

[0021] Then, derive the polar elevation using the polar radius:

[0022] , , ,

[0023] in, The satellite position components in the Earth-fixed coordinate system are expressed in meters (m). for The sign function, dimensionless. For symbolic functions, The distance from the Earth's center is in meters (m). The polar radius of the Earth is expressed in meters (m). This is the polar elevation, in meters (m). Polar projection longitude, unit: rad. It is the arctangent function in the fourth quadrant. The latitude is the polar projection latitude, in rad. The polar radius scaling factor is dimensionless. For initial latitude calculation, the unit is rad. This is the latitude threshold, in rad.

[0024] As a preferred embodiment of the satellite sub-satellite point calculation method based on the SGP4 orbit model described in this invention, wherein: in step S2, the SGP4 orbit model calculation step includes:

[0025] Analyzing the average number of orbital elements in TLE data;

[0026] An atmospheric drag term is introduced to dynamically compensate for the rate of change of the semi-major axis. This includes analyzing the drag coefficient term in the TLE data, calling the atmospheric density empirical model to output the density estimate of the satellite's current position, and inputting the drag coefficient, atmospheric density, and satellite relative velocity vector into the dynamic equation.

[0027] As a preferred embodiment of the satellite sub-satellite point calculation method based on the SGP4 orbit model described in this invention, wherein: in the SGP4 orbit model calculation step of step S2, dynamic compensation for the semi-major axis variation rate specifically includes:

[0028] Within each integration step Δt, a real-time drag model is constructed using the satellite aerodynamic ballistic coefficients:

[0029] ,

[0030] ,

[0031] ,

[0032] in, This is the ballistic coefficient, in units of... , The aerodynamic drag coefficient is dimensionless. This is the reference area for stress application, in units of... , Satellite mass, in kg. It is a relative velocity vector, with units of . , For the Earth's inertial velocity, , This is the Earth's rotational angular velocity vector, in units of . , This is the satellite position vector, in meters. Real-time atmospheric density, in units , Baseline density , For satellite altitude, As the reference height, The unit is height, in meters (m).

[0033] Define drag acceleration as:

[0034] , ,

[0035] in, The drag acceleration vector, in units of , For relative velocity scalar, ;

[0036] Convert the radial projection to the rate of change of the semi-major axis:

[0037] , ,

[0038] in, The rate of change of the semi-major axis due to resistance, in units of , This is the current semi-major axis, in meters. The Earth's gravitational constant, in units of 1000 ppm. , The intrinsic semi-major axis variation rate of SGP4 , This is the corrected total rate of change. ;when Time setting ;

[0039] In step S2, the SGP4 orbital model calculation also involves a higher-order derivation of the mean anterior angle, specifically including: calculating the mean angular velocity using the updated semi-major axis, with the following formula:

[0040] ,

[0041] ,

[0042] ,

[0043] in, For average motion, the unit is , The first derivative of the angle of approach is given by unit. , The second derivative of the angle of approach. ;

[0044] Update using Taylor expansion within the integration step Δt:

[0045] ,

[0046] in, The angle is the angle of approach, and the unit is rad. The integral step size is in seconds. The conversion between UTC and mechanical time is completed through the ΔAT table and epoch difference, without extrapolation truncation. ΔAT is the table of time difference parameters between International Atomic Time and UTC.

[0047] As a preferred embodiment of the satellite nadir point calculation method based on the SGP4 orbit model described in this invention, the atmospheric density empirical model in the dynamic compensation includes:

[0048] Use a standard atmospheric model database to obtain reference values ​​for atmospheric density corresponding to satellite altitude;

[0049] The reference value is scaled and adjusted by combining the real-time spatial environment index.

[0050] As a preferred embodiment of the satellite sub-satellite point calculation method based on the SGP4 orbit model described in this invention, it further includes:

[0051] Real-time acquisition of diurnal variation parameters from the International Earth Rotation Service bulletins;

[0052] When the network is interrupted, switch to linear extrapolation of the most recent 72 hours of historical parameters stored locally;

[0053] The Earth's rotation angular velocity is dynamically corrected based on the diurnal variation value, including converting the diurnal variation parameter into an angular velocity compensation coefficient and generating the corrected real-time angular velocity value through multiplication.

[0054] Greenwich Mean Time (GMT) is calculated using the corrected angular velocity.

[0055] As a preferred embodiment of the satellite sub-satellite point calculation method based on the SGP4 orbit model described in this invention, the real-time acquisition of the International Earth Rotation Service Bulletin includes:

[0056] Request the parameters corresponding to the target UTC time through the astronomical data API interface;

[0057] After verifying the validity of the parameters, they are written into the Earth's rotation angular velocity correction calculation module.

[0058] As a preferred embodiment of the satellite sub-satellite point calculation method based on the SGP4 orbit model described in this invention, it further includes:

[0059] When calculating the rate of change of the mean anomaly angle, a second-order time differential term is generated by combining the satellite orbital angular velocity and the rate of change of the semi-major axis.

[0060] Incremental interpolation is performed based on the time difference between the target UTC time and the TLE epoch.

[0061] As a preferred embodiment of the satellite nadir point calculation method based on the SGP4 orbit model described in this invention, all coordinate outputs are bound to UTC timestamps accurate to the millisecond level;

[0062] TLE data is acquired and cached in real time through the aerospace data platform interface.

[0063] As a preferred embodiment of the satellite nadir point calculation method based on the SGP4 orbit model described in this invention, the polar spherical projection operation and the WGS84 ellipsoid model conversion are executed mutually exclusively and are automatically switched through the latitude threshold judgment module.

[0064] The beneficial effects of this invention are as follows: This invention improves the reliability of high-latitude positioning; the polar spherical projector eliminates the longitude jump phenomenon in latitude >75°, ensuring continuous trajectory output in scenarios such as polar scientific expeditions and glacier monitoring; by directly mapping three-dimensional radial information to generate latitude, it avoids the iterative convergence problem of the ellipsoidal model, making the North and South Pole coordinates naturally symmetrical and the longitude domain continuous; the second-order differential compensation of the mean apogee angle combined with dynamic correction of the Earth's rotation angular velocity strictly aligns the UTC timestamp with the orbital position; it suppresses the cumulative error generated by traditional linear interpolation at perigee, meeting the millisecond-level timing synchronization requirements of 6G satellite communication, etc.

[0065] This invention enhances adaptability to complex perturbation environments. Atmospheric drag vectorization compensation captures real-time aerodynamic direction through relative velocity cross terms, and combined with continuous processing of density index model, significantly improves the orbital attenuation prediction accuracy of low-orbit satellites during geomagnetic storms. The energy conservation design of the semi-major axis change rate effectively suppresses position drift of long arc segment integrals.

[0066] This invention optimizes system robustness. The local extrapolation mechanism ensures uninterrupted operation of Earth rotation correction when the IERS data network is interrupted, while the automatic switching of latitude thresholds achieves seamless global coverage. While maintaining the SGP4 calculation framework, this invention achieves a leap in accuracy through preprocessing (TLE resolution) and postprocessing (coordinate transformation). Attached Figure Description

[0067] To more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings used in the following description of the embodiments will be briefly introduced. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0068] Figure 1 This is a flowchart illustrating the satellite sub-satellite point calculation method based on the SGP4 orbital model in Example 1. Detailed Implementation

[0069] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the specific embodiments of the present invention will be described in detail below with reference to the accompanying drawings.

[0070] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and those skilled in the art can make similar extensions without departing from the spirit of the invention. Therefore, the invention is not limited to the specific embodiments disclosed below.

[0071] Secondly, the term "one embodiment" or "embodiment" as used herein refers to a specific feature, structure, or characteristic that may be included in at least one implementation of the present invention. The phrase "in one embodiment" appearing in different places in this specification does not necessarily refer to the same embodiment, nor is it a single or selective embodiment that is mutually exclusive with other embodiments.

[0072] Example 1, referring to Figure 1 This embodiment provides a method for calculating the sub-satellite point based on the SGP4 orbit model, including:

[0073] Step S1: Obtain the satellite's TLE orbit data and target UTC time;

[0074] Step S2: Calculate the satellite's position in the Earth-fixed coordinate system based on the SGP4 orbit model;

[0075] In step S2, the SGP4 orbital model calculation steps include:

[0076] Analyzing the average number of orbital elements in TLE data;

[0077] The atmospheric drag term is introduced to dynamically compensate for the rate of change of the semi-major axis. This includes analyzing the drag coefficient term in the TLE data, calling the atmospheric density empirical model to output the density estimate of the satellite's current position, and inputting the drag coefficient, atmospheric density, and satellite relative velocity vector into the dynamic equation.

[0078] The rate of change of the mean anterior angle is calculated based on UTC time, and a second-order time differential term is generated.

[0079] In step S2, the SGP4 orbital model calculation step involves dynamically compensating for the rate of change of the semi-major axis, specifically including:

[0080] Within each integration step Δt, a real-time drag model is constructed using the satellite aerodynamic ballistic coefficients:

[0081] ,

[0082] ,

[0083] ,

[0084] in, This is the ballistic coefficient, in units of... , The aerodynamic drag coefficient is dimensionless. This is the reference area for stress application, in units of... , Satellite mass, in kg. It is a relative velocity vector, with units of . , For the Earth's inertial velocity, , This is the Earth's rotational angular velocity vector, in units of . , This is the satellite position vector, in meters. Real-time atmospheric density, in units , Baseline density , For satellite altitude, As the reference height, The unit is height, in meters (m).

[0085] Define drag acceleration as:

[0086] , ,

[0087] in, The drag acceleration vector, in units of , For relative velocity scalar, ;

[0088] Convert the radial projection to the rate of change of the semi-major axis:

[0089] , ,

[0090] in, The rate of change of the semi-major axis due to resistance, in units of , This is the current semi-major axis, in meters. The Earth's gravitational constant, in units of 1000 ppm. , The intrinsic semi-major axis variation rate of SGP4 , This is the corrected total rate of change. ;when Time setting ;

[0091] Specifically, the compensation captures the aerodynamic direction in real time with the relative velocity cross term, the density index model avoids density discrete jumps by processing the scale height continuously, the semi-major axis change rate directly quantifies energy loss through position acceleration dot product, the drag dissipation is fed back into the orbit size in real time, the spiral sinking error during multi-step integration in near-Earth orbit is suppressed, the threshold clipping freezes the drag component in the thin layer, maintains the numerical smoothness of the high orbit, and improves the overall consistency of energy conservation in the calculation of long arc segments.

[0092] In step S2, the SGP4 orbital model calculation also involves a higher-order derivation of the mean anterior angle, specifically including: calculating the mean angular velocity using the updated semi-major axis, with the following formula:

[0093] ,

[0094] ,

[0095] ,

[0096] in, For average motion, the unit is , The first derivative of the angle of approach is given by unit. , The second derivative of the angle of approach. ;

[0097] Update using Taylor expansion within the integration step Δt:

[0098] ,

[0099] in, The angle is the angle of approach, and the unit is rad. The integration step size is in seconds; the conversion between UTC and mechanical time is completed through the ΔAT table and epoch difference, without extrapolation truncation. ΔAT is the difference between International Atomic Time and UTC.

[0100] Specifically, the second derivative of the mean motion is directly coupled to the semi-major axis decay rate, capturing the higher-order perturbation of the orbit phase by the low-orbit drag fluctuation. The single-step Taylor expansion retains the second-order terms, which reduces the phase error of long-term propagation compared to single-order extrapolation. The time reference is unified to mechanical time, and the mean aperimeter angle is inverted at the sub-satellite point and maintains a strict one-to-one correspondence with the time mark, enhancing the consistency of attitude and orbit decoupling timing.

[0101] Empirical models of atmospheric density in dynamic compensation include:

[0102] Use a standard atmospheric model database to obtain reference values ​​for atmospheric density corresponding to satellite altitude;

[0103] The reference value is scaled and adjusted based on the real-time spatial environment index;

[0104] Step S3: When the absolute value of the satellite latitude is greater than or equal to the preset threshold, perform the polar spherical projection operation.

[0105] In step S3, the polar spherical projection operation includes:

[0106] Calculate the sign of the satellite's geocentric position vector component along the Z-axis;

[0107] Elevation is determined based on satellite geocentric distance and Earth's polar radius;

[0108] Output longitude based on the arctangent function of the X and Y coordinates in the Earth-fixed coordinate system;

[0109] The latitude value is generated based on the geometric relationship between elevation and polar radius, including dividing the elevation by a preset polar radius scaling factor, subtracting the scaled elevation ratio from π / 2, and then multiplying by the sign of the Z-axis component.

[0110] During the polar spherical projection operation in step S3, when the satellite's approximate latitude... When entering the polar projection branch, first obtain the geocentric coordinates:

[0111] , , ,

[0112] Then, derive the polar elevation using the polar radius:

[0113] , , ,

[0114] in, The satellite position components in the Earth-fixed coordinate system are expressed in meters (m). for The sign function, dimensionless. For symbolic functions, The distance from the Earth's center is in meters (m). The polar radius of the Earth is expressed in meters (m). This is the polar elevation, in meters (m). Polar projection longitude, unit: rad. It is the arctangent function in the fourth quadrant. The latitude is the polar projection latitude, in rad. The polar radius scaling factor is dimensionless. The value of 0.996 can be obtained by regression analysis of measured gravity anomalies in polar regions; the drag density threshold is determined according to the corresponding value of the MSIS-00 (Mass Spectrometer and Incoherent Scatter, MSIS-00) model at an altitude of 600 km. For initial latitude calculation, the unit is rad. Latitude threshold, in rad;

[0115] Specifically, this projection process no longer performs ellipsoid inverse calculation, but directly uses polar elevation to map the three-dimensional radial information into collinear latitude offset, avoiding the problem of slow convergence of ellipsoid iteration at high latitudes. The symmetrical sign term ensures the natural symmetry of the North and South pole coordinates, and the longitude still maintains the continuity of the arctangent domain, eliminating the 0° / 360° jump. The overall switching is smooth near the latitude threshold, eliminating boundary ring errors and maintaining coordinate monotonicity.

[0116] Step S4: When the absolute value of the satellite latitude is less than the preset threshold, latitude and longitude conversion is performed using the WGS84 ellipsoid model.

[0117] Step S5: Output the corrected latitude and longitude coordinates of the nadir point and the corresponding UTC time;

[0118] The method also includes:

[0119] Real-time acquisition of diurnal variation parameters from the International Earth Rotation Service bulletins;

[0120] When the network is interrupted, switch to linear extrapolation of the most recent 72 hours of historical parameters stored locally;

[0121] The Earth's rotation angular velocity is dynamically corrected based on the diurnal variation value, including converting the diurnal variation parameter into an angular velocity compensation coefficient and generating the corrected real-time angular velocity value through multiplication.

[0122] Greenwich Mean Time (GMT) was calculated using the corrected angular velocity;

[0123] Real-time access to the International Earth Rotation Service bulletins includes:

[0124] Request the parameters corresponding to the target UTC time through the astronomical data API interface;

[0125] After verifying the validity of the parameters, they are written into the Earth's rotation angular velocity correction calculation module.

[0126] The above methods also include:

[0127] When calculating the rate of change of the mean anomaly angle, a second-order time differential term is generated by combining the satellite orbital angular velocity and the rate of change of the semi-major axis.

[0128] Incremental interpolation is performed based on the time difference between the target UTC time and the TLE epoch;

[0129] All coordinate outputs are bound to UTC timestamps accurate to milliseconds;

[0130] TLE data is acquired and cached in real time through the aerospace data platform interface;

[0131] The polar spherical projection operation and the WGS84 ellipsoid model conversion are executed mutually exclusively, and are automatically switched through the latitude threshold judgment module.

[0132] Unless otherwise expressly defined, the various constants or thresholds involved in this invention (e.g.) =0.996、 threshold All of these parameters (etc.) are adjustable based on model fitting, experimental calibration, or engineering experience. The listed values ​​are merely preferred examples and do not constitute a limitation of the present invention. Those skilled in the art can determine them within a reasonable range according to specific task requirements.

[0133] The above parameters can be obtained by regression analysis using publicly available atmospheric models (such as the MSIS series), gravity field / ellipsoid parameter libraries, or mission simulation data.

[0134] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention, and all such modifications or substitutions should be covered within the scope of the claims of the present invention.

Claims

1. A method for calculating the sub-satellite point based on the SGP4 orbital model, characterized in that, Includes the following steps: Step S1: Obtain the satellite's TLE orbit data and target UTC time; Step S2: Calculate the satellite's position in the Earth-fixed coordinate system based on the SGP4 orbit model; Step S3: When the absolute value of the satellite latitude is greater than or equal to the preset threshold, perform the polar spherical projection operation. Step S4: When the absolute value of the satellite latitude is less than the preset threshold, latitude and longitude conversion is performed using the WGS84 ellipsoid model. Step S5: Output the corrected latitude and longitude coordinates of the nadir point and the corresponding UTC time.

2. The satellite nadir point calculation method based on the SGP4 orbit model as described in claim 1, characterized in that, In step S3, the polar spherical projection operation includes: Calculate the sign of the satellite's geocentric position vector component along the Z-axis; Elevation is determined based on satellite geocentric distance and Earth's polar radius; Output longitude based on the arctangent function of the X and Y coordinates in the Earth-fixed coordinate system; The latitude value is generated based on the geometric relationship between elevation and polar radius, including dividing the elevation by a preset polar radius scaling factor, and multiplying the scaled elevation ratio by π / 2 and then multiplying by the sign of the Z-axis component. During the polar spherical projection operation in step S3, when the satellite's approximate latitude... When entering the polar projection branch, first obtain the geocentric coordinates: ， ， ， Then, derive the polar elevation using the polar radius: ， ， ， in, The satellite position components in the Earth-fixed coordinate system are expressed in meters (m). for The sign function, dimensionless. For symbolic functions, The distance from the Earth's center is in meters (m). The polar radius of the Earth is expressed in meters (m). This is the polar elevation, in meters (m). Polar projection longitude, unit: rad. It is the arctangent function in the fourth quadrant. The latitude is the polar projection latitude, in rad. The polar radius scaling factor is dimensionless. For initial latitude calculation, the unit is rad. This is the latitude threshold, in rad.

3. The satellite nadir point calculation method based on the SGP4 orbit model as described in claim 1, characterized in that, In step S2, the SGP4 orbital model calculation steps include: Analyzing the average number of orbital elements in TLE data; The atmospheric drag term is introduced to dynamically compensate for the rate of change of the semi-major axis. This includes analyzing the drag coefficient term in the TLE data, calling the atmospheric density empirical model to output the density estimate of the satellite's current position, and inputting the drag coefficient, atmospheric density, and satellite relative velocity vector into the dynamic equation. The rate of change of the mean anterior angle is calculated based on UTC time, and a second-order time differential term is generated.

4. The method for calculating the satellite nadir point based on the SGP4 orbit model as described in claim 3, characterized in that, In step S2, the SGP4 orbital model calculation step involves dynamically compensating for the rate of change of the semi-major axis, specifically including: Within each integration step Δt, a real-time drag model is constructed using the satellite aerodynamic ballistic coefficients: ， ， ， in, This is the ballistic coefficient, in units of... , The aerodynamic drag coefficient is dimensionless. This is the reference area for stress application, in units of... , Satellite mass, in kg. It is a relative velocity vector, with units of . , For the Earth's inertial velocity, , This is the Earth's rotational angular velocity vector, in units of . , This is the satellite position vector, in meters. Real-time atmospheric density, in units , Baseline density , For satellite altitude, As the reference height, The unit is height, in meters (m). Define drag acceleration as: ， ， in, The drag acceleration vector, in units of , For relative velocity scalar, ; Convert the radial projection to the rate of change of the semi-major axis: ， ， in, The rate of change of the semi-major axis due to resistance, in units of , This is the current semi-major axis, in meters. The Earth's gravitational constant, in units of 1000 ppm. , The intrinsic semi-major axis variation rate of SGP4 , This is the corrected total rate of change. ;when Time setting ; In step S2, the SGP4 orbital model calculation also involves a higher-order derivation of the mean anterior angle, specifically including: calculating the mean angular velocity using the updated semi-major axis, with the following formula: ， ， ， in, For average motion, the unit is , The first derivative of the angle of approach is given by unit. , The second derivative of the angle of approach. ; Update using Taylor expansion within the integration step Δt: ， in, The angle is the angle of approach, and the unit is rad. The integral step size is in seconds. The conversion between UTC and mechanical time is completed through the ΔAT table and epoch difference, without extrapolation truncation. ΔAT is the table of time difference parameters between International Atomic Time and UTC.

5. The method for calculating the satellite nadir point based on the SGP4 orbit model as described in claim 3, characterized in that, The empirical model for atmospheric density in the dynamic compensation includes: Use a standard atmospheric model database to obtain reference values ​​for atmospheric density corresponding to satellite altitude; The reference value is scaled and adjusted by combining the real-time spatial environment index.

6. The satellite nadir point calculation method based on the SGP4 orbit model as described in claim 1, characterized in that, Also includes: Real-time acquisition of diurnal variation parameters from the International Earth Rotation Service bulletins; When the network is interrupted, switch to linear extrapolation of the most recent 72 hours of historical parameters stored locally; The Earth's rotation angular velocity is dynamically corrected based on the diurnal variation value, including converting the diurnal variation parameter into an angular velocity compensation coefficient and generating the corrected real-time angular velocity value through multiplication. Greenwich Mean Time (GMT) is calculated using the corrected angular velocity.

7. The method for calculating the satellite nadir point based on the SGP4 orbit model as described in claim 6, characterized in that, The real-time acquisition of the International Earth Rotation Service bulletins includes: Request the parameters corresponding to the target UTC time through the astronomical data API interface; After verifying the validity of the parameters, they are written into the Earth's rotation angular velocity correction calculation module.

8. The satellite nadir point calculation method based on the SGP4 orbit model as described in claim 1, characterized in that, Also includes: When calculating the rate of change of the mean anomaly angle, a second-order time differential term is generated by combining the satellite orbital angular velocity and the rate of change of the semi-major axis. Incremental interpolation is performed based on the time difference between the target UTC time and the TLE epoch.

9. The method for calculating the satellite nadir point based on the SGP4 orbit model as described in claim 1, characterized in that, All coordinate outputs are bound to UTC timestamps accurate to milliseconds; TLE data is acquired and cached in real time through the aerospace data platform interface.

10. A method for calculating the satellite nadir point based on the SGP4 orbit model as described in any one of claims 1 to 9, characterized in that, The polar spherical projection operation and the WGS84 ellipsoid model conversion are executed mutually exclusively and are automatically switched through the latitude threshold judgment module.

Images