A method for optimizing parameters of a binocular stereo vision camera

Abstract

The application discloses a kind of binocular stereo vision camera calibration parameter optimization method, corresponding geometric constraint is constructed to different situations, so that the result after three-dimensional reconstruction of calibration board is consistent with actual three-dimensional scene, i.e. the coordinate of reconstruction point should be consistent with actual space coordinate and the relative position between point and point is invariable, while the geometric characteristics of calibration board in space should be kept unchanged after calibration reconstruction, and all constraint terms are superimposed to form final optimization objective function, finally the objective function is optimized using Levenberg-Marquardt method, and the optimized camera calibration parameter is obtained.The camera calibration parameter optimization algorithm proposed in the application has higher accuracy compared to Zhang Zhengyou camera calibration method, and after calibration using the optimized parameter, the reprojection error and three-dimensional reconstruction measurement result error are greatly reduced, improving the camera calibration parameter precision.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of binocular stereo vision, and more particularly to a binocular stereo vision camera calibration parameter optimization method. BACKGROUND

[0002] With the rapid development of economy and technology, intelligent manufacturing technology is increasingly valued, among which, binocular stereo vision technology embodying artificial intelligence is particularly important. Three-dimensional reconstruction based on binocular vision can restore the three-dimensional information of a scene, and can effectively improve production efficiency and product quality, reduce operating costs and resource waste in intelligent production, traffic monitoring, machine navigation, aerospace, medical modeling, cultural relic restoration, and non-contact high-precision measurement.

[0003] Camera calibration plays a very important role in the field of binocular stereo vision, and is the primary step to realize three-dimensional scene reconstruction, and the accuracy of camera calibration directly affects the reconstruction result. Camera calibration is a process of solving camera internal and external parameters from known spatial position points and their image corresponding points. Among them, the camera calibration parameters to be optimized are as follows:

[0004] Left camera internal parameters:

[0005] Left camera distortion coefficient: D l =[k 1l k 2l p 1l p 2l k 3l ];

[0006] Right camera internal parameters:

[0007] Right camera distortion coefficient: D r =[k 1r k 2r p 1r p 2r k 3r ]:

[0008] Rotation parameters between left and right cameras: and translation vector: T = [t1 t2 t3] T .

[0009] Among them, f xl , f yl are the physical lengths of 1 pixel in the x and y directions on the left camera photosensitive plate, u 0l , v 0l are the horizontal and vertical pixel differences between the center pixel coordinates of the left image and the image origin pixel coordinates, k 1l , k2l , k 3l is a radial distortion coefficient of the left camera, p 1l , p 2l is a tangential distortion coefficient of the left camera; similarly, f xr , f yr are the physical lengths of one pixel in the x, y directions on the photosensitive plate of the right camera, u 0r , v 0r are the horizontal and vertical pixel numbers between the center pixel coordinates and the image origin pixel coordinates of the right image, k 1r , k 2r , k 3r is a radial distortion coefficient of the right camera, p 1r , p 2r is a tangential distortion coefficient of the right camera; r1,...,r9 are functions of the rotation angles of the x, y, z axes, and t1, t2, t3 are the translation dimensions in the x, y, z directions, respectively.

[0010] The coordinate conversion in the binocular camera calibration process is shown in Figure 1 The world coordinate system, the camera coordinate system, the image coordinate system and the pixel coordinate system are included in the calibration system, and in the calibration process, the calibration coordinates are converted to different coordinate systems through a transformation matrix composed of the camera internal and external parameters. Among them, the Zhang Zhengyou camera calibration method is the most commonly used camera calibration method, which has high robustness and accuracy, but this method only considers the 2D pixel constraint and introduces certain assumptions when estimating the parameters, so there is a certain error in the analytical parameters, and the parameter accuracy still has room for improvement. SUMMARY

[0011] The purpose of the present application is to overcome the shortcomings of the prior art, and to provide an optimization method for binocular stereo vision camera calibration parameters, which optimizes the calibrated camera calibration parameters under the condition of known calibrated camera internal and external parameters, and comprehensively considers various constraints, so as to improve the calibration accuracy of the camera calibration parameters and strengthen the accuracy and practicality of the camera calibration results.

[0012] To achieve the above-mentioned purpose of the application, the optimization method for binocular stereo vision camera calibration parameters comprises the following steps:

[0013] (1) Obtain b images of the calibration board, and obtain the initial values of the camera calibration parameters according to the camera calibration method based on the b images of the calibration board, and take the initial values as the camera calibration parameters to be optimized, wherein the number of corner points in each row is n and the number of corner points in each column is m in the corner point specification of the calibration board;

[0014] (2) Construct 2D pixel constraints

[0015] 2.1), for the i-th, i = 1, 2,..., m row and j-th, j = 1, 2,..., n column corner point of the q-th, q = 1, 2,..., b calibration board image q P ij whose world coordinate system spatial coordinates are q P ij = [ q X ij , q Y ij , q Z ij ] T , q X ij , q Y ij , q Z ij are the x, y, z axis coordinate values of the corner point q P ij , directly after corner point extraction on the q-th calibration board image, the obtained corner point q P ij has coordinate positions of on the left and right camera images respectively are the x, y axis coordinate values of the corner point , are the x, y axis coordinate values of the corner point ;

[0016] According to the pinhole camera model, the corner point q P ij is projected onto the left and right camera image planes respectively through the left and right camera calibration parameters, to obtain the left camera 2D projection point and the right camera projection point are the x, y axis coordinate values of the projection point , are the x, y axis coordinate values of the projection point ;

[0017] Taking the left camera coordinate system as the world coordinate system, the parameter conversion relationship of the corner point q P ij on the left camera is:

[0018]

[0019] The parameter conversion relationship on the right camera is:

[0020]

[0021] where wl and w r The scale factor is unknown.

[0022] 2.2) Compare the coordinates of all projected points with the coordinates of the corner points extracted directly from the image, calculate the reprojection error, and construct a 2D pixel constraint, i.e.:

[0023]

[0024] (3) Constructing 3D reconstruction constraints

[0025] 3.1) Through the external parameters of the left camera T lwc =[0 0 0] T The coordinate set Ω of the calibration plate corner points in the world coordinate system w ={ q P ij |[ q X ij , q Y ij , q Z ij ] T}, q=1,2,...,b, i=1,2,...m, j=1,2,...n, transformed to the left camera coordinate system, yields a set of three-dimensional coordinates of the ideal calibration board corner points in the left camera coordinate system. Points The x, y, and z axis coordinates are given, and the coordinate transformation formula is as follows:

[0026]

[0027] 3.2) Pixel coordinates of calibration board corner points extracted from the left and right images using camera distortion coefficients. Perform distortion correction, and denote the pixel coordinates of the ideal distortion-free left image as... And the pixel coordinates of the right image are Each pixel x and y coordinates Each pixel Given the x and y coordinates of the image, the radial distortion correction formula for the left image is as follows:

[0028]

[0029] These are the ideal corner points on the left image. The x and y coordinates after radial distortion;

[0030] The formula for radial distortion correction in the right image is as follows:

[0031]

[0032] respectively are the ideal corner points on the right image x, y axis coordinate values after radial distortion;

[0033] wherein, are the distances from the image pixel points to the image center point;

[0034] The tangential distortion correction formula of the left image is as follows:

[0035]

[0036] respectively are the corner points on the left image x, y axis coordinate values after tangential distortion;

[0037] The tangential distortion correction formula of the right image is as follows:

[0038]

[0039] respectively are the corner points on the right image x, y axis coordinate values after tangential distortion;

[0040] In combination with the radial distortion and tangential distortion correction formulas of the left and right cameras, there are the following relational expressions:

[0041]

[0042] The final ideal undistorted left image pixel coordinates and right image pixel coordinates

[0043] 3.3), the reconstructed three-dimensional point coordinates are calculated by the direct linear transformation method, and the projection matrices of the left and right cameras are denoted as P l and P r , and the expression formulas are as follows:

[0044]

[0045] wherein, p l1 , p l2 , and p l3 are row vectors of P l , and p r1 , p r2 , and p r3 are row vectors of P r ;

[0046] According to the direct linear transformation method, the following equation is obtained:

[0047]

[0048] Solving the above equation using singular value decomposition yields the reconstructed 3D points. Let the coordinate set of all reconstructed points be Ω={ q p ij |[ q x ij , q y ij , q z ij ] T}, q x ij , q y ij , q z ij Reconstruction points q p ij The x, y, and z axis coordinates;

[0049] 3.4) Calculate the current residual to obtain the constraints for constructing the 3D reconstruction:

[0050]

[0051] (4) Construct coplanar constraints for 3D points

[0052] 4.1) The normal vector of the plane containing the corner points of the ideal calibration plate in the q-th calibration plate image. Represented by the third column of the rotation matrix R of the camera extrinsic parameters, i.e., the vector

[0053] 4.2) The three-dimensional corner point set Ω reconstructed from the q-th calibration plate image q ={ q p i′ |[ q x i′ , q y i′ , q z i′ ] T In the context of}, denote the starting point from the corner point q p i′ Let i′ = 1, 2, ..., m × n, starting from the corner point. q p j′ The unit direction vector of the vector j′=1,2,...,m×n,j′≠i′ is Right now:

[0054]

[0055] 4.3) Transform the vector and unit direction vector The sum of the absolute values of the inner products, as an error value measuring the degree of perpendicularity between the surface vector and the normal vector of the surface, constitutes the coplanar constraint term:

[0056]

[0057] (5) Construction of the collinear constraint of 3D points

[0058] 5.1) For the q-th image of the calibration board, find the reconstructed three-dimensional corner point set Ω q ={ q p i′ | q x i′ , q y i′ , q z i′} T} in which the first end point and the tail end point of the r-th row of the calibration board in space are located q p rs and q p re , connect the first end point and the tail end point to get the collinear vector i.e.

[0059] 5.2) For the k-th intermediate point of the r-th row q p rk , similarly, connect the intermediate point q p rk and the first end point q p rs to get the second vector

[0060] 5.3) Then, calculate the projection point of the vector on the collinear vector of the r-th row , and the calculation formula is: where is the unit direction vector of the collinear vector of the row;

[0061] 5.4) Then, calculate the distance between the intermediate point q p rk and its projection point , that is, the distance between the intermediate point q p rk and the collinear vector of the r-th row, as one of the residual terms in the reconstructed collinear error, and the same processing is performed for all rows and columns and all points. Calculate the final collinear error term E line :

[0062]

[0063] where,q p k′c 、 respectively the k'th intermediate point of the c'th column of the q'th image of the calibration board in space and its projection point on the column vector of the column;

[0064] (6) Construct radial constraint of 3D point

[0065] 6.1) 3D point set Ω of the four corners of the q'th image of the calibration board q ={ q p i′ |[ q x i′ , q y i′ , q z i′ ] T} whose numbers are selected from left to right and from top to bottom q p1、 q p n 、 q p m×n 、 q p (m-1)×n+1 ;

[0066] 6.2) The two errors of the radial constraint are:

[0067]

[0068] where radius is the ideal radial length of the calibration board;

[0069] (7) Construct angle constraint of 3D point

[0070] 7.1) The four corners q p1、 q p n 、 q p m×n 、 q p (m-1)×n-1 selected in step 6.1) are respectively denoted as points A, B, C, D, and form a quadrilateral, and the direction vectors of the line segments AB, AD, BC, DC are respectively written as have:

[0071]

[0072] 7.2) Calculate the four angles of the quadrilateral, and calculate the residual error of each angle value with the right angle. For all calibration images, this angle error can be written as:

[0073]

[0074] (8), considering all the above 6 constraint terms, a new constraint is obtained, that is, the optimization objective function is:

[0075] min E total = E pix + E pnt + E pl + E line + E diag + E ang

[0076] The camera calibration parameters to be optimized are written as follows:

[0077]

[0078] The optimization objective function is simplified as:

[0079]

[0080] In the formula, θ represents the set of camera calibration parameters to be optimized {A l , A r , D l , D r , R, T}, and e is the set of square roots of all constraint functions

[0081] (9), the initial value of the camera calibration parameter variable is initialized as θ←θ0, the Levenberg-Marquardt method is used, the camera calibration parameters to be optimized are optimized according to the optimization objective function, and the optimal solution θ of the camera calibration parameters is obtained * , wherein θ0 is the camera calibration parameter after calibration.

[0082] The purpose of the application is achieved as follows:

[0083] To improve the accuracy of the camera calibration parameters, the result of the three-dimensional reconstruction of the calibration board should be consistent with the actual three-dimensional scene, that is, the coordinates of the reconstructed points should be consistent with the actual spatial coordinates, and the relative positions between the points should not change, and the geometric characteristics of the calibration board in space should remain unchanged after calibration and reconstruction. Considering the above situations, the optimization method of the binocular stereo vision camera calibration parameters constructs corresponding geometric constraints for different situations, and all the constraint terms are superimposed to form the final optimization objective function. Finally, the Levenberg-Marquardt method is used to optimize the objective function, and the optimized camera calibration parameters are obtained. The camera calibration parameter optimization algorithm proposed in the application has higher accuracy compared with the Zhang Zhengyou camera calibration method. After calibration using the optimized parameters, the reprojection error and the three-dimensional reconstruction measurement error are greatly reduced. BRIEF DESCRIPTION OF DRAWINGS

[0084] Figure 1is a schematic diagram of coordinate conversion in the binocular camera calibration process;

[0085] Figure 2 is a flow chart of a specific embodiment of the method for optimizing binocular stereo vision camera calibration parameters of the present application;

[0086] Figure 3 is a coplanar constraint schematic diagram;

[0087] Figure 4 is a collinear representation schematic diagram;

[0088] Figure 5 is a radial constraint schematic diagram;

[0089] Figure 6 is a binocular captured image: (a) left camera; (b) right camera;

[0090] Figure 7 is a curve of changes in partial internal parameters of the left and right cameras in the iteration process;

[0091] Figure 8 is a curve of changes in external parameters in the iteration process, since there are many parameters in the rotation matrix, the present application transforms it into a representation of quaternions for ease of display;

[0092] Figure 9 is a curve of changes in parameters in the optimization algorithm, wherein (a) is a curve of changes in the gain ratio ρ, (b) is a curve of changes in the correction factor ε, (c) is a curve of changes in the step adjustment factor v, and (d) is a curve of changes in the gradient value ||J(θ u )||;

[0093] Figure 10 is a curve of changes in the objective function in the iteration process;

[0094] Figure 11 is a comparison diagram of re-projection errors of each single calibration pose, wherein the horizontal coordinate represents the calibration pose, and the vertical coordinate represents the re-projection error, the lower curve is the re-projection error of the single calibration pose of the present application, and the upper curve is the re-projection error of the single calibration pose of the Zhang Zhengyou calibration method;

[0095] Figure 12 is the coordinate error of the re-projection point; (a)-(b) are results of the traditional method; (c)-(d) are results of the algorithm of the present application; (a), (c) are diagrams in which the horizontal coordinate axis represents the absolute value of the difference between the re-projection point and the ideal corner pixel coordinate in the horizontal coordinate, and the vertical coordinate represents the absolute value of the difference between the re-projection point and the ideal corner pixel coordinate in the vertical coordinate axis; (b), (d) are diagrams in which the horizontal coordinate represents the difference between the re-projection point coordinate value and the ideal coordinate, and the vertical coordinate is the number of re-projection points;

[0096] Figure 13is a point cloud instance reconstruction method. DETAILED DESCRIPTION

[0097] The specific embodiments of the present application are described below in conjunction with the accompanying drawings, so that those skilled in the art can better understand the present application. It should be particularly noted that in the following description, when the detailed description of the known functions and designs may obscure the main content of the present application, these descriptions will be omitted here.

[0098] Figure 1 is a flow chart of a specific embodiment of the binocular stereo vision camera calibration parameter optimization method of the present application.

[0099] In this embodiment, as shown in Figure 2 , the binocular stereo vision camera calibration parameter optimization method of the present application comprises the following steps:

[0100] Step S1: Obtain the calibration board image and the initial value of the camera calibration parameter after calibration

[0101] Obtain b calibration board images, and obtain the initial value of the camera calibration parameter according to the camera calibration method according to the b calibration board images, and take the camera calibration parameter as the camera calibration parameter to be optimized, wherein in the calibration board corner point specification, the number of corner points in each row is n, and the number of corner points in each column is m. In the specific implementation process, the camera calibration method can adopt Zhang Zhengyou camera calibration method to obtain the initial value of the camera calibration parameter.

[0102] Step S2: Construct 2D pixel constraints

[0103] Step S2.1: For the i-th, i = 1, 2,..., m row and the j-th, j = 1, 2,..., n column corner point of the q-th, q = 1, 2,..., b calibration board image q P ij , the world coordinate system space coordinates are q P ij = [ q X ij , q Y ij , q Z ij ] T , q X ij 、 q Y ij 、 q Z ij are the x, y, z axis coordinate values of the corner point q P ij , after directly extracting the corner point from the q-th calibration board image, the corner point q P ij on the left and right camera images the coordinate positions of the corner points are respectively the coordinate values of the x, y axes of the corner points are respectively the coordinate values of the x, y axes of the corner points

[0104] According to the pinhole camera model, the corner points q P ij are projected onto the left and right camera image planes by the left and right camera calibration parameters respectively, to obtain the left camera 2D projection point and the right camera projection point are respectively the coordinate values of the x, y axes of the projection points are respectively the coordinate values of the x, y axes of the projection points

[0105] Taking the left camera coordinate system as the world coordinate system, the parameter conversion relationship of the corner points q P ij on the left camera is:

[0106]

[0107] The parameter conversion relationship on the right camera is:

[0108]

[0109] where w l and w r are unknown scale factors.

[0110] Step S2.2: Compare the coordinates of all the projection points with the coordinates obtained by directly extracting the corner points from the image, calculate the re-projection error, and construct a 2D pixel constraint, i.e.:

[0111]

[0112] Step S3: Construct a 3D reconstruction constraint

[0113] Step S3.1: Through the left camera extrinsic parameter T lwc = [0 0 0] T , the coordinate point set Ω w = { q P ij |[ q X ij , q Y ij , q Z ij ​​​​] T}, q=1,2,...,b, i=1,2,...m, j=1,2,...n, transformed to the left camera coordinate system, yields a set of three-dimensional coordinates of the ideal calibration board corner points in the left camera coordinate system. Points The x, y, and z axis coordinates are given, and the coordinate transformation formula is as follows:

[0114]

[0115] Step S3.2: Extract the pixel coordinates of the calibration board corner points from the left and right images using the camera distortion coefficients. Perform distortion correction, and denote the pixel coordinates of the ideal distortion-free left image as... And the pixel coordinates of the right image are Each pixel x and y coordinates Each pixel Given the x and y coordinates of the image, the radial distortion correction formula for the left image is as follows:

[0116]

[0117] These are the ideal corner points on the left image. The x and y coordinates after radial distortion;

[0118] The formula for radial distortion correction in the right image is as follows:

[0119]

[0120] These are the ideal corner points on the right image. The x and y coordinates after radial distortion;

[0121] in, All are the distances from image pixels to the center point of the image;

[0122] The formula for correcting tangential distortion in the left image is as follows:

[0123]

[0124] The top corner of the left image The x and y coordinates after tangential distortion;

[0125] The formula for correcting tangential distortion in the right image is as follows:

[0126]

[0127] respectively, the upper corner point of the right image x, y axis coordinate values after tangential distortion correction;

[0128] Combining the radial distortion and tangential distortion correction formulas of the left and right cameras, there are the following relationships:

[0129]

[0130] Solving the above equation to get the final ideal undistorted left image pixel coordinates and right image pixel coordinates

[0131] Step S3.3: Calculate the reconstructed three-dimensional point coordinates by direct linear transformation method, and denote the projection matrices of the left and right cameras as P l and P r , respectively, and the expression formula is as follows:

[0132]

[0133] where p l1 , p l2 , and p l3 are row vectors of P l , and p r1 , p r2 , and p r3 are row vectors of P r .

[0134] According to the direct linear transformation method, the following equation is obtained:

[0135]

[0136] Solve the above equation using the singular value decomposition method to get the reconstructed three-dimensional points, and denote the coordinate set of all reconstructed points as Ω = { q p ij |[ q x ij , q y ij , q z ij ] T}, q x ij , q y ij , q z ij are the x, y, z axis coordinate values of the reconstructed point q p ij .

[0137] Step S3.4: Calculate the current residual to get the construction 3D reconstruction constraint:

[0138]

[0139] Step S4: Constructing coplanar constraint of 3D points

[0140] Step S4.1: Normal vector of the plane where the ideal corner point of the q-th image of the calibration board is located with the third column of the rotation matrix R of the camera extrinsic parameter, i.e. the vector

[0141] Step S4.2: 3D point set Ω reconstructed from the q-th image of the calibration board q ={ q p i′ |[ q x i′ , q y i′ , q z i′ ] T}where the unit direction vector of the vector from the corner point q p i′ , i' = 1, 2,..., m x n to the corner point q p j′ , j' = 1, 2,..., m x n, j' ≠ i' is i.e.

[0142]

[0143] Step S4.3: As shown in Figure 3 , in the ideal case, the normal vector of the plane and the vector between any two points on the plane are perpendicular, and the inner product should be equal to 0. In the case of error, the value will not be 0, so the absolute value of the inner product of the vector and the unit direction vector is summed up as an error value to measure the degree of perpendicularity of the vector on the plane to the normal vector of the plane, which constitutes the coplanar constraint term:

[0144]

[0145] Step S5: Constructing collinear constraint of 3D points

[0146] As shown in Figure 4 , the ideal 3D points of the reconstructed calibration board should keep the points in the same row on the same straight line, and the points in the same column also on the same straight line. This feature should not change after the rigid transformation and projection transformation and inverse transformation. However, due to the inaccuracy of the calibration parameters, this item will be introduced into the error, which will cause the points originally in the same column or the same row not to be strictly in the collinear constraint, as shown in Figure 4The left figure. The degree of reverse collinearity can represent the accuracy of the calibration parameters, and as a constraint for calibration parameter optimization. To represent the degree of collinearity of the reconstructed corner points of the calibration board that should be on the same line, the present application constructs a collinearity constraint model of 3D points as follows:

[0147] Step S5.1: For the qth calibration board image, find its reconstructed three-dimensional corner point set Ω q ={ q p i′ [ q x i′ , q y i′ , q z i′ ] T} in the first end point q p rs and the tail end point q p re of the rth row of the calibration board in space, connect the first and last end points to the collinear vector i.e.

[0148] Step S5.2: For the kth intermediate point q p rk of the rth row, similarly, connect the intermediate point q p rk and the first end point q p rs to get the second vector

[0149] Step S5.3: Then, calculate the projection point of the vector on the collinear vector of the rth row, the calculation formula is: where is the unit direction vector of the collinear vector of the row.

[0150] Step S5.4: Then, calculate the distance between the intermediate point q p rk and its projection point , that is, the distance between the intermediate point q p rk and the collinear vector of the rth row, as a residual term in the reconstructed collinearity error, and the same processing is performed for all rows and columns, and all points, to calculate the final collinearity error term E line :

[0151]

[0152] where, q p k′c , The k'th intermediate point of the c'th column of the q'th image of the calibration board in space and its projection point on the column-wise vector of the column.

[0153] Step S6: Constructing radial constraints of 3D points

[0154] The reconstructed corner points of the calibration board should also satisfy the radial constraints as shown in Figure 5 . The farthest distance between points of the reconstructed point cloud set should be consistent. For convenience of representation, the upper left corner of the calibration board is marked as A, the upper right corner as B, the lower right corner as C, and the lower left corner as D. The distance between AD and BC should be consistent with the ideal value.

[0155] Step S6.1: Three-dimensional corner point set Ω of the q'th image of the calibration board q ={ q p i′ | q x i′ , q y i′ , q z i′ ] T}, numbered from left to right and from top to bottom, four corner points q p1、 q p n 、 q p m×n 、 q p (m-1)×n+1 , i.e. Figure 4 , the upper left corner A, the upper right corner B, the lower right corner C, and the lower left corner D.

[0156] Step S6.2: Then the two items of error of the radial constraint are:

[0157]

[0158] where radius is the ideal radial length of the calibration board.

[0159] Step S7: Constructing angle constraints of 3D points

[0160] Step S7.1: The four corner points q p1、 q p n 、 q p m×n 、 q p (m-1)×n-1 selected in step S6.1 are respectively denoted as points A, B, C, and D, and form a quadrilateral. The direction vectors of line segments AB, AD, BC, and DC are respectively written as have:

[0161]

[0162] Step S7.2: Calculate the four angles of the quadrangle, and calculate the residual of each angle value with a right angle. For all the calibration images, the angle error can be written as:

[0163]

[0164] Step S8: Construct the optimization objective function

[0165] Considering all the above six constraint terms, a new constraint, i.e., the optimization objective function, is obtained as follows:

[0166] min E total = E pix + E pnt + E pl + E line + E diag + E ang

[0167] The camera calibration parameters to be optimized are written as follows:

[0168]

[0169] The optimization objective function is simplified as:

[0170]

[0171] In the formula, θ represents the set of the above-mentioned camera calibration parameters to be optimized {A l , A r , D l , D r , R, T}, and e is the set of square roots of all constraint term functions

[0172] Step S9: Take the calibrated camera calibration parameters as initial values, and use the Levenberg-Marquardt method to optimize according to the optimization objective function.

[0173] Initialize the initial value of the camera calibration parameter variable θ←θ0, and use the Levenberg-Marquardt method to optimize the camera calibration parameters to be optimized according to the optimization objective function, to obtain the optimal solution θ * of the camera calibration parameters, wherein θ0 is the calibrated camera calibration parameter.

[0174] Since the distortion model of the camera is nonlinear, after the partial derivative, the Levenberg-Marquardt (L-M) method can be used to optimize the objective function. According to the principle of the L-M method, the update rule of the camera calibration parameters in the optimization process has the following form:

[0175]

[0176] where u is the iteration number, ε u is the step size adaptive parameter for the u-th iteration, ε u I is a positive definite modification term, which is adjusted at each iteration to ensure the reduction of the residual. When ε u is large, the convergence is similar to the gradient descent method, the iteration number is reduced, θ u is quickly reduced to the vicinity of the minimum; when ε u is small or even close to 0, the step size is reduced, and the convergence is stable, similar to the Gauss-Newton method. The adjustment of ε u requires the introduction of a damping adjustment factor ρ u :

[0177]

[0178] ρ u represents the difference between the approximate model and the actual model function descent, the numerator is the actual model descent, and the denominator is the approximate model descent. The larger the value of ρ u , the closer to 1, which proves that the approximate model is reliable, and the convergence is close to the minimum, which can reduce the value of ε u , and perform similar Gauss-Newton convergence; if the value of ρ u is small, it means that it is still in the gradient descent stage, far from the extreme value, and the value of ε u can be increased; and when the value of ρ u is negative, it means that the convergence direction is opposite, and the value of ε u needs to be increased to make ρ u positive, correcting the convergence direction. J(θ u ) is the Jacobian matrix of the system, which can be written as:

[0179]

[0180] Let the initial estimate of the calibration parameter θ0∈R 30 , the optimization objective function term e={e1,...,e6}∈R 30 →R is a mapping from a 30-dimensional real space to a 1-dimensional real space, the iteration step size adaptive parameter u max , v in , τ, and the termination condition judgment threshold σ1, σ2, σ3. Then, the solution of the optimization algorithm system is obtained, i.e. θ * ∈R 30 calibration parameters, so that the objective function is minimized.

[0181] The detailed steps of the optimization algorithm are as follows:

[0182] Step S9.1: initialize the initial value of the calibration parameter variable θ <- θ0, the initial parameter of the adaptive iteration step size v <- v in and the iteration number u <- 0;

[0183] Step S9.2: to avoid dead loop, for the u-th iteration, when the amplitude J(θ u ) T e(θ u ) converges to a threshold value close to 0, it means that the optimization algorithm has reached a minimum value, and the termination condition flag η can be calculated according to the following formula:

[0184]

[0185]

[0186] When η = true, the iteration process of the algorithm will terminate;

[0187] Step S9.3: determine whether the iteration is terminated, the termination condition is: η = true or u ≥ u max , if true, jump to step S9.10; otherwise, go to the next step S9.4;

[0188] Step S9.4: according to the optimization formula, the LM method is a modification of the Gauss-Newton method, in order to ensure that the direction generated by each iteration is a descending direction, that is, to ensure that the Hessian matrix (J(θ u ) T J(θ new )+ε u I) is positive definite, a positive definite correction term factor ε u is introduced. At this time, the semi-positive definite matrix can be obtained from the product of the transpose and itself of the Jacobian matrix, and the correction factor ε u can be calculated by selecting the maximum value of the diagonal vector of the semi-positive definite matrix, the calculation process is as follows:

[0189]

[0190] The correction factor ε u can be assigned according to the input correction parameter τ:

[0191]

[0192] Step S9.5: calculate the update change value Δθ u of the calibration parameter in the current u-th iteration. According to the positive definite correction factor, the semi-positive definite matrix obtained in step S9.4 is corrected to a positive definite matrix to ensure that the overall solution is strictly minimum. And multiply its inverse matrix by the judgment matrix in step S9.2 to obtain the update change value, the process is represented by the formula:

[0193]

[0194] Where, Δθ u The dimension is 30, and the first element is Δf xl-u :

[0195]

[0196] Similarly, the updated values ​​of other elements can be calculated;

[0197] Step S9.6: Determine if the updated change value meets the iteration termination condition. When the updated change value is too small, it is desirable to terminate the iteration. Therefore, if max(Δf) xl-u ,Δf yl-u ,...,Δt 3-u )≤σ2·max(f xl-u ,f yl-u ,...,t 3-u If ), then s = true, proceed to step S9.10; otherwise, proceed to step S9.7.

[0198] Step S9.7: Update the calibration parameter θ based on the change value. new ←θ u +Δθ u And calculate the current gain ratio ρ u This parameter is used to determine the quality of the current step size and provides an indicator for adaptive step size updates. ρ u The calculation method is as follows:

[0199] Step S9.7.1: First, based on the judgment matrix in step S9.2, calculate:

[0200]

[0201] Step S9.7.2: Then, calculate the gain ratio:

[0202]

[0203] Step S9.8: From step S9.4, we know that J(θ) u ) T J(θ u () is positive semi-definite, in calculating (J(θ) u ) T J(θ u )) -1 During the process, J(θ) u ) T J(θ u) may be singular or ill-conditioned, which will lead to unstable increment of parameters and even divergence of algorithm. Therefore, the L-M method introduces a trust region on the basis of the Gauss-Newton method, which can limit the increment of parameters according to the value of gain ratio ρ u Update the trust region interval. When ρ u > 0, it indicates that the quality of the model is poor, and the trust region needs to be reduced, and the adaptive adjustment factor of the step size is set to the initial value v = v in ; otherwise, the search radius can be appropriately increased by the adaptive step size adjustment factor v, and the following is executed in turn: ε u+1 = ε u · v; v = v in · v;

[0204] Step S9.9: judge whether the termination condition is satisfied or not, when ||J(θ u ) T e(θ u )|| ∞ ≤ σ1 or e 2 (θ u+1 ) ≤ σ3, set the value of termination flag η to true, and enter step S9.10; otherwise, let the iteration number u = u + 1, and return to step S9.2 for new loop;

[0205] Step S9.10: output the calibration parameter value θ * = θ u at this time as the optimal solution.

[0206] Examples

[0207] In this example, the optimization of camera calibration parameters is performed for a convergent binocular stereo vision system. The adopted calibration board has a specification of 11 columns x 8 rows (i.e. m = 11, n = 8), and there are 88 inner corners on the board, and the distance between each inner corner point is 15.00 mm. The resolution of the convergent binocular camera used in this example is 4096 x 3000 pixels, and the pose of the calibration board is switched for a total of 10 times (b = 10) during the experiment, so there are a total of 20 images captured by the left camera I lk and the right camera I rk , as shown in Fig. 1. Figure 6

[0208] The initial camera calibration parameter estimation value, i.e. the camera calibration parameter result of the traditional calibration algorithm, is shown in Table 1.

[0209]

[0210] Table 1

[0211] ​Using the above calibration results to re-project the source data, the average re-projection error value is 1.61367 pixels.

[0212] According to the expression of each constraint term introduced in the application, the overall optimization objective function is established:

[0213]

[0214] The calibration results in Table 1 are input as initial parameter values θ0 into the calibration optimization algorithm, the maximum number of iterations of the optimization algorithm is set to u max = 50, the initial damping coefficient value is v in = 2, τ = 0.8, the termination threshold parameter gradient tolerance σ1 = 1.0e -10 , the parameter tolerance σ2 = 1.0e -8 , the function tolerance σ3 = 1.0e -6 , and the initial value is substituted into step (9) for cyclic iteration solution.

[0215] This example has a total of 27 iterations, combined with Figure 7 , Figure 8 and Figure 9 , it can be seen that in the first few iterations, the objective function value (as shown in Figure 7 and Figure 8 ) and the gradient value ||J(θ u )|| change greatly, so the trust region ε should also gradually increase, as shown in Figure 9 (a), the maximum value of ρ appears at this time, and the step size adaptive adjustment factor v is still the initial value. With the gradual iteration of the algorithm, the change value of the objective function gradually decreases, and the gradient value also gradually tends to 0, in this process, the trust region gradually decreases and approaches the optimal solution.

[0216] The output change curve of the objective function is shown in Figure 10 , it can be seen that the value of the objective function is finally stabilized near 20, which is due to the assumption that the objective function is 0 is only in the ideal case, in the actual situation, there is a certain error, which leads to the objective function cannot be accurately 0, the final convergence and stability of the value shows the feasibility of the algorithm. The optimized camera calibration parameters obtained by the application are shown in Table 2.

[0217]

[0218] Table 2

[0219] In order to verify the effectiveness of the binocular stereo vision camera calibration parameter optimization method of the application, the average re-projection error is calculated as 0.649198 pixels, as shown in Figure 11The re-projection error comparison chart of the re-projection error of the optimized result of the present application is shown, and the result is reduced by 59.77% compared with the average re-projection error of the traditional calibration algorithm.

[0220] In order to more intuitively compare the optimization results of the present application and the traditional Zhang calibration algorithm, the re-projection results of each calibration board corner point are compared with the ideal value, and the experimental results are shown in Figure 12 Figure 12 (a) It can be seen that the re-projection results of the traditional Zhang calibration algorithm wrap 95% of the points in a circle with a radius R, and the area of the coordinate difference set is S=πR 2 Figure 12 (c) In the present application, 95% of the re-projection results can be covered only by a sector with a radius r, and the area of the sector is s=0.25πr 2 <S, the coordinate error distribution of the output results of the present application is more concentrated than that of the traditional Zhang camera calibration parameter method, and the statistical histogram Figure 12 (b), 12(d) also intuitively verifies the above conclusion.

[0221] The accuracy of the camera calibration parameters needs to be further verified in three-dimensional reconstruction. In this example, the parameters of the Zhang calibration algorithm are used to reconstruct the three-dimensional corner points, and the average distance between the corner points is taken as the final reconstruction result, which is 14.9984 mm. Then the camera calibration parameters optimized by the present application are used to reconstruct the corner points, and the average corner point distance is 15.0006 mm. It is known that the inner focal point distance of the source calibration board is 15.0000 mm. Therefore, the reconstruction errors of the two camera calibration algorithms are -0.0016 mm and 0.0006 mm, respectively. Compared with the Zhang camera calibration algorithm, the reconstruction error of the binocular stereo vision camera calibration parameter optimization method of the present application is reduced by 62.50%.

[0222] In summary, compared with the traditional Zhang calibration method, the binocular stereo vision calibration parameter optimization algorithm proposed in the present application has higher accuracy, which is reflected in the reduction of the re-projection error by 59.77% under the experimental data, and the measurement result error of the three-dimensional reconstruction is reduced by 62.50%.

[0223] Although the above describes the specific embodiments of the present application for the purpose of facilitating the understanding of the present application by those skilled in the art, it should be clear that the present application is not limited to the scope of the specific embodiments, and for those skilled in the art, it is obvious that various changes are within the spirit and scope of the present application defined and determined by the appended claims, and all the inventions utilizing the concept of the present application are within the scope of protection.​​

Claims

1. A method for optimizing calibration parameters of a binocular stereo vision camera, characterized in that, include: (1) Obtain b calibration board images, and obtain the initial values ​​of camera calibration parameters based on the camera calibration method according to the b calibration board images, and use them as camera calibration parameters to be optimized. In the calibration board corner point specification, the number of corner points in each row is n, and the number of corner points in each column is m. (2) Constructing 2D pixel constraints 2.1) For the qth calibration plate image, the corner point in the i-th row and j-th column... q P ij Its world coordinate system spatial coordinates are q P ij =[ q X ij , q Y ij , q Z ij ] T , q X ij , q Y ij , q Z ij Corner points q P ij The x, y, and z axis coordinates are obtained by directly extracting corner points from the q-th calibration board image. q P ij Top corner of the left and right camera images The coordinates are respectively Corner points The x and y coordinates of the axis. Corner points The x and y coordinates of the , where, q=1,2,…,b, i=1,2,…m, j=1,2,…n; Based on the pinhole camera model, the corner points are determined by calibrating the left and right cameras respectively. q P ij Projected onto the image planes of the left and right cameras, the 2D projection point of the left camera is obtained. and right camera projection point They are the projection points The x and y coordinates of the axis. They are the projection points The x and y coordinates; If the left camera coordinate system is used as the world coordinate system, then the corner points q P ij The parameter transformation formula on the left camera is: The parameter conversion formula on the right camera is: Among them, w l and w r For an unknown scale factor, A l For the left camera's intrinsic parameters, A r R is the intrinsic parameter of the right camera, R is the rotation parameter between the left and right cameras, and T is the translation vector; 2.2) Compare the coordinates of all projected points with the coordinates of the corner points extracted directly from the image, calculate the reprojection error, and construct a 2D pixel constraint, i.e.: (3) Constructing 3D reconstruction constraints 3.1) Through the external parameters of the left camera T lwc =[0 0 0] T The coordinate set Ω of the calibration plate corner points in the world coordinate system w ={ q P ij |[ q X ij , q Y ij , q Z ij ] T }, q=1,2,…,b, i=1,2,…m, j=1,2,…n Transformed to the left camera coordinate system, a series of three-dimensional coordinate points of the ideal calibration board corner points in the left camera coordinate system are obtained. Points The x, y, and z axis coordinates are given, and the coordinate transformation formula is as follows: 3.2) Pixel coordinates of calibration board corner points extracted from the left and right images using camera distortion coefficients. Perform distortion correction, and denote the pixel coordinates of the ideal distortion-free left image as... And the pixel coordinates of the right image are Each pixel x and y coordinates Each pixel Given the x and y coordinates of the image, the radial distortion correction formula for the left image is as follows: These are the ideal corner points on the left image. The x and y coordinates after radial distortion, k 1l k 2l k 3l The radial distortion coefficient of the left camera; The formula for radial distortion correction in the right image is as follows: These are the ideal corner points on the right image. The x and y coordinates after radial distortion, k 1r k 2r k 3r The radial distortion coefficient of the right camera; in, All are the distances from image pixels to the center point of the image; The formula for correcting tangential distortion in the left image is as follows: The top corner of the left image The x and y coordinate values ​​after tangential distortion, p 1l p 2l p represents the tangential distortion coefficient of the left camera. 1r p 2r The tangential distortion coefficient of the right camera; The formula for correcting tangential distortion in the right image is as follows: in, The top corner of the right image The x and y coordinates after tangential distortion; Combining the formulas for correcting radial and tangential distortion of the left and right cameras, we have the following relationship: Solving the above equation yields the final ideal, distortion-free pixel coordinates of the left image. and right image pixel coordinates 3.3) The reconstructed 3D point coordinates are calculated using the direct linear transformation method. Let P be the projection matrices of the left and right cameras, respectively. l P r Its formula is as follows: Where, p l1 ,p l2 ,p l3 For P l row vector, p r1 ,p r2 ,p r3 For P r The row vector; According to the direct linear transformation method, we have the following equation: Solving the above equation using singular value decomposition yields the reconstructed 3D points. Let the coordinate set of all reconstructed points be Ω={ q p ij |[ q x ij , q y ij , q z ij ] T }, q x ij , q y ij , q z ij Reconstruction points q p ij The x, y, and z axis coordinates; 3.4) Calculate the current residual to obtain the constraints for constructing the 3D reconstruction: (4) Construct coplanar constraints for 3D points 4.1) The normal vector of the plane containing the corner points of the ideal calibration plate in the q-th calibration plate image. Represented by the third column of the rotation matrix R of the camera extrinsic parameters, i.e., the vector 4.2) The three-dimensional corner point set Ω reconstructed from the q-th calibration plate image q ={ q p i′ |[ q x i′ , q y i′ , q z i′ ] T In the context of}, denote the starting point from the corner point q p i′ Departure to Corner q p j′ The unit direction vector of the vector is Right now: in, i′=1,2,…,m×n, j′=1,2,…,m×n, j′≠i′; 4.3) Transform the vector and unit direction vector The sum of the absolute values ​​of the inner products, used as an error measure of the perpendicularity between a vector on the surface and the normal vector of the surface, constitutes the coplanarity constraint term: (5) Construct collinear constraints for 3D points 5.1) For the q-th calibration board image, find its reconstructed 3D corner point set Ω. q ={ q p i′ |[ q x i′ , q y i′ , q z i′ ] T The first endpoint of the calibration plate in the r-th row of space is located in} q p rs and tail endpoint q p re Connecting the first and last endpoints yields a collinear vector. Right now 5.2) For the k-th intermediate point of the r-th row q p rk Similarly, the intermediate point q p rk With the first endpoint q p rs The concatenation yields a second vector. 5.3) Next, calculate the vector. Collinear vectors in the r-th row The projection point on the surface is calculated using the following formula: in, It is the unit direction vector of the collinear vectors in that row; 5.4) Then, calculate the intermediate point. q p rk Its projection point The distance between, i.e., the midpoint q p rk The distance to the collinear vector in the r-th row is used as a residual term in the reconstructed collinearity error. This process is repeated for all rows, columns, and points to calculate the final collinearity error term E. line : in, q p k′c , These are the k′-th intermediate point of the c-th column in space of the q-th calibration board image and its projection point on the collinear vector of that column, respectively. (6) Constructing radial constraints for 3D points 6.1) The three-dimensional corner point set Ω of the q-th calibration plate image q ={ q p i′ |[ q x i′ , q y i′ ,qz i′ ] T The four corner points are selected from left to right and top to bottom. q p1、 q p n , q p m×n , q p (m-1)×n+1 ; 6.2) Then the two errors of the radial constraint are: Where radius is the radial length of the ideal calibration plate; (7) Constructing 3D point angle constraints 7.1) Select the four corner points from step 6.1). q p1、 q p n , q p m×n , q p (m-1)×n+1 Let A, B, C, and D be the points, forming a quadrilateral. The direction vectors of line segments AB, AD, BC, and DC are respectively written as... have: 7.2) Calculate the four angles of the quadrilateral and perform residual calculations between each angle value and the right angle. For all calibration images, this angle error can be written as: (8) Considering all six constraints mentioned above, a new constraint, i.e., the optimization objective function, is obtained as follows: min E total =E pix +E pnt +E pl +E line +E diag +E ang The camera calibration parameters to be optimized are written in the following form: in, f xl f yl These are the physical lengths of one pixel in the x and y directions on the left camera sensor, respectively. 0l v 0l f represents the difference in horizontal and vertical pixels between the center pixel coordinates and the origin pixel coordinates of the left image, respectively. xr f yr These are the physical lengths of one pixel in the x and y directions on the right camera sensor, respectively. 0r v 0r D represents the difference in horizontal and vertical pixels between the center pixel coordinates and the origin pixel coordinates of the right image, respectively. l D represents the distortion coefficient of the left camera. r R1, ..., R9 are the distortion coefficients of the right camera, and r1, ..., r9 are functions of the rotation angles of the x, y, and z axes, respectively. t1, t2, and t3 are the translation dimensions in the x, y, and z directions, respectively. The objective function can be simplified as follows: In the formula, θ represents the set of camera calibration parameters to be optimized {A}. l A r D l D r Let e ​​be the set of all square roots of the constraint term functions, R, T. (9) Initialize the initial values ​​of the camera calibration parameter variables θ←θ0. Using the Levenberg-Marquardt method, optimize the camera calibration parameters to be optimized according to the objective function to obtain the optimal solution θ of the camera calibration parameters. * Where θ0 is the camera calibration parameter after calibration.

Images