Satellite attitude guidance method and system for observing target star from same visual angle

Abstract

The application provides a satellite attitude guiding method and system for observing a target star in the same visual angle, comprising the following steps: calculating the position R m0 and velocity V m0 of the target star at the time of ephemeris data; performing numerical integration to obtain the position R mt0 and velocity V mt0 of the target star at the starting time of a task; calculating the orbit root number of the target star at the starting time of the task according to the position R mt0 and velocity V mt0 ; calculating the real-time position R m and velocity V m of the target star based on the orbit root number; and calculating the guiding attitude of the satellite according to the position R m and velocity V m . The application can improve the same visual angle observation tracking performance of the satellite and the target star by autonomous recursion of orbit information, and has practical significance for dynamic target tracking imaging tasks of the satellite by taking into account the on-board computing resources.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of satellite attitude guidance, and particularly relates to a satellite attitude guidance method and system for observing a target star from the same visual angle. BACKGROUND

[0002] For some advanced optical imaging technologies suitable for medium-high orbit satellites, when on-orbit verification is carried out, a low-orbit demonstration verification satellite is often used to verify the actual application of subsequent medium-high orbit satellites. This requires the low-orbit satellite to obtain observation data from the same visual angle as the medium-high orbit satellite, thereby involving an attitude guidance algorithm for observing a satellite and other target stars from the same visual angle. The low-orbit satellite has the characteristics of low launch cost, small launch difficulty, and high spatial resolution. Under the same angular resolution, the medium-high orbit satellite has a wider field of view.

[0003] In some satellite cross-radiometric calibration task scenarios, two satellites need to observe the same target at the same time, and the sensor of the satellite to be calibrated is accurately calibrated by the satellite with better calibration results. This process requires high-precision attitude guidance algorithm intervention on-orbit.

[0004] When a satellite observes a target star from the same visual angle, the motion trajectory of the target star needs to be obtained in advance. Due to the uncertainty of the task, it is difficult to track and point by looking up a large amount of target star position data on the ground, so the orbit information needs to be autonomously propagated during the task according to the orbit parameter information of the target star at a certain time. Existing methods give some orbit propagation calculation methods, but these methods are no longer applicable in the task of observing the target star from the same visual angle on the satellite in real time, and the main reasons are as follows:

[0005] Patent document CN104459732A discloses a desingularization method suitable for satellite orbit parameter propagation calculation. This scheme uses a glonass satellite position solution algorithm, but uses a numerical integration algorithm and relies on a glonass receiver.

[0006] Patent document CN103268407A discloses an orbit data interpolation method based on Lagrange interpolation and Kalman filtering. This scheme introduces an orbit data interpolation algorithm. Due to the variability of the task period, the interpolation algorithm only has a small orbit position deviation within a limited period. When the task changes beyond the original period, extrapolation based on the interpolation data will rapidly increase the orbit position deviation, affecting the tracking and pointing accuracy.

[0007] Patent document CN103995800A discloses a kind of on-board autonomous orbit extrapolation method suitable for circular orbit satellite, and the scheme introduces a kind of orbit extrapolation algorithm of circular orbit satellite, but the recursive model of this algorithm is relatively simple, not suitable for long time high-precision orbit extrapolation.

[0008] Patent document CN114002713A discloses a satellite orbit parameter extrapolation processing prediction system, and the scheme introduces a kind of satellite orbit parameter prediction method, which can reduce the amount of calculation by using normalization processing, but the amount of calculation cannot be met for the scene of real-time calculation of guidance attitude during the task.

[0009] Optical Precision Engineering, 2016, No. 10 discloses a low-orbit satellite orbit prediction algorithm based on orbit elements, and the scheme introduces an algorithm for predicting satellite orbit using elliptic curve, but the partial differential of coefficient needs to be calculated in the solving process, which is not suitable for on-board computer operation.

[0010] When satellite observes the target star with the same visual angle, the motion trajectory of the target star needs to be obtained in advance. Due to the uncertainty of the task, it is difficult to track and point by looking up the vast amount of target star position data on the ground, and the on-board needs to complete the autonomous extrapolation of orbit information during the task according to the orbit parameter information of the target star at a certain time.

[0011] As described above, the existing method gives some orbit extrapolation calculation methods, but these methods are no longer applicable in the task of real-time observation of the target star with the same visual angle on the satellite, and this problem needs to be solved urgently. SUMMARY

[0012] In view of the defects in the prior art, the purpose of the present application is to provide a satellite attitude guidance method and system for observing the target star with the same visual angle.

[0013] According to the satellite attitude guidance method for observing the target star with the same visual angle provided by the present application, the following steps are included:

[0014] Step S1: calculating the position R m0 and velocity V m0 of the target star at the ephemeris data time;

[0015] Step S2: numerically integrating to solve the position R mt0 and velocity V mt0 of the target star at the task start time;

[0016] Step S3: according to the position R mt0 and velocity V mt0 , the orbit plane element of the target star at the task start time is solved;

[0017] Step S4: based on the orbital elements, the real-time position R of the target star is calculated m and velocity V m ;

[0018] Step S5: according to the position R m and velocity V m , the guiding attitude of the satellite is calculated.

[0019] Preferably, in the step S1, the position R m0 is expressed as:

[0020] R m0 = Q·r p

[0021] wherein r p represents the position vector in the perifocal coordinate system, and Q represents the conversion matrix from the perifocal coordinate system to the inertial system;

[0022] The velocity V m0 is expressed as:

[0023] V m0 = Q·v p

[0024] The Q is expressed as:

[0025]

[0026] wherein the superscript T represents the matrix transpose symbol; R X is the X-axis rotation matrix function; R Z is the Z-axis rotation matrix function; Ω ms represents the right ascension of the ascending node; i ms represents the orbital inclination; ω ms represents the argument of perigee;

[0027] The R X , i.e. the X-axis rotation matrix function, is expressed as:

[0028]

[0029] The R Z , i.e. the Z-axis rotation matrix function, is expressed as:

[0030]

[0031] wherein α represents an angle;

[0032] The r p is expressed as:

[0033]

[0034] Among them, f ms Indicates the true nearest point angle; e ms The eccentricity is represented by μ; the gravitational constant of Earth is represented by h. ms Indicates the angular momentum of the target star;

[0035] The v p The mathematical expression is:

[0036]

[0037] Among them, v p Represents the velocity vector in the near-focal coordinate system;

[0038] The h ms The mathematical expression is:

[0039]

[0040] Among them, a ms Indicates the semi-major axis of the track;

[0041] In step S2, the numerical integration solution step includes:

[0042] Step A1: Calculate the integration step size from the start time of the task;

[0043] Step A2: Set the maximum single-step integration step size h max Then, the integration step size is decomposed.

[0044] Step A3: Recursively extrapolate the target star's orbit from the ephemeris data time to the mission start time, and calculate the target star's position R at the mission start time. mt0 With speed V mt0 ;

[0045] In step A1, the mathematical expression for the integration step size is:

[0046] h = t0 - t s

[0047] Where h represents the integration step size; t0 represents the task start time; t s Indicates the time of ephemeris data;

[0048] In step A2, the mathematical expression for the decomposition process is:

[0049] h = N·h max +h0

[0050] Where N represents a natural number, and 0 ≤ h < h max h0 represents the initial integration step size;

[0051] In the step A3, the N+1 times value calculation method is called to obtain the position R of the target star at the starting time of the task mt0 With the velocity V mt0 ;

[0052] The numerical calculation method is mathematically expressed as:

[0053]

[0054] Wherein, f represents a function; y0 represents an integral initial value; the symbol Δ represents an increment symbol; the differential equation group f(p) is a 6×1 column vector, p is the item in the brackets among k1-k4, and p is set as

[0055] The mathematical expression of f(p) is:

[0056]

[0057] Wherein, R e is the radius of the Earth equator, and J2 is the 2nd order banding harmonic coefficient of the non-spherical gravity perturbation term of the Earth;

[0058] In the step A3, y0 is the integral initial value used for the first time, and the mathematical expression is:

[0059] y0=[R m0 ; V m0 ]

[0060] In the first N times, the single-step maximum integral step h max is used, that is, Δ=h max , and the output y1 of each time is used as the integral initial value y0 of the next time; in the last time, Δ=h0 is set, and the output y1 is calculated to obtain the position R mt0 and the velocity V mt0 of the target star at the starting time of the task.

[0061] The mathematical expression of R mt0 is:

[0062] R mt0 =[y1(1); y1(2); y1(3)]

[0063] The mathematical expression of the velocity V mt0 is:

[0064] V mt0 =[y1(4); y1(5); y1(6)]

[0065] Wherein, y1(1), y1(2), y1(3), y1(4), y1(5) and y1(6) represent six components of the output y1.

[0066] Preferably, in the step S3, the position R mt0 with the speed V mt0 , the orbital instantaneous root number and the amplitude angle information are obtained; according to the orbital instantaneous root number and the amplitude angle information, the orbital plane root number of the target star at the mission starting moment is obtained;

[0067] The orbital instantaneous root number comprises e', a', i', Ω' and ω';

[0068] The mathematical expression of the e' is:

[0069] e' =‖E‖

[0070] wherein e' represents the eccentricity at the mission starting moment of the target star, and is the orbital instantaneous root number; E represents the eccentricity vector at the mission starting moment of the target star;

[0071] The mathematical expression of the E is:

[0072]

[0073] wherein E represents the eccentricity vector at the mission starting moment of the target star; the symbol ||| represents the modulus operator; H represents the angular momentum vector at the mission starting moment of the target star;

[0074] The mathematical expression of the H is:

[0075] H = R mt0 × V mt0

[0076] The mathematical expression of the a' is:

[0077]

[0078] wherein a' represents the orbit semi-major axis at the mission starting moment of the target star;

[0079] The mathematical expression of the i' is:

[0080]

[0081] wherein i' represents the orbit inclination at the mission starting moment of the target star;

[0082] The mathematical expression of the Ω' is:

[0083]

[0084] wherein Ω' represents the right ascension of the ascending node at the mission starting moment of the target star;

[0085] The mathematical expression of the ω' is:

[0086]

[0087] wherein ω' represents the argument of perigee at the starting time of the target star mission;

[0088] The argument information comprises: u', φ' and λ';

[0089]

[0090] wherein u' represents the argument of true latitude at the starting time of the target star mission;

[0091]

[0092] wherein φ' represents the argument of declination at the starting time of the target star mission;

[0093] λ' = φ' - e' sin (φ' - ω')

[0094] wherein λ' represents the argument of latitude at the starting time of the target star mission.

[0095] Preferably, in the step S3, the orbital elements comprise: a mt0 , i mt0 , Ω mt0 , ω mt0 , e mt0 and M mt0 .

[0096] The mathematical expression of the a mt0 is:

[0097]

[0098] wherein a mt0 represents the semi-major axis at the starting time of the target star mission; x5 represents the intermediate variable five; and x6 represents the intermediate variable six.

[0099] The mathematical expression of the i mt0 is:

[0100] i mt0 = i' + (1 / 6) x2Ω1 sin i'

[0101] wherein i mt0 represents the inclination at the starting time of the target star mission; and x2 represents the intermediate variable two.

[0102] The mathematical expression of the Ω mt0 is:

[0103] Ω mt0 = Ω' - Ω1 (x1 - x3 / 6)

[0104] wherein Ωmt0 Ω1 represents the right ascension of the ascending node at the beginning of the target star mission; x1 is an intermediate variable one; x3 is an intermediate variable three;

[0105] The mathematical expression of the Ω1 is:

[0106]

[0107] The mathematical expression of the ω mt0 is:

[0108] ω mt0 = arctan(-η / ξ)

[0109] The mathematical expression of the ω mt0 is:

[0110] The mathematical expression of the η is:

[0111] η = η' - (-x 10 -(x4+x7)ξ' + x8η')

[0112] The mathematical expression of the η is: 10 x10 is an intermediate variable ten; x4 is an intermediate variable four; x7 is an intermediate variable seven; x8 is an intermediate variable eight;

[0113] The mathematical expression of the e mt0 is:

[0114]

[0115] The mathematical expression of the e mt0 is:

[0116] The mathematical expression of the ξ is:

[0117] ξ = ξ' - (x9 + (x4+x7)η' + x8ξ')

[0118] x8 is an intermediate variable eight; x9 is an intermediate variable nine;

[0119] The mathematical expression of the M mt0 is:

[0120] M mt0 = λ-ω mt0

[0121] The mathematical expression of the M mt0 is:

[0122] The mathematical expression of the λ is:

[0123] λ = λ' - (x4+ε(x7e'2 +x9η′+x 10 ξ′))

[0124] where x9 is an intermediate variable nine.

[0125] Preferably, in said step S4, said R m is mathematically expressed as:

[0126]

[0127] where E mt is the eccentric anomaly;

[0128] said V m is mathematically expressed as:

[0129]

[0130] said P m is mathematically expressed as:

[0131]

[0132] where P m is an intermediate variable;

[0133] said Q m is mathematically expressed as:

[0134]

[0135] where Q m is an intermediate variable;

[0136] In said step S5, the guiding attitude of said satellite, i.e. the roll angle the pitch angle θ and the yaw angle ψ of said satellite in the inertial frame are mathematically expressed as:

[0137]

[0138] θ = arcsin (i z (1))

[0139]

[0140] where i x , i y , i z are mathematically expressed as:

[0141]

[0142] i x = i y x i z

[0143] wherein, i x , i y , i z respectively represent the direction vectors of the three axes of the satellite body in the inertial system; R s and V s are respectively the position vector and the velocity vector of the satellite in the J2000 inertial coordinate system.

[0144] According to the satellite attitude guiding system for observing the target star in the same visual angle provided by the application, the system comprises:

[0145] Module M1: calculating the position R m0 and the velocity V m0 of the target star at the ephemeris data moment;

[0146] Module M2: numerically integrating to solve the position R mt0 and the velocity V mt0 of the target star at the task starting moment;

[0147] Module M3: according to the position R mt0 and the velocity V mt0 , solving the orbital plane root number of the target star at the task starting moment;

[0148] Module M4: based on the orbital plane root number, solving the real-time position R m and the velocity V m of the target star;

[0149] Module M5: according to the position R m and the velocity V m , calculating the guiding attitude of the satellite.

[0150] Preferably, in the module M1, the mathematical expression of the position R m0 is:

[0151] R m0 = Q·r p

[0152] wherein, r p represents the position vector under the perifocal coordinate system, and Q represents the conversion matrix from the perifocal coordinate system to the inertial system;

[0153] The mathematical expression of the velocity V m0 is:

[0154] V m0 = Q·v p

[0155] The mathematical expression of the Q is:

[0156]

[0157] wherein the superscript T is the matrix transpose symbol; R X is the X-axis rotation matrix function; R Z is the Z-axis rotation matrix function; Ω ms denotes the right ascension of the ascending node; i ms denotes the orbital inclination; ω ms denotes the argument of perigee;

[0158] The R X , i.e. the mathematical expression of the X-axis rotation matrix function, is:

[0159]

[0160] The R Z , i.e. the mathematical expression of the Z-axis rotation matrix function, is:

[0161]

[0162] wherein α denotes an angle;

[0163] The mathematical expression of the r p is:

[0164]

[0165] wherein f ms denotes the true anomaly; e ms denotes the eccentricity; μ denotes the Earth gravitational constant; h ms denotes the angular momentum of the target star;

[0166] The mathematical expression of the v p is:

[0167]

[0168] wherein v p denotes the velocity vector in the perifocal coordinate system;

[0169] The mathematical expression of the h ms is:

[0170]

[0171] wherein a ms denotes the semi-major axis of the orbit;

[0172] In the module M2, the numerical integration solving step comprises:

[0173] Module A1: calculating the integration step from the task start time;

[0174] Module A2: setting single-step maximum integral step size h max , and then decomposing the integral step size;

[0175] Module A3: recursively calculating the target star orbit from the ephemeris data time to the mission start time, and obtaining the position R mt0 and velocity V mt0 of the target star at the mission start time;

[0176] In the module A1, the mathematical expression of the integral step size is:

[0177] h = t0-t s

[0178] wherein h represents the integral step size; t0 represents the mission start time; t s represents the ephemeris data time;

[0179] In the module A2, the mathematical expression of the decomposition processing is:

[0180] h = N·h max +h0

[0181] wherein N represents a natural number, and 0≤h0 max ; h0 represents the initial integral step size;

[0182] In the module A3, the N+1 numerical calculation method is called to obtain the position R mt0 and velocity V mt0 of the target star at the mission start time;

[0183] The numerical calculation method has the mathematical expression:

[0184]

[0185] wherein f represents a function; y0 represents an integral initial value; the symbol Δ represents an increment symbol; the differential equation group f(p) is a 6×1 column vector, p is an item in the brackets in k1-k4, and let

[0186] The mathematical expression of the f(p) is:

[0187]

[0188] wherein R e is the Earth equatorial radius, and J2 is the 2nd order zonal harmonic coefficient of the Earth non-spherical gravitational perturbation term;

[0189] In the module A3, y0 is the integral initial value used for the first time, and the mathematical expression is:

[0190] y0 = [Rm0 ; V m0 ]

[0191] The first N times, use single-step maximum integral step h max , that is, let Δ = h max , and take each output y1 as the integral initial value y0 of the next time; the last time let Δ = h0, calculate the output y1 to obtain the position R of the target star at the starting time of the task mt0 and the velocity V mt0 ;

[0192] The mathematical expression of the R mt0 is:

[0193] R mt0 = [y1(1); y1(2); y1(3)]

[0194] The mathematical expression of the velocity V mt0 is:

[0195] V mt0 = [y1(4); y1(5); y1(6)]

[0196] Wherein, y1(1), y1(2), y1(3), y1(4), y1(5) and y1(6) represent six components of the output y1.

[0197] Preferably, in the module M3, according to the position R mt0 and the velocity V mt0 , the orbital instant root number and the amplitude angle information are obtained; according to the orbital instant root number and the amplitude angle information, the orbital plane root number of the target star at the starting time of the task is obtained;

[0198] The orbital instant root number includes e', a', i', Ω' and ω';

[0199] The mathematical expression of the e' is:

[0200] e' = ‖E‖

[0201] Wherein, e' represents the eccentricity of the target star at the starting time of the task, which is the orbital instant root number; E represents the eccentricity vector of the target star at the starting time of the task;

[0202] The mathematical expression of the E is:

[0203]

[0204] Wherein, E represents the eccentricity vector of the target star at the starting time of the task; the symbol |||| represents the modulus operator; H represents the angular momentum vector of the target star at the starting time of the task;

[0205] The mathematical expression of the H is:

[0206] H = R mt0 x V mt0

[0207] The mathematical expression of a' is:

[0208]

[0209] Wherein, a' represents the orbit semi-major axis at the starting moment of the target star task;

[0210] The mathematical expression of i' is:

[0211]

[0212] Wherein, i' represents the orbit inclination at the starting moment of the target star task;

[0213] The mathematical expression of Ω' is:

[0214]

[0215] Wherein, Ω' represents the longitude of the ascending node at the starting moment of the target star task;

[0216] The mathematical expression of ω' is:

[0217]

[0218] Wherein, ω' represents the argument of perigee at the starting moment of the target star task;

[0219] The argument information includes: u', φ' and λ';

[0220]

[0221] Wherein, u' represents the true latitude argument at the starting moment of the target star task;

[0222]

[0223] Wherein, φ' represents the argument of declination at the starting moment of the target star task;

[0224] λ' = φ' - e' sin (φ' - ω')

[0225] Wherein, λ' represents the argument of latitude at the starting moment of the target star task.

[0226] Preferably, in the module M3, the orbit plane elements include: a mt0 , i mt0 , Ω mt0 , ω mt0 , e mt0 and Mmt0 ;

[0227] the a mt0 Mathematical expression is:

[0228]

[0229] Wherein, a mt0 Indicate the target star mission starting time of the orbit semi-major axis; x5 indicates the intermediate variable five; x6 is the intermediate variable six;

[0230] the i mt0 Mathematical expression is:

[0231] i mt0 = i'+ (1 / 6) x2Ω1sin i'

[0232] Wherein, i mt0 Indicate the target star mission starting time of the orbit inclination; x2 is the intermediate variable two:

[0233] the Ω mt0 Mathematical expression is:

[0234] Ω mt0 = Ω'- Ω1(x1-x3 / 6)

[0235] Wherein, Ω mt0 Indicate the target star mission starting time of the ascending node right ascension; x1 is the intermediate variable one; x3 is the intermediate variable three;

[0236] the mathematical expression of Ω1 is:

[0237]

[0238] the ω mt0 Mathematical expression is:

[0239] ω mt0 = arctan (-η / ξ)

[0240] Wherein, ω mt0 Indicate the target star mission starting time of the perigee amplitude angle;

[0241] the mathematical expression of η is:

[0242] η = η'- (-x 10 -(x4+x7)ξ'+x8η')

[0243] Wherein, x 10 Is the intermediate variable ten; x4 is the intermediate variable four; x7 is the intermediate variable seven; x8 is the intermediate variable eight;

[0244] the e mt0The mathematical expression of R is:

[0245]

[0246] Wherein, e mt0 represents the eccentricity of the target star mission starting time;

[0247] The mathematical expression of ξ is:

[0248] ξ = ξ' - (x9 + (x4 + x7)η' + x8ξ')

[0249] Wherein, x8 is the intermediate variable eight; x9 is the intermediate variable nine;

[0250] The mathematical expression of M mt0 is:

[0251] M mt0 = λ - ω mt0

[0252] Wherein, M mt0 represents the mean anomaly of the target star mission starting time;

[0253] The mathematical expression of λ is:

[0254] λ = λ' - (x4 + ε(x7e' 2 +x9η' + x 10 ξ')

[0255] Wherein, x9 is the intermediate variable nine.

[0256] Preferably, in the module M4, the mathematical expression of R m is:

[0257]

[0258] Wherein, E mt is the eccentric anomaly;

[0259] The mathematical expression of V m is:

[0260]

[0261] The mathematical expression of P m is:

[0262]

[0263] Wherein, P m is the intermediate variable;

[0264] The mathematical expression of Q m is:

[0265]

[0266] wherein Q m is an intermediate variable;

[0267] In the module M5, the guiding attitude of the satellite, i.e. the roll angle of the satellite in inertia

[0268]

[0269] = arcsin (i z (1))

[0270]

[0271] wherein i x , i y , i z are the mathematical expressions for

[0272]

[0273] i x = i y x i z

[0274] wherein i x , i y , i z denote the direction vectors of the three axes of the satellite body in the inertial system; R s and V s are the projections of the position vector and the velocity vector of the satellite in the J2000 inertial coordinate system, respectively.

[0275] Compared with the prior art, the application has the following beneficial effects:

[0276] 1. The orbit of the target star is recursively propagated to the start time of the observation task by means of the numerical integration method, and the orbit plane root is obtained by means of the plane instantaneous root conversion.

[0277] 2. The position and the velocity of the target star at the current time are quickly calculated by means of the analytical method, and finally the three-axis guiding attitude for the same angle observation of the satellite and the target star is solved in combination with the real-time orbit of the star.

[0278] 3. The orbit information is recursively propagated autonomously, so that the same angle observation tracking performance of the satellite and the target star can be improved, and the calculation resources on the satellite are taken into account, so that the application has practical significance for the dynamic target tracking imaging task of the satellite. BRIEF DESCRIPTION OF DRAWINGS

[0279] Other features, objects, and advantages of the application will become more apparent from the following detailed description of non-limiting embodiments thereof, when read in conjunction with the accompanying drawings:

[0280] Figure 1 The flowchart provided by the present application. DETAILED DESCRIPTION

[0281] The application will be described in detail below with specific embodiments. The following examples will help those skilled in the art to further understand the application, but do not limit the application in any form. It should be noted that for those skilled in the art, without departing from the concept of the present application, a number of changes and improvements can be made. These are within the scope of the present application.

[0282] The present application provides a satellite attitude guidance method for observing the target star in the same visual angle, which recursively propagates the orbit of the target star to the starting time of the observation task by a numerical integration method, obtains the orbital flat root number by using the flat instantaneous root conversion, and then uses an analytical method to quickly calculate the position and velocity of the target star at the current time. Finally, combined with the real-time orbit of the star, the three-axis guidance attitude of the satellite for observing the target star in the same visual angle is solved. The present application can improve the tracking performance of the satellite and the target star in the same visual angle observation, and takes into account the on-board computing resources, and has practical significance for the satellite dynamic target tracking imaging task.

[0283] In other words, the present application can improve the tracking performance of the satellite and the target star in the same visual angle observation when the satellite maintains the same visual angle observation attitude of the target star, and takes into account the on-board computing resources, and has practical significance for the satellite dynamic target tracking imaging task.

[0284] As shown in the accompanying drawings, Figure 1 According to the attitude guidance method for observing the target star in the same visual angle provided by the present application, comprising:

[0285] Step S1: calculating the position R m0 and velocity V m0 of the target star at the ephemeris data time;

[0286] Step S2: numerically integrating to solve the position R mt0 and velocity V mt0 of the target star at the starting time of the task;

[0287] Step S3: calculating the orbital flat root number of the target star at the starting time of the task;

[0288] Step S4: analytically solving the position R m and velocity V m of the target star at the current time;

[0289] Step S5: calculating the three-axis guidance attitude θ and ψ.

[0290] This invention is based on the target star at time t in the ephemeris data. s The instantaneous root of orbit in the J2000.0 inertial frame, the projections of the satellite's position and velocity vectors into the J2000 inertial coordinate system, and other parameters are used to calculate the three-axis attitude in the inertial frame when the satellite and the target satellite are observed from the same angle. The instantaneous root of orbit includes: the semi-major axis a. ms eccentricity e ms Track inclination angle i ms Right ascension of ascending node Ω ms Argument of perigee ω ms True nearest angle f ms .

[0291] First, it is necessary to use numerical integration methods to recursively calculate the target star's orbit back to the mission's start time.

[0292] (1) Calculate the position R of the target star at the time of the ephemeris data. m0 and speed V m0

[0293] The target star is recorded at time t in the ephemeris data from the ground. s Calculate the instantaneous root of the orbit in the J2000.0 inertial frame, and calculate the position R in the inertial frame. m0 and speed V m0 The formula is as follows:

[0294] Calculate the angular momentum h of the target star ms The transformation matrix Q from the near-focal coordinate system to the inertial system;

[0295]

[0296] Where μ represents the Earth's gravitational constant;

[0297]

[0298] Where the superscript T is the matrix transpose symbol; the rotation matrix function R with α as the independent variable... X (α), R Z The expressions for (α) are as follows:

[0299]

[0300] Then calculate the position vector r in the near-focal coordinate system. p The velocity vector v in the near-focal coordinate system p α represents the angle;

[0301]

[0302] The position R of the system is obtained m0and velocity V m0 ;

[0303] R m0 = Q · r p

[0304] V m0 = Q · v p

[0305] numerical integration to solve the position R mt0 and velocity V mt0 of the target star at the mission start time t0.

[0306] The integral step h is calculated from the mission start time t0, and the set maximum integral step h max The integral step is decomposed h = N · h max + h0, where N is a non-negative integer and 0 ≤ h0 < h max The following numerical calculation method is called N+1 times to recursively propagate the target star orbit from the ephemeris time to the mission start time, and the first time uses the integral initial value y0 = [R m0 ; V m0 ].

[0307] Specifically, the first N times use the maximum integral step, that is, let Δ = h max , and the output y1 each time is used as the integral initial value y0 next time. The last time is Δ = h0, and the output y1 is calculated to obtain the position R mt0 = [y1(1); y1(2); y1(3)] and velocity V mt0 = [y1(4); y1(5); y1(6)] of the target star at the mission start time t0. The y1(1), y1(2), y1(3), y1(4), y1(5), y1(6) respectively represent the six components of the output y1;

[0308] The mathematical expression of the integral step h is:

[0309] h = t0-t s

[0310] Where h represents the integral step; t0 represents the mission start time; t s represents the ephemeris time;

[0311] In the step S2, the following numerical calculation method is used

[0312]

[0313] f represents a function; y0 represents an integral initial value; the symbol Δ represents an increment symbol; the differential equation group f(p) is a 6x1 column vector, p is an item in the brackets among k1-k4, and let The simplified approximate expression of f(p) is

[0314]

[0315] wherein R e is the equatorial radius of the earth, and J2 is the 2nd order zonal harmonic coefficient of the non-spherical gravitational perturbation term of the earth.

[0316] Secondly, the inertial system position and velocity of the target star at the current time are solved by using an analytical method.

[0317] In the step S3, the orbit instantaneous root number is calculated according to the position R mt0 and the velocity V mt0 of the target star at the starting time of the task;

[0318] H=R mt0 x V mt0

[0319] wherein H represents the angular momentum vector of the target star at the starting time of the task;

[0320]

[0321] wherein E represents the eccentricity vector of the target star at the starting time of the task; the symbol |||| represents a modulus operator;

[0322] e'=‖E‖

[0323] wherein e' represents the eccentricity of the target star at the starting time of the task, and is the instantaneous root number;

[0324]

[0325] wherein a' represents the semi-major axis of the orbit of the target star at the starting time of the task, and is the instantaneous root number;

[0326]

[0327] wherein i' represents the orbit inclination of the target star at the starting time of the task, and is the instantaneous root number;

[0328]

[0329] wherein Ω' represents the longitude of the ascending node of the target star at the starting time of the task (the instantaneous root number);

[0330]

[0331] wherein ω' represents the argument of perigee of the target star at the starting time of the task (the instantaneous root number);

[0332]

[0333] wherein u' represents the true latitude amplitude angle at the starting moment of the target star task;

[0334]

[0335] wherein φ' represents the deflection latitude amplitude angle at the starting moment of the target star task;

[0336] λ' = φ' - e' sin (φ' - ω')

[0337] wherein λ' represents the plane latitude amplitude angle at the starting moment of the target star task;

[0338] In the step S3, (2) the orbit plane root number is calculated from the orbit instantaneous root number, and the algorithm is simplified as follows:

[0339]

[0340] ξ' = e' cos ω'

[0341] η' = -e' sin ω'

[0342] δ = a' (1 - e ′2 )

[0343] γ = a' (1 - ξ' cos φ' + η' sin φ')

[0344]

[0345] Specifically, according to the above, ε, ξ', η', δ, γ, Ω1 and ω1 are all intermediate variables, and have no real meaning;

[0346] x1 = u' - λ' + ξ' sin u' + η' cos u'

[0347] x2 = ξ' (3 cos u' + cos3u' ) + η' (3 sin u' - sin3u' ) + 3 cos2u' + ε (2 - ε) (ξ' 2 - η' 2 )

[0348] x3 = ξ' (3 sin u' + sin3u' ) - η' (3 cos u' - cos3u' ) + 3 sin2u' - 2ε (2 - ε) ξ' η'

[0349] x4 = - (5 / 6) Ω1 x3 cos i' + ω1 (x1 - 0.5 x3)

[0350]

[0351] x7 = -x5(ξ'sinu' + η'cosu') - 2x6ε(1 - ε) 2 ξ'η'

[0352] x8 = 0.5x5(3 - ε(2 - ε)e' 2 ) + 3x6(δ / γ)cos2u'

[0353]

[0354] wherein x1~x 10 are all intermediate variables without real meaning; specifically, x1~x 10 correspond to intermediate variable one~intermediate variable ten;

[0355]

[0356] wherein a mt0 represents the orbital semi-major axis at the starting moment of the target star mission; (dimensionless) ;

[0357] i mt0 = i' + (1 / 6)x2Ω1sini'

[0358] wherein i mt0 represents the orbital inclination at the starting moment of the target star mission; (dimensionless) ;

[0359] Ω mt0 = Ω' - Ω1(x1 - x3 / 6)

[0360] wherein Ω mt0 represents the right ascension of the ascending node at the starting moment of the target star mission; (dimensionless) ;

[0361] ξ = ξ' - (x9 + (x4 + x7)η' + x8ξ')

[0362] η = η' - (-x 10 -(x4 + x7)ξ' + x8η')

[0363] λ = λ' - (x4 + ε(x7e' 2 +x9η' + x 10 ξ')

[0364] wherein λ, ξ and η are all intermediate variables without real meaning;

[0365] ω mt0 = arctan(-η / ξ)

[0366] wherein ω mt0 represents the argument of perigee at the starting moment of the target star mission; (dimensionless) ;

[0367]

[0368] wherein e mt0 represents the eccentricity of the target star at the starting time of the task; (the root number)

[0369] M mt0 = λ-ω mt0

[0370] wherein M mt0 represents the mean anomaly of the target star at the starting time of the task; (the root number)

[0371] (3) Analytical method is used to solve the position R m and the velocity V m of the target star at the current time.

[0372] The auxiliary vector is calculated

[0373]

[0374] wherein P m is an intermediate variable without practical meaning.

[0375]

[0376] wherein Q m is an intermediate variable without practical meaning.

[0377] Up to now, the above steps only need to be calculated once in one observation task, and the calculation steps of each control period only need to perform the following steps, thereby saving the calculation resources.

[0378] In the step S4, the recursive time Δt is calculated from the current time t, and the mean anomaly is analytically recursively calculated Then, the mean anomaly M mt is converted into the eccentric anomaly E mt , and the calculation method of the conversion of the mean anomaly into the eccentric anomaly is a conventional algorithm, which will not be described herein.

[0379] The mathematical expression of the recursive time Δt is as follows:

[0380] Δt = t-t0

[0381] The mathematical expression of the mean anomaly M mt is as follows:

[0382]

[0383] By using the orbit root number of the target star, the formula of the analytical algorithm is used to calculate the position R m and the velocity V m of the target star at the current time, i.e. the position R m and the velocity Vm The formula is as follows:

[0384]

[0385] Finally, the three-axis attitude of the satellite in the inertial system when observing the target star with the same visual angle is calculated according to the position and velocity information of the target star and the position and velocity information of the satellite.

[0386] (1) The expected attitude rotation matrix of the satellite when observing the target star with the same visual angle is calculated, the visual axis (+Z axis) of the satellite camera is reversely directed to the target star, and the rotation direction of the visual axis during the task is parallel to the +Y axis as a constraint condition, vector i x y z The calculation method is as follows:

[0387]

[0388] x y z

[0389] Wherein, R s and V s are the projections of the position vector and the velocity vector of the satellite in the J2000 inertial coordinate system respectively.

[0390] In the step S5, (2) the three-axis attitude in the inertial system is calculated, taking the Euler angle of 1-2-3 rotation sequence as an example, that is, the roll angle of the satellite when observing the target star with the same visual angle in the inertial system The pitch angle θ and the yaw angle ψ are calculated according to the following formula:

[0391]

[0392] θ = arcsin (i z (1))

[0393]

[0394] Wherein, the mathematical expressions of i x , i y , i z are as follows:

[0395]

[0396] x y z

[0397] Wherein, i x , i y , i z ​​​​​​​​respectively represent the direction vectors of the three axes of the satellite body in the inertial system; R s and V s respectively represent the projection of the position vector and the velocity vector of the satellite in the J2000 inertial coordinate system.

[0398] The application also provides a satellite attitude guiding system for observing the target star in the same visual angle, which can be realized by executing the process steps of the satellite attitude guiding method for observing the target star in the same visual angle, i.e., the satellite attitude guiding method for observing the target star in the same visual angle can be understood as the preferred implementation of the satellite attitude guiding system for observing the target star in the same visual angle by those skilled in the art.

[0399] According to the application, a satellite attitude guiding system for observing a target star in the same visual angle is provided, which comprises:

[0400] Module M1: calculating the position R m0 and the velocity V m0 of the target star at the ephemeris data moment;

[0401] Module M2: numerically integrating to solve the position R mt0 and the velocity V mt0 of the target star at the task starting moment;

[0402] Module M3: according to the position R mt0 and the velocity V mt0 , the orbital root number of the target star at the task starting moment is solved;

[0403] Module M4: based on the orbital root number, the real-time position R m and the velocity V m of the target star are solved;

[0404] Module M5: according to the position R m and the velocity V m , the guiding attitude of the satellite is calculated.

[0405] Those skilled in the art know that, in addition to implementing the system provided by the present application and each device, module and unit thereof in the form of pure computer readable program code, the system provided by the present application and each device, module and unit thereof can also be implemented in the form of logic gates, switches, application specific integrated circuits, programmable logic controllers and embedded microcontrollers, etc. by logically programming the method steps to achieve the same functions. Therefore, the system provided by the present application and each device, module and unit thereof can be considered as a hardware component, and the devices, modules and units included therein for achieving various functions can also be considered as structures within the hardware component; the devices, modules and units for achieving various functions can also be considered as both software modules implementing methods and structures within hardware components.

[0406] The specific embodiments of the present application are described above. It needs to be understood that the present application is not limited to the specific embodiments described above, and various changes or modifications can be made by those skilled in the art within the scope of the claims, which does not affect the essential content of the present application. The embodiments of the present application and the features in the embodiments can be combined with each other in any manner without conflict.

Claims

1. A satellite attitude guidance method for observing a target star with the same visual angle, characterized by, The method comprises the following steps: Step S1: Calculate the position R of the target star at the ephemeris data time m0 and velocity V m0 ; Step S2: Numerically integrate to solve the position R of the target star at the start of the mission mt0 and velocity V mt0 ; Step S3: determining the position R of the target star at the start of the mission based on the velocity V mt0 and the time T mt0 , the orbital plane of the target star at the start of the mission is determined. Step S4: Based on the orbital square root number, calculate the real-time position R of the target star. m and speed V m ; Step S5: determining the position R of the satellite according to the position P of the ground station and the distance D between the ground station and the satellite m and the velocity V m calculating the guiding attitude of the satellite In said step S3, the position R mt0 with the velocity V mt0 , the orbital moment and the amplitude information are determined; According to the orbital instantaneous element and the argument information, the orbital plane element of the target star at the starting time of the task is obtained; The orbital instantaneous element comprises e', a', i', Ω' and ω'; The mathematical expression of e' is: e'=‖E‖ Wherein, e' represents the eccentricity of the target star at the starting time of the task, and E represents the eccentricity vector of the target star at the starting time of the task; The mathematical expression of E is: Wherein, E represents the eccentricity vector of the target star at the starting time of the task; the symbol || || represents the modulus operator; and H represents the angular momentum vector of the target star at the starting time of the task; The mathematical expression of H is: H = R mt0 x V mt0 The mathematical expression of a' is: Wherein, a' represents the semi-major axis of the target star at the starting time of the task; The mathematical expression of i' is: Wherein, i' represents the inclination of the target star at the starting time of the task; The mathematical expression of Ω' is: Wherein, Ω' represents the right ascension of the ascending node of the target star at the starting time of the task; The mathematical expression of ω' is: Wherein, ω' represents the argument of perigee of the target star at the starting time of the task; The argument information comprises u', φ' and λ'; Wherein, u' represents the argument of true latitude of the target star at the starting time of the task; Wherein, φ' represents the argument of eccentric latitude of the target star at the starting time of the task; λ'=φ'-e'sin(φ'-ω') Wherein, λ' represents the argument of plane latitude of the target star at the starting time of the task; In the step S3, the track flat number comprises: a mt0 , i mt0 , Ω mt0 , ω mt0 , e mt0 and M mt0 ; The a mt0 The mathematical expression is: wherein a mt0 x5 represents an intermediate variable five; x6 is an intermediate variable six; The i mt0 The mathematical expression is: i mt0 = i' + (1 / 6)x2Ω1sini' where i mt0 x2is an intermediate variable two: The Ω mt0 The mathematical expression is: Ω mt0 = Ω' - Ω1(x1 - x3 / 6) where Ω mt0 denotes the right ascension of the ascending node at the beginning of the mission; x1 is an intermediate variable one; x3 is an intermediate variable three; The mathematical expression of Ω1 is: The omega mt0 The mathematical expression is: ω mt0 = arctan(-η / ξ) where ω mt0 denotes the argument of perigee at the beginning of the mission. The mathematical expression of η is: η = η' - (-x 10 - (x4+x7)ξ' + x8η' wherein x 10 is an intermediate variable ten; x4 is an intermediate variable four; x7 is an intermediate variable seven; x8 is an intermediate variable eight; The e mt0 The mathematical expression is: where e mt0 eccentricity of the target star mission start time; The mathematical expression of ξ is: ξ=ζ'-(x9+(x4+x7)η'+x8ξ') Wherein, x8 is the eighth intermediate variable; and x9 is the ninth intermediate variable; The M mt0 The mathematical expression is: M mt0 = λ - ω mt0 wherein M mt0 denotes the mean anomaly at the beginning of the mission; The mathematical expression of λ is: λ = λ' - (x4+ ε(x7e' + x9η' + x1 1 ζ' ) ). 2 + x9η' + x 10 ζ' ) ). ζ' ) ). Wherein, x9 is the ninth intermediate variable; In said step S4, said R m The mathematical expression is: wherein E mt is the eccentric anomaly; The V m The mathematical expression is: The P m The mathematical expression is: P = P + P m is an intermediate variable; The Q m The mathematical expression is: wherein Q m is an intermediate variable; In said step S5, the satellite's guiding attitude, i.e. the satellite's roll angle in inertial The mathematical expressions of the pitch angle θ and the yaw angle ψ are respectively: θ = arcsin(i z (1)) wherein i x , i y , i z The mathematical expression is: i x = i y x i z where i x , i y , i z denote the direction vectors of the satellite body three-axes in the inertial frame; R s and V s are the projections of the satellite position and velocity vectors in the J2000 inertial coordinate system, respectively.

2. The method according to claim 1, wherein the target star is a star in the constellation of Ursa Major. In said step S1, the position R m0 The mathematical expression of R is: R m0 = Q - r p where r p represents the position vector in the near-focus coordinate system, and Q represents the conversion matrix from the near-focus coordinate system to the inertial system; The speed V m0 The mathematical expression is: V m0 = Q - v p The mathematical expression of Q is: where the superscript T denotes the matrix transpose; P X is the X-axis rotation matrix function; R Z is the Z-axis rotation matrix function; Ω ms denotes the right ascension of the ascending node; i ms denotes the inclination of the orbit; ω ms denotes the argument of the perigee; The R X That is, the mathematical expression of the X-axis rotation matrix function is: The R Z The mathematical expression of the Z-axis rotation matrix function is: Wherein, α represents an angle; The r p The mathematical expression is: where f ms denotes the true anomaly; e ms denotes the eccentricity; μ denotes the Earth gravitational constant; h ms denotes the target star angular momentum; The v p The mathematical expression is: where v p denotes the velocity vector in the near-focus coordinate system; The h ms The mathematical expression is: wherein a ms denotes the semi-major axis of the orbit; In the step S2, the numerical integral solving method comprises the following steps: Step A1: calculating the integral step length from the starting time of the task; Step A2: Set the single step maximum integration step size h max and further decompose the integration step size; Step A3: Propagate the target star orbit from the ephemeris data time to the mission start time, calculate the position of the target star at the mission start time R mt0 with the velocity V mt0 ; In the step A1, the mathematical expression of the integral step length is: h = t0 - t s wherein h denotes the integration step; t0denotes the task start time; t s denotes the time of the ephemeris data; In the step A2, the mathematical expression of the decomposition processing is: h = N - h max + h0 where N denotes a natural number, and 0 < h0< h max ; h0denotes an initial integration step size; In the step A3, the N+1 times value calculation method is called to obtain the position R of the target star at the beginning of the task mt0 With the speed V mt0 ; The numerical calculation method is mathematically expressed as: wherein f denotes a function; y0 denotes an initial value of integration; the symbol Δ denotes an increment symbol; the differential equation set f(p) is a 6 x 1 column vector, p is a term in the brackets among k1 to k4, and let The mathematical expression of f(p) is: wherein R e is the Earth equatorial radius, and J2is the second order zonal coefficient of the Earth's non-spherical gravitational perturbation term. In the step A3, y0, that is, the first used integral initial value, is mathematically expressed as: y0 = [R m0 ; V m0 ] The first time, use the single-step maximum integration step size h max That is, let Δ = h max And output y1 each time as the next integration initial value y0; the last time, let Δ = h0, and calculate the output y1 to obtain the position R of the target star at the start time of the task mt0 And the velocity V mt0 ; The R mt0 The mathematical expression is: R mt0 = [y1(1); y1(2); y1(3)] The speed V mt0 The mathematical expression is: V mt0 = [y1(4); y1(5); y1(6)] Wherein, y1(1), y1(2), y1(3), y1(4), y1(5) and y1(6) represent six components of the output y1.

3. A satellite attitude guidance system for observing a target star in the same visual field, characterized by, The method comprises the following steps: Module M1 : Compute the position R of the target star at the time of the ephemeris data m0 and the velocity V m0 ; Module M2: Numerical integration to solve the position R of the target star at the beginning of the mission mt0 and the velocity V mt0 ; Module M3: determining the position R of the target star at the start of the mission mt0 with the speed V mt0 , the orbital plane of the target star at the start of the mission is determined Module M4: Based on the stated orbital square root number, calculate the real-time position R of the target star. m and speed V m ; Module M5: calculating the guiding attitude of said satellite as a function of said position R m with the speed V m calculating the guiding attitude of said satellite; In said module M3, from said position R mt0 with the speed V mt0 , the orbital instantaneous number and amplitude angle information are found; According to the orbital instantaneous element and the argument information, the orbital plane element of the target star at the starting time of the task is obtained; The orbital instantaneous element comprises e', a', i', Ω' and ω'; The mathematical expression of e' is: e'=‖E‖ Wherein, e' represents the eccentricity of the target star at the starting time of the task, and E represents the eccentricity vector of the target star at the starting time of the task; The mathematical expression of E is: Wherein, E represents the eccentricity vector of the target star at the starting time of the task; the symbol || || represents the modulus operator; and H represents the angular momentum vector of the target star at the starting time of the task; The mathematical expression of H is: H = R mt0 x V mt0 The mathematical expression of a' is: Wherein, a' represents the orbit semi-major axis of the target star task starting time; The mathematical expression of i' is: Wherein, i' represents the orbit inclination of the target star task starting time; The mathematical expression of Ω' is: Wherein, Ω' represents the ascending node right ascension of the target star task starting time; The mathematical expression of ω' is: Wherein, ω' represents the perigee amplitude angle of the target star task starting time; The amplitude angle information includes: u', φ' and λ'; Wherein, u' represents the true latitude amplitude angle of the target star task starting time; Wherein, φ' represents the deflection latitude amplitude angle of the target star task starting time; λ' = φ' - e' sin (φ' - ω') Wherein, λ' represents the plane latitude amplitude angle of the target star task starting time; In said module M3, said track flat number comprises: a mt0 , i mt0 , Ω mt0 , ω mt0 , e mt0 and M mt0 ; The a mt0 The mathematical expression is: wherein a mt0 x5 represents an intermediate variable five; x6 is an intermediate variable six; The i mt0 The mathematical expression is: i mt0 = i' + (1 / 6)x2Ω1sini' where i mt0 x2is an intermediate variable two: The Ω mt0 The mathematical expression is: Ω mt0 = Ω' - Ω1(x1 - x3 / 6) wherein Ω mt0 denotes the right ascension of the ascending node at the beginning of the mission; x1 is an intermediate variable one; x3 is an intermediate variable three; The mathematical expression of Ω1 is: The omega mt0 The mathematical expression is: ω mt0 = arctan(-η / ξ) where ω mt0 denotes the argument of perigee at the beginning of the mission. The mathematical expression of η is: η = η' - (-x 10 - (x4+x7)ξ' + x8η' wherein x 10 is an intermediate variable ten; x4 is an intermediate variable four; x7 is an intermediate variable seven; x8 is an intermediate variable eight; The e mt0 The mathematical expression is: where e mt0 eccentricity of the target star mission start time; The mathematical expression of ξ is: ξ = ξ' - (x9 + (x4 + x7) η' + x8 ξ') Wherein, x8 is the intermediate variable eight; x9 is the intermediate variable nine; The M mt0 The mathematical expression is: M mt0 = λ - ω mt0 wherein M mt0 denotes the mean anomaly at the beginning of the mission; The mathematical expression of λ is: λ = λ' - (x4+ ε(x7e' + x9η' + x 2 + x 10 ξ' )) Wherein, x9 is the intermediate variable nine; In said module M4, said R m The mathematical expression is: wherein E mt is the eccentric anomaly; The V m The mathematical expression is: The P m The mathematical expression is: P = P + P m is an intermediate variable; The Q m The mathematical expression is: wherein Q m is an intermediate variable; In said module M5, the satellite's guided attitude, i.e. the satellite's roll angle in inertial The mathematical expressions of the pitch angle θ and the yaw angle ψ are respectively: θ = arcsin(i z (1)) wherein i x , i y , i z The mathematical expression is: i x = i y x i z where i x , i y , i z denote the direction vectors of the satellite body three-axes in the inertial frame; R s and V s are the projections of the satellite position and velocity vectors in the J2000 inertial coordinate system, respectively.

4. The satellite attitude guidance system for observing the target star with the same visual angle as claimed in claim 3, wherein In said module M1, said position R m0 is mathematically expressed by: R m0 = Q - r p where r p represents the position vector in the near-focus coordinate system, and Q represents the conversion matrix from the near-focus coordinate system to the inertial system; The speed V m0 The mathematical expression is: V m0 = Q - v p The mathematical expression of Q is: where the superscript T is the matrix transpose symbol; R X is the X-axis rotation matrix function; R Z is the Z-axis rotation matrix function; Ω ms denotes the right ascension of the ascending node; i ms denotes the orbital inclination; ω ms denotes the argument of perigee; The R X That is, the mathematical expression of the X-axis rotation matrix function is: The R Z The mathematical expression of the Z-axis rotation matrix function is: Wherein, α represents the angle; The r p The mathematical expression is: where f ms denotes the true anomaly; e ms denotes the eccentricity; μ denotes the Earth gravitational constant; h ms denotes the target star angular momentum; The v p The mathematical expression is: where v p denotes the velocity vector in the near-focus coordinate system; The h ms The mathematical expression is: wherein a ms denotes the semi-major axis of the orbit; In the module M2, the numerical integral solving method step includes: Module A1: calculating the integral step length from the task starting time; Module A2: Setting the single-step maximum integration step size h max and further decomposing the integration step size; Module A3: Propagate the target star's orbit from the ephemeris data time to the mission start time, calculate the target star's position R at the mission start time mt0 with the velocity V mt0 ; In the module A1, the mathematical expression of the integral step length is: h = t0 - t s wherein h denotes the integration step; t0denotes the task start time; t s denotes the time of the ephemeris data; In the module A2, the mathematical expression of the decomposition processing is: h = N - h max + h0 where N denotes a natural number, and 0 < h0< h max ; h0denotes an initial integration step size; In the module A3, the N+1 times value calculation method is called to obtain the position R of the target star at the beginning of the task mt0 With the speed V mt0 ; The numerical calculation method, mathematical expression is: wherein f denotes a function; y0 denotes an initial value of integration; the symbol Δ denotes an increment symbol; the differential equation set f(p) is a 6 x 1 column vector, p is a term in the brackets among k1 to k4, and let The mathematical expression of f (p) is: wherein R e is the Earth equatorial radius, and J2is the second order zonal coefficient of the Earth's non-spherical gravitational perturbation term. In the module A3, y0 is the first time using the integral initial value, the mathematical expression is: y0 = [R m0 ; V m0 ] The first time, use the single-step maximum integration step size h max That is, let Δ = h max And output y1 each time as the next integration initial value y0; the last time, let Δ = h0, and calculate the output y1 to obtain the position R of the target star at the start time of the task mt0 And the velocity V mt0 ; The R mt0 The mathematical expression is: R mt0 = [yl(l); yl(2); yl(3)] The speed V mt0 The mathematical expression is: V mt0 = [yl(4); yl(5); yl(6)] Wherein, y1 (1), y1 (2), y1 (3), y1 (4), y1 (5) and y1 (6) represent six components of the output y1.

Images