Camera calibration is a fundamental issue in visual technology, a key step in linking 2D image information with 3D spatial information, and a necessary condition for 3D reconstruction. The checkerboard target is not only easy to make, but also provides rich and easily detectable feature points in the image. Therefore, this article uses a checkerboard target with clear corner features for calibration. The parameters that need to be calibrated include the camera intrinsic matrix K , the distortion coefficients k 1 , k 2 , and the extrinsic parameters R h C T , t h C T . Firstly, randomly place a flat checkerboard target in the public field of view, and ensure that the size and position of the target can cover most of the camera’s field of view, in order to collect sufficient calibration data. Capture calibration images at different positions from the same perspective by moving the target. Afterwards, we use the camera calibration method proposed by Zhang [21] to obtain the initial estimates of K , k 1 , k 2 and R h C T , t h C T . The LM algorithm [22] was used to perform nonlinear optimization on K , k 1 , k 2 and R h C T , t h C T through Equation 7, minimizing the reprojection error of all feature points. L = ∑ h = 1 H ∑ j = 1 J e h j − e ^ K , k 1 , k 2 , R h C T , t h C T , E j 2 (7) where E j denotes the three-dimensional coordinate of the j t h feature point in the calibration image, e ^ K , k 1 , k 2 , R h C T , t h C T , E j denotes the projection of E j in the h t h calibration image, and e h j denotes the pixel coordinate of the j t h feature point in the h t h calibration image.

The edge contour of CDSP comprises circular contours of end faces and straight contours. The method presented in this paper completes 3D reconstruction based on the constraint relationship between the contour lines and the axis of CDSP. Therefore, this method does not require detecting the circular contours of end faces of CDSP and can achieve 3D reconstruction for CDSP of any length. In this paper, we use the subpixel edge detection method proposed by Trujillo-Pino [23] to obtain the pixel coordinates u ˘ , v ˘ of contour feature points in CDSP image. This method estimates the subpixel position of the edge by considering the partial area effects around the edge. Due to the distortion of the camera lens, the distortion coefficients k 1 , k 2 obtained from the camera calibration module, and the coordinate u 0 , v 0 of the principal point in the intrinsic matrix K , are used to correct the distortion of the obtained contour feature points by Equation 8. u ˘ = u + u − u 0 k 1 x 2 + y 2 + k 2 x 2 + y 2 2 v ˘ = v + v − v 0 k 1 x 2 + y 2 + k 2 x 2 + y 2 2 (8) where x , y denotes the undistorted normalized image coordinates, and u , v denotes the undistorted pixel coordinates.

Due to the apparent contours of CDSP in the image are two straight lines, we use the RANSAC algorithm [24] to estimate one of the straight line models by randomly selecting sample points from the contour feature points. Then, the distance from other points to this model is calculated to determine whether they belong to the model’s feature points. After selecting the points that comply, we obtain the feature points set u g ′ , v g ′ g ∈ 1 , … , G of the straight line model, where u g ′ , v g ′ denotes the pixel coordinates of the feature points in this model. Among the remaining feature points, the RANSAC algorithm [24] is used again to sift the feature points of another straight line model and obtain a feature points set u m ″ , v m ″ m ∈ 1 , … , M , where u m ″ , v m ″ denotes the pixel coordinates of the feature points in the other straight line model. For CDSP images with complex backgrounds, we can also obtain the contour feature point sets of CDSP through manual filtering methods. Finally, based on the least squares method, the feature point sets in the two straight line models are fitted separately. The equations for fitting lines l 1 and l 2 are expressed in Equation 9. l 1 : a ′ u g ′ + b ′ v g ′ + c ′ = 0 l 2 : a ″ u m ″ + b ′ v m ″ + c ″ = 0 (9) where a ′ , b ′ , c ′ denote the coefficients of the equation of the fitting straight line l 1 , and a ″ , b ″ , c ″ denote the coefficients of the equation of the fitting straight line l 2 . The system of equations, which is constructed using the homogeneous coordinates of the feature points of the l 1 model as the coefficient matrix, is expressed as, u 1 ′ v 1 ′ 1 ⋮ ⋮ ⋮ u G ′ v G ′ 1 a ′ b ′ c ′ = 0 (10)

Similarly, the system of equations, which is constructed using the homogeneous coordinates of the feature points of the l 2 model as the coefficient matrix, is expressed as, u 1 ″ v 1 ″ 1 ⋮ ⋮ ⋮ u M ″ v M ″ 1 a ″ b ″ c ″ = 0 (11)

For the above two systems of equations, the eigenvector corresponding to the smallest eigenvalue of the coefficient matrix represents the least squares solution of the equations. Alternatively, we can perform SVD decomposition on the coefficient matrix and then take the vector in the right singular matrix that corresponds to the smallest singular value as the optimal solution. After solving from Equation 10 and Equation 11, we can obtain two contour lines l i in the CDSP image.

In this section, based on the geometric constraints provided by the perspective projection imaging model of CDSP, the process of solving the axis 3D pose using the contour lines of CDSP in the image as input is derived. According to the contour lines l i obtained from the contours extraction module, the back-projection planes π i can be expressed in Equation 12. π i = K T l i (12) where T denotes vector transpose. Due to n i ⊥ L , the direction vector v L of L can be obtained through Equation 1. In addition, since n m ⊥ O P C → , n m ⊥ L , and O P C → ⊥ L , the direction vector v O P C of O P C → can be obtained by Equation 3. Based on the v L obtained from Equation 1 and the n m obtained from Equation 2, v O P C obtained from Equation 3 can be expressed as, v O P C = n 1 × n 2 n 1 × n 2 × n 1 n 1 + n 2 n 2 (13)

Since α is the angle between the planes π 1 and π 2 , and also the angle between the normal vectors n 1 and n 2 , the angle α obtained by Equation 4 is, α = cos − 1 n 1 · n 2 n 1 · n 2 (14)

Besides, since the radius r of the cross section circle C is known and the angle α is obtained from Equation 14, the distance d O P C between point O and point P C can be obtained through Equation 5. can be expressed as, d O P C = r sin cos − 1 n 1 · n 2 n 1 · n 2 / 2 (15)

According to the distance d O P C obtained from Equation 15 and the v O P C obtained from Equation 13, the three-dimensional coordinate of point P C in the world coordinate system obtained from Equation 6 can be expressed as, P C = r sin cos − 1 n 1 · n 2 n 1 · n 2 / 2 · n 1 × n 2 n 1 × n 2 × n 1 n 1 + n 2 n 2 n 1 × n 2 n 1 × n 2 × n 1 n 1 + n 2 n 2 (16) where n 1 = K T l 1 , n 2 = K T l 2 .

Based on the above derivation, under the premise that the radius r is known, taking the contour lines in the CDSP image as input, the 3D pose of the CDSP axis can be determined by the direction vector v L obtained from Equation 1 and the 3D coordinate of point P C obtained from Equation 16.

In this part, the proposed method is tested on simulated data. The cell size is set as ∆ u , ∆ v = 0.00345 , 0.00345 mm , the focal length is set as f = 12 mm , the image size is set as N u , N v = 2448 , 2048 pixel, and the coordinate of principal point is set as u 0 , v 0 = N u / 2 , N v / 2 , then the camera intrinsic matrix can be expressed in Equation 22. K = f / ∆ u 0 N u / 2 0 f / ∆ v N v / 2 0 0 1 (22)

In addition, the radius and length of CDSP are set as 8 mm and 500 mm respectively. The experimental scene is shown in Figure 3. We establish the world coordinate system O W − X W Y W Z W in the scene, and place a CDSP in the plane Z W = 0 mm . We denote the CDSP with its center located at the origin O W and its axis parallel to the X W axis as p 1 , and the CDSP obtained by translating p 1 along the positive direction of Y W as p 2 , as shown in Figure 3A. Similarly, we denote the CDSP with its center located at the origin O W and its axis parallel to the Y W axis as p 3 , and the CDSP obtained by translating p 3 along the positive direction of X W as p 4 , as illustrated in Figure 3B. The optical center O C is located at the world coordinate system − 150 , − 200 , − 300 mm and points to O W .

Based on the scene of the above simulation experiment, The images of p 1 , p 2 , p 3 and p 4 generated from simulated data are shown in Figure 4.

In order to verify the correctness of the proposed method, the projected contour lines of CDSP are obtained from the simulated data through the cylindrical target perspective projection model [25]. Then, we use the proposed method to achieve 3D reconstruction of four pipes p 1 , p 2 , p 3 and p 4 . The final experimental results are shown in Figure 5. Figure 5 displays the 3D pipes reconstructed using the proposed method, as well as the true 3D pipes generated using simulated data. It can be seen that the reconstructed 3D pipes completely overlap with the true 3D pipes.

In order to consider the impact of radius measurement errors on reconstruction, we set the radius measurement errors to ±   0.02 mm , ±   0.05 mm , and ±   0.10 mm respectively, and then conducted experiments with radius measurement values containing errors. After completing the reconstruction, we calculated the RMSE using Equation 17, and the calculation results are shown in Table 1. It can be seen from the table that when the radius measurement errors are ±0.02mm, ±0.05mm, and ±0.10mm, the average RMSEs are 0.944mm, 2.361mm, and 4.721mm, respectively. As the radius measurement error increases, the accuracy of reconstruction gradually decreases.

In addition, in order to approach the actual experimental conditions more closely, we add gaussian random noise with a mean of 0 and a noise level of 0.1 to the contour feature points, and then conduct a comparison experiment in terms of accuracy and speed using the Doignon [19] method, the Cheng [20] method, and the proposed method respectively. In the accuracy comparison experiment, each method was tested 10 times and the RMSE was calculated using Equation 17. The calculation results are shown in Figure 6. Figure 6 displays the RMSE for 10 reconstructions of each pipe using three different methods, as well as the average RMSE over the 10 times. The average RMSE for reconstructing the four pipes using the proposed method is 0.013mm, while the average RMSE for the Doignon [19] method is 0.057mm, and the average RMSE for the Cheng [20] method is 0.033 mm. Compared with the other two methods, the RMSE of the proposed method is the lowest, which also indicates that the proposed method is less sensitive to noise. This is because the proposed method adopts simple straight line fitting model, which has low sensitivity to noise. The Cheng [20] method affects the matching accuracy of the coupled point pairs due to noise, which in turn affects the accuracy of the reconstruction results. The Doignon [19] method uses curve model for fitting, and its parameter estimation is easily affected by noise, which in turn affects the accuracy of pose parameters. Compared with the other two methods, the proposed method is faster and less sensitive to noise. In terms of applicability, Cheng [20] method can achieve 3D reconstruction of constant-diameter pipelines, while the proposed method and Doignon [19] method are only applicable to 3D reconstruction of constant-diameter straight pipelines. This is a potential drawback of the method proposed in this article and also an area for improvement in our future work.

Afterwards, we calculate the time taken by the three methods to complete the reconstruction of pipes p 1 , p 2 , p 3 and p 4 respectively, and the results are shown in Table 2. Table 2 displays the time taken to complete the reconstruction of each pipe with three methods separately, and take the average of the time taken to reconstruct four pipes as the time-consuming of each method. The time-consuming of the proposed method is 9.984 × 10 − 4 s , while the Doignon [19] method and Cheng [20] method are 2.451 s and 2.367 s respectively. The time-consuming for reconstructing each pipe by the proposed method is less than that of the other two methods, and there is a significant increase in reconstruction speed. This is because the straight line fitting used in the proposed method is mainly based on solving linear equations, and the computational complexity is relatively small. The Doignon [19] method uses curve fitting with equations containing multiple coefficients, and the solving process involves nonlinear operations. The Cheng [20] method requires finding the coupled point pairs for each cross section circle, making the reconstruction process cumbersome. Through the comparison of accuracy and speed experiments, we have verified the correctness of the proposed method and can achieve rapid 3D reconstruction of CDSP.

To test the noise resistance of the proposed method in this article, we added Gaussian noise with an average value of 0 and a standard deviation range of 0 to 3 pixel (step size of 0.2 pixel) to the feature point coordinates of the synthesized image, and conducted 10 independent experiments for each noise level. Calculate the RMSE for different noise levels after completing the reconstruction. The RMSE for each noise level is shown in Figure 7. As can be seen from Figure 7, the accuracy gradually decreases with the increase of noise level, and it has a linear relationship with the noise level. When the accuracy requirement reaches 0.2mm, 1.5 pixel of noise can be tolerated. When the accuracy requirement reaches 0.4mm, 3 pixel of noise can be tolerated.

The above simulation experiment has verified the correctness of the proposed method. In order to further verify the feasibility of the proposed method, this section will conduct real experiment based on real data. The camera used in the experiment is Daheng MER-503-23 GM-P camera with resolution of 2448 × 2048 pixels, equipped with HN-P-1628-6M-C2/3 lens with focal length of 16 mm . We set up a real experimental scene that is consistent with the scene of simulation experiment as shown in Figure 8A. When the CDSP is placed at two different positions in the scene, they are denoted as p ∼ 1 and p ∼ 3 respectively. By translating p ∼ 1 and p ∼ 3 separately, we obtain p ∼ 2 and p ∼ 4 . Since there is no true value in the real experiment, we use the distance of pipe translation as the true distance d in Equation 17. Measured by a vernier caliper with an accuracy of 0.02 mm , the diameter and length of CDSP are 16 mm and 300 mm respectively, the distance between the two support columns of p ∼ 1 and p ∼ 2 is 112.67 mm , the distance between the two support columns of p ∼ 3 and p ∼ 4 is 112.68 mm , and the diameter of the support column is 12.68 mm . Based on this, we can calculate that the translation distances for pipes p ∼ 1 and p ∼ 3 are 99.99 mm and 100.00 mm respectively.

Based on the scene of the above real experiment, we first calibrate the camera used in the experiment. The calibration parameters obtained using the method in the camera calibration module are shown in Table 3.

Then, we collect images of p ∼ 1 , p ∼ 2 , p ∼ 3 and p ∼ 4 respectively, and the collected images are shown in Figure 8B. After feature extraction from the collected images using the method in the contours extraction module, the contour lines of p ∼ 1 , p ∼ 2 , p ∼ 3 and p ∼ 4 in the image are obtained. Subsequently, taking the obtained contour lines as input, the pose of each pipe axis is calculated using the method in the pipeline 3D reconstruction module. The contour lines and projected axes are shown in Figure 8C. Finally, the reconstructed 3D effects of p ∼ 1 , p ∼ 2 , p ∼ 3 and p ∼ 4 based on the collected images are shown in Figure 9. Figure 9 displays the reconstructed 3D pipes and their axes from the above images.

In addition, based on the scene of the above real experiment, we collect five groups of images for p ∼ 1 and p ∼ 2 , and another five groups of images for p ∼ 3 and p ∼ 4 . Afterwards, we use the proposed method, the Doignon [19] method, and Cheng [20] method to perform 3D reconstruction on these 10 groups of images. The reconstruction error results of one group of p ∼ 1 , p ∼ 2 and one group of p ∼ 3 , p ∼ 4 are shown in Figure 10. In Figure 10, the x-axis represents the point number, and the y-axis represents the error of each point. Figure 10A displays the distance deviation of each point on the reconstructed axis of p ∼ 1 to the reconstructed axis of p ∼ 2 , and Figure 10B displays the distance deviation of each point on the reconstructed axis of p ∼ 3 to the reconstructed axis of p ∼ 4 . Meanwhile, the reconstruction time-consuming for the set of p ∼ 1 , p ∼ 2 and the set of p ∼ 3 , p ∼ 4 are shown in Table 4.

After completing the 3D reconstruction of 10 groups of images using three different methods, we calculate the RMSE of each group using Equation 17. Meanwhile, we calculate the time taken to complete the reconstruction for each pipe and used the average of the time taken to complete the reconstruction for the two pipes as the time-consuming for each group of experiments. In the experimental results shown in Figure 11, the x-axis represents the group number, and the y-axis shows the RMSE of the reconstructed pipe for each group and the average time-consuming to complete the pipe reconstruction for each group. The average RMSE for the 10 groups using the proposed method is 0.165 mm , and the average time-consuming is 1.917 × 10 − 3 s . The average RMSE of the Doignon [19] method is 0.314 mm , with the average time-consuming of 2.034 s . The average RMSE of the Cheng [20] method is 0.180 mm , with the average time-consuming of 2.071 s . The accuracy of the proposed method is similar to Cheng [20] method and slightly higher than Doignon [19] method, and the time-consuming of the proposed method is the least. In the real experiment, the three-dimensional reconstruction accuracy of the pipe axis is affected by many factors such as camera calibration accuracy, lens distortion, edge contour extraction accuracy, image noise, etc., resulting in the improvement of accuracy is not as obvious as that of the simulation experiment. The obvious advantage of the proposed method compared to the three methods is the improvement in speed. The comparison results of accuracy and speed are consistent with the simulation experiment, proving that proposed method can effectively achieve rapid 3D reconstruction of CDSP.

To discuss the robustness of the proposed method under different lighting conditions and backgrounds, we collected CDSP images under different lighting conditions and backgrounds for experiments. For more complex scenes, we use manual filtering to extract contour features, and then use the proposed method for reconstruction. The collected images and their reconstructed effects are shown in Figure 12. Figure 12A shows the collected CDSP images, and Figure 12B shows the 3D reconstruction effect of CDSP. From the experimental results, it can be seen that the method proposed in this paper can effectively achieve 3D reconstruction of CDSP under different lighting conditions and backgrounds, provided that the contour features can be accurately extracted.