Method for calculating visible time period between satellite and any area on ground

Abstract

The invention discloses a method for calculating a visible time period between a satellite and any area of the ground, relates to the technical field of visible time period calculation, and aims to solve the problems of orbit perturbation errors, low intensive transit efficiency and coverage blind areas in satellite visible time period calculation. Through three technical means of dynamic resistance correction, space-time clustering optimization and inter-satellite link cooperation, a resistance compensation mechanism fuses real-time change of atmospheric density and a satellite motion direction, an accumulative error of orbit prediction is corrected, and it is ensured that starting and ending moments of a visible window accurately correspond to an actual transit event; the space-time clustering algorithm performs grouping processing on dense transit satellites, so that repeated calculation of satellite coverage in the same time-space domain is avoided, and the large-scale constellation scheduling response time is remarkably shortened; the inter-satellite link cooperation mechanism breaks through the geometric limitation of single-satellite direct connection, indirectly covers the traditional blind area through the relay satellite, and improves the service continuity of high-latitude and complex terrain areas.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of visible time period calculation, in particular to a satellite and ground any area visible time period calculation method. BACKGROUND

[0002] The current satellite and ground area visible time period calculation mainly serves real-time scenarios such as remote sensing monitoring, emergency communication, etc.; with the large-scale deployment of super low orbit satellites, its wide coverage and low transmission delay advantages are obvious, but the orbit characteristics introduce new challenges: such satellites are affected by the high atmosphere resistance perturbation, and the orbit root number changes nonlinearly; the existing method drives the SGP4 model based on TLE data, iterates the satellite position through fixed step, and judges whether the ground area elevation angle exceeds the threshold, and finally outputs the visible time window in UTC format.

[0003] Some schemes introduce an atmospheric resistance compensation factor, such as the NRLMSISE-00 model, to correct the SGP4 output position, and increase the calculation step size during the satellite overpass period; however, the resistance model relies on static atmospheric parameters and cannot respond to instantaneous density changes caused by sudden space weather events such as solar storms; at the same time, the fixed step strategy still has a response lag when the satellite enters and exits the dense atmosphere, resulting in a deviation in the determination of the start and end points of the visible window; some systems use Kalman filtering to predict the orbit, but the calculation load is significantly increased.

[0004] To balance accuracy and efficiency, some schemes propose a dynamic step size adjustment mechanism: when the satellite height drops to a critical value, such as 350km, the iteration step size is automatically shortened, this method alleviates the error accumulation under normal working conditions, but the response to the sudden change of orbit parameters is insufficient; when the satellite encounters sudden atmospheric density fluctuations, the position prediction error will still be transmitted to the time window calculation layer, causing the start and end time of the visible period to shift; in addition, the dynamic step size algorithm needs to preset the height threshold, which is difficult to adapt to the maneuvering characteristics of different satellite configurations. SUMMARY

[0005] In view of the above existing problems, the present application is proposed.

[0006] The present application provides a satellite and ground any area visible time period calculation method to solve the problem of orbit mutation of super low orbit satellites caused by atmospheric resistance, and the deviation and low efficiency of existing methods caused by static compensation and fixed step iteration in visible window calculation.

[0007] To solve the above technical problems, the present application provides the following technical scheme:

[0008] The satellite and ground any area visible time period calculation method provided by the embodiment of the present application comprises,

[0009] Step S1, obtaining TLE data of satellites and longitude and latitude boundary of target ground area;

[0010] Step S2, calculating orbit position sequence of satellites in preset time range based on SGP4 orbit model;

[0011] Step S3, dynamically correcting instantaneous position of satellites according to orbit position sequence and atmospheric resistance perturbation compensation factor;

[0012] Step S4, grouping transit satellites by using space-time clustering algorithm based on corrected satellite position;

[0013] Step S5, calculating elevation angle variation curve of target ground area according to satellite group; if the representative satellite is invisible to the target area, searching for a satellite having effective inter-satellite link with the representative satellite and verifying whether the satellite is visible to the target area and the representative satellite at the same time;

[0014] Step S6, outputting UTC format visible time window according to continuous time period in which the elevation angle exceeds a set threshold.

[0015] As a preferred scheme of the satellite and ground any area visible time period calculation method, in step S3, the atmospheric resistance perturbation compensation factor is generated by the following way:

[0016] NRLMSISE-00 atmospheric density model is called to calculate resistance proportion coefficient according to current orbit height and speed of the satellite, and the proportion coefficient is embedded into the correction operation example of SGP4 position output result.

[0017] As a preferred scheme of the satellite and ground any area visible time period calculation method, in step S3, the generation way of the atmospheric resistance perturbation compensation factor comprises:

[0018] The instantaneous atmospheric density is obtained based on NRLMSISE-00:

[0019] ,

[0020] Wherein, represents the outer atmospheric density of the i-th iteration, unit kg , is the orbit height, unit km, , , are the sub-point latitude and longitude respectively, is the current time, is the daily 10.7cm solar flux index sfu, is the daily planetary activity index, represents the NRLMSISE-00 density evaluation function;

[0021] The relative velocity after the correction of the earth rotation is calculated by the formula:

[0022] ,

[0023] Where, is the relative velocity scalar, unit , is the satellite geocentric velocity vector, unit , is the earth rotation angular velocity vector , is the satellite geocentric position vector, unit km, represents the vector cross product, represents the vector module;

[0024] The area mass ratio is:

[0025] ,

[0026] Where, is the area mass ratio, unit , is the satellite force reference area, unit , is the satellite mass, unit kg;

[0027] The aerodynamic drag coefficient is:

[0028] ,

[0029] Where, is the instantaneous aerodynamic drag coefficient, is the reference drag coefficient, is the high altitude correction amplitude, is the natural exponential function, is the atmospheric scale height, unit km;

[0030] The drag proportionality coefficient is:

[0031] ,

[0032] Where, is the th iteration drag proportionality coefficient, unit ;

[0033] Position correction embedded SGP4:

[0034] ,

[0035] Where, is the corrected satellite position vector, is a step size.

[0036] As a preferred scheme of the method for calculating the time period during which any area on the ground is visible to a satellite, the spatio-temporal clustering algorithm in step S4 specifically comprises:

[0037] The spatio-temporal neighborhood is defined as a satellite geocentric angular difference of ≤5° and a transit time difference of ≤60 seconds. A density-based clustering algorithm is used to group satellites satisfying the neighborhood condition into the same group, and a satellite at the geometric center position of each group is selected as a representative satellite.

[0038] As a preferred scheme of the method for calculating the time period during which any area on the ground is visible to a satellite, the inter-satellite link cooperative coverage judgment is performed after step S5:

[0039] When the representative satellite is not visible to the target area, a neighboring satellite having an inter-satellite link with the representative satellite is searched, and it is verified whether the neighboring satellite and the representative satellite are simultaneously visible to the target area. If they are simultaneously visible, the visible time period of the neighboring satellite is included in the output window.

[0040] The judgment criterion for simultaneous visibility is:

[0041] The instantaneous elevation angle of the neighboring satellite to the target area is >5°, the line connecting the neighboring satellite and the representative satellite is not blocked by the Earth, and the sum of the propagation delays of the two links is < the maximum allowed delay.

[0042] As a preferred scheme of the method for calculating the time period during which any area on the ground is visible to a satellite, the effectiveness of the inter-satellite link needs to meet the following two conditions simultaneously:

[0043] The geometric distance between the satellites is less than an orbit height-related threshold, and the link has a continuous visible arc segment outside the Earth's shadow zone.

[0044] As a preferred scheme of the method for calculating the time period during which any area on the ground is visible to a satellite, the dynamic adjustment method of the elevation angle threshold in step S6 is:

[0045] If the target area is located within the polar circle, i.e., the latitude is ≥65°, the threshold is lowered to 3°, and if the satellite orbit height is lower than 350 km, the threshold is raised to 7°.

[0046] As a preferred scheme of the method for calculating the time period during which any area on the ground is visible to a satellite, the orbit position calculation in step S2 adopts a variable step size iteration strategy:

[0047] When the satellite orbit height is lowered to below 400 km, the iteration step size is shortened to 1 / 10 of the original step size.

[0048] When the satellite elevation angle rate of change exceeds 0.5° per second, an automatic step halving mechanism is triggered.

[0049] As a preferred scheme of the satellite and ground any area visible time period calculation method, in step S2, the variable step iteration strategy is used to calculate the orbit position, including:

[0050] a) initialization, formula:

[0051] ,

[0052] Where, is the satellite geocentric position vector at the epoch time, unit km, is the satellite geocentric velocity vector at the same time, unit , is the starting time of the epoch, unit s, is the two-line element set;

[0053] b) single step and double half step propagation, propagation formula:

[0054] ,

[0055] Where, is the satellite position vector obtained by single step propagation, unit km, is the corresponding velocity vector, unit , is the starting time of the th iteration, unit s, is the current step, unit s,

[0056] ,

[0057] Where, is the satellite position vector obtained by double half step propagation, unit km, is the corresponding velocity vector, unit ;

[0058] c) local truncation error is:

[0059] ,

[0060] Where, is the th local truncation error estimate;

[0061] d) define height, elevation angle and rate of change respectively:

[0062] ,

[0063] Where, is the orbital height, is the module length, is the average radius of the earth,

[0064] ,

[0065] where, is the elevation angle of the satellite to the target area, is the geocentric coordinate of the geometric center of the target area, in km, is the unit vector of the target center,

[0066] ,

[0067] where, is the rate of change of the elevation angle, in , is the last iteration elevation angle, is the last step length, in s;

[0068] e) Step length update is performed in the following way:

[0069] ,

[0070] where, is the error proportion coefficient, is the error tolerance, is the constant to prevent zero denominator km,

[0071] ,

[0072] where, is the error-driven scaling factor,

[0073] ,

[0074] ,

[0075] where, is the height-driven scaling factor, is the height step length scaling coefficient, is the height threshold, is the elevation angle-driven scaling factor, is the elevation angle step length scaling coefficient, is the rate of change of the elevation angle threshold;

[0076] ,

[0077] where, is the next iteration step length;

[0078] f) Termination with storage: when decreases to a set minimum value or exceeds a preset time range, iteration is stopped, and the output wherein, is the sequence time.

[0079] As a preferred scheme of the satellite and ground any area visible time period calculation method of the application, wherein: it further contains a visibility conflict resolution mechanism:

[0080] Detecting the coverage overlap of multiple satellites to the same area in the same period, automatically assigning communication frequency bands according to the satellite payload type and area priority, and marking the main satellite identifier in the output time window.

[0081] The application has the advantages that: the application aims at the orbit perturbation error, low efficiency of dense transit and coverage blind area in satellite visible period calculation, and solves the problems through three technical means of dynamic drag correction, time-space clustering optimization and inter-satellite link cooperation: the drag compensation mechanism fuses the real-time change of atmospheric density and the motion direction of the satellite, corrects the cumulative error of orbit prediction, and ensures that the visible window start and end time accurately corresponds to the actual transit event; the time-space clustering algorithm processes the dense transit satellites in groups, avoids repeated calculation of satellite coverage in the same time-space domain, and significantly shortens the large-scale constellation scheduling response time; the inter-satellite link cooperation mechanism breaks through the geometric limitation of single-satellite direct connection, indirectly covers the traditional blind area through relay satellites, and improves the service continuity of high-latitude and complex terrain areas; the variable step iteration strategy dynamically adjusts the calculation granularity according to the orbit height and elevation angle change rate, and balances the precision and efficiency independently, and adapts to different task scene requirements. BRIEF DESCRIPTION OF DRAWINGS

[0082] In order to more clearly illustrate the technical solutions of the embodiments of the application, the following will briefly introduce the drawings needed to be used in the embodiment description. Obviously, the drawings in the following description are only some embodiments of the application, and other drawings can also be obtained by those skilled in the art without creative labor.

[0083] Figure 1 It is the flowchart of the satellite and ground any area visible time period calculation method in embodiment 1. DETAILED DESCRIPTION

[0084] In order to make the above-mentioned purposes, features and advantages of the application more obvious and easy to understand, the specific embodiments of the application will be described in detail below with reference to the drawings of the specification.

[0085] In the following description, numerous specific details are set forth in order to provide a thorough understanding of the present application. However, it will be apparent to one skilled in the art that the present application can be practiced without the specific details set forth in this description. In other instances, well-known methods, procedures, components, and circuits have not been described in detail so as not to unnecessarily obscure aspects of the present application.

[0086] It is also noted that, as used herein, "one embodiment" or "an embodiment" means

[0087] Embodiment 1 provides a method for calculating the time period during which a satellite is visible from any region on the ground, comprising:

[0088] Step S1, obtaining TLE data of the satellite and longitude and latitude boundaries of the target ground region;

[0089] Step S2, calculating a sequence of orbital positions of the satellite within a preset time range based on the SGP4 orbit model;

[0090] The orbit position calculation of Step S2 adopts a variable step iteration strategy:

[0091] When the satellite orbit height falls below 400 km, the iteration step is shortened to 1 / 10 of the original step;

[0092] When the satellite elevation angle change rate exceeds 0.5° per second, the step halving mechanism is automatically triggered;

[0093] In Step S2, the variable step iteration strategy is used to calculate the orbit position, comprising:

[0094] a) initialization, formula:

[0095] ,

[0096] wherein, is the satellite geocentric position vector at the epoch time, unit km, is the satellite geocentric velocity vector at the same time, unit , is the starting time of the epoch, unit s, is the two-line element set;

[0097] b) single-step and double-half-step propagation, propagation formula:

[0098] ,

[0099] wherein, is the satellite position vector in km obtained by single-step propagation, is the corresponding velocity vector in , is the initial time of the th iteration in s, is the current step size in s,

[0100] ,

[0101] where, is the satellite position vector in km obtained by double-half-step propagation, is the corresponding velocity vector in ;

[0102] c) The local truncation error is:

[0103] ,

[0104] where, is the th local truncation error estimate;

[0105] d) The altitude, elevation, and rate of change are defined as:

[0106] ,

[0107] where, is the th orbital altitude, denotes the modulus, is the mean radius of the Earth,

[0108] ,

[0109] where, is the th satellite elevation to the target area, is the geocentric coordinate of the geometric center of the target area in km, is the unit vector of the target center,

[0110] ,

[0111] where, is the rate of change of elevation in , is the elevation of the last iteration, is the step size of the last step in s;

[0112] e) The step size is updated as:

[0113] ,

[0114] in, This is the error proportionality coefficient. For error tolerance, To prevent zero denominator constant km,

[0115] ,

[0116] in, For error-driven scaling factor,

[0117] ,

[0118] ,

[0119] in, The scaling factor is driven by height. This is the height step scaling factor. km is the altitude threshold. The scaling factor is driven by the elevation angle. This is the scaling factor for the elevation angle step size. The threshold for the rate of change of elevation angle;

[0120] ,

[0121] in, This is the step size for the next iteration.

[0122] f) Termination and Storage: When Decrease to the set minimum value or Iteration stops if the preset time range is exceeded, and the output is... ,in, For sequence time points;

[0123] Specifically, this iterative strategy immediately reduces the step size in sections with altitudes below 400km or drastic elevation changes, allowing for a full characterization of trajectory details caused by Earth's non-spherical perturbations and atmospheric drag. In high-orbit and gently sloping elevation sections, the step size is widened, significantly reducing computational calls. The double half-step method provides a self-consistent error estimate with a proportionality coefficient. Together with threshold pruning, it suppresses step size oscillations, and the truncation error is always limited to within the tolerance; the whole mechanism can dynamically balance time resolution and computing resources in global coverage missions without manual intervention;

[0124] Step S3: Dynamically correct the satellite's instantaneous position based on the orbital position sequence and atmospheric drag perturbation compensation factor;

[0125] The atmospheric drag perturbation compensation factor in step S3 is generated in the following way:

[0126] Call the NRLMSISE-00 atmospheric density model;

[0127] The drag ratio coefficient is calculated in real time based on the satellite's current orbital altitude and speed.

[0128] Example of correcting the output result by embedding the scaling factor into the SGP4 position;

[0129] In step S3, the atmospheric drag perturbation compensation factor is generated in the following ways:

[0130] Instantaneous atmospheric density was obtained based on NRLMSISE-00:

[0131] ,

[0132] in, Indicates the first The outer atmosphere density of the next iteration, in kg. , The altitude is the orbital height, in kilometers. , These are the latitude and longitude of the nadir point, respectively. For the current moment, The solar flux index (SFU) for that day was 10.7 cm. This is the planetary activity index for the day. This represents the density evaluation function for NRLMSISE-00;

[0133] The formula for calculating the relative velocity after Earth's rotation correction is:

[0134] ,

[0135] in, For relative velocity scalar, units , This is the satellite's geocentric velocity vector, in units of... , The vector of Earth's rotational angular velocity , The satellite's geocentric position is indicated by the unit "km". Represents the cross product of vectors. Represents the magnitude of a vector;

[0136] The area-to-mass ratio is:

[0137] ,

[0138] in, The area-to-mass ratio is expressed in units of... , The force reference area of ​​the satellite, in units of , Satellite mass, unit: kg;

[0139] The aerodynamic drag coefficient is:

[0140] ,

[0141] in, This is the instantaneous aerodynamic drag coefficient. For reference drag coefficient, For high-altitude correction amplitude, It is a natural exponential function. Atmospheric altitude, measured in km;

[0142] The drag ratio coefficient is:

[0143] ,

[0144] in, For the first The resistance ratio coefficient for each iteration, in units ;

[0145] Position correction embedding SGP4:

[0146] ,

[0147] in, This is the corrected satellite position vector. Step size;

[0148] Specifically, the density is obtained in real time by NRLMSISE-00 and coupled with space weather indicators, enabling timely capture of abrupt changes in thinness in high-latitude auroral regions and during the solar expansion period; the relative velocity obtained after rotation correction avoids the omission of geostationary components, making the drag proportionality coefficient completely conform to the aerodynamic model; the exponential decay expansion of the drag coefficient directly correlates scale height with altitude variables, with a continuous transition in the 120km to 500km range, without the need for piecewise constants; the correction term adopts a velocity square dependence form, seamlessly spliced ​​with the SGP4 native spherical harmonic perturbation solution, and embedded only with a first-order Taylor approximation to ensure numerical stability;

[0149] Step S4: Based on the corrected satellite positions, a spatiotemporal clustering algorithm is used to group the transiting satellites; transiting satellites refer to satellites whose nadir points intersect with the target ground area within a preset time range;

[0150] The spatiotemporal clustering algorithm in step S4 specifically includes:

[0151] The spatiotemporal neighborhood is defined as a satellite geocentric angle difference ≤ 5° and a transit time difference ≤ 60 seconds;

[0152] A density-based clustering algorithm is used to group satellites that meet the neighborhood condition into the same group;

[0153] The satellite at the geometric center of each group was selected as the representative satellite;

[0154] Step S5: Calculate the elevation angle variation curve of the target ground area according to the satellite group; if the representative satellite is not visible to the target area, search for satellites with effective inter-satellite links with it, and verify whether the satellite can simultaneously see the target area and the representative satellite.

[0155] After step S5, perform inter-satellite link cooperative coverage determination:

[0156] When the representative satellite is not visible to the target area, search for neighboring satellites that have inter-satellite links with the representative satellite.

[0157] Verify whether adjacent satellites can simultaneously see the target area and the representative satellite.

[0158] If both are visible simultaneously, the visible time periods of the adjacent satellites will be included in the output window;

[0159] The visible criteria for judgment are:

[0160] The instantaneous elevation angle of adjacent satellites to the target area is >5°.

[0161] The line connecting adjacent satellites and the representative satellite is not obstructed by the Earth.

[0162] The sum of the propagation delays of the two links (satellite-region, satellite-satellite) is less than the maximum permissible delay;

[0163] The effectiveness of inter-satellite links requires the following conditions to be met simultaneously:

[0164] The geometric distance between satellites is less than the orbital altitude-related threshold.

[0165] The continuous visible arc of the link outside the Earth's obstruction zone;

[0166] Step S6: Output the visible time window in UTC format based on the continuous time period during which the elevation angle exceeds the set threshold;

[0167] This method also includes a visibility conflict resolution mechanism:

[0168] Detecting overlapping coverage of the same area by multiple satellites within the same time period;

[0169] Communication frequency bands are automatically allocated based on satellite payload type and regional priority.

[0170] The primary satellite identifier is displayed in the output time window;

[0171] The dynamic adjustment method for the elevation angle setting threshold in step S6 is as follows:

[0172] If the target area is located within the polar circle, i.e., latitude ≥ 65°, the threshold is reduced to 3°;

[0173] If the satellite's orbital altitude is below 350 kilometers, the threshold is raised to 7°.

[0174] 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 a time period during which a satellite and a ground are visible to each other, characterized in that, The method comprises the steps of: Step S1, acquiring TLE data of satellites and longitude and latitude boundary of a target ground area; Step S2, calculating orbit position sequence of the satellites in a preset time range based on an SGP4 orbit model; Step S3, dynamically correcting an instantaneous position of the satellites according to the orbit position sequence and an atmospheric resistance perturbation compensation factor; Step S4, grouping the passing satellites by using a space-time clustering algorithm based on the corrected satellite positions; Step S5, calculating an elevation angle change curve of the target ground area according to the satellite groups; If the representative satellite is invisible to the target area, a satellite having an effective inter-satellite link with the representative satellite is searched, and it is verified whether the satellite is simultaneously visible to the target area and the representative satellite; Step S6, outputting a UTC format visible time window according to a continuous time period in which the elevation angle exceeds a set threshold.

2. The method of claim 1, wherein the time period is calculated by the following equation: ###0001### where T is the time period, N is the number of satellites, and Tj is the time period of the i-th satellite. The atmospheric resistance perturbation compensation factor in step S3 is generated by the following method: ​ An NRLMSISE-00 atmospheric density model is called, a resistance proportion coefficient is calculated in real time according to a current orbit height and speed of the satellite, and the proportion coefficient is embedded in a correction operation example of an SGP4 position output result.

3. The method for calculating the visible time period of any region on the ground and satellite as described in claim 2, characterized in that, In step S3, the generation method of the atmospheric resistance perturbation compensation factor comprises: Instantaneous atmospheric density is acquired based on the NRLMSISE-00: ， wherein, denotes the outer atmosphere density at the ith iteration, in kg / m3 , is the orbit altitude, in km, , are the subpoint latitude and longitude, respectively, is the current time, is the daily 10.7 cm solar flux index, in sfu, is the daily planetary activity index, denotes the NRLMSISE-00 density evaluation function; A relative speed after a rotation correction of the earth is calculated, and a formula is: ， where is the relative velocity scalar, in km / s , is the satellite geocentric velocity vector, in km / s , is the Earth rotation angular velocity vector , is the satellite geocentric position vector, in km denotes the vector cross product denotes the vector modulus An area mass ratio is: ， wherein, is the area mass ratio, unit is , is the satellite force reference area, unit is , is the satellite mass, unit kg; An aerodynamic resistance coefficient is: ， wherein, is the instantaneous aerodynamic drag coefficient, is the reference drag coefficient, is the high-altitude correction amplitude, is the natural exponential function, is the atmospheric scale height in km; A resistance proportion coefficient is: ， wherein is the first is the second iteration resistance proportionality coefficient, unit Position correction is embedded in the SGP4: ， wherein, is the modified satellite position vector, is the step size.

4. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 1, characterized in that, The space-time clustering algorithm in step S4 specifically comprises: A space-time adjacent domain is defined as a satellite geocentric angle difference ≤ 5° and a passing time difference ≤ 60 seconds, a density-based clustering algorithm is used to group the satellites meeting the adjacent domain condition into the same group, and a satellite having a geometric center position in each group is selected as a representative satellite.

5. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 1, characterized in that, The inter-satellite link cooperative coverage judgment is performed after step S5: When the representative satellite is invisible to the target area, an adjacent satellite having an inter-satellite link with the representative satellite is searched, it is verified whether the adjacent satellite is simultaneously visible to the target area and the representative satellite, if the adjacent satellite is simultaneously visible, a visible time period of the adjacent satellite is included in the output window; The judgment standard for the simultaneous visibility is: An instantaneous elevation angle of the adjacent satellite to the target area > 5°, a line connecting the adjacent satellite and the representative satellite is not blocked by the earth, and a total propagation time delay of the two links < a maximum allowable time delay.

6. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 5, characterized in that, The effectiveness of the inter-satellite link needs to meet the following conditions: A geometric distance between the satellites is less than an orbit height related threshold, and a continuous visible arc segment of the link is outside an earth blocking area.

7. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 1, characterized in that, In step S6, a dynamic adjustment method of the elevation angle set threshold is as follows: If the target area is located within the polar circle, that is, the latitude ≥ 65°, the threshold is reduced to 3°, and if the satellite orbit height is lower than 350 km, the threshold is increased to 7°.

8. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 1, characterized in that, In step S2, a variable step iteration strategy is used to calculate the orbit position: When the satellite orbit height is reduced to below 400 km, the iteration step is shortened to 1 / 10 of the original step; When the satellite elevation angle change rate exceeds 0.5° per second, a step halving mechanism is automatically triggered.

9. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 8, characterized in that, In step S2, the variable step iteration strategy is used to calculate the orbit position, comprising: a) initialization, a formula is: ， wherein, is the satellite geocentric position vector at the epoch time, in km, is the satellite geocentric velocity vector at the same time, in km / s , is the epoch start time, in s, is the set of two-line elements; b) single-step and double-half-step propagation, a propagation formula is: ， wherein, is the satellite position vector obtained from the single step propagation, in km, is the corresponding velocity vector, in km / s, , is the time of the start of the th iteration, in s, is the current step size, in s, ， wherein, is the satellite position vector in km obtained from the double half-step propagation, is the corresponding velocity vector in km / s, ; c) local truncation error is: ， wherein is the first local truncation error estimate; d) the height, elevation angle and change rate are defined as: ， wherein is the orbital height, denotes the module length, is the average radius of the earth, ， wherein, is the number of satellites, is the elevation angle of the nth satellite to the target area, is the geocentric coordinate of the geometric center of the target area, in km, is the unit vector of the target center, ， wherein is the rate of change of elevation, in units of , is the elevation of the last iteration, is the step size of the last iteration, in units of s; e) performing step length update in the following way: ， wherein, is an error proportionality coefficient, is an error tolerance, is a zero-prevention denominator constant km, ， wherein is an error drive scaling factor, ， ， wherein, is a height drive scaling factor, is a height step scaling coefficient, is a height threshold, is an elevation drive scaling factor, is an elevation step scaling coefficient, is an elevation rate of change threshold; ， wherein is the next iteration step; f) Termination with storage: when falling below a set minimum or exceeding a preset time range, the iteration is stopped and the output wherein, is the sequence time.

10. The method for calculating the visible time period of a satellite and any area on the ground as described in claim 1, characterized in that, Also contains the visibility conflict resolution mechanism: Detect the coverage overlap of multiple satellites in the same period to the same area, automatically assign communication frequency bands according to satellite load type and regional priority, and mark the main satellite identifier in the output time window.

Images