Adaptive state noise compensation maneuvering spacecraft orbit determination method

Abstract

The invention discloses a maneuvering spacecraft orbit determination method for adaptive state noise compensation, and belongs to the technical field of aerospace. The implementation method comprises the following steps: establishing a general dynamical model and an angle-only observation model for unknown mutation of the random nonlinear dynamical system; predicting a system state and a covariance matrix based on the dynamic model; calculating a system measurement residual error and a covariance matrix at the tk + 1 moment, and representing the matching degree of a prediction result of the system based on the kinetic model and the maneuver-free model through an index function lambda k + 1; comparing the index function at the tk + 1 moment with a mean value of the index functions from t1 to tk, when the difference between an evaluation index and the mean value is greater than a threshold value, judging that the state of the system is abruptly changed, and performing state noise compensation on a system prediction result; and calculating orbit determination algorithm gain based on a system state noise compensation result, weighting and outputting a tk + 1 system state and a covariance estimation value based on a system prediction result # imgabs0 # noise compensation result # imgabs1 # and a measurement result zk + 1, and iteratively solving a state estimation value at each moment so as to realize high-precision orbit determination of the maneuvering spacecraft.

Description

Technical Field

[0001] The invention relates to a method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation, and is particularly suitable for tracking a maneuvering spacecraft that uses pulse thrust for maneuvering, and belongs to the field of aerospace technology. Background Art

[0002] There are more than 1,000 maneuverable spacecraft in orbit, and the number of such spacecraft is still increasing. This has triggered a strong demand for situational awareness of space targets, which requires the ability to accurately predict the trajectory of maneuverable spacecraft, so as to perform tasks such as intention identification and function prediction of target spacecraft.

[0003] At present, many spacecraft maneuvering tracking algorithms have been developed to estimate the position and velocity of maneuvering spacecraft. The traditional interactive multiple model (IMM) algorithm has been successfully applied to the tracking problem of maneuvering spacecraft. However, the traditional IMM algorithm mixes the spacecraft state models under different motion states based on artificial experience given values, and the tracking and orbit determination accuracy of the spacecraft is low. The improved IMM algorithm assisted by prior knowledge establishes a state mutation probability model, which makes the state mutation depend on the prior information of the system. Although this method effectively improves the tracking and orbit determination accuracy of the system, the covariance matrix and scaling coefficient in the probability model are given based on empirical values, which makes the system tracking accuracy low near the state mutation area, and for most maneuvering orbit change problems, the prior knowledge of the system is difficult to obtain.

[0004] Based on this, this patent proposes a maneuverable spacecraft tracking method with adaptive state noise compensation. First, the state prediction is performed based on the system dynamics equation; then, a model matching evaluation method based on the system residual is proposed, and based on this result, it is judged whether the system has a state mutation; subsequently, for the case where the system has a state mutation, the system state noise is solved by inversion based on the system model matching, and the state noise compensation is performed on the prediction result based on the system dynamics equation; finally, the final weighted result is output in combination with the measurement result, and the above process is iterated to achieve high-precision orbit determination of the maneuverable spacecraft.

[0005] Among the developed maneuvering spacecraft tracking methods, the technique [1] (see: [1] LEE S, LEE J, HWANGI. Maneuvering Spacecraft Tracking via State-Dependent Adaptive Estimation [J / OL]. Journal of Guidance, Control, and Dynamics, 2016, 39 (9): 2034-2043.) proposed an improved IMM algorithm assisted by prior knowledge, established a probability model to describe the sudden change of the system state, and then integrated this state-based probability model with the IMM algorithm, replacing the model transition probability given by artificial experience in the traditional IMM algorithm that mixes the spacecraft state models under different motion states. However, this method relies on the prior information of the system, and the covariance matrix and scaling coefficient of the probability model are based on given empirical values. Therefore, the actual tracking accuracy is low near the sudden change of the system state. Summary of the invention

[0006] The technical problem to be solved by the method for orbit determination of maneuvering spacecraft with adaptive state noise compensation disclosed in the present invention is to track and determine the orbit of a spacecraft that uses pulse thrust for maneuvering, so that the target spacecraft can still achieve high-precision tracking and orbit determination when the prior maneuvering information is unknown. The present invention has the following advantages: (1) high tracking accuracy and good convergence; (2) compared with the traditional IMM algorithm and the parallel operation of two different state models in its improved algorithm, the method performs state noise compensation based on the evaluation index, and has high tracking and orbit determination efficiency; (3) wide application range, even when the prior maneuvering information of the target spacecraft is unknown, it can still achieve high-precision tracking and orbit determination.

[0007] The objective of the present invention is achieved through the following technical solutions:

[0008] The invention discloses a method for orbit determination of a maneuverable spacecraft with adaptive state noise compensation, which establishes a general dynamics model and an angle-only observation model for unknown mutations of a random nonlinear dynamics system; predicts the system state and covariance matrix based on the dynamics model; calculates t k+1 The system measures the residual and covariance matrix at each moment and uses the indicator function Λ k+1 Characterizes the matching degree between the prediction results based on the dynamic model and the non-maneuvering model of the system; t k+1 The index function at time t1 to t k When the difference between the evaluation index and the mean is greater than the threshold value, it is determined that the system has a sudden change in state, and the state noise compensation is performed on the system prediction result; based on the result of the system state noise compensation, the orbit determination algorithm gain is calculated, and based on the system prediction result Noise compensation results And the measurement results z k+1 Weighted output t k+1 The system state and covariance estimates are iteratively solved to achieve high-precision orbit determination for maneuverable spacecraft.

[0009] The present invention discloses a method for orbit determination of a maneuverable spacecraft with adaptive state noise compensation, comprising the following steps:

[0010] Step 1: Establish a general dynamic model and an angle-only observation model for unknown mutations of random nonlinear dynamic systems;

[0011] Step 1.1 establish a general dynamic model for unknown mutations of random nonlinear dynamic systems;

[0012] The state variable differential equation of the space target is expressed as:

[0013]

[0014] In the formula, Table system state vector; r k represents the position vector of the spatial target at time k, [x k y k z k ] T is its scalar form; v k represents the velocity vector of the space target at time k, is its scalar form;

[0015] Establish a general dynamics model for unknown mutations in random nonlinear dynamic systems:

[0016] x k+1 =F(x k +δ t,k Gu)+Bw k (2)

[0017] Among them, δ t,k represents the Kronecker symbol. When t≠k, δ t,k =0, when t = k, δ t,k =1; is a constant matrix; represents the unknown input vector of the system; is a constant matrix; represents white noise with an expected value of 0; δ t,k Gu represents the random maneuvers that occur in the system; F(x k +δ t,k Gu) represents the state recursion based on the state variable differential equation (1) of the space target. The specific process is to transform the system state x at time k into k +δt,k Gu is used as the initial value, and the numerical integrator is used to integrate equation (1) to obtain the recursive result x based on dynamics at the k+1th moment. k+1 ;

[0018] Step 1.2: Establish an angle measurement observation model only;

[0019] The angle information of the maneuvering target relative to the observation satellite is measured as:

[0020] z k+1 =[α k+1 β k+1 ] T (3)

[0021] In the formula, z k+1 Indicates the quantity measurement, α k+1 and β k+1 are the pitch angle α of the maneuvering target relative to the observation satellite k+1 and azimuth β k+1 ;

[0022] The observation model is expressed as:

[0023] z k+1 =h(x k+1 )+ν k+1 (4)

[0024] Where:

[0025]

[0026] [x k+1 y k+1 z k+1 ] T is the position vector of the maneuvering target, [X k+1 Y k+1 Z k+1 ] T is the position vector of the observed star; ν k+1 is two-dimensional zero-expectation white noise, R k+1 is the measurement noise covariance matrix, σ z is the standard deviation of the measured value; h(x k+1 ) and ρ k+1 is an intermediate variable.

[0027] Step 2: Predict the system state and covariance matrix based on the dynamic model;

[0028] Predict the system status:

[0029]

[0030] In the formula, represents the predicted value of the state at time k+1; represents the estimated value of the state at time k; It represents the state prediction based on the state variable differential equation (1) of the space target. The specific process is to convert the state estimate value at time k into As the initial value, the numerical integrator is used to integrate equation (1) to obtain the state prediction value at the k+1th moment:

[0031] Make predictions on the system covariance matrix:

[0032]

[0033] in, Represents the predicted value of the covariance matrix at time k+1; represents the estimated value of the covariance matrix at time k; Φ k represents the system state transfer matrix from time k to time k+1, and its differential equation is expressed as:

[0034]

[0035] In the formula, represents the Jacobian matrix, which is specifically expressed as:

[0036]

[0037] Step 3: Calculate t k+1 The system measures the residual and covariance matrix at each moment and uses the indicator function Λ k+1 Characterizes the matching degree between the prediction results based on the dynamic model and the non-maneuvering model of the system; t k+1 The index function at time t1 to t k When the difference between the evaluation index and the mean value is greater than the threshold value, it is considered that the system has a sudden change in state, and the state noise compensation is performed on the system prediction result;

[0038] Step 3.1: Calculate t k+1 The time system measurement residual γ k+1 and the covariance matrix S k+1 , and through the indicator function Λ k+1 Characterize the degree of match between the prediction results of the system based on the dynamic model and the non-maneuvering model;

[0039] System measurement residual γ k+1 and the covariance matrix S k+1 It is expressed as:

[0040]

[0041] Among them, the matrix H k+1 It is expressed as:

[0042]

[0043] t k+1 Time index function Λ k+1 It is expressed as:

[0044]

[0045] Among them, tr(S k+1 ) represents the trace of the measurement residual covariance matrix;

[0046] Step 3.2: Obtain t through equation (15) k+1 The indicator function of the moment; calculate t1 to t k The mean of the indicator function; when the difference between the indicator function and the mean is greater than the threshold value, the system state changes suddenly, and the system predicts the covariance matrix For state noise compensation;

[0047] From time t1 to t k The mean of the indicator function at time is expressed as:

[0048] vale k =(vale k-1 ×(k-1)+Λ k ) / k (16)

[0049] Among them, the initial value of the indicator function mean is set to 0;

[0050] State Noise Compensation Q k+1 It is expressed as:

[0051]

[0052] Where η is the threshold value for judging whether a mutation occurs in the system; D k+1 Indicates the state noise compensation of the position when the system state suddenly changes; V k+1 It represents the state noise compensation of the position when the system state suddenly changes;

[0053] State noise compensation for position when the system state suddenly changes k+1 It is expressed as:

[0054]

[0055] In the formula, the intermediate variable dr k+1 It is expressed as:

[0056]

[0057] In the formula, the intermediate variable Ω k+1 It can be expressed as:

[0058]

[0059] State transition matrix Φ k+1 The scalar form of is:

[0060]

[0061] State noise compensation for speed when the system state suddenly changes V k+1 It is expressed as:

[0062]

[0063] In the formula, the intermediate variable dv1 k+1 ,dv2 k+1 and dv3 k+1 It is expressed as:

[0064]

[0065] In the formula, △v represents the estimated value of position deviation, which is specifically expressed as:

[0066]

[0067] Compensate the system prediction results for state noise:

[0068]

[0069] Step 4: Calculate the orbit determination algorithm gain based on the system state noise compensation result; Noise compensation results And the measurement results z k+1 Weighted output t k+1 The system state and covariance estimates are brought into step 2 to iteratively solve the estimated values ​​of the system state variables at each moment, and high-precision orbit determination of the maneuverable spacecraft is achieved based on the estimated values ​​of the system state variables at each moment.

[0070] Step 4.1 Orbit determination algorithm gain:

[0071]

[0072] Step 4.2 Based on the system prediction results Noise compensation results And the measurement results z k+1 Weighted output t k+1 System state and covariance estimates:

[0073]

[0074] Will Substitute into formula (8), and change P k+1 Substitute into equation (9), iteratively solve to obtain the estimated value of the system state variable at each moment, set the observation stop time as the orbit determination termination time, and realize high-precision orbit determination of the maneuverable spacecraft.

[0075] Beneficial effects:

[0076] 1. Compared with the traditional IMM maneuverable spacecraft tracking method, the present invention discloses a maneuverable spacecraft tracking method with adaptive state noise compensation. First, the system state model in the non-maneuvering state is evaluated by measuring the residual error. Then, based on the evaluation result, a step function is constructed to fuse the system state noise to compensate for the state mutation in the maneuvering state. This method can improve the tracking and orbit determination accuracy of the system while effectively identifying the state mutation in the maneuvering state of the spacecraft.

[0077] 2. Compared with the traditional IMM algorithm and its improved algorithm that perform parallel operations on two different state models, the present invention discloses a maneuverable spacecraft tracking method with adaptive state noise compensation, which only compensates for state noise based on evaluation indicators, thereby improving the computational efficiency of maneuverable spacecraft tracking.

[0078] 3. The present invention discloses a method for tracking a maneuverable spacecraft with adaptive state noise compensation. The method performs state noise compensation based on the residual evaluation results of the system state model measured in a non-maneuvering state. The method does not rely on the prior information of the system state mutation. It can still achieve high-precision tracking and orbit determination when the prior maneuver information of the target spacecraft is unknown.

[0079] 4. Compared with the traditional IMM algorithm and its improved algorithm, the present invention discloses a method for tracking a maneuverable spacecraft with adaptive state noise compensation. It solves the system state noise by inverting according to the matching degree of the system model, which can more accurately fit the actual state mutation of the system and improve the orbit determination accuracy of the spacecraft during maneuvering. BRIEF DESCRIPTION OF THE DRAWINGS

[0080] Figure 1 The present invention discloses a flow chart of a maneuvering spacecraft tracking method with adaptive state noise compensation.

[0081] Figure 2 This is the distribution diagram of observed constellations under the CRTBP model in Example 1.

[0082] Figure 3 This is the target spacecraft maneuvering at the apogee to maintain orbit diagram in Example 1.

[0083] Figure 4 It is the root mean square error result of the position obtained by 100 Monte Carlo simulations in Example 1.

[0084] Figure 5 This is the root mean square error result of the 100 Monte Carlo simulations obtained in Example 1. DETAILED DESCRIPTION

[0085] In order to better illustrate the purpose and advantages of the present invention, a simulation analysis of the apogee maneuver maintenance scenario of the NRHO orbit in the cis-lunar space is performed below, and the present invention is explained in detail.

[0086] The target spacecraft is tracked under the CRTBP model, and the observation satellite orbit is set to the L1 point Halo orbit (orbital amplitude is 23200km, orbital initial phase is 240deg) and the L2 point plane Lyapunov orbit (orbital amplitude is 5000km, orbital initial phase is 120deg). According to this setting, the observation scene is as follows Figure 2 The target spacecraft orbit is set to NRHO. Since the spacecraft speed is fast near the perigee, and considering the low fuel consumption, maneuvers are applied near the apogee. The nominal orbit of the target orbit is obtained by the differential correction method. The initial orbit data is shown in Table 1. According to the above settings, the nominal orbit of the target spacecraft is obtained as follows Figure 3 The observation interval is selected as 0.1 hour, and the target is observed continuously for 6 hours.

[0087] Table 1 Initial data for maintaining apogee orbit

[0088]

[0089] Embodiment 1:

[0090] like Figure 1 As shown, the present embodiment discloses a method for tracking a maneuverable spacecraft with adaptive state noise compensation, and the specific implementation steps are as follows:

[0091] Step 1: Establish a general dynamic model and an angle-only observation model for unknown mutations of random nonlinear dynamic systems;

[0092] Step 1.1 Establish a general dynamic model for unknown mutations of random nonlinear dynamic systems;

[0093] The state variable differential equation of the space target is expressed as:

[0094]

[0095] In the formula, Table system state vector; r k represents the position vector of the spatial target at time k, [x k y k z k ] T is its scalar form; vk represents the velocity vector of the space target at time k, is its scalar form;

[0096] Establish a general dynamics model for unknown mutations in random nonlinear dynamic systems:

[0097] x k+1 =F(x k +δ t,k Gu)+Bw k (30)

[0098] Among them, δ t,k represents the Kronecker symbol. When t≠k, δ t,k =0, when t = k, δ t,k =1; is a constant matrix. In order to make the setting of G change only the speed but not the position, it is specifically set to G = [O 3×3 ;I 3×3 ]; represents the unknown input vector of the system; is a constant matrix. In this embodiment, the system state noise from other sources is not considered, and is specifically set to B=0 6×6 ; represents white noise with an expected value of 0; δ t,k Gu represents the random maneuvers that occur in the system; F(x k +δ t,k Gu) represents the state recursion based on the state variable differential equation (29) of the space target. The specific process is to transform the system state x at time k into k +δ t,k Gu is used as the initial value, and the numerical integrator is used to integrate equation (29) to obtain the recursive result x based on dynamics at the k+1th moment: k+1 , select ode45 integrator in MATLAB as the preferred numerical integrator;

[0099] Step 1.2: Establish an angle measurement observation model only;

[0100] The angle information of the maneuvering target relative to the observation satellite is measured as:

[0101] z k+1 =[α k+1 β k+1 ] T (31)

[0102] In the formula, z k+1 Indicates the quantity measurement, α k+1 and β k+1 are the pitch angle α of the maneuvering target relative to the observation satellite k+1 and azimuth βk+1 ;

[0103] The observation model is expressed as:

[0104] z k+1 =h(x k+1 )+ν k+1 (32)

[0105] Where:

[0106]

[0107] [x k+1 y k+1 z k+1 ] T is the position vector of the maneuvering target, [X k+1 Y k+1 Z k+1 ] T is the position vector of the observed star; ν k+1 is two-dimensional zero-expectation white noise, R k+1 is the measurement noise covariance matrix, σ z is the standard deviation of the measurement error, specifically set to σ z =10 -6 rad;h(x k+1 ) and ρ k+1 is an intermediate variable.

[0108] Step 2: Predict the system state and covariance matrix based on the dynamic model;

[0109] Predict the system status:

[0110]

[0111] In the formula, represents the predicted value of the state at time k+1; represents the estimated value of the state at time k; It represents the state prediction based on the state variable differential equation (29) of the space target. The specific process is to convert the state estimate value at time k into As the initial value, the numerical integrator is used to integrate equation (29) to obtain the state prediction value at the k+1th moment

[0112] Make predictions on the system covariance matrix:

[0113]

[0114] in, Represents the predicted value of the covariance matrix at time k+1; represents the estimated value of the covariance matrix at time k; Φk It represents the state transfer matrix of the system from time k to time k+1, and its differential equation is expressed as follows: the initial prior covariance matrix of the system is set to a diagonal matrix, the 1-σ position standard deviation is set to 100km, and the speed standard deviation is set to 0.1m / s.

[0115]

[0116] In the formula, represents the Jacobian matrix, which is specifically expressed as:

[0117]

[0118] Step 3: Calculate t k+1 The system measures the residual and covariance matrix at each moment and uses the indicator function Λ k+1 Characterizes the matching degree between the prediction results based on the dynamic model and the non-maneuvering model of the system; t k+1 The index function at time t1 to t k When the difference between the evaluation index and the mean value is greater than the threshold value, it is considered that the system has a sudden change in state, and the state noise compensation is performed on the system prediction result;

[0119] Step 3.1: Calculate t k+1 The time system measurement residual γ k+1 and the covariance matrix S k+1 , and through the indicator function Λ k+1 Characterize the degree of match between the prediction results of the system based on the dynamic model and the non-maneuvering model;

[0120] System measurement residual γ k+1 and the covariance matrix S k+1 It is expressed as:

[0121]

[0122] Among them, the matrix H k+1 It is expressed as:

[0123]

[0124] t k+1 Time index function Λ k+1 It is expressed as:

[0125]

[0126] Among them, tr(S k+1 ) represents the trace of the measurement residual covariance matrix;

[0127] Step 3.2: Obtain t through equation (43) k+1 The indicator function of the moment; calculate t1 to tk The mean of the indicator function; when the difference between the indicator function and the mean is greater than the threshold value, the system state changes suddenly, and the system predicts the covariance matrix For state noise compensation;

[0128] From time t1 to t k The mean of the indicator function at time is expressed as:

[0129] vale k =(vale k-1 ×(k-1)+Λ k ) / k (44)

[0130] Among them, the initial value of the indicator function mean is set to 0;

[0131] State Noise Compensation Q k+1 It is expressed as:

[0132]

[0133] Wherein, η is the threshold value for judging whether a sudden change occurs in the system. In this embodiment, the threshold value is set to η=20; D k+1 Indicates the state noise compensation of the position when the system state suddenly changes; V k+1 It represents the state noise compensation of the position when the system state suddenly changes;

[0134] State noise compensation for position when the system state suddenly changes k+1 It is expressed as:

[0135]

[0136] In the formula, the intermediate variable dr k+1 It is expressed as:

[0137]

[0138] In the formula, the intermediate variable Ω k+1 It can be expressed as:

[0139]

[0140] State transition matrix Φ k The scalar form of is:

[0141]

[0142] State noise compensation for speed when the system state suddenly changes V k+1 It is expressed as:

[0143]

[0144] In the formula, the intermediate variable dv1 k+1 ,dv2 k+1 and dv3 k+1 It is expressed as:

[0145]

[0146] In the formula, △v represents the estimated value of position deviation, which is specifically expressed as:

[0147]

[0148] Compensate the system prediction results for state noise:

[0149]

[0150] Step 4: Calculate the orbit determination algorithm gain based on the system state noise compensation result; Noise compensation results And the measurement results z k+1 Weighted output t k+1 The system state and covariance estimation values ​​are brought into step 2 for iterative solution to achieve high-precision orbit determination for maneuverable spacecraft.

[0151] Step 4.1 Orbit determination algorithm gain:

[0152]

[0153] Step 4.2 Based on the system prediction results Noise compensation results And the measurement results z k+1 Weighted output t k+1 System state and covariance estimates:

[0154]

[0155] Will Substitute into formula (8), and change P k+1 Substitute into equation (9) and iteratively solve the estimated value of the system state variable at each moment. In this embodiment, the end moment is set to k=60, thereby achieving high-precision orbit determination for the maneuverable spacecraft.

[0156] According to the above settings, the steps proposed in the present invention are used, and the results are simulated 100 times by Monte Carlo simulation to obtain statistical results. Figure 4 represents the position root mean square error of 100 Monte Carlo simulation results, Figure 5Indicates the velocity root mean square error of 100 Monte Carlo simulation results. In the figure, the red solid line represents the traditional IMM algorithm for maneuvering spacecraft orbit determination, the green dotted line represents the improved IMM algorithm assisted by prior knowledge, and the blue dotted line shows the orbit determination result using the method of this patent. This method has few limitations and higher orbit determination accuracy. It can realize the state estimation of maneuvering targets in various scenarios in space, and then accurately predict the trajectory of maneuvering spacecraft, so as to perform tasks such as intention discrimination and function prediction of the target spacecraft.

[0157] The specific description above further illustrates in detail the purpose, technical solutions and beneficial effects of the invention. It should be understood that the above is only a specific embodiment of the present invention, which is used to explain the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation, characterized in that: The steps include: Step 1: Establish a general dynamic model and an angle-only observation model for unknown mutations of random nonlinear dynamic systems; Step 2: t based on the kinetic model k+1 The system state and covariance matrix are used for prediction at each moment; Step 3: Calculate t k+1 The system measures the residual and covariance matrix at each moment and uses the indicator function Λ k+1 Characterizes the matching degree between the prediction results based on the dynamic model and the non-maneuvering model of the system; t k+1 The index function at time t1 to t k When the difference between the evaluation index and the mean value is greater than the threshold value, it is determined that the system has a sudden change in state, and the state noise compensation is performed on the system prediction result; Step 4: Calculate the orbit determination algorithm gain based on the system state noise compensation result; Noise compensation results And the measurement results z k+1 Weighted output t k+1 The system state and covariance estimates are brought into step 2 to iteratively solve the system state estimate at each moment, and high-precision orbit determination of the maneuverable spacecraft is achieved based on the state estimate at each moment.

2. The method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation according to claim 1, characterized in that: The implementation method of step one is: Step 1.1 establish a general dynamic model for unknown mutations of random nonlinear dynamic systems; The state variable differential equation of the space target is expressed as: In the formula, Table system state vector; r k represents the position vector of the spatial target at time k, [x k y k z k ] T is its scalar form; v k represents the velocity vector of the space target at time k, is its scalar form; Establish a general dynamics model for unknown mutations in random nonlinear dynamic systems: x k+1 =F(x k +d t,k Gu)+Bw k (2) Among them, δ t,k represents the Kronecker symbol. When t≠k, δ t,k =0, when t = k, δ t,k =1; is a constant matrix; represents the unknown input vector of the system; is a constant matrix; represents white noise with an expected value of 0; δ t,k Gu represents the random maneuvers that occur in the system; F(x k +δ t,k Gu) represents the state recursion based on the state variable differential equation (1) of the space target. The specific process is to transform the system state x at time k into k +δ t,k Gu is used as the initial value, and the numerical integrator is used to integrate equation (1) to obtain the recursive result x based on dynamics at the k+1th moment. k+1 ; Step 1.2: Establish an angle measurement observation model only; The angle information of the maneuvering target relative to the observation satellite is measured as: z k+1 ＝[α k+1 β k+1 ] T (3) In the formula, z k+1 Indicates the quantity measurement, α k+1 and β k+1 are the pitch angle α of the maneuvering target relative to the observation satellite k+1 and azimuth β k+1 ; The observation model is expressed as: z k+1 =h(x k+1 )+ν k+1 (4) Where: [x k+1 y k+1 z k+1 ] T is the position vector of the maneuvering target, [X k+1 Y k+1 Z k+1 ] T is the position vector of the observed star; ν k+1 is a two-dimensional zero-expectation white noise, Rk +1 is the measurement noise covariance matrix, σz is the standard deviation of the measurement error; h(xk +1 ) and ρ k+1 is an intermediate variable.

3. The method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation according to claim 2, characterized in that: The implementation method of step 2 is: Predict the system status: In the formula, represents the predicted value of the state at time k+1; represents the estimated value of the state at time k; It represents the state prediction based on the state variable differential equation (1) of the space target. The specific process is to convert the state estimate value at time k into As the initial value, the numerical integrator is used to integrate equation (1) to obtain the state prediction value at the k+1th moment: Make predictions on the system covariance matrix: in, Represents the predicted value of the covariance matrix at time k+1; represents the estimated value of the covariance matrix at time k; Φ k represents the system state transfer matrix from time k to time k+1, and its differential equation is expressed as: In the formula, Represents the Jacobian matrix, which is specifically expressed as:

4. The method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation according to claim 3, characterized in that: The implementation method of step three is: Step 3.1: Calculate t k+1 The time system measurement residual γ k+1 and the covariance matrix S k+1 , and through the indicator function Λ k+1 Characterize the degree of match between the prediction results of the system based on the dynamic model and the non-maneuvering model; System measurement residual γ k+1 and the covariance matrix S k+1 It is expressed as: Among them, the matrix H k+1 It is expressed as: t k+1 Time index function Λ k+1 It is expressed as: Among them, tr(S k+1 ) represents the trace of the measurement residual covariance matrix; Step 3.2: Obtain t through equation (15) k+1 The indicator function of the moment; calculate t1 to t k The mean of the indicator function; when the difference between the indicator function and the mean is greater than the threshold value, the system state changes suddenly, and the system predicts the covariance matrix For state noise compensation; From time t1 to t k The mean of the indicator function at time is expressed as: vale k =(vale k-1 ×(k-1)+Λ k ) / k (16) Among them, the initial value of the indicator function mean is set to 0; State Noise Compensation Q k+1 It is expressed as: Where η is the threshold value for judging whether a mutation occurs in the system; D k+1 Indicates the state noise compensation of the position when the system state suddenly changes; V k+1 It represents the state noise compensation of the position when the system state suddenly changes; State noise compensation for position when the system state suddenly changes k+1 It is expressed as: In the formula, the intermediate variable dr k+1 It is expressed as: In the formula, the intermediate variable Ω k+1 It can be expressed as: State transition matrix Φ k+1 The scalar form of is: State noise compensation for speed when the system state suddenly changes V k+1 It is expressed as: In the formula, the intermediate variable dv1 k+1 ,dv2 k+1 and dv3 k+1 It is expressed as: In the formula, △v represents the estimated value of position deviation, which is specifically expressed as: Compensate the system prediction results for state noise:

5. The method for orbit determination of a maneuvering spacecraft with adaptive state noise compensation according to claim 4, characterized in that: Step 4 is implemented as follows: Step 4.1 Orbit determination algorithm gain: Step 4.2 Based on the system prediction results Noise compensation results And the measurement results z k+1 Weighted output t k+1 System state and covariance estimates: Will Substitute into formula (8), and change P k+1 Substitute into equation (9), iteratively solve to obtain the estimated value of the system state at each moment, set the observation stop time as the orbit determination termination time, and realize high-precision orbit determination of the maneuverable spacecraft based on the estimated value of the system state at each moment.

Images