A method for tracking non-cooperative pulse-maneuvering targets in space.

Abstract

The application discloses a space non-cooperative impulsive maneuvering target tracking method and relates to the technical field of aerospace. According to target ranging and speed information obtained by a ground station, a target state prediction value is obtained through a nonlinear state equation propagation, and a measurement prediction value is generated in combination with a measurement model; a Huber robust cost function is used to replace a mean square error cost function of a traditional filter, a performance index function is reconstructed, a corrected measurement noise covariance matrix is solved by minimizing the function, then the target impulsive maneuvering is regarded as an unknown input of a nonlinear system, an impulsive estimation value is obtained by using an unbiased minimum input estimation method, a detection index is constructed by calculating Mahalanobis distance of the impulsive estimation error, a detection threshold is determined by using an error upper limit, and whether the impulsive maneuvering occurs is judged; finally, the impulsive estimation value and error covariance are compensated to a filter state prediction process, a Kalman gain matrix is calculated, a target state estimation value and a covariance estimation value are updated in combination with measurement information, and the filtering is completed.

Description

Technical Field

[0001] This invention relates to the field of aerospace technology, and more specifically to a method for tracking non-cooperative pulse-maneuvering targets in space. Background Technology

[0002] As space activity continues to increase, the space environment is becoming increasingly crowded and competitive, leading to a corresponding increase in the demand for space situational awareness (SSA). Tracking non-cooperative maneuvering targets in space is fundamental to SSA, and accurately determining and predicting the orbits of space targets is necessary to ensure the safety of satellites in orbit.

[0003] One of the key challenges in tracking non-cooperative maneuvering targets in space is determining and predicting the trajectory of a maneuvering spacecraft solely based on sensor-collected observation data in the absence of prior information about the maneuver. Traditional tracking filters do not compensate for spacecraft maneuvers, resulting in a mismatch between nominal dynamics and actual maneuvers. Consequently, measurements following unknown maneuvers are rejected by the filter, leading to degraded tracking performance and even divergence. Furthermore, traditional Kalman filters are derived based on the Gaussian noise assumption; however, measurements obtained from sensors often deviate from this assumption. In actual tracking, measurement noise contaminated by outliers follows a heavy-tailed non-Gaussian distribution. Accurately tracking non-cooperative maneuvering spacecraft becomes a significant challenge when prior information about the target maneuver is lacking and measurements are contaminated by non-Gaussian noise.

[0004] Currently, research on tracking non-cooperative maneuvering spacecraft has been conducted both domestically and internationally from different perspectives. Existing literature (Jiang, Y., Ma, P., and Baoyin, H., “Residual-Normalized Strong TrackingFilter for Tracking a Noncooperative Maneuvering Spacecraft,” Journal of Guidance, Control, and Dynamics, Vol. 42, No. 10, 2019, pp. 2304–2309) uses a residual-normalized strong tracking filter to detect and compensate for pulse maneuvers, achieving real-time tracking of maneuvering spacecraft. However, this method is only applicable when the measurement noise follows a Gaussian distribution. When the measurement noise is contaminated and exhibits a non-Gaussian distribution, frequent false detections of maneuvers occur, leading to a decrease in filtering performance. Existing literature (Qiang, Q., Lin, B., Liu, Y., Lin, X., and Wang, S., “Robust UKF Orbit Determination Method with Time-Varying Forgetting Factor for Angle / Range-Based Integrated Navigation System,” Chinese Journal of Aeronautics, Vol. 37, No. 11, 2024, pp. 420–434) proposes an adaptive robust unscented Kalman filter based on a time-varying forgetting factor, achieving accurate tracking of spacecraft under non-Gaussian noise. However, this method can only suppress the measurement model uncertainty caused by non-Gaussian measurement noise. When the target spacecraft maneuvers, the filter diverges rapidly and loses its tracking capability. Existing literature (Wang, Y., Sun, S., and Li, L. “Adaptively Robust Unscented Kalman Filter for Tracking a Maneuvering Vehicle.” Journal of Guidance, Control, and Dynamics, Vol. 37, No. 5, 2014, pp. 1696–1701) proposes an adaptive robust unscented Kalman filter, employing a strong tracking filter (STF) structure to reduce the impact of errors in the dynamics and measurement models on tracking accuracy. However, this method cannot estimate the maneuver vector, and when the target maneuver amplitude is large, redundancy in robustness occurs, causing a decrease in filter tracking performance. Summary of the Invention

[0005] To address this, the present invention provides a method for tracking non-cooperative pulse maneuvering targets in space, thereby solving the problems of false detection of maneuvers, degraded tracking performance, or even divergence of the filter under non-Gaussian noise conditions in existing technologies.

[0006] To achieve the above objectives, the present invention provides the following technical solution: a method for tracking a non-cooperative pulse maneuvering target in space, comprising the following steps:

[0007] Step 1: Based on the ranging and velocity information of the non-cooperative target satellite obtained from the ground station, the current state prediction value of the target is obtained through nonlinear state equation propagation, and the measurement prediction value is calculated based on the measurement model;

[0008] Step 2: Replace the mean square error cost function of the traditional filter with the Huber robust cost function and reconstruct the index function;

[0009] Step 3: Obtain the modified measurement noise covariance matrix by minimizing the reconstructed robustness index function;

[0010] Step 4: Treat the target pulse maneuver as an unknown input to a nonlinear system, and calculate the maneuver estimate using the unbiased minimum input estimation method;

[0011] Step 5: Using the Mahalanobis distance of the maneuver estimation error as the computer motion detection index, the upper bound of the computer motion estimation error is used to obtain the maneuver detection threshold. The maneuver estimate is compared with the detection threshold to determine whether the maneuver has occurred.

[0012] Step 6: Compensate the maneuver estimate and the maneuver estimate error covariance into the state prediction process of the filter to obtain the state prediction value and state covariance prediction value containing the maneuver.

[0013] Step 7: Calculate the Kalman gain matrix at the current time, update the target's state estimate and state covariance estimate based on the measurement information, and repeat the process time-by-time until the filtering ends.

[0014] Preferably, step one specifically includes:

[0015] Step 1.1: Ground stations collect measurements of non-cooperative pulse maneuvering targets in space using the North-Sky-East coordinate system. ,in The measurement value at time k includes the distance measurement. and speed measurement ; Let K be the measurement noise of the ground station at time k. The measurement noise is zero-mean Gaussian white noise, which satisfies the following conditions: ,in, To measure the noise covariance matrix;

[0016] Step 1.2, the state-space model of the nonlinear system is represented as follows:

[0017] ;

[0018] in This indicates the spacecraft's status, including its position in the J2000 coordinate system. and speed Where x, y, and z represent the positions of the three axes. , , Indicates the three-axis velocity. The nominal dynamic model of the target spacecraft. For unknown maneuver acceleration, The input coefficient matrix for acceleration. Let k be the state vector at time k. For measurement model;

[0019] Step 1.3: Based on the state equation, calculate the predicted state value and the predicted state covariance value at the current time based on the state estimate of the system at the previous time step, as shown in the following formula:

[0020] ;

[0021] in This is the state estimate at time k-1. The state prediction value at time k. Let k be the predicted state covariance value at time k. This represents the numerical integral of the nonlinear dynamic model, used to propagate the state to the next time step. This represents the state transition matrix from time k-1 to time k. This represents the estimated state covariance at time k-1. This represents the process noise covariance matrix.

[0022] Preferably, step two specifically includes: considering the non-Gaussian nature of the measurement noise, designing the filter's performance index function using the Huber function. For measurement outliers, the Huber function is equivalent to the l1 norm of the residual, possessing the ability to suppress measurement outliers. The standardized state error is denoted as... and standardized measurement error ,in Representing the measurement mapping matrix, the performance index function is constructed as follows:

[0023] ;

[0024] in, The state estimate at time k. for The i-th component, for The j-th component, The Huber function is defined as follows:

[0025] ;

[0026] in, The independent variable of the function, This is an adjustable parameter.

[0027] Preferably, step three specifically includes:

[0028] Minimize the filter performance index function constructed in step two to obtain the modified measurement noise covariance. The solution process is as follows:

[0029] ;

[0030] in, , Represents Huber function The derivative of the function is defined. ,matrix The above equation can be rewritten in matrix form:

[0031] ;

[0032] Huber weight function The specific expression is:

[0033] ;

[0034] Will and Substituting the values, we get the following formula:

[0035] ;

[0036] The modified measurement noise covariance matrix is ​​obtained:

[0037] .

[0038] Preferably, step four specifically includes: treating the target spacecraft maneuver as an unknown input to a nonlinear system, and calculating the maneuver estimate using the unbiased minimum input estimation method. The specific calculation process is as follows:

[0039] After linearizing the nonlinear system, we obtain the following linear discrete-time system:

[0040] ;

[0041] ;

[0042] in The coefficient matrix, For motor acceleration, For process noise, For measuring the mapping matrix;

[0043] The basic framework of unbiased minimum input estimation is as follows:

[0044] ;

[0045] ;

[0046] ;

[0047] ;

[0048] in For the estimated value of maneuver, The state vector after maneuver compensation. and This is the gain matrix;

[0049] By using the modified measurement noise covariance matrix obtained in step three Substituting the values ​​into the calculation, we obtain the estimated error covariance matrix:

[0050] ;

[0051] The gain matrix is ​​calculated. and maneuver estimates :

[0052] ;

[0053] .

[0054] Preferably, step five specifically includes: using the Mahalanobis distance of the maneuver estimation error to calculate the maneuver index, calculating the upper bound of the Mahalanobis distance during the process where no maneuver occurred and using it as the maneuver detection threshold, comparing the maneuver estimation value with the detection threshold to determine whether a maneuver has occurred, the specific process is as follows:

[0055] Based on Mahalanobis distance, the motion detection index is defined as follows:

[0056] ;

[0057] If the spacecraft does not perform a maneuver at time k-1, i.e., the actual maneuver value Then the above formula can be rewritten as:

[0058] ;

[0059] in Indicates the error in maneuver estimation. For motor detection, represents the Mahalanobis distance of the acceleration estimation error. Expressing expectations, The covariance, representing the acceleration estimation error, is expressed as follows:

[0060] ;

[0061] Analyzing the statistical characteristics of the target's acceleration estimation error during the spacecraft's non-maneuvering phase reveals that its probability distribution approximately follows a Gaussian distribution. Therefore, the upper bound of the acceleration estimation error during the non-maneuvering phase is determined by the following formula:

[0062] ;

[0063] in, , and These represent the standard deviations of the acceleration estimation error along the x, y, and z axes, respectively.

[0064] To simplify the threshold analysis, the covariance of the acceleration estimation error is approximated as a diagonal matrix, assuming that the components of the estimation error are independent of each other:

[0065] ;

[0066] Based on the above analysis, the approximate motion detection threshold can be calculated as follows:

[0067] ;

[0068] When the maneuver detection index is less than the detection threshold, it is considered that the maneuver has not occurred, and the estimated maneuver value is... and its error covariance Set to zero; conversely, when the maneuver detection index is greater than the detection threshold, a maneuver is considered to have occurred, and the estimated maneuver value is... and its error covariance The following formula is used to calculate:

[0069] ;

[0070] .

[0071] Preferably, step six specifically includes: adjusting the maneuver estimate... and the covariance of the estimation error During the state prediction process of the compensation filter, the state prediction value including the maneuver and the state covariance prediction value are obtained, as shown in the formula:

[0072]

[0073] .

[0074] Preferably, step seven specifically includes: calculating the Kalman gain matrix at the current time, and updating the target's state estimate and state covariance estimate at the next time step based on the measurement information, using the following formula:

[0075] ;

[0076] ;

[0077] ;

[0078] in Represents the identity matrix. Represents the filter gain matrix;

[0079] For each step in the filtering process, it is necessary to calculate the maneuver estimate and maneuver detection index corresponding to the current step, and determine whether the maneuver detection index of the current step is greater than the detection threshold, i.e. whether a pulse maneuver has occurred. The maneuver estimate and maneuver estimate error covariance are then compensated into the filter, and the state estimate and state covariance estimate are updated until the filtering process ends.

[0080] This invention is based on an input estimation robust filter. Through the Huber robust performance index function, it corrects the measurement noise covariance, suppresses the influence of measurement outliers and heavy-tailed measurement noise, and utilizes an unbiased minimum input estimation method to detect and estimate impulsive maneuvers. It also compensates for the mismatch between the nominal dynamic model and the actual maneuver caused by unknown maneuvers, effectively solving the problems of decreased tracking performance, filter divergence, and false detection of maneuvers in existing Kalman filters under non-Gaussian noise conditions. The Huber robust cost function improves the filter's ability to suppress measurement model uncertainties, and the unbiased minimum input estimation method enables rapid and accurate detection and estimation of maneuvers. Simultaneously, the estimated maneuver is compensated into the filtering process, reducing the mismatch between the nominal and actual dynamic models caused by impulsive maneuvers and enhancing the filter's convergence capability. The state estimate and state covariance estimate of the target at the next moment are obtained by updating the Kalman gain matrix and measurement information at the current moment, realizing real-time tracking of impulsive maneuvering spacecraft under non-Gaussian conditions. The overall approach is novel and has strong innovation and engineering application prospects. Attached Figure Description

[0081] Figure 1 This is a flowchart of a method for tracking a non-cooperative pulse maneuvering target in space, provided by the present invention.

[0082] Figure 2 This diagram illustrates a comparison of the position estimation errors of the ARKFIE tracking method of this invention with existing adaptive robust unscented Kalman filters ARUKF and TFFARUKF based on time-varying forgetting factors during the tracking process.

[0083] Figure 3 This diagram illustrates a comparison of the velocity estimation errors of the ARKFIE tracking method described in this invention with existing adaptive robust unscented Kalman filters (ARUKF) and TFFARUKF based on time-varying forgetting factors during the tracking process.

[0084] Figure 4 This is a schematic diagram illustrating the changes in motion detection indicators obtained by the ARKFIE tracking method described in this invention during the tracking process.

[0085] Figure 5 This is a schematic diagram illustrating the changes in the maneuver estimates obtained by the ARKFIE tracking method described in this invention during the tracking process. Detailed Implementation

[0086] The following specific embodiments illustrate the implementation of the present invention. Those skilled in the art can easily understand other advantages and effects of the present invention from the content disclosed in this specification. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0087] This invention proposes a method for tracking non-cooperative pulse maneuvering targets in space (i.e., ARKFIE), comprising the following steps:

[0088] Step 1: Based on the ranging and velocity information of the non-cooperative target satellite obtained from the ground station, the current state prediction value of the target is obtained through nonlinear state equation propagation, and the measurement prediction value is calculated based on the measurement model;

[0089] Step 2: Replace the mean square error cost function of the traditional filter with the Huber robust cost function and reconstruct the index function;

[0090] Step 3: Obtain the modified measurement noise covariance matrix by minimizing the reconstructed robustness index function;

[0091] Step 4: Treat the target pulse maneuver as an unknown input to a nonlinear system, and calculate the maneuver estimate using the unbiased minimum input estimation method;

[0092] Step 5: Using the Mahalanobis distance of the maneuver estimation error as the computer motion detection index, the upper bound of the computer motion estimation error is used to obtain the maneuver detection threshold. The maneuver estimate is compared with the detection threshold to determine whether the maneuver has occurred.

[0093] Step 6: Compensate the maneuver estimate and the maneuver estimate error covariance into the state prediction process of the filter to obtain the state prediction value and state covariance prediction value containing the maneuver.

[0094] Step 7: Calculate the Kalman gain matrix at the current time, update the target's state estimate and state covariance estimate based on the measurement information, and repeat the process time-by-time until the filtering ends.

[0095] The following will provide a detailed analysis and explanation of each step:

[0096] Step 1: Establish a simulation scenario, which requires determining the pre-maneuver orbit, pulse maneuver time, pulse maneuver vector, and post-maneuver orbit of the non-cooperative target satellite. Ground stations collect measurements of the non-cooperative pulse maneuvering target in the North-Sky-East coordinate system, specifically including... ,in The measurement value at time k includes the distance measurement. and speed measurement . Let K be the measurement noise of the ground station at time k. The measurement noise is zero-mean Gaussian white noise, which satisfies the following conditions: ,in To measure the noise covariance matrix.

[0097] The state-space model of a nonlinear system is expressed as follows:

[0098] ;

[0099] in This indicates the spacecraft's status, including its position in the J2000 coordinate system. and speed Where x, y, and z represent the positions of the three axes, and v x v y v z Indicates the three-axis velocity. The nominal dynamic model of the target spacecraft. For unknown maneuver acceleration, The input coefficient matrix for acceleration. Let k be the state vector at time k. For measurement model.

[0100] Based on the state equation, the predicted state value and the predicted state covariance value at the current moment can be calculated from the state estimate of the system at the previous moment, as shown in the following formula:

[0101] ;

[0102] in This is the state estimate at time k-1. The state prediction value at time k. Let k be the predicted state covariance value at time k. This represents the numerical integral of the nonlinear dynamic model, used to propagate the state to the next time step. This represents the state transition matrix from time k-1 to time k. This represents the estimated state covariance at time k-1. This represents the process noise covariance matrix.

[0103] Step 2: Considering the non-Gaussian nature of measurement noise, design the filter's performance index function using the Huber function. For measurement outliers, the Huber function is equivalent to the l1 norm of the residual, exhibiting suppression capability for measurement outliers. For non-measurement outliers, the Huber function is equivalent to the traditional l2 norm of the residual. Let the standardized state error be... and standardized measurement error ,in Representing the measurement mapping matrix, the performance index function can be constructed as follows:

[0104] ;

[0105] in, The state estimate at time k. for The i-th component, for The j-th component, The Huber function is defined as follows:

[0106] ;

[0107] in The independent variable of the function, This is an adjustable parameter.

[0108] Step 3: Minimize the filter robustness index function constructed in Step 2 to obtain the modified measurement noise covariance. The solution process is as follows:

[0109] ;

[0110] in , Represents Huber function The derivative of the function is defined. ,matrix The above formula can be rewritten in matrix form:

[0111] ;

[0112] Huber weight function The specific expression is:

[0113] ;

[0114] Will and Substituting the values, we get the following formula:

[0115] ;

[0116] The modified measurement noise covariance matrix can be obtained as follows:

[0117] ;

[0118] Step 4: After linearizing the nonlinear system, we can obtain the following linear discrete-time system:

[0119] ;

[0120] ;

[0121] in The coefficient matrix, For motor acceleration, For process noise, For measuring the mapping matrix.

[0122] The basic framework of unbiased minimum input estimation is as follows:

[0123] ;

[0124] ;

[0125] ;

[0126] ;

[0127] in For the estimated value of maneuver, The state vector after maneuver compensation. and This is the gain matrix.

[0128] By using the modified measurement noise covariance matrix obtained in step three Substituting the values ​​into the calculation, we obtain the estimated error covariance matrix:

[0129] ;

[0130] The gain matrix is ​​calculated. and maneuver estimates :

[0131] ;

[0132] ;

[0133] Step 5: Based on Mahalanobis distance, the following motion detection index can be defined:

[0134] ;

[0135] If the spacecraft does not perform a maneuver at time k-1, i.e., the actual maneuver value Then the above formula can be rewritten as:

[0136] ;

[0137] in Indicates the error in maneuver estimation. For motor control, the Mahalanobis distance represents the acceleration estimation error. Expressing expectations, The covariance, representing the acceleration estimation error, is expressed as follows:

[0138] ;

[0139] Analyzing the statistical characteristics of the target's acceleration estimation error during the spacecraft's non-maneuvering phase reveals that its probability distribution approximately follows a Gaussian distribution. Therefore, the upper bound of the acceleration estimation error during the non-maneuvering phase can be determined by the following formula:

[0140] ;

[0141] in, , and These represent the standard deviations of the acceleration estimation error along the x, y, and z axes, respectively.

[0142] To simplify the threshold analysis, the covariance of the acceleration estimation error can be approximated as a diagonal matrix, assuming that the components of the estimation error are independent of each other:

[0143] ;

[0144] Based on the above analysis, the motion detection threshold can be approximately calculated:

[0145] ;

[0146] When the maneuver detection index is less than the detection threshold, it is considered that the maneuver has not occurred, and the estimated maneuver value is... and its error covariance Set to zero. Conversely, when the maneuver detection index is greater than the detection threshold, a maneuver is considered to have occurred, and the estimated maneuver value is... and its error covariance The following formula is used to calculate:

[0147] ;

[0148] ;

[0149] Step 6: Estimate the maneuver value and the covariance of the estimation error During the state prediction process of the compensation filter, the state prediction value including the maneuver and the state covariance prediction value are obtained, as shown in the formula:

[0150] ;

[0151] ;

[0152] Step 7: Calculate the Kalman gain matrix at the current time step, and update the target's state estimate and state covariance estimate at the next time step based on the measurement information. The formula is as follows:

[0153] ;

[0154] ;

[0155] ;

[0156] in Represents the identity matrix. Represents the filter gain matrix;

[0157] For each step in the filtering process, it is necessary to calculate the maneuver estimate and maneuver detection index corresponding to the current step, and determine whether the maneuver detection index of the current step is greater than the detection threshold, i.e. whether a pulse maneuver has occurred. The maneuver estimate and maneuver estimate error covariance are then compensated into the filter, and the state estimate and state covariance estimate are updated until the filtering process ends.

[0158] The feasibility of this invention is illustrated by the following examples.

[0159] The following calculation conditions and technical parameters are set:

[0160] 1. The initial position and velocity of the non-cooperative target satellite in the J2000 coordinate system are [-5726510.59, -3472068.72, 1156512.92] m and [-1594.7645, -150.6295, -7498.2957] m / s, respectively;

[0161] 2. The initial position estimation error of the target is [1000, 1000, 1000] m, and the initial velocity estimation error is [1, 1, 1] m / s;

[0162] 3. A total of three ground stations were set up, located in the geodetic coordinate system at longitude 38.6°, latitude 87.2°, altitude 1100.0 m; longitude 40.1°, latitude 98.2°, altitude 1200.0 m; and longitude 35.4°, latitude 97.6°, altitude 1000.0 m.

[0163] 4. The standard deviation of the distance measurement noise at the ground station is 3 m, the standard deviation of the speed measurement noise is 0.01 m / s, and the measurement interval is 1 s;

[0164] 5. The measurement noise conforms to a non-Gaussian noise distribution. The probability density function of non-Gaussian noise is defined as follows:

[0165] ;

[0166] in and These are the standard deviations of independent Gaussian distributions. Set as , Here, represents the disturbance parameter, indicating the degree of noise pollution. .

[0167] 6. The tracking process began at 00:19:00 on December 5, 2024, and lasted for 1000 seconds. At 500 seconds, the target spacecraft generated a pulse maneuver with an amplitude of 100 m / s along the velocity direction. A measurement outlier was observed at 500 seconds. .

[0168] Figure 2 , Figure 3 , Figure 4 , Figure 5 Simulation results for the present invention are presented respectively.

[0169] Figure 2 and Figure 3The target position estimation error and velocity estimation error obtained by the tracking method for a non-cooperative pulse maneuvering target in space described in this invention are presented respectively. It can be seen that, compared with existing adaptive robust unscented Kalman filters (ARUKF) and adaptive robust unscented Kalman filters based on time-varying forgetting factors (TFFARUKF), the tracking method of this invention has higher estimation accuracy. Compared with existing methods, the tracking method disclosed in this invention reduces the root mean square error of position by approximately 99.5% and the root mean square error of velocity by approximately 74.2% throughout the tracking process, demonstrating the superiority of the method described in this invention.

[0170] Figure 4 The changes in maneuver detection indicators obtained by the tracking method described in this invention are presented. It can be seen that after a maneuver occurs, the proposed filter detects the maneuver promptly and accurately, and there are no false detections of maneuver during periods when no maneuver occurs.

[0171] Figure 5 The changes in the maneuver estimates obtained using the tracking method described in this invention are presented. It can be seen that when a maneuver occurs, the unbiased minimum input estimation method described in this invention accurately estimates the maneuver vector, demonstrating that the method described in this invention can effectively estimate the direction and magnitude of the maneuver.

[0172] Therefore, by using the method proposed in this invention, tracking of space non-cooperative pulse maneuvering spacecraft under non-Gaussian noise conditions can be achieved using only ground station measurement data, and the tracking accuracy is significantly improved compared to existing methods. No false maneuvering detections occur in the non-maneuvering segment, and maneuvering can be detected and maneuvering vectors estimated simultaneously.

[0173] Although the present invention has been described in detail above with general descriptions and specific embodiments, modifications or improvements can be made to it, which will be obvious to those skilled in the art. Therefore, all such modifications or improvements made without departing from the spirit of the present invention fall within the scope of protection claimed by the present invention.

Claims

1. A method of tracking a spatially non-cooperative impulsive maneuvering target, characterized in that: The steps comprise the following: Step one: according to the ground station, the ranging and velocity information of the non-cooperative target satellite is obtained, the state prediction value of the target at the current time is obtained through the propagation of the nonlinear state equation, and the measurement prediction value is calculated according to the measurement model; Step two: according to the Huber robust cost function, the mean square error cost function of the traditional filter is replaced, and the index function is reconstructed; Step three: the modified measurement noise covariance matrix is obtained by minimizing the reconstructed robust performance index function; Step four: the target pulse maneuver is regarded as an unknown input of the nonlinear system, and the maneuver estimation value is calculated by the unbiased minimum input estimation method; Step five: the maneuver detection index is calculated by using the Mahalanobis distance of the maneuver estimation error, the maneuver detection threshold is obtained by calculating the upper bound of the maneuver estimation error, and the maneuver estimation value is compared with the detection threshold to judge whether the maneuver occurs; Step six: the maneuver estimation value and the maneuver estimation error covariance are compensated to the state prediction process of the filter to obtain the state prediction value and the state covariance prediction value containing the maneuver; Step seven: the Kalman gain matrix at the current time is calculated, and the state estimation value and the state covariance estimation value of the target are updated according to the measurement information, and the time is recursively propagated until the filtering is finished.

2. The method of claim 1, wherein: The step one specifically comprises: Step 1.

1. The ground station collects the measurements of the space non-cooperative impulsive maneuvering target in the North-Earth-East coordinate system where is the measurement at time k, which includes the range and the velocity ; is the measurement noise of the ground station at time k, which is zero-mean Gaussian white noise and satisfies where is the measurement noise covariance matrix; Step 1.2, the state space model of the nonlinear system is represented as follows: ； wherein denotes the spacecraft state, including its position in the J2000 coordinate system and velocity where x, y, z denote the three-axis position, , , denotes the three-axis velocity, is a nominal dynamic model of the target spacecraft, is an unknown maneuvering acceleration, is an input coefficient matrix for the acceleration, is the state vector at time k, is a measurement model; Step 1.3, based on the state equation, the state prediction value and the state covariance prediction value at the current time are calculated according to the state estimation of the system at the previous time, as shown in the following formula: ； wherein is a state estimate at time k - 1, is a state prediction at time k, is a state covariance prediction at time k, denotes a numerical integration of the nonlinear dynamics model for propagating the state to the next time instant, denotes a state transition matrix from time k - 1 to time k, denotes a state covariance estimate at time k - 1, denotes a process noise covariance matrix.

3. The method of claim 2, wherein: The step two specifically comprises: considering the non-Gaussian characteristic of the measurement noise, designing a performance index function of the filter by a Huber function, the Huber function is equivalent to the l1 norm of the residual for the measurement outliers, has the suppression ability for the measurement outliers, and recording the normalized state error as and the normalized measurement error as wherein The measurement mapping matrix is denoted as H, and the performance index function is constructed as follows: ； wherein, is the state estimation value at time k, is is the i-th component of is is the j-th component of is the Huber function, defined as follows: ； wherein is a function argument, is an adjustable parameter.

4. The method of claim 3, wherein: The step three specifically comprises: The modified measurement noise covariance is obtained by minimizing the filter performance index function constructed in step two, and the solving process is as follows: ； where , denotes the derivative of the Huber function , the function , the matrix , the above equation is rewritten in matrix form: ； Huber weight function The specific expression for Huber weight function is: ； Substituting and yields the following equation: ； The modified measurement noise covariance matrix is obtained: 。 5. The method of claim 4, wherein: The step four specifically comprises: the target spacecraft maneuver is regarded as an unknown input of the nonlinear system, and the maneuver estimation value is calculated by the unbiased minimum input estimation method, and the specific calculation process is as follows: After linearizing the nonlinear system, the following linear discrete-time system is obtained: ； ； wherein is a coefficient matrix, is a motor acceleration, is a process noise; The basic framework of the unbiased minimum input estimation is as follows: ； ； ； ； wherein is a maneuver estimate, is a maneuver-compensated state vector, and is a gain matrix; The modified measurement noise covariance matrix obtained in step three is substituted into the equation Substituting into the calculation, the estimated error covariance matrix is obtained: ； The gain matrix is calculated and the motor estimate : ； 。 6. The method of claim 5, wherein: The step five specifically comprises: the maneuver index is calculated by using the Mahalanobis distance of the maneuver estimation error, the upper bound of the Mahalanobis distance in the non-manipulation process is calculated and taken as the maneuver detection threshold, and the maneuver estimation value is compared with the detection threshold to judge whether the maneuver occurs, and the specific process is as follows: Based on the Mahalanobis distance, the maneuver detection index is defined as follows: ； If no maneuver occurs at time k-1, i.e. the true value of the maneuver is zero, then the above equation is rewritten as: ； wherein denotes the maneuvering estimation error, is a maneuver detection indicator, denoting the Mahalanobis distance of the acceleration estimation error, denotes the expectation, denotes the covariance of the acceleration estimation error, which is expressed as follows: ； In the non-manipulation stage of the spacecraft, the statistical characteristics of the acceleration estimation error of the target are analyzed, and it is found that the probability distribution approximately obeys the Gaussian distribution, therefore, the upper bound of the acceleration estimation error in the non-manipulation stage is determined by the following formula: ； wherein, , and respectively represent the standard deviation of the acceleration estimation error along the x, y and z axis directions; In order to simplify the threshold analysis, the covariance of the acceleration estimation error is approximated as a diagonal matrix, and it is assumed that each component of the estimation error is independent of each other: ； Through the above analysis, the maneuver detection threshold is approximately calculated as follows: ； When the maneuver detection index is less than the detection threshold, no maneuver is considered to have occurred, and the maneuver estimate and its error covariance is set to zero; conversely, when the maneuver detection index is greater than the detection threshold, a maneuver is considered to have occurred, and the maneuver estimate and its error covariance is calculated by the following equation: ； 。 7. The method of claim 6, wherein: The step six specifically includes: compensating the maneuver estimation value and the maneuver estimation error covariance to a state prediction process of the filter, to obtain a state prediction value and a state covariance prediction value containing the maneuver, and the formula is: ； 。 8. The method of claim 7, wherein: The step seven specifically comprises: the Kalman gain matrix at the current time is calculated, and the state estimation value and the state covariance estimation value of the target at the next time are updated according to the measurement information, and the formula is as follows: ； ； ； wherein denotes the identity matrix, denotes a filter gain matrix; The maneuver estimation value and the maneuver detection index corresponding to each step in the filtering process are calculated, and it is judged whether the maneuver detection index of the current step is greater than the detection threshold, that is, whether the pulse maneuver occurs, and the maneuver estimation value and the maneuver estimation error covariance are compensated to the filter, and the state estimation value and the state covariance estimation value are updated until the filtering process ends.

Images