An In-Orbit Stereo Navigation Camera Self-Calibration Method for Planetary Rovers with Multiple Constraints

: In order to complete the high-precision calibration of the planetary rover navigation camera using limited initial data in-orbit, we proposed a joint adjustment model with additional multiple constraints. Speciﬁcally, a base model was ﬁrst established based on the bundle adjustment model, second-order radial and tangential distortion parameters. Then, combining the constraints of collinearity, coplanarity, known distance and relative pose invariance, a joint adjustment model was constructed to realize the in orbit self-calibration of the navigation camera. Given the problem of directionality in line extraction of the solar panel due to large differences in the gradient amplitude, an adaptive brightness-weighted line extraction method was proposed. Lastly, the Levenberg-Marquardt algorithm for nonlinear least squares was used to obtain the optimal results. To verify the proposed method, ﬁeld experiments and in-orbit experiments were carried out. The results suggested that the proposed method was more accurate than the self-calibration bundle adjustment method, CAHVOR method (a camera model used in machine vision for three-dimensional measurements), and vanishing points method. The average error for the ﬂag of China and the optical solar reﬂector was only 1 mm and 0.7 mm, respectively. In addition, the proposed method has been implemented in China’s deep space exploration missions.


Introduction
With the continuous development of aerospace technologies, many countries have carried out extensive exploration of extraterrestrial bodies, such as to the moon and Mars [1][2][3]. To successfully conduct these tasks, such as planetary rover navigation, high-precision mapping of the detection area, the determination of interior orientation (IO) elements, exterior orientation (EO) elements, and distortion parameters of the cameras are of great significance [4,5]. Because of the influences of environmental factors, such as temperature and pressure, and engine vibration during landing, the camera parameters can be inevitable changed. In 2004, China officially launched the lunar exploration project and named it the "Chang'e Project". Chang'e-3 is a lunar probe launched in the second phase of China's lunar exploration project. It consists of a lander and a patrol ("Yutu" lunar rover). The "Yutu" rover patrolled the lunar surface at a speed of 200 m per hour and a rhythm of about 7 m per "step", and carried out various scientific exploration tasks, including lunar surface morphology and geological structure, lunar surface material composition, available resources, earth plasma layer, etc. The planetary rover is equipped with a pair of stereo navigation cameras for navigation and positioning [3]. Taking Chang'e-3 as an example, the landing process can be divided into the main deceleration stage (15 km), the approaching stage (2 km), the hovering stage (100 m), the obstacle avoidance stage (100 m-30 m), the slow-falling stage (30 m-3 m), and the free-falling stage (<3 m). The This study proposed a self-calibration bundle adjustment method with multiple feature constraints for in-orbit camera calibration. We summarize our contributions as follows: (1) A combined self-calibration bundle adjustment method with multiple feature constraints is constructed, which makes full use of the limited conditions in orbit, can complete the calibration under the condition of a small amount of observation data, and effectively improves the calibration accuracy in-orbit.
(2) An adaptive intensity weighted line extraction method is proposed, which ensures the completeness of line detection in different directions and effectively improves the extraction effect.
(3) The redundant observations with different constraints are analyzed, and the relationship between the redundant observations and the standard deviation of the model adjustment result is obtained, which provides a reference for the optimization of other bundle adjustment models.
(4) The proposed method can provide a reference for in-orbit camera calibration of planetary rovers in the deep space exploration missions of all countries.
The remaining sections of this paper are as follows: Section 2 describes the overall process and each part of the proposed method. Section 3 presents a series of experiments on different experimental conditions, and the results are compared. Section 4 provides the conclusions.

Stereo Vision System and Coordinate Systems of Planetary Rovers
This section describes the structure of the planetary rover's stereo vision system and the parameters of each component, and then different coordinate systems and their transformation relationship are explained.

Navcam System of Planetary Rovers
Currently, the Navcam system is one of the key components for navigation and positioning of planetary rovers. Taking China's Chang'e-3 rover as an example, its Navcam system is composed of a pair of navigation cameras, a pair of panoramic cameras, a pair of obstacle avoidance cameras, a mast system and an inertial measurement unit (IMU), as shown in Figure 1. The navigation cameras use fixed-focus lenses with a medium focal length and are installed on the camera gimbal, which is fixed to the top of the planetary rover through the mast system. The height of the navigation cameras is approximately 1500 mm above the ground [32]. The cameras can rotate 360 • in the horizontal plane, and −70 • to +80 • around the pitch axis. The mast system consists of three joints to enable roll, yaw, and pitch movement of the navigation cameras and the panoramic cameras [33]. In the remote-control mode, the planetary rover maps the surrounding landscape within 120 • after moving a fixed distance, and the data are sent back to the ground to generate a digital orthophoto map (DOM), digital elevation model (DEM), and other cartographic products for path planning for the next navigation point [34]. The navigation cameras on the Chang'e-3 rover are parallel binocular stereo and can be considered pinhole imaging cameras, where the two main optical axes The navigation cameras on the Chang'e-3 rover are parallel binocular stereo cameras and can be considered pinhole imaging cameras, where the two main optical axes are parallel to each other, the horizontal axes are collinear, and cameras are placed with a certain baseline length, as shown in Figure 2. B lr is the photogrammetric baseline vector, and the length is the distance between the optical centers of the left and right cameras.

Posture Model of the Stereo Navigation Camera
The navigation cameras on the Chang'e-3 rover are parallel binocular stereo cameras and can be considered pinhole imaging cameras, where the two main optical axes are parallel to each other, the horizontal axes are collinear, and cameras are placed with a certain baseline length, as shown in Figure 2. Blr is the photogrammetric baseline vector, and the length is the distance between the optical centers of the left and right cameras. Assume that the rotation matrices between the image coordinate systems of the left and right cameras and the ground coordinate system of the control points are R l and R r . The center of the left and right images are S l and S r . Taking the image coordinate system of the right camera as a reference, then, the relative rotation of the left camera is R lr . Based on the baseline vector between the left and right navigation cameras, the transformation relationship between the rotation matrices of the stereo navigation cameras is

Bundle Adjustment Model for Binocular Cameras
The bundle adjustment model is a classic high-precision adjustment model used in the field of photogrammetry. In this study, the bundle adjustment model was used as the basic model for combined adjustment [35]. According to the definitions of the image plane coordinate system and the photogrammetric coordinate system, the classic collinear equation can be established using the stereo navigation images: Assume that the rotation matrices between the image coordinate systems of the left and right cameras and the ground coordinate system of the control points are R l and R r . The center of the left and right images are S l and S r . Taking the image coordinate system of the right camera as a reference, then, the relative rotation of the left camera is R lr . Based on the baseline vector between the left and right navigation cameras, the transformation relationship between the rotation matrices of the stereo navigation cameras is According to Equation (1), the transformation relationship between the EO elements can be obtained:

Bundle Adjustment Model for Binocular Cameras
The bundle adjustment model is a classic high-precision adjustment model used in the field of photogrammetry. In this study, the bundle adjustment model was used as the basic model for combined adjustment [35]. According to the definitions of the image plane coordinate system and the photogrammetric coordinate system, the classic collinear equation can be established using the stereo navigation images: where (x, y) are the coordinates of an image point in the image plane, (X, Y, Z) are the corresponding coordinates under the photogrammetric coordinate system, (X S , Y S , Z S ) are the translation components of the EO elements, and (a1, b1, c1), . . . , (a3, b3, c3) are the elements of the rotation matrix R. (x 0 , y 0 ) are the IO elements, and (∆x, ∆y) are the coordinate differences of image points in the left and right cameras: where k 1 , k 2 are the first-order and second-order radial distortion parameters, respectively, and p 1 , p 2 are the first-order and second-order tangential distortion parameters, respectively.
Through a Taylor series expansion of Equation (4), the linearized bundle adjustment error model for the binocular camera is obtained: where V = [v x , v y ] T is the correction of the image point residuals, and t 1 = [∆X l S , ∆Y l S , ∆Z l S , ∆ϕ l , ∆ω l , ∆κ l ] T and t 2 = [∆X r S , ∆Y r S , ∆Z r S , ∆ϕ r , ∆ω r , ∆κ r ] T are the corrections of the EO elements in the left and right images, respectively. X 1 = ∆x l 0 , ∆y l 0 , ∆ f l T and X 2 = ∆x r 0 , ∆y r 0 , ∆ f r T are the corrections of the IO elements in the left and right images, respectively. X 3 = [∆X, ∆Y, ∆Z] T is the correction of the photogrammetric coordinate system, D and F are the coefficient matrices corresponding to X 2 and X 3 , respectively, L = [l x , l y ] T is the coordinate vector of the image point, and (A, B, C) is the corresponding coefficient matrix.

Adaptive Line Extraction and Constraint Model
Lens distortion is one of the key factors in the calibration of navigation cameras. The straight line feature is a typical pattern seen in the image of the solar panel of the lander and the planetary rover, and it is not subject to the influence of vibrations. Thus, using line features as calibration constraints can ensure the robustness of the calibration algorithm and improve the accuracy of calibration parameters. In the current line extraction methods, the overall image is used to determine the threshold. However, the images containing solar panels can be affected by the lunar illumination such that the line features of solar panels in different directions are disrupted, as shown in Figure 3. This disruption can prevent the extraction of lines with poor gradients. Moreover, if manual extraction is used, the extraction result might deviate from the true position of the line. Therefore, an adaptive brightness-weighted line extraction method was proposed.

Adaptive Line Extraction and Constraint Model
Lens distortion is one of the key factors in the cal straight line feature is a typical pattern seen in the im and the planetary rover, and it is not subject to the i line features as calibration constraints can ensure the rithm and improve the accuracy of calibration param methods, the overall image is used to determine the th taining solar panels can be affected by the lunar illumi solar panels in different directions are disrupted, as s can prevent the extraction of lines with poor gradients used, the extraction result might deviate from the tru adaptive brightness-weighted line extraction method During extraction, the extracted lines are optimi extracted edges. Specifically, based on the gradient te dient magnitude and direction of the pixels are cal During extraction, the extracted lines are optimized and straightened based on the extracted edges. Specifically, based on the gradient template shown in Figure 4, the gradient magnitude and direction of the pixels are calculated, and the average gradient direction of the 3 × 3 neighborhood pixels is used to replace the original gradient of the pixel. By traversing the entire image, the gradient magnitude and direction of all pixels are obtained. Figure 5 shows the histograms of gradient directions of four navigation images. direction of the 3 × 3 neighborhood pixels is used to replace the original gradient of the pixel. By traversing the entire image, the gradient magnitude and direction of all pixels are obtained. Figure 5 shows the histograms of gradient directions of four navigation images.   Figure 5 and the statistical results of many experiments show that the distribution of the gradient directions of the rover navigation images in the Chinese lunar exploration project satisfies a Gaussian distribution. Therefore, the images are divided into different gradient intervals according to the confidence interval (CI) of a Gaussian distribution, and the pixel gradient threshold of each interval is determined [36][37][38]. According to the properties of a Gaussian distribution, 99.7% of the area is within the range of three standard deviations (SDs) 3δ . The subintervals of the gradient direction distribution are established according to the SD, as shown in Figure 6. In total, there are eight subintervals, i.e.
The red vertical line is the expected value of the gradient direction of the entire image.  Remote Sens. 2022, 13, x FOR PEER REVIEW 7 of 26 direction of the 3 × 3 neighborhood pixels is used to replace the original gradient of the pixel. By traversing the entire image, the gradient magnitude and direction of all pixels are obtained. Figure 5 shows the histograms of gradient directions of four navigation images.   Figure 5 and the statistical results of many experiments show that the distribution of the gradient directions of the rover navigation images in the Chinese lunar exploration project satisfies a Gaussian distribution. Therefore, the images are divided into different gradient intervals according to the confidence interval (CI) of a Gaussian distribution, and the pixel gradient threshold of each interval is determined [36][37][38]. According to the properties of a Gaussian distribution, 99.7% of the area is within the range of three standard deviations (SDs) 3δ . The subintervals of the gradient direction distribution are established according to the SD, as shown in Figure 6. In total, there are eight subintervals, i.e., The red vertical line is the expected value of the gradient direction of the entire image.   Figure 5 and the statistical results of many experiments show that the distribution of the gradient directions of the rover navigation images in the Chinese lunar exploration project satisfies a Gaussian distribution. Therefore, the images are divided into different gradient intervals according to the confidence interval (CI) of a Gaussian distribution, and the pixel gradient threshold of each interval is determined [36][37][38]. According to the properties of a Gaussian distribution, 99.7% of the area is within the range of three standard deviations (SDs) 3δ. The subintervals of the gradient direction distribution are established according to the SD, as shown in Figure 6 Figure 5 and the statistical results of many experiments show the gradient directions of the rover navigation images in the Chi project satisfies a Gaussian distribution. Therefore, the images are gradient intervals according to the confidence interval (CI) of a Gau the pixel gradient threshold of each interval is determined [36][37][38]. erties of a Gaussian distribution, 99.7% of the area is within the r deviations (SDs) 3δ . The subintervals of the gradient direction lished according to the SD, as shown in Figure 6. In total, there are [2 ,3 ] δ δ an tical line is the expected value of the gradient direction of the entir  To minimize the influence of illumination on the image gradient, the normalized brightness is used as a basis for weight determination, and the dual thresholds for edge extraction in different regions are also determined. The calculation of the weight and upper threshold is shown in Equation (6), and the lower threshold is 0.4 times the upper threshold. Lastly, edge extraction in the dark part of the image is optimized.
where P i is the weight at (x, y), m is the number of subintervals, Th m is the upper threshold for the m-th interval, and n is the number of pixels in the image. After obtaining the edge extraction results, the edge pixels are constructed into an edge chain. Specifically, according to the distribution of the gradient directions of the edge chain elements, pixelwise tracking is used to perform the initial straight-line detection. The process is as follows: (1) A histogram of the gradient directions of the edge pixels is established (Figure 7), and the corresponding seed point is set in each edge chain along the peak gradient direction ( Figure 8a). To minimize the influence of illumination on the image gradient, the normalized brightness is used as a basis for weight determination, and the dual thresholds for edge extraction in different regions are also determined. The calculation of the weight and upper threshold is shown in Equation (6), and the lower threshold is 0.4 times the upper threshold. Lastly, edge extraction in the dark part of the image is optimized.
where i P is the weight at ( , ) x y , m is the number of subintervals, m Th is the upper threshold for the m-th interval, and n is the number of pixels in the image. After obtaining the edge extraction results, the edge pixels are constructed into an edge chain. Specifically, according to the distribution of the gradient directions of the edge chain elements, pixelwise tracking is used to perform the initial straight-line detection. The process is as follows: (1) A histogram of the gradient directions of the edge pixels is established (Figure 7), and the corresponding seed point is set in each edge chain along the peak gradient direction ( Figure 8a). (2) Starting from the seed point, the pixels on both sides are tracked. If the gradient direction difference in the connection line between the seed point and the pixel is less than the empirical threshold of 22.5° [39], the tracking continues. Otherwise, the tracking is stopped, and a line index is created.
(3) All edge chains in the image are traversed for line detection. The edge growth process is shown in Figure 8.  To minimize the influence of illumination on the image gradient, the normalized brightness is used as a basis for weight determination, and the dual thresholds for edge extraction in different regions are also determined. The calculation of the weight and upper threshold is shown in Equation (6), and the lower threshold is 0.4 times the upper threshold. Lastly, edge extraction in the dark part of the image is optimized.
where i P is the weight at ( , ) x y , m is the number of subintervals, m Th is the upper threshold for the m-th interval, and n is the number of pixels in the image. After obtaining the edge extraction results, the edge pixels are constructed into an edge chain. Specifically, according to the distribution of the gradient directions of the edge chain elements, pixelwise tracking is used to perform the initial straight-line detection. The process is as follows: (1) A histogram of the gradient directions of the edge pixels is established (Figure 7), and the corresponding seed point is set in each edge chain along the peak gradient direction ( Figure 8a). (2) Starting from the seed point, the pixels on both sides are tracked. If the gradient direction difference in the connection line between the seed point and the pixel is less than the empirical threshold of 22.5° [39], the tracking continues. Otherwise, the tracking is stopped, and a line index is created.
(3) All edge chains in the image are traversed for line detection. The edge growth process is shown in Figure 8.  (2) Starting from the seed point, the pixels on both sides are tracked. If the gradient direction difference in the connection line between the seed point and the pixel is less than the empirical threshold of 22.5 • [39], the tracking continues. Otherwise, the tracking is stopped, and a line index is created.
(3) All edge chains in the image are traversed for line detection. The edge growth process is shown in Figure 8.
To obtain continuous straight lines, the detected line chains need to be processed via connection and fitting. First, the obtained initial straight lines are linearly fitted to obtain the slope. Then, the initial straight lines are connected based on the constraints, such as the distance between straight lines and the consistency between the slope of the line connecting the two centers of gravity (CG). Lastly, length inspection is performed to obtain the final straight lines.
(1) Fitting and connection of initial straight lines. Based on the principle of least squares, the slopeâ of the lines is calculated using Equation (7).
where n is the number of elements in the line, (x i , y i ) are the pixel coordinates, and (x, y) are the coordinates of the CG of the line. Then, the direction of the line is calculated based on the slope. To increase the calculation speed, the direction difference between the lines is taken as the initial criterion for line connection. If the direction difference between two lines is less than 22.5 • , they can be considered the same line. Otherwise, they are not the same line.
Assume that the initial line chain includes the yellow elements in Figure 9, and the line chain to be connected includes the blue elements. After connecting the center of the two chains, if the distance between every pixel on the chain and the CG line is less than one pixel, then the two lines can be connected. In addition, the blank pixels between the two lines are evaluated. If the gradient direction difference between these pixels and the two lines is less than 22.5 • , then the two lines are considered the same line. Lastly, the line chain to be connected and the blank pixels are merged into the initial line chain. All initial line chains are traversed in this manner until all connections are completed.
the distance between straight lines and the consistency between the sl necting the two centers of gravity (CG). Lastly, length inspection is p the final straight lines.
(1) Fitting and connection of initial straight lines. Based on th squares, the slope â of the lines is calculated using Equation (7).
where n is the number of elements in the line, ( , i i x y ) are the pixel coo ) are the coordinates of the CG of the line. Then, the direction of the line on the slope. To increase the calculation speed, the direction differenc is taken as the initial criterion for line connection. If the direction diffe lines is less than 22.5°, they can be considered the same line. Otherwi same line. Assume that the initial line chain includes the yellow elements i line chain to be connected includes the blue elements. After connectin two chains, if the distance between every pixel on the chain and the one pixel, then the two lines can be connected. In addition, the blank two lines are evaluated. If the gradient direction difference between t two lines is less than 22.5°, then the two lines are considered the same chain to be connected and the blank pixels are merged into the initial l line chains are traversed in this manner until all connections are comp NFA l k ≤ and the length is the shortest, k = l, i.e., NF . Then, the following inequality can be obtained: (2) Deletion of short line chains. Based on the Helmholtz principle, the number of false alarms (NFA) is introduced [39]. Assuming that A is a line chain with a length of l, in an N × N image, at least k pixels are on the line chain. Then, the NFA is defined as where N 4 is the total number of all possible line chains, and p is the direction accuracy of a line. Since the threshold used above is 22.5 • , p = 0.125. If NFA(l, k) > 1, the probability that the line belongs to the background is high, and it can be considered a false line. If NFA(l, k) ≤ 1, the probability that the line belongs to the foreground is high, and it can be considered a straight line with high robustness. When NFA(l, k) ≤ 1 and the length is the shortest, k = l, i.e., NFA(l, l) = N 4 * p l ≤ 1. Then, the following inequality can be obtained: Then, the minimal line length L min is obtained through Equation (9). If the length of an initial line is greater than L min , the line is retained; otherwise the line is discarded. After all the lines are checked for length, the final line extraction result is obtained.
The feature lines on the solar panels of the lander and the rover were extracted using the method above. Then, based on the collinear nature of the points on a line, the lines can be used as an adjustment constraint to improve the calibration accuracy. Assuming that the points are on the same line, the ground coordinates of these three points should satisfy the requirement that the slope of the line connecting any two points is the same: represent the three-dimensional coordinates of the start and end points of the line, respectively. The error model of collinear constraint was constructed according to Equation (10): where X 3 is the correction of a point in the photogrammetric coordinate system, G 1 is the corresponding coefficient matrix, and L 1 is the residual matrix. When using collinear conditions to construct an adjustment model, we do not use all points on the line to achieve constraints, but only use the starting point, midpoint and end point.

Models of Distance Constraint, Coplanar Constraint and Relative Pose Constraint
In the traditional bundle adjustment model, the control points and the pixel coordinates are considered observations, which lead to many unknowns and requires more control points or more images, thus causing the rank reduction in the coefficient matrix of the normal equation [40]. Therefore, this study jointly used non-photogrammetric and photogrammetric observations or conditions to carry out the adjustment. In addition to the line feature constraints described in Section 2.3, the other non-photogrammetric observations used as feature constraints include the following three types: (1) Distance constraint. In the calibration of the Chang'e-3 planetary rover's navigation cameras, a high-precision total station measurement system can be used to accurately measure the corners of the national flag and the solar panels and other feature points, which are then used as the virtual control points for in-orbit calibration. Since the dimensions of the flag and solar panels do not change with vibration, the distances between these virtual control points should remain the same as the laboratory measurements. Hence, the distance between two points can be used as a constraint to improve the calibration accuracy. In addition, according to the adjustment theory, the error equation must have two basic conditions, i.e., a coordinate system and a length reference, while the adjustment based on the collinear equation only has a coordinate reference. Thus, the actual length between the control points should be used in combination with the collinear equation during adjustment [41]. The distance between two virtual control points is calculated in Equation (12): where (X m , Y m , Z m ) and (X n , Y n , Z n ) are the coordinates of the control points in the photogrammetric coordinate system, and D mn is the high-precision laboratory-measured distance between the two control points, and ∆D mn is the residual between the distance obtained from the observation and the true distance.
The coordinates of the virtual control points are considered virtual observations with errors, and by performing a Taylor series expansion on Equation (12), the linearized error equation can be obtained: where V 3 is the correction of the distance residual, X 3 is the correction of the coordinates of the virtual control point, F 3 is the corresponding coefficient matrix, and L 3 is a constant term obtained by substituting the virtual control point's coordinates before correction into Equation (12).
(2) Planar conditions. On the solar panels of the lander and the planetary rover are many coplanar feature points. Therefore, based on the coplanar property, the coordinates of virtual control points or connection points can be substituted into the coplanar equation to realize coordinate adjustment, thereby improving the accuracy. Assuming that the points are on the same plane, the coplanar equation of the four points is as follows: This expression can be rewritten in a function format: After linearization of Equation (15), the following equation is obtained: where X 3 is the correction of the 3D coordinates, G 2 is the corresponding coefficient matrix, and L 2 is the corresponding residual.
(3) Constraint on the relative position and orientation of the cameras. In every photographic station, the relative position and orientation between the left and right navigation cameras are the same. Assume that two image pairs are used for adjustment, and both satisfy Equations (1) and (2). Then, based on the orthogonal property of the rotation matrix, the following equations can be obtained: where the subscripts i and j are the number of the stereo image pair. After linearization of Equation (17), the error equation of the relative position and orientation constraint of the binocular stereo cameras is obtained: where (t i 1 , t i 2 , t j 1 , t j 2 ) represents the left and right EO elements of the i-th and j-th image pair, (A l i , B r i , A l j , B r j ) represents the first-order partial derivative of the EO elements of the i-th and j-th image pair, and L ij is the residual error obtained after substituting the result of the previous iteration.

Final Calibration Model and Weighting of Observations
More radial and tangential distortion parameters do not necessarily mean a higher accuracy, because too many nonlinear parameters can reduce the robustness of the solution [42]. Therefore, the final calibration parameters used in this study included the IO Remote Sens. 2022, 14, 402 12 of 25 elements of the left and right cameras (u 0 , v 0 , f ), the radial distortion parameters k 1 , k 2 , and the tangential distortion parameters p 1 , p 2 . To improve the overall accuracy of the solution, the coordinates of the virtual control points and the connection points (X, Y, Z) are considered observations with errors and adjusted. Based on Equation (3), the above parameters, the constraint Equations (11), (13), (16) and (18), and the weight matrix P of the observations, the final indirect adjustment model with multiple constraints is where (19) can be simplified as according to the principle of least squares, the parameters can be solved: where N 1 = H T PH, N 2 = MN −1 1 M T , and W = H T PL. Then, the LM algorithm is used to optimize the parameters. The accuracy of the parameters is estimated by Equation (22).
where σ 0 is the unit weight RMSE, n is the number of observations used in the adjustment, r is the rank of the coefficient matrix of the unknown, and Q xx = N −1 is a hybrid matrix. The above adjustment model requires the weight matrix P to be determined. The virtual control points on the lander and planetary rover are measured consistently, yet the distance between the corresponding image point and the camera is quite different. Therefore, equal weight cannot be used, and the weights need to be determined separately to realize optimal solutions of parameters.
The weights of virtual control points and connection points used in the adjustment cannot be directly determined, yet the corresponding image points can reflect the accuracy. Therefore, the weight determination of these points can be converted into the weight determination of image points. In the field of close-range photogrammetry, the resolution of an image point depends on the distance between the target and the camera. Under the same conditions, the smaller the distance means the higher image resolution and better accuracy of subsequent data processing, and vice versa. In the experiment on the ground, the distance can be measured with a laser rangefinder, yet the distance is difficult to obtain under in-orbit conditions. The approximate distance (depth) can be obtained through image matching and binocular vision, as well as the parameters obtained before launch. Therefore, the depth can be used as a basis for determining the weight of the observation value of an image point, which is calculated as follows.
where Z is the depth value, i.e., the vertical distance from the control point to the optical center of the camera, T is the baseline length between the photography centers, f is the focal length, and x l , x r are the column coordinates of homologous points in the left and right images, respectively. According to the principle that high-resolution image points have large weights, the normalized weight of image point i in different images is where P i is the final weight, Z min is the minimum depth, and Z i is the depth value of the current image point.

Materials and Methods
To test the accuracy of the proposed feature-based self-calibration model, various types of experiments were carried out in a general laboratory, the simulated experimental field built by the China Academy of Space Technology (CAST), and under in-orbit conditions using the Chang'e-3 rover. In the general laboratory test, the calibration results of highprecision calibration plates were used as the initial values. In the lunar test site experiment, high-precision checkpoints were used to verify the calibration accuracy. For the in-orbit test, the dimensions of the national flag, solar panels, and other components were used.

Experiments in General Laboratory
The Zhang calibration method was used to perform an accurate calibration. A highprecision (0.1 mm) standard calibration plate with 15 × 11 grids and an edge length of 50 mm and a simulated planetary rover ( Figure 10) equipped with stereo navigation cameras were used. In this study, 13 image pairs were obtained at a fixed focal length. Then, the stereo camera was calibrated using MATLAB calibration toolbox [43]. The calibration process is mainly composed of three steps, i.e., corner detection, calibration of IO and EO elements and distortion parameters, and optimal solution calculation via reprojection.
According to the principle that high-resolution image points have larg normalized weight of image point i in different images is min min min where i P is the final weight, min Z is the minimum depth, and i Z is the d the current image point.

Materials and Methods
To test the accuracy of the proposed feature-based self-calibration m types of experiments were carried out in a general laboratory, the simulated field built by the China Academy of Space Technology (CAST), and under tions using the Chang'e-3 rover. In the general laboratory test, the calibra high-precision calibration plates were used as the initial values. In the lun periment, high-precision checkpoints were used to verify the calibration acc in-orbit test, the dimensions of the national flag, solar panels, and other com used.

Experiments in General Laboratory
The Zhang calibration method was used to perform an accurate calibr precision (0.1 mm) standard calibration plate with 15 × 11 grids and an edg mm and a simulated planetary rover ( Figure 10) equipped with stereo navig were used. In this study, 13 image pairs were obtained at a fixed focal len stereo camera was calibrated using MATLAB calibration toolbox [43]. The c cess is mainly composed of three steps, i.e., corner detection, calibration elements and distortion parameters, and optimal solution calculation via re We analyzed the calibration results by re-projection error. The statisti be concluded that the average error of reprojection obtained by using 13 p is 0.3 pixels. The reprojection error of Camera 1 was slightly larger than tha indicating that the quality of Camera 1 was lower than that of Camera 2. reprojection errors of the 13th image pair were above 0.4 pixels. To obtain bration results, the 13th image pair was discarded. After that, the reprojecti left and right cameras was reduced to 0.2 pixels. The IO elements and disto ters obtained during calibration are shown in Table 1. We analyzed the calibration results by re-projection error. The statistical results can be concluded that the average error of reprojection obtained by using 13 pairs of images is 0.3 pixels. The reprojection error of Camera 1 was slightly larger than that of Camera 2, indicating that the quality of Camera 1 was lower than that of Camera 2. Moreover, the reprojection errors of the 13th image pair were above 0.4 pixels. To obtain accurate calibration results, the 13th image pair was discarded. After that, the reprojection error of the left and right cameras was reduced to 0.2 pixels. The IO elements and distortion parameters obtained during calibration are shown in Table 1. After the calibration, the camera parameters could be used as the initial values of the proposed combined adjustment method. The method mainly includes four steps: (1) extraction of the coordinates of the virtual control points and the matching points, (2) determination of the initial EO elements, (3) self-calibration bundle adjustment with multiple constraints, and (4) accuracy assessment.
Specifically, the image point coordinate extraction software was used to extract the coordinates of the virtual control points with sub-pixel precision. Then the proposed method in Section 2.3 was used to extract the straight line of the solar panel area. There were 10 pairs of virtual control points in the lander images and 13 pairs in the solar panel images. The virtual control points are shown in Figure 11. Since it is impossible to simulate the actual lighting environment of solar panels, the comparison experiment of line extraction is implemented in Section 3.3. After the calibration, the camera parameters could be used as the initial values of the proposed combined adjustment method. The method mainly includes four steps: (1) extraction of the coordinates of the virtual control points and the matching points, (2) determination of the initial EO elements, (3) self-calibration bundle adjustment with multiple constraints, and (4) accuracy assessment.
Specifically, the image point coordinate extraction software was used to extract the coordinates of the virtual control points with sub-pixel precision. Then the proposed method in Section 2.3 was used to extract the straight line of the solar panel area. There were 10 pairs of virtual control points in the lander images and 13 pairs in the solar panel images. The virtual control points are shown in Figure 11. Since it is impossible to simulate the actual lighting environment of solar panels, the comparison experiment of line extraction is implemented in Section 3.3.  Before adjustment, the initial value of EO must be obtained. In this study, the DLT algorithm was used to calculate the initial values of the EO elements. The error equations of the left and right images were constructed separately. Based on the above initial data, the proposed self-calibration bundle adjustment model with multiple constraints was used to carry out parameters calibration, and the LM algorithm was used for adjustment. Since the first image pair did not contain line features, only distance and coplanar constraint conditions were used. The model converged after 10 iterations, and the unit weight RMSE was 2.0 mm. The average reprojection error of the images was 2.7 pixels by the left and right camera. The second image pair contained line features, so collinear constraint was included in the calibration. The adjustment was then carried out based on the two pairs of images, and RMSE was reduced to 1.0 mm. The average reprojection error of the images was 1.1 pixels by the left camera and 0.9 pixels by the right camera. The average reprojection error indicated that the accuracy of the result was higher when the collinear constraint was included (using two pairs of images). Thus, the use of collinear constraints can effectively improve the accuracy of the calibration parameters.
The errors of calibration results when using one pair and two pairs of stereo images are shown in Table 2.  Table 2 shows the RMSE of each parameter after the adjustment of the two pairs of images, and it was substantially lower than that of the first image pair. Among the key parameters, the relative accuracy of the camera focal length f increased by 43.6% on average, the coordinates of the principal points u 0 and v 0 increased by 44.6% and 49.7%, respectively, and k 1 , k 2 , p 1 , and p 2 increased by 38.4%, 34.1%, 56.6%, and 50.1%, respectively. Thus, the addition of a collinear constraint significantly improved the calibration results. In addition, since there were more constraints for calibration using two pairs of images, this approach is beneficial to the calibration of the parameters.
To verify the accuracy of the calibration results, the proposed method was compared with the classic self-calibration bundle adjustment method, CAHVOR, and Vanishing Points. The results of these methods were compared with the coordinates obtained by the total station measurement system. Table 3 shows the errors at 20 checkpoints obtained by different methods.  Table 3 shows that the proposed method had a higher accuracy than the other three methods, with an average error of less than 1 mm. Compared with the other three methods, the accuracy of our method is improved by 25.0%, 40.0%, and 43.8% respectively. The error of the bundle adjustment method was slightly larger than that of the proposed method and the errors of CAHVOR and vanishing points were much larger.
The above results show that the proposed method can effectively realize camera calibration using multiple pairs of images. The constraints increased the robustness of the calibration results and the accuracy of the parameters. Thus, it can perform the calculation and inspection of the IO elements and lens distortion in the stereo camera under the in-orbit condition.

Experiments in Simulated Lunar Environment
To simulate the various working conditions of the planetary rover in a real deep space environment, the CAST has established a large indoor test field with parallel light arrays, volcanic ash-simulated lunar soil, rocks, craters, and various types of planetary vehicle prototypes. Figure 12 shows the indoor test field and the prototype rover.
CAHVOR [30] 1.2 1.6 1 Vanishing Points [17] 1.4 2.0 1 Table 3 shows that the proposed method had a higher accuracy than the other methods, with an average error of less than 1 mm. Compared with the other three ods, the accuracy of our method is improved by 25.0%, 40.0%, and 43.8% respectivel error of the bundle adjustment method was slightly larger than that of the pro method and the errors of CAHVOR and vanishing points were much larger.
The above results show that the proposed method can effectively realize came ibration using multiple pairs of images. The constraints increased the robustness calibration results and the accuracy of the parameters. Thus, it can perform the calcu and inspection of the IO elements and lens distortion in the stereo camera under t orbit condition.

Experiments in Simulated Lunar Environment
To simulate the various working conditions of the planetary rover in a real space environment, the CAST has established a large indoor test field with paralle arrays, volcanic ash-simulated lunar soil, rocks, craters, and various types of plan vehicle prototypes. Figure 12 shows the indoor test field and the prototype rover. The prototype rover took eight pairs of images at different photographic stati the indoor test field. Sixty-four checkpoints were arranged in the field of view of the navigation cameras to verify the calibration accuracy. The coordinates of the check were measured with a high-precision total station measurement system (Leica TS5 ure 13a), and the angle measurement accuracy was 0.5 seconds. Figure 13b show moving trajectory of the planetary rover and the photographic stations. Figure 14 s the 64 checkpoints for stereo cameras. The prototype rover took eight pairs of images at different photographic stations of the indoor test field. Sixty-four checkpoints were arranged in the field of view of the stereo navigation cameras to verify the calibration accuracy. The coordinates of the checkpoints were measured with a high-precision total station measurement system (Leica TS50, Figure 13a), and the angle measurement accuracy was 0.5 s. Figure 13b shows the moving trajectory of the planetary rover and the photographic stations. Figure 14 shows the 64 checkpoints for stereo cameras.    The stereo navigation cameras of the Chang'e-3 prototype were used to colle pairs of stereo images for calibration. Then, the proposed model was used (withou collinear constraint since the prototype was not equipped with solar panels) to verif IO elements and the radial and tangential distortion parameters. The calibrated par ters of the stereo navigation cameras are shown in Table 4. The stereo navigation cameras of the Chang'e-3 prototype were used to collect 27 pairs of stereo images for calibration. Then, the proposed model was used (without the collinear constraint since the prototype was not equipped with solar panels) to verify the IO elements and the radial and tangential distortion parameters. The calibrated parameters of the stereo navigation cameras are shown in Table 4. According to the calibration results, it can be concluded that the differences in the calibration parameters of the proposed method and the calibration plate were between 2.8 pixels and 6.2 pixels. Specifically, the differences in the focal length f of the left and right cameras were 6.2 pixels and 4.6 pixels, respectively, the differences in the coordinates of principal points u 0 and v 0 were 2.8 pixels, 5.0 pixels, 3.1 pixels, and 4.3 pixels, respectively, and the results of k 1 , k 2 , p 1 , and p 2 were close.
Since the calibration method based on the calibration plate and the proposed method are different types of models, it is difficult to directly compare the accuracy from the differences. Therefore, checkpoints were used to verify the accuracy of the results. The space intersection appropriate for multi-images was used to calculate the coordinates of 64 checkpoints in the object coordinate system. Taking the coordinates of the checkpoints measured by the total station measurement system as the true values, the calibration results of the four methods were compared, and the results are shown in Table 5. From the statistical results of checkpoint errors, it can be concluded that the proposed method with multiple constraints was significantly more accurate than the self-calibration bundle adjustment model, which was more accurate than the other two methods. Compared with the other three methods, the average accuracy increased by 31.58%, 49.72%, and 45.50%, the maximum error decreased by 26.46%, 59.11%, and 43.41%, and the RMS decreased by 36.44%, 45.09%, and 43.86%. Overall, the proposed method was 42.27%, 42.99%, and 42.35% better than the other methods. The test results indicated that the constraints improved the robustness and accuracy of the adjustment model, and the proposed multiple constraints were effective.
Four pairs of checkpoints in Figure 15 were used for verification. The distance was calculated based on their coordinates obtained by the adjustment. The accuracy of the distance is shown in Table 6.
intersection appropriate for multi-images was used to calculate the coordinates of 6 checkpoints in the object coordinate system. Taking the coordinates of the checkpoint measured by the total station measurement system as the true values, the calibration re sults of the four methods were compared, and the results are shown in Table 5. From the statistical results of checkpoint errors, it can be concluded that the propose method with multiple constraints was significantly more accurate than the self-calibratio bundle adjustment model, which was more accurate than the other two methods. Com pared with the other three methods, the average accuracy increased by 31.58%, 49.72% and 45.50%, the maximum error decreased by 26.46%, 59.11%, and 43.41%, and the RM decreased by 36.44%, 45.09%, and 43.86%. Overall, the proposed method was 42.27% 42.99%, and 42.35% better than the other methods. The test results indicated that the con straints improved the robustness and accuracy of the adjustment model, and the propose multiple constraints were effective.
Four pairs of checkpoints in Figure 15 were used for verification. The distance wa calculated based on their coordinates obtained by the adjustment. The accuracy of th distance is shown in Table 6.    According to the statistical results of the distance error, it can be concluded that the distance deviation between checkpoints was between 0.9 mm and 1.8 mm, with an average deviation of 1.4 mm. Thus, the camera parameters and checkpoint coordinates calculated by the proposed method have high accuracy, demonstrating the reliability of the proposed method.

In-Orbit Calibration and Analysis
After the planetary rover separated from the lander, it adjusted its moving direction and took sequential images at five photographic stations, i.e., photographic stations X, A, B, C and D (Figure 16). Ten pairs of stereo images from the top of the lander, photographic stations C and D are used to verify the effectiveness of the proposed method under inorbit conditions.

In-Orbit Calibration and Analysis
After the planetary rover separated from the lander, it adjusted its moving direction and took sequential images at five photographic stations, i.e., photographic stations X, A, B, C and D ( Figure 16). Ten pairs of stereo images from the top of the lander, photographic stations C and D are used to verify the effectiveness of the proposed method under inorbit conditions. The in-orbit calibration process of the navigation cameras can be divided into four steps: (1) extraction of the image point coordinates of the virtual control points and the matching connection points, (2) calculation of the initial EO elements, (3) self-calibration bundle adjustment with multiple constraints, and (4) accuracy assessment. Since no checkpoints with known coordinates were on the lunar surface, the accuracy of the camera parameters was evaluated based on the dimensions of the devices on the rover and the positioning results of the planetary rover.
(1) Manual extraction of the coordinates of the virtual control points. The sub-pixel coordinate extraction of the connecting points is carried out with the software specially designed for the Chang'e-3 mission. A total of 9 pairs of connection points were extracted from the images taken at station D and were taken as the control points. Another 8 pairs of feature points were extracted as distance constraints from the images of the solar panel taken at station C. Figure 17 shows the distribution of virtual control points on the solar panel of the lander and planetary rover whose coordinates are measured on the ground. (2) Automatic matching of the connection points. The parameters to be solved included the position and orientation of the camera and the virtual control points were not The in-orbit calibration process of the navigation cameras can be divided into four steps: (1) extraction of the image point coordinates of the virtual control points and the matching connection points, (2) calculation of the initial EO elements, (3) self-calibration bundle adjustment with multiple constraints, and (4) accuracy assessment. Since no checkpoints with known coordinates were on the lunar surface, the accuracy of the camera parameters was evaluated based on the dimensions of the devices on the rover and the positioning results of the planetary rover.
(1) Manual extraction of the coordinates of the virtual control points. The sub-pixel coordinate extraction of the connecting points is carried out with the software specially designed for the Chang'e-3 mission. A total of 9 pairs of connection points were extracted from the images taken at station D and were taken as the control points. Another 8 pairs of feature points were extracted as distance constraints from the images of the solar panel taken at station C. Figure 17 shows the distribution of virtual control points on the solar panel of the lander and planetary rover whose coordinates are measured on the ground. proposed method.

In-Orbit Calibration and Analysis
After the planetary rover separated from the lander, it adjusted its moving direction and took sequential images at five photographic stations, i.e., photographic stations X, A, B, C and D ( Figure 16). Ten pairs of stereo images from the top of the lander, photographic stations C and D are used to verify the effectiveness of the proposed method under inorbit conditions. The in-orbit calibration process of the navigation cameras can be divided into four steps: (1) extraction of the image point coordinates of the virtual control points and the matching connection points, (2) calculation of the initial EO elements, (3) self-calibration bundle adjustment with multiple constraints, and (4) accuracy assessment. Since no checkpoints with known coordinates were on the lunar surface, the accuracy of the camera parameters was evaluated based on the dimensions of the devices on the rover and the positioning results of the planetary rover.
(1) Manual extraction of the coordinates of the virtual control points. The sub-pixel coordinate extraction of the connecting points is carried out with the software specially designed for the Chang'e-3 mission. A total of 9 pairs of connection points were extracted from the images taken at station D and were taken as the control points. Another 8 pairs of feature points were extracted as distance constraints from the images of the solar panel taken at station C. Figure 17 shows the distribution of virtual control points on the solar panel of the lander and planetary rover whose coordinates are measured on the ground. (2) Automatic matching of the connection points. The parameters to be solved included the position and orientation of the camera and the virtual control points were not (2) Automatic matching of the connection points. The parameters to be solved included the position and orientation of the camera and the virtual control points were not uniformly distributed in the image; therefore, to improve the accuracy of the model solution, the stereo images were matched to build a set of connection points. Figure 18 shows the matching results of two pairs of stereo images. uniformly distributed in the image; therefore, to improve the accuracy of the model solution, the stereo images were matched to build a set of connection points. Figure 18 shows the matching results of two pairs of stereo images.  Figure 19 shows the result of line extraction using the proposed method and the classic LSD (line segment detector) method.  Figure 19 shows the result of line extraction using the proposed method and the classic LSD (line segment detector) method. Figure 18. Matching results of connection points in images taken at adjacent stereo images (Th same color represents the corresponding matching points).
(3) Extraction of line features on the solar panel. The adaptive line extraction metho proposed in Section 2.3 was used to extract the line features in the images of the sola panel. Figure 19 shows the result of line extraction using the proposed method and th classic LSD (line segment detector) method. The results show that the proposed method extracted more lines. Especially, the line with smaller gradients were also extracted. The main reason was that the threshold of th proposed method was adaptively determined according to the gradient direction of th line, which effectively overcame the problem that the global unified threshold cannot b used for extraction because of the large brightness difference in the solar panel area in th stereo image due to the color and illumination.
According to the extracted coordinates of the image point, the coordinates of the con nection point, and the initial camera parameters obtained before launching, the initial EO elements were solved by the DLT algorithm. Then, the self-calibration bundle adjustmen method with multiple constraints was used for 3 image pairs of the top of the lander an 7 image pairs from each of photographic stations C and D. The IO elements and distortio parameters and the corresponding RMSEs obtained by the combined adjustment ar shown in Table 7. The results show that the proposed method extracted more lines. Especially, the lines with smaller gradients were also extracted. The main reason was that the threshold of the proposed method was adaptively determined according to the gradient direction of the line, which effectively overcame the problem that the global unified threshold cannot be used for extraction because of the large brightness difference in the solar panel area in the stereo image due to the color and illumination.
According to the extracted coordinates of the image point, the coordinates of the connection point, and the initial camera parameters obtained before launching, the initial EO elements were solved by the DLT algorithm. Then, the self-calibration bundle adjustment method with multiple constraints was used for 3 image pairs of the top of the lander and 7 image pairs from each of photographic stations C and D. The IO elements and distortion parameters and the corresponding RMSEs obtained by the combined adjustment are shown in Table 7.  Table 7 shows that the proposed method obtained the camera's IO elements and distortion parameters with high accuracy. The average RMSE of focal length f and coordinates of the principal points u 0 and v 0 was only 0.1 pixel, 0.2 pixels, and 0.2 pixels, respectively. Thus, the proposed method was feasible under in-orbit conditions and can provide technical references for in-orbit calibration of planetary rovers in deep space exploration missions. Furthermore, the influence of different constraints on the calibration was analyzed. Combining the above observations and different constraints, 7 different models were established. Since the same points are used for collinear and coplanar constraints, both are used as constraints at the same time. Finally, the median errors of various models are obtained after adjustment, and the details are shown in Table 8. From the above experimental results, it can be concluded that BA can be solved without constraints, but the standard deviation is larger than other models with additional constraints, indicating that the proposed constraints are beneficial to the adjustment.
Compared with the BA, the standard deviation of the various models with different constraints is reduced by 8%, 61%, 63%, 62%, 75%, and 76%, respectively, indicating that the accuracy of the proposed method is significantly improved. Since there are more redundant observations than other models, CCC and RC have a greater effect on the improvement of calibration accuracy. we found the sums of squared residuals of BA with different additional constraints are pretty close. The larger the redundancy number, the smaller the estimation precision σ 0 of calibration parameters. So the redundancy number of BA with different additional constraints is larger than that of BA. Thus, the calibration parameters estimation precision of BA with different additional constraints is better, and the adoption of various models proposed in this paper as geometric constraints is more rational.
In addition, the accuracy of the parameters was indirectly verified using the size of the national flag. Figure 20 shows the manually extracted corner points of the national flag, and Table 9 shows the design size and calculated size of the national flag and the differences.  According to Table 9, the calculation errors were not more than 2 mm, indicating that using the IO elements and distortion parameters calculated by the proposed method to solve the space coordinates of other points can achieve good accuracy.
Based on the line extraction results of the solar panel from images taken by stereo cameras at different photographic stations and angles, the coordinates of the intersection points between lines, and the calculated parameters, the lengths of multiple OSRs were  According to Table 9, the calculation errors were not more than 2 mm, indicating that using the IO elements and distortion parameters calculated by the proposed method to solve the space coordinates of other points can achieve good accuracy.
Based on the line extraction results of the solar panel from images taken by stereo cameras at different photographic stations and angles, the coordinates of the intersection points between lines, and the calculated parameters, the lengths of multiple OSRs were calculated and compared with the true value. Figure 21 shows the endpoints of the OSRs after line extraction. The length and width of the OSR obtained by the 3D coordinates calculated by the binocular vision measurement method were compared with the design values, as shown in Table 10. According to Table 9, the calculation errors were not more than 2 mm, indicating that using the IO elements and distortion parameters calculated by the proposed method to solve the space coordinates of other points can achieve good accuracy.
Based on the line extraction results of the solar panel from images taken by stereo cameras at different photographic stations and angles, the coordinates of the intersection points between lines, and the calculated parameters, the lengths of multiple OSRs were calculated and compared with the true value. Figure 21 shows the endpoints of the OSRs after line extraction. The length and width of the OSR obtained by the 3D coordinates calculated by the binocular vision measurement method were compared with the design values, as shown in Table 10.    Table 10 shows that the calculated OSR dimensions from images of different angles and distances were accurate. The average length of calculated by the parameters obtained by the proposed method was 31.2 mm, and the average width was 39.8 mm, with an average deviation from the designed value of 0.4 mm and 1.0 mm, respectively. Thus, the   Table 10 shows that the calculated OSR dimensions from images of different angles and distances were accurate. The average length of calculated by the parameters obtained by the proposed method was 31.2 mm, and the average width was 39.8 mm, with an average deviation from the designed value of 0.4 mm and 1.0 mm, respectively. Thus, the camera parameters obtained by the proposed method have good accuracy. In particular, the lens distortion parameters have good reliability and accuracy.

No. Calculated (mm) Designed (mm) Deviation (mm) No. Calculated (mm) Designed (mm) Deviation (mm)
In summary, the classic bundle adjustment model has good accuracy. However, in the absence of constraints, the point position may still deviate in the light direction during the adjustment process. With the addition of distance constraints and coplanar constraints, the spatial geometric structure of the observations is more stable and reliable, thereby reducing the point deviation during the solution, and ultimately improving the accuracy of camera parameter calibration. In addition, the relative pose constraints of the left and right navigation cameras can also realize the spatial geometric constraints of the EO, making the model during the adjustment more stable and the solution more accurate.

Conclusions
Accurate IO elements of the stereo camera are crucial for the in-orbit navigation and positioning of planetary rovers. High-precision in-orbit parameter calibration can provide a reference for path planning of the planetary rover and play a key role in ensuring the safety of the planetary rover and improving the efficiency of scientific missions. In deep space exploration missions, the variables that can be used for in-orbit camera calibration are very limited. Therefore, an in-orbit stereo camera self-calibration method of planetary rovers with multiple constraints was proposed. Based on the bundle method of the binocular stereo camera, four types of constraints were included to construct the self-calibration bundle adjustment model with multiple constraints, i.e., the constraint of the distance between feature points on the lander and the planetary rover, the collinear constraint of the connection points on the solar panel, the coplanar constraint of feature points on the solar panel and national flag, and the relative position and orientation constraint of the left and right cameras. Given the problem of unsatisfactory distribution of line gradient features in the solar panel images, an adaptive line extraction method based on brightness weighting was proposed to ensure the accuracy of line extraction in different directions. Through various types of experiments, the following conclusions were obtained: (1) In the simulations in a general laboratory, two pairs of stereo images were used in the adjustment. The results showed that using the extracted linear features as additional constraints can effectively improve the accuracy of the obtained camera parameters, and the reprojection error was reduced to 1.0 pixels. Overall, the accuracy of various parameters was improved by 34.1-56.6%.
(2) In the experiments at the simulated experimental field built by CAST, 64 checkpoints were used to test the calibration accuracy. The results showed that the average error of the proposed method was only 2.1 mm, which was 31.58%, 49.72%, and 45.50% better than the classic self-calibration bundle adjustment method, CAHVOR, and Vanishing Points, respectively. Four pairs of checkpoints were selected for distance verification, and the average distance error was only 1.4 mm.
(3) In the in-orbit calibration test of the Chang'e-3 planetary rover, the average RMSE of the focal length and the principal point coordinates was only 0.2 pixels for the left camera and right camera. The maximum error in flag dimensions was 2.0 mm and the minimum error was 0.2 mm. The maximum error in the length of the OSR of the solar panel was 1.5 mm, and the minimum error was 0.1 mm.
Based on the above results, the proposed method can obtain high-precision IO elements and lens distortion parameters of the stereo cameras under in-orbit conditions, which provide a reference for in-orbit camera calibration of planetary rovers in the deep space exploration missions of all countries.