Method for calculating time window of agile satellite earth observation

Abstract

The application discloses a kind of agile satellite earth observation time window calculation method, belongs to satellite earth observation technical field.The present application is in view of the method for calculating agile satellite earth observation time window by orbit recursion in the prior art, it is difficult to realize the problem of long time mission planning in orbit due to large amount of calculation.Six numbers of orbit and flight time are considered according to the orbit of earth observation satellite, and orbit recursion is carried out considering J2 perturbation, the position vector is converted into geographical longitude and latitude considering the influence of earth rotation, and the analytical expression of subsatellite point trajectory about time is established;According to the expression of satellite subsatellite point trajectory, the time when the latitude of ground target is equal is screened out, and whether the relative distance constraint is satisfied at this time is judged, so as to estimate the time when the observation condition is met;The observation attitude angle is calculated according to the relative position vector of satellite pointing target, the maximum roll angle and pitch angle constraint are considered, and the over-the-top time and observation time window are iteratively solved near the estimated time.The application is used to calculate satellite earth observation time window.

Description

TECHNICAL FIELD

[0001] The present application relates to a method for calculating a time window for earth observation of an agile satellite, and belongs to the technical field of satellite earth observation. BACKGROUND

[0002] With the increasing number of on-orbit earth observation satellites, how to realize autonomous mission planning of the satellite becomes a key problem to be solved urgently to reduce the cost of artificial control and management on the ground. In order to meet the needs of users, the earth observation satellite needs to complete the observation tasks of multiple ground targets within a specified time. Due to the limited performance of the on-board computer, the on-orbit satellite needs an efficient and simple method to calculate the observation time window of the ground target, and to plan the mission accordingly.

[0003] The earth observation satellite is equipped with a high-resolution payload, which can be directed to the target for observation by adjusting the attitude when flying over the space above the ground target. Compared with traditional non-agile satellites, agile satellites can quickly perform three-axis maneuvers, greatly increasing the observation time window of the target and planning the observation task time independently. The traditional method for calculating the earth observation time window is to complete the orbit recursion within a specified time and traverse all time nodes to see if the observation constraints are met. Due to the large amount of calculation, it is difficult to realize long-term mission planning on orbit. SUMMARY

[0004] In view of the problem that the existing method for calculating the earth observation time window of an agile satellite by orbit recursion is difficult to realize long-term mission planning on orbit due to the large amount of calculation, the present application provides a method for calculating the earth observation time window of an agile satellite.

[0005] The method for calculating the earth observation time window of an agile satellite according to the present application comprises:

[0006] According to the initial values of the orbit six elements of the earth observation satellite and the initial time of flight, the orbit is recursively calculated considering the J2 perturbation to obtain the orbit six elements of the satellite at time t; and the position vector in the earth-centered inertial coordinate system and the velocity vector in the earth-centered inertial coordinate system of the earth observation satellite at time t are calculated based on the orbit six elements of the satellite at time t;

[0007] The position vector in the earth-centered fixed coordinate system of the earth observation satellite is calculated from the position vector in the earth-centered inertial coordinate system; and the star spot represented by the geographic longitude and the geographic latitude is obtained by converting the position vector in the earth-centered fixed coordinate system, to obtain an analytical expression of the star spot trajectory with respect to time;

[0008] Based on the analytical expression of the star spot trajectory with respect to time, all time instants at which the star spot has the same latitude as the ground target are screened out; and the sines of the distances between the star spot and the ground target corresponding to the screened time instants are calculated; and the time instants corresponding to the sines of the distances less than a set distance threshold are taken as estimated time instants at which the visible time window exists.

[0009] The three-axis unit vector of the earth observation satellite in the orbit coordinate system is calculated from the position vector in the earth-centered inertial coordinate system and the velocity vector in the earth-centered inertial coordinate system of the earth observation satellite at the estimated moment, and the projection of the relative position vector between the sub-satellite point and the ground target in the orbit coordinate system at the estimated moment is calculated; the expressions of the target roll angle and the target pitch angle at the estimated moment are obtained in combination with the attitude angle of the earth observation satellite; the over-the-top moment is obtained by iteratively using a numerical method on the target pitch angle expression, and the target roll angle is determined based on the over-the-top moment; the over-the-top moment corresponding to the maximum roll angle constraint of the target roll angle is taken as the selected over-the-top moment of the observation time window;

[0010] The observation time window corresponding to the selected over-the-top moment is iteratively calculated based on the start and end time constraints of the observation time window, and is taken as the final agile satellite earth observation time window.

[0011] According to the agile satellite earth observation time window calculation method of the application, the satellite orbit six-element number calculation method at the moment t is:

[0012] The satellite orbit six-element number is expressed as [a e i Omega omega M], wherein a is the orbit semi-major axis, e is the eccentricity, i is the inclination, Omega is the ascending node right ascension, omega is the pericenter angle distance, and M is the perigee argument;

[0013] The J2 perturbation term of the non-spherical gravity is considered in the orbit recursion, the moment t0 is taken as the flight initial moment, and the satellite orbit six-element number at the moment t is calculated based on the linearized J2 perturbation equation:

[0014]

[0015] In the formula, a0 is the orbit semi-major axis initial value at the moment t0, e0 is the eccentricity initial value at the moment t0, i0 is the inclination initial value at the moment t0, Omega0 is the ascending node right ascension initial value at the moment t0, omega0 is the pericenter angle distance initial value at the moment t0, M0 is the perigee argument initial value at the moment t0, J2 is the perturbation coefficient, n is the average orbit angular velocity, μ is the earth gravity constant; is the earth equatorial radius, p is the semi-latus rectum, p=a0(1-e 2 ).

[0016] According to the agile satellite earth observation time window calculation method of the application, the position vector of the earth observation satellite in the earth-centered inertial coordinate system at the moment t is expressed as r ECI :

[0017]

[0018] In the formula, f is the true anomaly;

[0019] The velocity vector in the geocentric inertial coordinate system is represented as v ECI :

[0020]

[0021] According to the agile satellite earth observation time window calculation method of the present application, the position vector of the earth observation satellite in the geocentric fixed coordinate system is represented as r ECEF :

[0022] r ECEF = R Z (α Gt )r ECI ,

[0023] wherein R Z is the rotation matrix of the earth observation satellite about the Z axis in the geocentric inertial coordinate system, α Gt is the Greenwich mean sidereal time, α G0 is the initial Greenwich mean sidereal time, and ω is the earth rotation angular velocity.

[0024] According to the agile satellite earth observation time window calculation method of the present application, the method for obtaining the analytical expression of the sub-satellite point trajectory with respect to time is as follows:

[0025] According to r ECEF = [r x , r y , r z ] T the geocentric longitude λ c and the geocentric latitude

[0026] λ c = arctan2(r y , r x ),

[0027] wherein r x is the X axis coordinate of the geocentric fixed coordinate system, r y is the Y axis coordinate of the geocentric fixed coordinate system, and r z is the Z axis coordinate of the geocentric fixed coordinate system;

[0028] The geographic longitude λ d is:

[0029]

[0030] wherein u is the latitude amplitude angle: u = ω + f;

[0031] According to the latitude amplitude angle u, the geocentric latitude is represented as:

[0032]

[0033] the geographic latitude is:

[0034]

[0035] wherein is the first eccentricity of the earth ellipsoid;

[0036] the star below point represented by geographic longitude λ d and geographic latitude at all times obtains an analytical expression of the star below point trajectory with respect to time.

[0037] The agile satellite earth observation time window calculation method according to the present application, the method for screening all times when the star below point is equal to the latitude of the ground target is:

[0038] When the star below point is equal to the latitude of the ground target, the latitude amplitude angle u satisfies:

[0039]

[0040] wherein N is the number of satellite flight orbit circles;

[0041] For a near-circular orbit: f similar to M; then the latitude amplitude angle u and t are in a linear relationship, and all times when the star below point is equal to the latitude of the ground target are calculated by the linearized J2 perturbation equation.

[0042] The agile satellite earth observation time window calculation method according to the present application, the determination method for the estimation time when the visible time window exists is:

[0043] Calculate the secant distance d:

[0044]

[0045] wherein is the difference between the current geographic latitude of the star below point and the ground target, is the current geographic latitude of the star below point, is the current geographic latitude of the ground target; and Δλ is the difference between the current geographic longitude of the star below point and the ground target;

[0046] If the secant distance d is less than a set distance threshold, the current time is taken as the estimation time.

[0047] The present application has the following advantages: the method of the present application uses analytical approximate model instead of orbit integral, uses constraint screening instead of traversal judgment, and has extremely small amount of calculation, and can be applied to long-time on-orbit satellite earth observation time window calculation, and provides sufficient preparation for autonomous task planning. The method of the present application describes the subsatellite point trajectory as a time-dependent analytical expression based on the linearized J2 perturbation model and considering the influence of the earth rotation. The time points with observation conditions are screened according to the latitude and longitude difference and the relative distance, and the observation time window is obtained by iterative calculation in the vicinity thereof. BRIEF DESCRIPTION OF DRAWINGS

[0048] Figure 1 is a flow chart of the agile satellite earth observation time window calculation method of the present application;

[0049] Figure 2 is a schematic diagram of the observation time window;

[0050] Figure 3 is a simulation diagram of the observation time window in a specific embodiment. DETAILED DESCRIPTION

[0051] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative work fall within the scope of protection of the present application.

[0052] It should be noted that the embodiments in the present application and the features in the embodiments can be combined with each other without conflict.

[0053] The present application will be further described below with reference to the drawings and specific embodiments, but is not limited by the present application.

[0054] DETAILED DESCRIPTION Figure 1 and Figure 2 The present application provides an agile satellite earth observation time window calculation method, which comprises the following steps:

[0055] According to the initial value of the orbit six elements of the earth observation satellite and the initial time of flight, the orbit is recursively calculated considering the J2 perturbation to obtain the orbit six elements of the satellite at time t; and the position vector in the earth-centered inertial coordinate system and the velocity vector in the earth-centered inertial coordinate system of the earth observation satellite at time t are calculated based on the orbit six elements of the satellite at time t;

[0056] The position vector in the earth-fixed coordinate system of the satellite is calculated from the position vector in the earth-centered inertial coordinate system; the position vector is converted into geographical longitude and latitude by considering the influence of the earth rotation, i.e. the position vector in the earth-fixed coordinate system is converted into the sub-satellite point represented by the geographical longitude and latitude, to obtain an analytical expression of the sub-satellite point trajectory with respect to time;

[0057] Based on the analytical expression of the sub-satellite point trajectory with respect to time, all time instants at which the sub-satellite point has the same latitude as the ground target are screened out; the corresponding sub-satellite point and the ground target are calculated based on the screened time instants to obtain the subtangent distance; the time instant corresponding to the subtangent distance less than a set distance threshold is taken as an estimated time instant at which a visible time window exists; the estimated time instant is a time instant at which the observation condition is met;

[0058] The three-axis unit vectors of the satellite in the orbit coordinate system are calculated from the position vector and the velocity vector in the earth-centered inertial coordinate system of the satellite at the estimated time instant, and the projection of the relative position vector between the sub-satellite point and the ground target at the estimated time instant in the orbit coordinate system is obtained; the expressions of the target roll angle and the target pitch angle are obtained in combination with the attitude angle of the satellite; the over-the-top time instant is obtained by iteratively using a numerical method for the expression of the target pitch angle, and the target roll angle is determined based on the over-the-top time instant; the over-the-top time instant corresponding to the target roll angle satisfying the maximum roll angle constraint is taken as a selected over-the-top time instant at which an observation time window exists;

[0059] The observation time window corresponding to the selected over-the-top time instant is iteratively calculated based on the start and end time instant constraints of the observation time window, and is taken as the final agile satellite observation time window on the ground.

[0060] Further, the calculation method of the satellite orbit six elements at the time instant t is as follows:

[0061] The satellite orbit six elements are expressed as [a e i Ω ω M], which are used to describe the orbit state, wherein a is the orbit semi-major axis, e is the eccentricity, i is the inclination, Ω is the ascending node right ascension, ω is the pericenter argument, and M is the mean anomaly;

[0062] In the orbit recursion, the influence of the non-spherical gravity J2 perturbation term is considered, t0 is taken as the initial time instant of flight, and the satellite orbit six elements at the time instant t are calculated based on the linearized J2 perturbation equation (Gauss perturbation equation):

[0063]

[0064] In the formula, a0 is the initial value of the orbit semi-major axis at the time instant t0, e0 is the initial value of the eccentricity at the time instant t0, i0 is the initial value of the inclination at the time instant t0, Ω0 is the initial value of the ascending node right ascension at the time instant t0, ω0 is the initial value of the pericenter argument at the time instant t0, M0 is the initial value of the mean anomaly at the time instant t0, and J2 is the perturbation coefficient, J2 = 1.082627 × 10-3 , n is the average orbital angular velocity, μ is the Earth's gravitational constant, μ=3.986×10 14 m 3 / s -2 ; is the Earth's equatorial radius, p is the semi-diameter, p=a0(1-e 2 ).

[0065] In this embodiment, the position vector of the Earth observation satellite at time t in the Earth Centered Inertial Coordinate System (ECI) is expressed as r ECI :

[0066]

[0067] Where f is the true anomaly, which is equivalent to the mean anomaly;

[0068] The velocity vector in the geocentric inertial coordinate system is represented by v ECI :

[0069]

[0070] If only the effect of the Earth's rotation is considered, the position vector of the satellite in the Earth-centered fixed coordinate system (ECEF) can be calculated from the ECI position vector. The position vector of the Earth observation satellite in the Earth-centered fixed coordinate system is expressed as r ECEF :

[0071] r ECEF =R Z (α Gt )r ECI ,

[0072] Where R Z is the rotation matrix of the Earth observation satellite around the Z axis in the Earth-centered inertial coordinate system, α Gt is Greenwich Mean Sidereal Time (GMST), α G0 is the initial Greenwich mean sidereal time, which can be calculated from Coordinated Universal Time (UTC); is the angular velocity of the Earth's rotation.

[0073] Furthermore, the analytical expression of the subsatellite point trajectory with respect to time is obtained as follows:

[0074] According to r ECEF =[r x ,r y ,r z ] T Calculate the geocentric longitude λ c and geocentric latitude

[0075] c = arctan2(r y , r x ),

[0076] where r x is the X-axis coordinate of the geocentric fixed coordinate system, r y is the Y-axis coordinate of the geocentric fixed coordinate system, and r z is the Z-axis coordinate of the geocentric fixed coordinate system.

[0077] The geographic longitude λ d can be expressed by the orbital elements as follows:

[0078]

[0079] where u is the latitude amplitude angle: u = ω + f.

[0080] According to the latitude amplitude angle u, the geocentric latitude φ can be expressed as:

[0081]

[0082] The geographic latitude φ is:

[0083]

[0084] where e is the first eccentricity of the Earth ellipsoid, and

[0085] In summary, the geographic longitude λ d and the geographic latitude φ of the satellite at any time, i.e., the subspace point trajectory, can be obtained from the initial orbital state and the flight time. The analytical expression of the subspace point trajectory with respect to time is obtained from the subspace point represented by the geographic longitude λ d and the geographic latitude φ .

[0086] Furthermore, the method for screening all the times when the subspace point is equal to the latitude of the ground target is as follows:

[0087] First, select the time when the subspace point trajectory of the satellite is close to the latitude of the ground target. When the subspace point is equal to the latitude of the ground target, the latitude amplitude angle u satisfies:

[0088]

[0089] where N is the number of the satellite flight orbit;

[0090] For a near-circular orbit: f ≈ M; then the latitude amplitude angle u and t are linearly related, and all the times when the subspace point is equal to the latitude of the ground target are calculated from the linearized J2 perturbation equation.

[0091] The method for determining the estimation time with the visible time window is as follows:

[0092] It is determined whether the distance between the subspace point and the ground target is less than a set threshold at the estimation time. The distance d is calculated by the spherical triangle formula:

[0093]

[0094] In the formula, Δφ is the difference between the current geographic latitude of the subspace point and the ground target, is the current geographic latitude of the subspace point, is the current geographic latitude of the ground target, and Δλ is the difference between the current geographic longitude of the subspace point and the ground target.

[0095] If the distance d is less than the set distance threshold, the current time is taken as the estimation time, that is, it is considered that the visible time window exists near the current time.

[0096] In the embodiment, the satellite attitude is defined in the orbit coordinate system (VVLH), and the three-axis unit vectors are represented as [i j k] T which can be obtained from the position and velocity vectors in the ECI system; in the formula, i is the X-axis unit vector of the orbit coordinate system, j is the Y-axis unit vector of the orbit coordinate system, and k is the Z-axis unit vector of the orbit coordinate system:

[0097] i = j × k,

[0098] The relative position vector δr of the earth observation satellite pointing to the ground target in the ECI coordinate system is: ECI

[0099] δr ECI = r T-ECI - r ECI ,

[0100] In the formula, r T-ECI is the position vector of the ground target in the ECI coordinate system;

[0101] The projection of the relative position vector of the subspace point and the ground target at the estimation time in the orbit coordinate system is δr VVLH :

[0102] δr VVLH = [i j k] T · δr ECI .

[0103] ​​The position vector of the ground target in the Earth-Centered Inertial (ECI) coordinate system needs to be converted from the Earth-Centered Earth-Fixed (ECEF) coordinate system. For the subsequent accurate calculation of the observation time window, the conversion process needs to fully consider the effects of orbit perturbation, precession, nutation, polar motion, etc., rather than using the above analytical model.

[0104] Further, the determination method of the selected over-the-top moment is:

[0105] The attitude angle of the Earth observation satellite includes roll angle φ, pitch angle θ and yaw angle ψ. Assuming ψ = 0, the attitude rotation matrix R of the Earth observation satellite in the orbit coordinate system is yx = R y (θ)R x (φ), R x is the rotation matrix of the Earth observation satellite around the X-axis in the orbit coordinate system, and R y is the rotation matrix of the Earth observation satellite around the Y-axis in the orbit coordinate system.

[0106] Assuming that the initial attitude of the Earth observation satellite payload points to the Z-axis of the orbit coordinate system, then:

[0107]

[0108] The expression of the estimated moment target roll angle and target pitch angle is obtained:

[0109]

[0110] In the formula, δr x is the X-axis coordinate of the projection δr VVLH at the estimated moment, δr y is the Y-axis coordinate of the projection δr VVLH at the estimated moment, and δr z is the Z-axis coordinate of the projection δr VVLH at the estimated moment.

[0111] The over-the-top moment is defined as the moment when the pitch angle θ = 0, at which time only roll maneuver can be performed to observe the ground target. In the neighborhood of the estimated moment, the over-the-top moment t over is obtained by using a numerical method to iteratively obtain the expression of the pitch angle θ; the over-the-top moment t over is substituted into the expression of the target roll angle to calculate the target roll angle; and the numerical method can be the secant method.

[0112] If the target roll angle satisfies the maximum roll angle constraint, the current over-the-top moment t over is selected as the selected over-the-top moment with the observation time window.

[0113] The method for calculating the Earth observation time window of the agile satellite is:

[0114] The start time tstart and the end time t end Satisfy the maximum pitch angle constraint of the Earth observation satellite:

[0115] θ(t start )=θ max ,θ(t end )=-θ max ,

[0116] θ max is the maximum pitch angle of the Earth observation satellite;

[0117] Set the starting time t start and the end time t end and the overhead moment t over The time interval is Δt:

[0118]

[0119] In t over ±Δt neighborhood, the final starting time t is determined iteratively using a numerical method for the pitch angle θ start and the final termination time t end , the final [t start ,t end ] as the final agile satellite earth observation time window. The numerical method can be the secant method.

[0120] Example:

[0121] Assume that the initial orbital numbers of the agile earth observation satellite at UTC time 2024 / 7 / 1 4:00:00 are The maximum pitch and roll angles are both 45°, the geographical latitude and longitude of the ground target are (116.388°E, 39.9289°N), and the total mission duration is 3 days.

[0122] After analysis and calculation, there are 6 earth observation time windows in the total mission duration, and the start and end times and durations are shown in Table 1. The simulation results of the observation time windows are shown in Table 1. Figure 3 As shown, Figure 3 The trajectory of the subsatellite point is represented by the blue dotted line, and the trajectory within the observation time window is represented by the green solid line. They are all distributed near the ground target represented by the star. Figure 3 It can be seen that this method can correctly calculate the earth observation time window.

[0123] Table 1 Calculation results of observation time window

[0124] Time window Start time (UTC) End time (UTC) Duration (seconds) 1 2024 / 7 / 1 8:41:22 2024 / 7 / 1 8:44:49 208 2 2024 / 7 / 2 2:06:21 2024 / 7 / 2 2:09:14 173 3 2024 / 7 / 2 8:51:41 2024 / 7 / 2 8:55:18 218 4 2024 / 7 / 3 2:16:48 2024 / 7 / 3 2:19:45 176 5 2024 / 7 / 3 7:21:22 2024 / 7 / 3 7:24:41 198 6 2024 / 7 / 4 0:46:20 2024 / 7 / 4 0:49:25 185

[0125] While the application has been described with reference to particular embodiments thereof, it is to be understood that these embodiments are merely illustrative of the principles and applications of the present application. It will be apparent to those skilled in the art that numerous modifications can be made within the scope of the present application as defined by the appended claims. It is intended that all such modification fall within the spirit and scope of the present application. It will be understood that the features described in connection with one embodiment can be used in connection with another embodiment.

Claims

1. A method for calculating the time window of an agile satellite earth observation, characterized in that: include: According to the initial values ​​of the six orbital elements of the Earth observation satellite and the initial flight time, orbit recursion is performed considering the J2 perturbation to obtain the six orbital elements of the satellite at time t; based on the six orbital elements of the satellite at time t, the position vector and velocity vector of the Earth observation satellite in the Earth-centered inertial coordinate system at time t are calculated; The position vector of the Earth observation satellite in the geocentric fixed coordinate system is calculated from the position vector in the geocentric inertial coordinate system; the position vector in the geocentric fixed coordinate system is converted into a sub-satellite point represented by geographic longitude and geographic latitude, and an analytical expression of the sub-satellite point trajectory with respect to time is obtained; Based on the analytical expression of the sub-satellite point trajectory with respect to time, all moments when the sub-satellite point and the ground target have the same latitude are selected. The haversine distance between the sub-satellite point and the ground target is calculated based on all selected moments. The moment when the haversine distance is less than the set distance threshold is used as the estimated moment when the visible time window exists. The three-axis unit vectors of the Earth observation satellite in the orbital coordinate system are calculated from the position vector and velocity vector of the Earth observation satellite in the Earth-centered inertial coordinate system at the estimated time. The projection of the relative position vector between the sub-satellite point and the ground target in the orbital coordinate system at the estimated time is then calculated. The expressions of the target roll angle and target pitch angle at the estimated time are obtained by combining the attitude angle of the Earth observation satellite. The target pitch angle expression is iterated by numerical method to obtain the over-the-top moment, and then the target roll angle is determined based on the over-the-top moment. The over-the-top moment corresponding to when the target roll angle satisfies the maximum roll angle constraint is taken as the selected over-the-top moment in the observation time window; The observation time window corresponding to the selected overpass time is obtained by iterative calculation based on the start and end time constraints of the observation time window, and is used as the final agile satellite earth observation time window; The method for determining the estimated time when there is a visible time window is: Calculate the haversine distance d: In the formula is the current geographic latitude difference between the sub-satellite point and the ground target, is the current geographic latitude of the subsatellite point, is the current geographical latitude of the ground target; Δλ is the current geographic longitude difference between the sub-satellite point and the ground target; is the equatorial radius of the Earth; If the haversine distance d is less than the set distance threshold, the current time is used as the estimated time.

2. The method for calculating the time window for agile satellite earth observation according to claim 1, characterized in that: The calculation method of the six satellite orbit numbers at time t is: The six satellite orbit numbers are expressed as [ae iΩωM], where a is the semi-major axis of the orbit, e is the eccentricity, i is the inclination, Ω is the right ascension of the ascending node, ω is the angular distance of the periapsis, and M is the mean anomaly. In the orbit recursion, the influence of the non-spherical gravity J2 perturbation term is considered, and the time t0 is taken as the initial time of the flight. The six satellite orbital elements at time t are calculated based on the linearized J2 perturbation equation: Where a0 is the initial value of the orbital semi-major axis at time t0, e0 is the initial value of the eccentricity at time t0, i0 is the initial value of the inclination at time t0, Ω0 is the initial value of the right ascension of the ascending node at time t0, ω0 is the initial value of the angular distance of the periapsis at time t0, and M0 is the initial value of the mean anomaly at time t0; J2 is the perturbation coefficient, n is the mean orbital angular velocity, μ is the earth's gravitational constant; p is the semi-diameter, p=a0(1-e 2 ).

3. The method for calculating the time window for agile satellite earth observation according to claim 2, characterized in that: The position vector of the Earth observation satellite in the geocentric inertial coordinate system at time t is expressed as r ECI : Where f is the true anomaly; The velocity vector in the geocentric inertial coordinate system is represented by v ECI :

4. The method for calculating the time window for agile satellite earth observation according to claim 3, characterized in that: The position vector of the Earth observation satellite in the geocentric fixed coordinate system is expressed as r ECEF : r ECEF ＝R Z (α Gt )r ECI ， Where R Z is the rotation matrix of the Earth observation satellite around the Z axis in the Earth-centered inertial coordinate system, α Gt is Greenwich Mean Sidereal Time, α G0 is the initial Greenwich mean sidereal time, is the angular velocity of the Earth's rotation.

5. The method for calculating the time window for agile satellite earth observation according to claim 4, characterized in that: The analytical expression of the subsatellite point trajectory with respect to time is obtained as follows: According to r ECEF =[r x ,r y ,r z ] T Calculate the geocentric longitude λ c and geocentric latitude Where r x is the X-axis coordinate of the geocentric fixed coordinate system, r y is the Y-axis coordinate of the geocentric fixed coordinate system, r z is the Z-axis coordinate of the geocentric fixed coordinate system; Geographic longitude λ d for: Where u is the latitude angle: u = ω + f; According to the latitude argument u, the geocentric latitude Expressed as: Geographic latitude for: In the formula is the first eccentricity of the Earth ellipsoid; By the geographical longitude λ at all times d and geographic latitude The analytical expression of the sub-satellite point trajectory with respect to time is obtained by denoting the sub-satellite point.

6. The method for calculating the time window for agile satellite earth observation according to claim 5, characterized in that: The method to filter out all the moments when the sub-satellite point is equal to the ground target latitude is: When the latitude of the sub-satellite point is equal to that of the ground target, the latitude argument u satisfies: Where N is the number of satellite orbits; For near-circular orbits: f≈M; then the latitude argument u is linearly related to t, and the linearized J2 perturbation equation is used to calculate all the times when the sub-satellite point and the ground target latitude are equal.

7. The method for calculating the time window for agile satellite earth observation according to claim 6, characterized in that: The three-axis unit vector is expressed as [ijk] T , where i is the X-axis unit vector of the orbital coordinate system, j is the Y-axis unit vector of the orbital coordinate system, and k is the Z-axis unit vector of the orbital coordinate system: The relative position vector δr of the Earth observation satellite pointing to the ground target in the Earth-centered inertial coordinate system ECI for: δr ECI =r T-ECI -r ECI , Where r T-ECI is the position vector of the ground target in the geocentric inertial coordinate system; The projection of the relative position vector between the sub-satellite point and the ground target in the orbital coordinate system at the estimated time is δr VVLH : δr VVLH =[i j k] T ·δr ECI 。 8. The method for calculating the time window for agile satellite earth observation according to claim 7, characterized in that: The method for determining the time of passing the top is: The attitude angle of the Earth observation satellite includes the roll angle φ, pitch angle θ and yaw angle ψ. Assume ψ = 0, and the attitude rotation matrix R in the orbital coordinate system is yx =R y (θ)R x (φ), R x is the rotation matrix of the Earth observation satellite around the X axis in the orbital coordinate system, R y is the rotation matrix of the Earth observation satellite around the Y axis in the orbital coordinate system; Assuming that the initial attitude of the Earth observation satellite payload points to the Z axis of the orbital coordinate system, then: The expressions of target roll angle and target pitch angle at the estimation moment are obtained: Where δr x is the estimated time projection δr VVLH X-axis coordinate, δr y is the estimated time projection δr VVLH Y-axis coordinate, δr z is the estimated time projection δr VVLH The Z-axis coordinate of The over-the-top moment is defined as the moment when the pitch angle θ = 0. In the neighborhood of the estimated time, the expression of the pitch angle θ is iterated by a numerical method to obtain the over-the-top moment t over ; Set the top time t over Substitute the target roll angle expression into the target roll angle and calculate the target roll angle; If the target roll angle satisfies the maximum roll angle constraint, the current over-the-top moment t over The selected passing moment as the existence observation time window.

9. The method for calculating the time window for agile satellite earth observation according to claim 8, characterized in that: The method for calculating the time window of agile satellite earth observation is: Let the starting time of the observation time window be t start and the end time t end Satisfy the maximum pitch angle constraint of the Earth observation satellite: θ(t start )=θ max ,θ(t end )=-θ max , θ max is the maximum pitch angle of the Earth observation satellite; Set the starting time t start and the end time t end and the overhead moment t over The time interval is Δt: In t over ±Δt neighborhood, the final starting time t is determined iteratively using a numerical method for the pitch angle θ start and the final termination time t end , the final [t start ,t end ] as the ultimate agile satellite earth observation time window.

Images