Method for estimating relative attitude of failed spacecraft under constraint of computing resources

Abstract

The invention relates to the technical field of spacecraft navigation, and discloses a method for estimating the relative attitude of a failed spacecraft under the constraint of computing resources, and the method comprises the steps: carrying out the preprocessing of measurement data, obtaining the observation data of the relative attitude of a camera based on an ellipse recognition algorithm, and obtaining the measurement data of the relative attitude of a laser radar based on an iterative nearest point algorithm. Pseudo measurement factors are designed for the camera and the laser radar, and variational integral prediction factors are designed for attitude increments at adjacent moments of the failed spacecraft; factor graph optimization based on a moving window is carried out, and the relative attitude of the failed spacecraft is estimated; the momentum conservation characteristics of the failed spacecraft can be mined through a designed variational integral predictive factor, the constraint relation of attitude increments at adjacent moments is constructed, and the problem that the estimation precision of an existing scheme depends on high-frequency measurement input is solved; and low-frequency measurement data of the camera and the laser radar are fused through a built failure spacecraft full-state factor graph model, so that the relative pose estimation of the non-cooperative target under the low computing power load is realized.

Description

Technical Field

[0001] The present invention relates to the technical field of spacecraft navigation, and more particularly to a method for estimating the relative attitude of a failed spacecraft under computational resource constraints. Background Art

[0002] The proximity operation of a failed spacecraft is an urgent need to ensure the healthy development of the space industry. Since the failed target has significant non-cooperative characteristics, the use of multiple vision sensors for fusion observation can significantly improve the detection and evaluation capabilities, but there are still many difficult challenges in its observation and estimation. On the one hand, the computing performance of the on-board platform is insufficient to meet the high-frequency sampling of vision sensors, resulting in serious limitations in attitude estimation. On the other hand, the current multi-source fusion method based on Bayesian filtering cannot handle the missing measurement sources, and often ignores the spacecraft dynamics identification, resulting in the need for additional time registration when processing the asynchronous update of multi-measurement information, and cannot solve the problem of estimation divergence caused by missing measurement information. Therefore, it has great theoretical and engineering value to design a multi-vision sensor fusion scheme for the relative attitude estimation of a failed spacecraft under limited computing performance conditions.

[0003] In the published document with the publication number CN111678512B, a method for determining the attitude of a satellite by combining a star sensor and a gyroscope based on a factor graph is disclosed. An incremental smoothing optimization method based on a factor graph is designed for the attitude determination by fusing the star sensor and the gyroscope, realizing the high-precision attitude determination of the spacecraft itself. However, the method in this published document only considers its own attitude determination and cannot handle the relative attitude determination of a failed spacecraft lacking information interaction. In the published document with the publication number CN114114311A, a method for measuring the relative pose of a non-cooperative spacecraft based on multi-source information fusion is disclosed. A relative pose measurement method is designed for the fusion of a binocular camera and a lidar to realize the determination of the relative pose information of a non-cooperative spacecraft. However, the estimation accuracy of the method in this published document is seriously affected by observation noise and is extremely dependent on the high-frequency input of visual measurement information. Facing the complex and changeable space lighting environment and the severely limited on-board computing power resources, the actual engineering implementation accuracy is not high.

[0004] Most of the existing research methods for the problem of multi-sensor fusion estimation of spacecraft are based on the cooperative assumption of complete information interaction, and cannot handle the fusion estimation of multiple vision sensors under low computing power loads, affecting the execution of the proximity operation task of a failed spacecraft. Summary of the Invention

[0005] In order to overcome the above-mentioned defects of the prior art, the present invention provides a method for estimating the relative attitude of a failed spacecraft under computational resource constraints to solve the problems existing in the above-mentioned background art.

[0006] The present invention provides the following technical solution: A method for estimating the relative attitude of a failed spacecraft under computational resource constraints, comprising the following steps:

[0007] Step S01: Perform preprocessing of measurement data, obtain camera relative attitude observation data based on an ellipse recognition algorithm, and obtain lidar relative attitude measurement data based on the iterative closest point algorithm;

[0008] Step S02: Design pseudo-measurement factors for the camera and lidar, and design a variational integration prediction factor for the attitude increment of the failed spacecraft at adjacent moments;

[0009] Step S03: Perform factor graph optimization based on a moving window to estimate the relative attitude of the failed spacecraft.

[0010] Preferably, the ellipse recognition algorithm includes arc segment extraction, ellipse detection, and clustering fusion; first perform Gaussian filtering on the original image once, then use the Canny operator to extract the image edges, divide the boundary points into positive and negative categories according to the phase of the boundary point gradient, detect the 8-neighborhood connectivity of the boundary points with the same gradient positive and negative and connect them into arc segments for extraction;

[0011] Design two indicators to screen the extracted arc segments: Indicator 1: The total length of the arc segment ; Indicator 2: The length of the short side of the arc segment boundary box , combined with the concavity and convexity of the alternative arc segments, the quadrant (Ⅰ, Ⅱ, Ⅲ, Ⅳ) to which it belongs can be divided;

[0012] Based on the extracted arc segments, first combine the arc segments in a counterclockwise order of quadrants in pairs, expressed as: ((Ⅰ, Ⅱ), (Ⅱ, Ⅲ), (Ⅲ, Ⅳ), (Ⅳ, Ⅰ)), and eliminate the combinations that do not meet the position constraints; list the possible combinations (Ⅰ, Ⅱ, Ⅲ), (Ⅱ, Ⅲ, Ⅰ), (Ⅲ, Ⅳ, Ⅱ), (Ⅳ, Ⅰ, Ⅲ), and eliminate the combinations that do not meet the concentricity constraints.

[0013] Preferably, the acquisition method of the camera relative attitude measurement data is:

[0014] Perform clustering fusion through two-stage screening;

[0015] The first stage: fuse the ellipses with overly close center positions;

[0016] The second stage: check the contour integrity and exclude the ellipses below the set threshold;

[0017] The contour integrity is the ratio of the arc segments detected in the image to the estimated ellipse edge. The distance between the centers of the ellipses to be fused is used as the clustering index, and the ellipse with the highest contour integrity is retained. By finally identifying the five-element parameters of the ellipse, the rotation matrix for affine transformation into a perfect circle is calculated, which is the relative attitude measurement data obtained from the camera image.

[0018] Preferably, the method for obtaining the relative attitude measurement data of the lidar is as follows:

[0019] The attitude measurement data is obtained based on the iterative closest point algorithm; based on two point sets with a rotational correspondence, the following optimization loss function is constructed: , where R represents the relative attitude, E(R) represents the loss function of the relative attitude, represents the relative attitude input parameter for which the loss function obtains the minimum value; U is the known point set of the three-dimensional model of the failed spacecraft, V is the lidar measurement point set, F represents the Frobenius norm of the matrix, "s.t." indicates that the variables in the objective function on the left satisfy the constraints on the right, and I 3 is the third-order identity matrix; , where n represents the number of feature points, m = 1, 2, 3,..., n; u m is the m-th three-dimensional feature point in the known three-dimensional model of the failed spacecraft, and v m is the m-th three-dimensional feature point measured by the lidar;

[0020] The rotation matrix that minimizes the loss function after iterative optimization is the relative attitude measurement data obtained by the lidar.

[0021] Preferably, in step S02, the pseudo-measurements of the camera and the lidar are expressed as:

[0022] ; ;

[0023] where , respectively represent the observations of the camera and the lidar at time k; , respectively represent the noises of the two vision sensors; , respectively represent the observation functions of the two sensors; x k represents the motion state of the failed spacecraft at time k;

[0024] The corresponding pseudo-measurement odometry factor is expressed as:

[0025] ; ;

[0026] Among them, represents the cost function; f cam represents the pseudo-measurement odometry factor corresponding to the camera, and f LiDAR represents the pseudo-measurement odometry factor corresponding to the lidar.

[0027] Preferably, if the measurement data satisfies the Gaussian distribution, the cost function of the odometry factor is expressed as:

[0028] ;

[0029] Among them, is the error function, expressed as the error between the sensor measurement value and the state prediction value, that is ; is the Mahalanobis distance of the error function, is the measurement noise variance matrix of the corresponding sensor, represents the mileage.

[0030] Preferably, in the step S02, the variational integration prediction factor is obtained as follows:

[0031] Apply variational integration to the rigid body rotation Euler equation, and is the rotation increment at the current moment and

[0032] ;

[0033] Among them, represents the vector part of the attitude quaternion at time k; the rotation increment is defined as ; represents quaternion multiplication, is 's conjugate quaternion; represents the moment of inertia of the failed spacecraft; represents the cross product matrix of three-dimensional vectors;

[0034] Let , then:

[0035] ;

[0036] Design the variational integration prediction factor as:

[0037] .

[0038] Preferably, in the step S03, the factor acquisition process in the factor graph is as follows:

[0039] Solve the maximum a posteriori probability distribution of the state quantity according to the observation quantity:

[0040] ;

[0041] wherein, represents the posterior probability of the state quantity determined according to the observation quantity; represents the joint probability of the observation quantity and the state quantity; represents the prior probability of the state quantity before there is no measurement quantity information; represents the probability of the occurrence of the measurement quantity, which is also called the marginal probability; represents the probability of the measurement quantity when the state quantity has occurred, which is also called the likelihood function;

[0042] ; ;

[0043] wherein, represents the state quantity when it has occurred, the measurement quantity the likelihood probability of occurrence; represents the state from the previous moment changing to the state at the next moment the transition probability; represents the prior probability of the initial state ;

[0044] ; wherein represents the state quantity when it has occurred, the measurement quantity the likelihood probability of occurrence, represents the measurement quantity other than the IMU in the whole system; z j represents the measurement information obtained through the camera or lidar and processed as the odometer; therefore, the posterior probability of the state inferred from the observation quantity is:

[0045] ;

[0046] (1);

[0047] wherein, represents the IMU measurement information at the moment; each term in the above formula corresponds to a factor in the factor graph.

[0048] Preferably, in the step S03, the specific method for optimizing the factor graph is:

[0049] After the factor graph construction is completed, the optimization of the entire graph is carried out; the states and constraints at all times are transformed into a least squares problem, and the final optimized estimation result is obtained through iterative solution by numerical optimization methods; to avoid the huge computational burden caused by the synchronous optimization of all states, the sliding window method is adopted to limit the optimization to time to time, while the states from the initial time to time are regarded as prior factors and added to the graph;

[0050] Taking the negative logarithm of the posterior probability in Equation (1), the optimal state solution is:

[0051] ;

[0052] wherein, represents the optimal state that minimizes the cost function.

[0053] Preferably, in step S03, the relative attitude of the failed spacecraft is estimated by obtaining the minimum value of the following formula, and the formula is expressed as:

[0054] .

[0055] The technical effects and advantages of the present invention:

[0056] (1) By providing step S01, the present invention is beneficial to effectively improving the robustness of visual attitude recognition through the established method based on arc segment detection, through three stages of arc segment extraction, ellipse detection and clustering fusion.

[0057] (2) By providing step S02, the present invention is beneficial to mining the momentum conservation characteristics of the failed spacecraft by designing the variational integral prediction factor, constructing the constraint relationship of the attitude increment at adjacent times, and solving the problem that the estimation accuracy of the existing scheme depends on high-frequency measurement input.

[0058] (3) By providing step S03, the present invention is beneficial to realizing the relative pose estimation of non-cooperative targets under low computing power load by fusing the low-frequency measurement data of cameras and lidars through the established full-state factor graph model of the failed spacecraft. BRIEF DESCRIPTION OF THE DRAWINGS

[0059] Figure 1 is the flow chart of the method for estimating the relative attitude of a failed spacecraft under the computing resource constraint of the present invention.

[0060] Figure 2 is the flow chart of the ellipse recognition algorithm of the present invention.

[0061] Figure 3 is the full-state factor graph model of the relative navigation of the failed spacecraft of the present invention.

[0062] Figure 4 It is a curve graph of attitude estimation error. Specific implementation manners

[0063] Next, in combination with the accompanying drawings in the present invention, the technical solutions in the present invention will be clearly and completely described. Additionally, the forms of each structure described in the following implementation manners are merely examples. A method for relative attitude estimation of a failed spacecraft under computational resource constraints involved in the present invention is not limited to the structures described in the following implementation manners. All other implementation manners obtained by those of ordinary skill in the art without creative efforts fall within the scope of protection of the present invention.

[0064] As Figure 1 shown, the present invention provides a method for relative attitude estimation of a failed spacecraft under computational resource constraints, including the following steps:

[0065] Step S01: Perform measurement data preprocessing, obtain camera relative attitude observation data based on an ellipse recognition algorithm, and obtain lidar relative attitude measurement data based on an iterative closest point algorithm; the measurement data is relative attitude measurement data;

[0066] Take pictures of circular features on the surface of the failed spacecraft by a vision camera carried on the servicing spacecraft, detect ellipse parameters from a single-frame image, and calculate relative attitude measurement data during the process of affine-transforming the ellipse into a perfect circle; the circular features include but are not limited to thruster nozzles and docking rings, etc., and the ellipse parameters include but are not limited to centroid pixel coordinates, semi-major axis, semi-minor axis, and ellipse inclination angle; the relative attitude measurement data includes camera relative attitude measurement data and lidar relative attitude measurement data;

[0067] Step S02: Design pseudo-measurement factors for the camera and lidar, and design a variational integration prediction factor for the attitude increment of the failed spacecraft at adjacent moments; design odometer factors for the camera image information and lidar point cloud information respectively, and then add the pseudo-measurement odometer factors to the full-state factor graph;

[0068] Step S03: Define all states from the initial moment to the kth moment, and define the camera, lidar, and variational integration predictor as observations of the system state; perform factor graph optimization based on a moving window to effectively reduce the computational resource requirements and achieve accurate estimation of the relative attitude of the failed spacecraft during proximity operations.

[0069] In this embodiment, it should be specifically noted that the acquisition method of the ellipse parameters is specifically as follows:

[0070] As Figure 2As shown, an ellipse recognition algorithm is adopted, including arc segment extraction, ellipse detection, and clustering fusion. First, perform a Gaussian filter on the original image, then use the Canny operator to extract the image edges. According to the phase of the gradient of the boundary points, they can be divided into two categories: positive and negative. Then, detect the 8-neighborhood connectivity of the boundary points with the same gradient positivity and negativity and connect them into arc segments for extraction. Design two indicators to screen the extracted arc segments: Indicator 1: The total length of the arc segment , this indicator is used to filter out false detections caused by noise; Indicator 2: The length of the short side of the bounding box of the arc segment , this indicator ensures that the arc segment has enough curvature to form an ellipse rather than a straight line segment. Combining with the concavity and convexity of the candidate arc segments, the quadrant (Ⅰ, Ⅱ, Ⅲ, Ⅳ) to which it belongs can be divided;

[0071] Based on the extracted arc segments, first combine the arc segments in a pairwise manner in the counterclockwise order of the quadrants, expressed as: ((Ⅰ, Ⅱ), (Ⅱ, Ⅲ), (Ⅲ, Ⅳ), (Ⅳ, Ⅰ)), and eliminate the combinations that do not meet the position constraints; The position constraint is expressed as the constraint situation of the arc segment position. Taking the (Ⅰ, Ⅱ) combination as an example, the right end of arc segment Ⅱ needs to be to the left of the left end of arc segment Ⅰ; Taking the (Ⅱ, Ⅲ) combination as an example, the lower end of arc segment Ⅱ needs to be above the upper end of arc segment Ⅲ. On this basis, list the possible combinations (Ⅰ, Ⅱ, Ⅲ), (Ⅱ, Ⅲ, Ⅰ), (Ⅲ, Ⅳ, Ⅱ), (Ⅳ, Ⅰ, Ⅲ), and eliminate the combinations that do not meet the concentricity constraints; The concentricity constraint utilizes the property that the midpoint locus of the parallel chords of an ellipse passes through the center of the circle. Two arc segments can determine a center of the circle. Therefore, the position error of the two centers of the circle determined by three arc segments needs to be within a certain range, and then estimate the five parameters of this ellipse.

[0072] In this embodiment, it should be specifically noted that the acquisition method of the relative pose measurement data of the camera is as follows:

[0073] Perform clustering fusion through two-stage screening;

[0074] The first stage: fuse the ellipses with too close center positions;

[0075] The second stage, check the contour integrity and exclude the ellipses below the set threshold;

[0076] The contour integrity is the proportion of the arc segments detected in the image to the estimated ellipse edge. The lower this proportion, the higher the probability of false matching. Use the distance between the centers of the ellipses to be fused as the clustering index and retain the ellipse with the highest contour integrity. By finally identifying the five-element parameters of the ellipse, calculate the rotation matrix for affine transformation into a perfect circle, which is the relative pose measurement data obtained from the image.

[0077] In this embodiment, it should be specifically noted that the acquisition method of the relative pose measurement data of the lidar is as follows:

[0078] The point cloud data of the failed spacecraft is obtained by the lidar carried on the service spacecraft, and the attitude measurement data is obtained based on the iterative closest point algorithm; based on two point sets with a rotational correspondence relationship, the following optimization loss function is constructed: , where R represents the relative attitude, and E(R) represents the loss function of the relative attitude, represents the relative attitude input parameter for obtaining the minimum value of the loss function; U is the known point set of the three-dimensional model of the failed spacecraft, V is the point set measured by the lidar, F represents the Frobenius norm of the matrix, "s.t." means that the variables in the objective function on its left satisfy the constraints on the right, and I 3 is the third-order identity matrix; , where n represents the number of feature points, and m = 1, 2, 3,..., n; u m is the m-th three-dimensional feature point in the known three-dimensional model of the failed spacecraft, and v m is the m-th three-dimensional feature point measured by the lidar;

[0079] The rotation matrix that minimizes the loss function after iterative optimization is the relative attitude measurement data obtained by the lidar.

[0080] In this embodiment, it should be specifically noted that the pseudo-measurements of the camera and the lidar in step S02 are expressed as:

[0081] ; ;

[0082] Among them, , respectively represent the observations of the camera and the lidar at time k; , respectively represent the noises of the two vision sensors; , respectively represent the observation functions of the two sensors; x k represents the motion state of the failed spacecraft at time k;

[0083] The corresponding pseudo-measurement odometry factor is expressed as:

[0084] ; ;

[0085] Among them, represents the cost function; f cam represents the pseudo-measurement odometry factor corresponding to the camera, and f LiDAR represents the pseudo-measurement odometry factor corresponding to the lidar; if the measurement data all satisfy the Gaussian distribution, the cost function of the odometry factor is expressed as:

[0086] ;

[0087] Among them, is the error function, which is expressed as the error between the sensor measurement value and the state prediction value, that is ; is the Mahalanobis distance of the error function, is the measurement noise variance matrix of the corresponding sensor, represents mileage.

[0088] In this embodiment, it should be specifically noted that in step S03, considering the relative navigation and positioning problem of the failed spacecraft, all states from the initial time to time k are defined as The camera, lidar, and variational integration predictor are all defined as observations of the system state ;

[0089] The specific method of variational integration prediction is as follows:

[0090] Considering that the failed spacecraft is freely tumbling in space, without considering external torque and disturbance torque, it is in angular momentum conservation and momentum conservation during the proximity operation of the servicing spacecraft; applying variational integration to the rigid body rotation Euler equation, the rotation increment at the current time and the rotation increment at the next time Construct an equality constraint:

[0091] ;

[0092] Among them, represents the vector part of the attitude quaternion at time k; the rotation increment is defined as , represents quaternion multiplication, is 's conjugate quaternion; represents the moment of inertia of the failed spacecraft; represents the cross product matrix of three-dimensional vectors;

[0093] Let , then:

[0094] ;

[0095] Design the variational integration prediction factor as:

[0096] ;

[0097] When using the numerical optimization method to solve the variable that minimizes the error function, the Jacobian matrix still needs to be calculated. To further reduce the computational amount, the Jacobian matrix is given in an analytical form:

[0098] 。

[0099] In this embodiment, it should be specifically noted that the state estimation problem is actually to solve the maximum a posteriori probability distribution of the state quantity according to the observed quantity:

[0100] ;

[0101] Among them, represents the posterior probability of the state quantity determined according to the observed quantity; represents the joint probability of the observed quantity and the state quantity; represents the prior probability of the state quantity before there is no measurement information; represents the probability of the occurrence of the measurement quantity, which is also called the marginal probability; represents the probability of the measurement quantity in the case where the state quantity has occurred, which is also called the likelihood function;

[0102] ; ;

[0103] Among them, represents the state quantity in the case where it has occurred, the likelihood probability of the occurrence of the measurement quantity ; represents the transition probability from the state at the previous moment to the state at the next moment; represents the prior probability of the initial state ;

[0104] ; Among them represents the state quantity in the case where it has occurred, the likelihood probability of the occurrence of the measurement quantity ; represents the measurement quantity other than the IMU in the entire system; z j represents the measurement information obtained by the camera or lidar and processed into an odometer; Therefore, the posterior probability of the state inferred from the observed quantity is:

[0105] ;

[0106] (1);

[0107] Among them, represents the IMU measurement information at the Figure 3 moment; Each term in the above formula corresponds to Figure 3As shown; after the factor graph construction is completed, further optimization of the entire graph is carried out, which is also called smoothing; the states and constraints at all times are transformed into a least squares problem, and iterative solutions are obtained through numerical optimization methods such as the Gauss-Newton method and the LM method to finally obtain the optimized estimation result; to avoid the huge computational burden caused by the synchronous optimization of all states, a sliding window method is used to limit the optimization to time to time, and the states from the initial time to time are regarded as prior factors and added to the graph;

[0108] Taking the negative logarithm of the posterior probability in Equation (1), the optimal state solution is:

[0109] ;

[0110] where represents the optimal state that minimizes the cost function; the relative attitude estimation problem of the finally failed spacecraft is described as the following minimum value solving problem: ;

[0111] Finally, a relative attitude estimation experiment is carried out on the ground dual-manipulator hardware-in-the-loop platform. The camera and lidar are placed on the servicing spacecraft, and relative measurement data of the target spacecraft are collected during the proximity operation. In addition, eight Vicon infrared cameras are installed on the periphery of the experimental scene, and together with the retroreflective target balls placed on the surface of the spacecraft in advance, accurate ground truth of the relative attitude is obtained. The initial values of the experimental scene are set as follows: target attitude quaternion [-0.0041, -0.7067, -0.0053, 0.7075], target relative position [7.03316, 0.03439, 1.90147] m, relative angular velocity [0.5, 1.5, 1] * 10-3 rad / s, relative linear velocity [-1, 1, -1] * 10-2 m / s. The moment of inertia of the target spacecraft is set as follows:

[0112] ;

[0113] The ellipse of the recognized circular nozzle is converted into a perfect circle through affine transformation, so as to obtain the attitude measurement value of the camera. The image processing time for each frame is about 250 milliseconds, so the camera measurement frequency is 4 Hz. Although the acquisition frequency of repeated scanning can reach tens of thousands of Hz, considering the computational efficiency of ICP, the actual LiDAR measurement frequency is 0.5 Hz.

[0114] By fusing the measurement data of the camera and LiDAR, the attitude estimation error is as Figure 4As shown in the figure. The present invention uses the classical multiplicative extended Kalman filter method (MEKF) in the field of spacecraft attitude estimation for comparison to highlight the accuracy of the estimation method proposed by the present invention. From Figure 4 It can be seen that the average attitude estimation error of the proposed method is 0.2264 degrees. Compared with the steady-state estimation error of MEKF which is 1.0829 degrees, the attitude estimation accuracy has been significantly improved. In addition, all filtering methods require a certain convergence time to ensure better estimation performance. The attitude estimation error of MEKF stabilizes at about 1 degree after 80 seconds. In contrast, the proposed factor graph-based fusion estimation method adopts global optimization, does not require a convergence time, and maintains a stable and low estimation error throughout the proximity operation estimation process.

[0115] Finally: The above are only the preferred embodiments of the present invention and are not used to limit the present invention. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present invention shall be included within the protection scope of the present invention.

[0116] The above is only the specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any person skilled in the art can easily think of changes or replacements within the technical scope disclosed in the present application, and all of them should be covered by the protection scope of the present application. Therefore, the protection scope of the present application shall be subject to the protection scope of the claims.

Claims

1. A method for estimating relative attitude of a failed spacecraft under computing resource constraints, characterized by: The method comprises the following steps: Step S01: performing measurement data preprocessing, obtaining camera relative attitude observation data based on an ellipse recognition algorithm, and obtaining lidar relative attitude measurement data based on an iterative closest point algorithm; Step S02: designing pseudo-measurement factors for the camera and lidar, and designing variational integral prediction factors for the attitude increments of the failed spacecraft at adjacent moments; Step S03: performing factor graph optimization based on a moving window to estimate the relative attitude of the failed spacecraft.

2. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 1, characterized in that: The ellipse recognition algorithm includes arc segment extraction, ellipse detection and cluster fusion; firstly, the original image is subjected to a Gaussian filter, and then the image edge is extracted using the Canny operator, and the boundary point gradient is divided into positive and negative categories according to the phase, and the connectivity of the eight neighborhoods of the boundary points with the same positive and negative gradients is detected and connected into arc segments for extraction; two indicators are designed to screen the extracted arc segments: Indicator 1: total length of the arc segment ; Indicator 2: The length of the short side of the arc segment bounding box , combined with the concavity and convexity of the candidate arc segments, the quadrants to which they belong can be divided (Ⅰ, Ⅱ, Ⅲ, Ⅳ); based on the extracted arc segments, first combine the arc segments in pairs in the counterclockwise order of the quadrants, expressed as: ((Ⅰ, Ⅱ), (Ⅱ, Ⅲ), (Ⅲ, Ⅳ), (Ⅳ, Ⅰ)), and eliminate the combinations that do not meet the position constraints; list the possible Combine (Ⅰ, Ⅱ, Ⅲ), (Ⅱ, Ⅲ, Ⅰ), (Ⅲ, Ⅳ, Ⅱ), (Ⅳ, Ⅰ, Ⅲ), and eliminate the combinations that do not satisfy the concentricity constraint.

3. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 2, characterized in that: The camera relative posture measurement data is obtained by: clustering fusion through two-stage screening; the first stage: fusing ellipses whose center positions are too close; In the second stage, the contour completeness is checked and ellipses below the set threshold are excluded; the contour completeness is the ratio of the arc segments detected in the image to the edge of the estimated ellipse. The distance between the centers of the ellipses to be fused is used as the clustering indicator, and the ellipse with the highest contour completeness is retained; by finally identifying the five element parameters of the ellipse, its affine transformation into the rotation matrix of a perfect circle is calculated, which is the relative posture measurement data obtained by the camera image.

4. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 3 is characterized by: The laser radar relative attitude measurement data is obtained by: obtaining attitude measurement data based on an iterative closest point algorithm; based on two point sets with a rotation correspondence relationship, constructing an optimization loss function as follows: , where R represents the relative posture, E(R) represents the loss function of the relative posture, represents the relative attitude input parameter for obtaining the minimum value of the loss function; U is the known point set of the three-dimensional model of the failed spacecraft, V is the laser radar measurement point set, F represents the Frobenius norm of the matrix, "st" indicates that the variables in the objective function on the left side satisfy the constraints on the right side, and I3 is the third-order unit matrix; , where n represents the number of feature points, m=1, 2, 3, ..., n; u m is the mth 3D feature point in the known 3D model of the failed spacecraft, v m is the mth three-dimensional feature point measured by the laser radar; the rotation matrix that minimizes the loss function after iterative optimization is the relative attitude measurement data obtained by the laser radar.

5. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 4, characterized in that: In step S02, the camera and lidar pseudo measurements are expressed as: ; ;in, , Respectively represent the observations of the camera and lidar at time k; , Represent the noise of two visual sensors respectively; , Respectively represent the observation functions of the two sensors; x k represents the motion state of the failed spacecraft at time k; the corresponding pseudo-measurement odometer factor is expressed as: ; ;in, represents the cost function; f cam represents the pseudo-measurement odometer factor corresponding to the camera, f LiDAR Indicates the pseudo-measurement odometry factor corresponding to the LiDAR.

6. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 5, characterized in that: If the measurement data satisfies the Gaussian distribution, the cost function of the odometer factor is expressed as: ;in, is the error function, which is expressed as the error between the sensor measurement value and the state prediction value, that is, ; is the Mahalanobis distance of the error function, is the measurement noise variance matrix of the corresponding sensor, Indicates mileage.

7. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 6, characterized in that: In step S02, the variational integral prediction factor is obtained by applying variational integral to the rigid body rotation Euler equation, which is the rotation increment at the current moment. and the next moment rotation increment Construct an equality constraint: ;in, Represents the vector part of the quaternion attitude at time k; rotation increment Defined as , represents quaternion multiplication, for The conjugate quaternion of ; represents the moment of inertia of the failed spacecraft; Represents the cross product matrix of three-dimensional vectors; let ,but: ; The designed variational integral prediction factor is: 。 8. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 7, characterized in that: In step S03, the process of obtaining factors in the factor graph is: solving the maximum a posteriori probability distribution of the state quantity according to the observed quantity: ;in, It indicates the posterior probability of determining the state quantity according to the observed quantity; Represents the joint probability of observation and state; It represents the prior probability of the state quantity before there is any measurement information; It represents the probability of the measured quantity occurring, also known as marginal probability; It indicates the probability of the measured quantity when the state quantity has occurred, also known as the likelihood function; ; ;in, Indicates the state Measurement of the amount of the situation that has already occurred Likelihood of occurrence; Indicates the state from the previous moment Change to the next state The transition probability of Indicates the initial state The prior probability of ;in Indicates the state Measurement of the amount of the situation that has already occurred The likelihood of occurrence, Represents the measurement quantity of the entire system except IMU; j represents the measurement information obtained by the camera or lidar, which is processed as an odometer; therefore, the posterior probability of the state inferred from the observation is: ； (1); among which, express IMU measurement information at the moment; each term in the above formula corresponds to a factor in the factor graph.

9. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 8, characterized in that: In step S03, the specific method of factor graph optimization is as follows: after completing the factor graph construction, optimize the entire graph; convert the states and constraints at all times into a least squares problem, and iterate through the numerical optimization method to solve the final optimization estimation result; in order to avoid the huge computational burden caused by the simultaneous optimization of all states, the sliding window method is used to limit the optimization to Time to time, and the initial time to The state at the moment is regarded as a priori factor added to the figure; the posterior probability of equation (1) is taken as the negative logarithm, and the optimal state is solved as: ;in, Represents the optimal state that minimizes the cost function.

10. The method for estimating relative attitude of a failed spacecraft under computing resource constraints according to claim 9, characterized in that: In step S03, the relative attitude of the failed spacecraft is estimated by solving the minimum value of the following formula, which is expressed as: 。

Images