A camera calibration method suitable for a palletizing robot

Abstract

The application discloses a camera calibration method suitable for a palletizing robot, and relates to the technical field of camera calibration.The method comprises the following steps: selecting a calibration board, wherein the calibration board is composed of two concentric circles and six auxiliary straight lines passing through the center of the concentric circles and uniformly distributed along the circumferences of the concentric circles; acquiring an image of the calibration board, and pre-processing the image to extract an ellipse and a straight line sub-pixel point set; combining an axial distance separation strategy and a concentric constraint to complete concentric ellipse fitting; solving the auxiliary straight line through angle clustering and gradient amplitude weighted SVD; extracting the intersection of the auxiliary straight line and the concentric ellipse as a feature point, and establishing a mapping relationship between the feature point in a pixel coordinate system and a point in a world coordinate system; and completing camera calibration by using a perspective projection model according to the matching pairs of the feature points in the pixel coordinate system and the world coordinate system.The application effectively alleviates the problem of insufficient discrimination accuracy of edge confusion points caused by the change of the axial distance due to perspective distortion, and effectively improves the calibration accuracy and stability.

Description

Technical Field

[0001] This invention relates to the field of camera calibration technology, and in particular to a camera calibration method suitable for palletizing robots. Background Technology

[0002] In industrial palletizing operations, cameras responsible for acquiring visual information are typically mounted above or around the palletizing equipment to identify targets within the palletizing area. Camera calibration enables the palletizing equipment to accurately identify the target's position in the world coordinate system, providing precise positional information for subsequent operations such as grasping, handling, and stacking. Therefore, camera calibration is a crucial prerequisite and essential guarantee for ensuring the accurate operation of palletizing equipment.

[0003] Based on whether a target is used during the calibration process, existing calibration methods can be divided into two categories. One category is the targetless method. During palletizing operations, due to factors such as changes in lighting, irregular material shapes, and cluttered backgrounds, targetless calibration is difficult to implement and leads to unstable calibration results and insufficient accuracy.

[0004] Another type is target-based calibration methods, with two-dimensional targets being the most widely used. Two-dimensional targets are typically square or circular in shape. Square targets, using a checkerboard pattern, rely on corner point extraction, which can easily lead to false detections or loss of targets in industrial environments with uneven lighting, dust, or blurred images, thus affecting calibration accuracy to some extent. Compared to square targets, circular targets are more suitable for camera calibration in complex environments, offering advantages in contour feature extraction, sensitivity to plane flatness, and anti-interference capabilities. Therefore, they are more suitable for camera calibration of palletizing equipment.

[0005] Existing circular target calibration methods are generally divided into two categories: one is a circular matrix composed of multiple independent circles, and the other is a target composed of multiple nested concentric circles. These methods typically rely on extracting the center points of the circles within the target as feature points for calibration. However, the following problems exist during calibration:

[0006] (1) When extracting the center point of the circle as a feature point, the center point has an unavoidable positioning error due to camera perspective projection and lens distortion, which limits the calibration accuracy.

[0007] (2) The axial distance variation of the segmented ellipse was ignored during ellipse fitting, which resulted in limited classification accuracy of edge confusion points, low ellipse fitting quality, and affected the accuracy of subsequent calibration.

[0008] (3) When optimizing the ellipse fitting parameters, the concentric geometric constraints are ignored, which affects the accuracy of the parameters and restricts the improvement of the fitting quality.

[0009] Existing methods for fitting ellipses to the inner and outer circles of a concentric target require optimizing the center positions, rotation angles, and major and minor axis parameters of the inner and outer circles. These methods solve for the ellipse parameters of the inner and outer circles independently, without utilizing the concentric geometric constraints of the circles. This approach not only leads to decreased accuracy in parameter calculation and affects the quality of the ellipse fitting, but also deviates from the original intention of the concentric circle design in the target, thus impacting subsequent calibration accuracy. Summary of the Invention

[0010] To overcome the aforementioned problems in the prior art, this invention proposes a camera calibration method suitable for palletizing robots.

[0011] The technical solution adopted by this invention to solve its technical problem is: a camera calibration method suitable for palletizing robots, comprising the following steps:

[0012] Step 1, Calibration plate selection: The calibration plate consists of two concentric circles and six auxiliary straight lines passing through the center of the circles, evenly distributed along the circumference of the concentric circles.

[0013] Step 2: Obtain the calibration board image and preprocess the image to extract the sub-pixel sets of ellipses and lines;

[0014] Step 3: Combine the axial distance separation strategy with concentric constraints to complete the concentric ellipse fitting;

[0015] Step 4: Solve for the auxiliary line using angle clustering and gradient magnitude-weighted SVD;

[0016] Step 5: Extract the intersection points of the auxiliary line and the concentric ellipse as feature points, and establish the mapping relationship between feature points in the pixel coordinate system and points in the world coordinate system;

[0017] Step 6: Based on the feature point matching pairs in the pixel coordinate system and the world coordinate system, complete the camera calibration using the perspective projection model.

[0018] In the aforementioned camera calibration method for palletizing robots, step 1 involves auxiliary lines evenly distributed at angles of 0°, 30°, 60°, 90°, 120°, and 150°. All lines pass through the center of the circle, and each line intersects the inner and outer circles at two points. The auxiliary lines and concentric circles form 24 feature points. The origin of the world coordinate system is... o w Located at the center of concentric circles, x w The axis to the left is positive. y w The axis pointing downwards is positive. z w The axis perpendicular to the plane of the calibration plate and pointing outwards toward the camera is positive.

[0019] The above-mentioned camera calibration method for palletizing robots, specifically step 2, involves:

[0020] Step 2.1: After converting the calibration board image into a grayscale image, Gaussian filtering and edge detection are performed to obtain the edge map.

[0021] Step 2.2: Extract the edge point set from the edge map obtained in Step 2.1;

[0022] Step 2.3: Perform gradient-based sub-pixel refinement on the edge points in the edge point set obtained in Step 2.2;

[0023] Step 2.4: Separate the elliptical and linear point sets from the middle edge points of the sub-pixel refined edge point set obtained in Step 2.3.

[0024] In the aforementioned camera calibration method for palletizing robots, step 3 specifically comprises:

[0025] Step 3.1: Based on the elliptical sub-pixel point set obtained in Step 2, the segmented ellipse is fitted using a rotated ellipse geometric model. The rotated ellipse geometric model is iteratively solved using a nonlinear least squares criterion to obtain the geometric parameters of the segmented ellipse. ,in( , () represents the coordinates of the center of the dividing ellipse. , These represent the lengths of the semi-major axis and the semi-minor axis, respectively. The rotation angle of the ellipse is the difference between the local coordinate system of the ellipse and the pixel coordinate system. x The angle between the axes;

[0026] Step 3.2, based on the coordinates of the center of the segmented ellipse and ellipse rotation angle The edge points in the sub-pixel point set of the ellipse are transformed from pixel coordinates to the local coordinate system of the ellipse, the axial distance is calculated, and the point set of the ellipse inside and outside the calibration plate is divided according to the axial distance.

[0027] Step 3.3: Set the initial values ​​of the center and rotation angle of the segmented ellipse. Using the initial values ​​of the center and rotation angle of the inner and outer ellipses, the point sets of the inner and outer ellipses obtained in step 3.2 are fitted with ellipses using the geometric model of the rotated ellipse to obtain the initial semi-axis length of the outer ellipse. Length of the initial semi-axis of the inner ellipse An ellipse fitting optimization model is established based on the normalized geometric residuals of ellipse fitting, which establishes concentric geometric constraints.

[0028] The aforementioned camera calibration method for palletizing robots specifically uses the following formula for the rotating ellipse geometric model:

[0029]

[0030] Among them, optimization variables For the geometric parameters of the ellipse, ( x , y ) represents the pixel coordinates of the sub-pixel point concentration point of the ellipse.

[0031] The aforementioned camera calibration method for palletizing robots, specifically the ellipse fitting optimization model based on the ellipse fitting normalized geometric residuals to establish concentric geometric constraints in step 3.3, is as follows:

[0032] The normalized geometric residual of the ellipse fitting is:

[0033]

[0034] in, For an elliptical sub-pixel set, For point The axial distance to the ellipse;

[0035] Ellipse fitting optimization model that establishes concentric geometric constraints using normalized residuals. The joint parameters to be optimized are ,in, The center coordinates are shared by the inner and outer ellipses. θ For the shared rotation angle, and These represent the lengths of the semi-major axis and semi-minor axis of the outer ellipse and the inner ellipse, respectively.

[0036] The initial parameters for model optimization are set to The LM algorithm is used to perform nonlinear least squares optimization on the ellipse fitting optimization model, and the parameters are updated iteratively. ε J Until convergence, the final parameters of the inner and outer ellipses are obtained:

[0037] .

[0038] The above-mentioned camera calibration method for palletizing robots, specifically step 4, involves:

[0039] Step 4.1, for the points in the sub-pixel set of the straight line obtained in Step 2 ( Find the polar angle of this point relative to the center of the ellipse. :

[0040]

[0041] in, , For the sub-pixel set of straight lines i Each point is relative to the center of the ellipse. The coordinate difference The polar angle is measured in degrees, with a value ranging from -180° to 180°. Mapping to [-0°, 360°) yields Based on the concentration point of sub-pixel points of the straight line and The angular differences are assigned to the corresponding straight lines and points. The corresponding line number is This yields the edge point sets of different straight lines. ,in k =1,2,...,6, N k The number of points in each group;

[0042] Step 4.2, regarding the first k edge point set of a straight line Calculate the weight of each point one by one. w k,i :

[0043]

[0044] in, For the first k The first straight line i Points The gradient magnitude;

[0045] No. k Weighted centroid of a set of points on a straight line C k ( x , y )for:

[0046]

[0047] Regarding the first k Construct a weighted decentralized coordinate matrix using straight lines. Q k :

[0048]

[0049] in, For the first k The first straight line i The weighted centroid of each point relative to this line. C k ( x , y The decentralized coordinate vector of ).

[0050] For matrix Q k Perform singular value decomposition and take the right singular vector corresponding to the minimum singular value. As the unit normal vector of the fitted line, the final... k The straight lines are: The parameters satisfy: .

[0051] In the aforementioned camera calibration method for palletizing robots, step 5 specifically comprises:

[0052] Step 5.1, let the radius of the inner circle of the calibration plate in the world coordinate system be... The outer radius is The target is marked with six evenly distributed auxiliary straight lines, whose orientation angles satisfy... θ k =30•( k -1),( k =1,2,...,6), the world coordinates of any feature point on the annulus are represented as:

[0053]

[0054] In the formula, the radius parameter r Desirable or These correspond to feature points on the inner and outer circles, respectively. Following the rule of sorting all feature points from the inner circle to the outer circle, and from smallest to largest direction angle within the same annulus, an ordered world feature point set p is constructed. w This provides a unified index benchmark for subsequent feature matching;

[0055] Step 5.2: The feature points in the pixel coordinate system are obtained by solving the equations of the ellipse and the line simultaneously. The polar angle of the feature points in the pixel coordinate system relative to the center of the ellipse is calculated, and the polar angle is mapped to the interval [-0°, 360°). The feature points are then transformed from the pixel coordinate system to the local elliptical coordinate system, and the axial distance is calculated. Belongs to the outer circle, The points belong to the inner circle; following the rule of first the inner circle, then the outer circle, and then sorting the points within the same annulus by direction angle from smallest to largest, the feature points in the pixel coordinate system are indexed hierarchically. First, the inner circle point set and the outer circle point set are divided, and then the points within the same annulus point set are sorted by polar angle from smallest to largest to obtain the feature point set. p pix ;

[0056] Step 5.3: Sort the set of pixel sequences. p pix World Feature Point Set p w Map feature points one-to-one by index to construct feature point matching pairs. , ).

[0057] In the aforementioned camera calibration method for palletizing robots, step 6 specifically comprises:

[0058] Based on the extrinsic parameters, the feature points are matched to the {( , World feature points within )} Transformed to the camera coordinate system, after normalization and second-order radial distortion mapping, and combined with the parameter vector to be optimized. φ =[ f x , f y , c x , c y , k 1 , k 2 ] T Obtain the predicted pixel coordinates corresponding to the world point. Construct the least squares objective function for reprojection error and solve for the optimal parameters:

[0059]

[0060] in, Let i be the coordinates of the measured pixel. For the current parameter vector The calculated predicted pixel coordinates;

[0061] The LM algorithm is used to iteratively minimize the reprojection error, and the optimal parameters are output after the iteration converges. That is, the optimal camera intrinsic radial distortion parameters.

[0062] The beneficial effects of the present invention are: (1) The present invention designs a novel feature point extraction method based on the intersection of a circle and an auxiliary line. This method does not rely on the traditional circle center as a feature point, thus avoiding the problem of decreased calibration accuracy caused by eccentricity error;

[0063] (2) This invention proposes a method for dividing elliptical edge confusion points based on axial distance. This method effectively alleviates the problem of insufficient accuracy in identifying edge confusion points caused by changes in axial distance due to perspective distortion.

[0064] (3) This invention constructs an ellipse fitting method based on concentric geometric constraints. This method utilizes a shared center and inclination angle to optimize the parameters of the inner and outer circles during the ellipse fitting process, thereby improving the ellipse fitting quality.

[0065] (4) A gradient magnitude-weighted SVD line fitting method is proposed. This method improves the accuracy and robustness of line fitting by assigning higher weights to high gradient points, thereby suppressing the influence of low gradient points.

[0066] (5) The calibration scheme proposed in this invention, when applied to a palletizing robot system, can effectively resist industrial interference factors such as dust, changes in lighting, and material obstruction, significantly improving the overall calibration accuracy and exhibiting strong environmental robustness. The calibration board has a simple structure, is easy to deploy, and is highly adaptable to various palletizing operation scenarios.

[0067] In summary, this invention effectively improves calibration accuracy and stability; ablation experiments confirm that the three modules of axial distance constraint, concentric ellipse constraint, and gradient magnitude-weighted SVD linear fitting work together to improve calibration accuracy. Attached Figure Description

[0068] Figure 1 This is a schematic diagram of the process of this invention;

[0069] Figure 2 This is a schematic diagram of the calibration plate of the present invention;

[0070] Figure 3 This is a comparison of the impact of the number of calibration images on the intrinsic parameter calibration error using different methods in the embodiments of the present invention; wherein, (a) is the impact of the number of calibration images on the intrinsic parameter calibration error. f x (a) Effect of calibration error; (b) Effect of the number of calibration images on the intrinsic parameters f y The effect of calibration error; (c) the effect of the number of calibration images on the intrinsic parameters C x The effect of calibration error; (d) is the effect of the number of calibration images on the intrinsic parameters. C y Influence of calibration error;

[0071] Figure 4 This is a comparison of the impact of noise levels on the intrinsic parameter calibration error of different methods in the embodiments of the present invention; wherein, (a) is the effect of noise level on focal length parameter f x (a) Effect of calibration error; (b) Effect of noise level on internal parameters f y Calibration error influence; (c) Noise level effect on internal parameters C x Calibration error influence; (d) is the effect of noise level on internal parameters C y Calibration error has an impact. Detailed Implementation

[0072] To enable those skilled in the art to better understand the technical solution of the present invention, the present invention will be described in detail below with reference to the accompanying drawings and specific embodiments.

[0073] This embodiment discloses a camera calibration method suitable for palletizing robots, which specifically includes the following steps:

[0074] First, there's the calibration plate design. The circular calibration plate designed in this embodiment is as follows: Figure 2 As shown. This calibration plate consists of two concentric circles and six straight lines passing through the centers of the circles. The origin of the world coordinate system. o w Located at the center of the double circle, x w The axis to the left is positive. y w The axis pointing downwards is positive. z w The axis is perpendicular to the plane of the calibration plate and points outward toward the camera. Six auxiliary straight lines are designed, evenly distributed at angles of 0°, 30°, 60°, 90°, 120°, and 150° respectively. All lines pass through the center of the circle, and each line intersects the inner and outer circles at two points, thus forming a total of 24 geometric feature points.

[0075] The overall calibration process is as follows: Figure 1 As shown, the acquired image is first filtered, edge detected, and sub-pixel sets of ellipses and lines are extracted. Then, high-precision ellipse fitting is completed by combining axial distance separation and concentric constraints. Next, auxiliary lines are solved by angle clustering and gradient magnitude weighted SVD, and the intersection of the line and ellipse is extracted as feature points. Finally, the pixel-world coordinate matching relationship is constructed, and the minimum reprojection error is optimized based on LM to solve the camera intra-camera distortion coefficient.

[0076] During camera calibration, the calibration board is first photographed to obtain... Figure 2 middle o pix -uv Concentric circles in the coordinate system. Then, the intersection points of the solved circles and auxiliary lines are used as the calibrated feature points, thereby establishing the pixel coordinate system. o pix -uv Feature points under ) With the world coordinate system ( o w -x w y w z w ) click next The mapping relationship between them. Finally, based on the feature point matching pairs in the pixel coordinate system and the world coordinate system ( , Camera calibration is completed using the perspective projection model formula (1).

[0077] (1)

[0078] in, s As a scale factor, K For the camera intrinsic parameter matrix, R and t These represent the rotation matrix and translation, respectively.

[0079] Unlike existing methods, the concentric circle calibration plate designed in this embodiment no longer uses the center of the circle as the feature point during calibration. Instead, it uses the intersection of the circle and an auxiliary straight line as the feature point, effectively avoiding the impact of eccentricity error on calibration accuracy. The entire calibration solution process can be broken down into five execution stages, implemented step by step:

[0080] (1) Obtaining the set of points of the line and the ellipse in the calibration plate

[0081] (2) Fitting concentric ellipses in the calibration plate

[0082] (3) Auxiliary linear fitting in calibration plate

[0083] (4) Feature point matching for hierarchical indexes

[0084] (5) Camera intrinsic parameter calibration solution, the specific process is as follows: Figure 2 As shown.

[0085] I. Obtaining the points of the straight line and ellipse on the calibration plate

[0086] (1) First, convert the calibration plate image into a grayscale image. I ( x , y ), and then to I ( x , y Gaussian filtering is performed to obtain a smooth grayscale image. I g ( x , y ).

[0087] (1.1) The Gaussian kernel used is shown in formula (2), where i , j This represents the relative offset of the Gaussian kernel center. σ The standard deviation of the Gaussian function is used to control the smoothness of the filter.

[0088] (2)

[0089] (1.2) Using formula (3) G (i , j Gaussian kernel normalization processing

[0090] (3)

[0091] in R w The window radius of the Gaussian kernel is given by the standard deviation. σ Sure.

[0092] (1.3) Compare the normalized Gaussian kernel with the original grayscale image I ( x , y Perform a two-dimensional convolution operation to obtain a smooth image. I g ( x , y The convolution process is shown in formula (4):

[0093] (4)

[0094] The Canny edge detection algorithm is used to process smooth grayscale images. I g ( x , y Edge extraction is performed to obtain the constituent points of ellipses and auxiliary lines in the image.

[0095] (2.1) Calculate the smoothed grayscale image I g ( x , y arbitrary pixel coordinates in ) x , y gradient magnitude M ( x , y ) and gradient direction θ ( x , y ) 。

[0096] (2.1.1) For I g ( x , y arbitrary pixel coordinates in ) x , y Construct a 3x3 local window matrix centered on this pixel. M x,y The matrix is ​​composed of the gray values ​​of the neighborhood of the center pixel, which are offset by 1 unit above, below, left, and right, arranged in rows and columns. If the coordinates after the offset exceed the image boundary, the gray values ​​at the corresponding positions are filled with zeros.

[0097] (2.1.2) Pixel ( x , y ) horizontal gradient G x with vertical gradient G y Composed of Sobel convolution kernels and local windows M x,y The expression for two-dimensional convolution is shown in formula (5):

[0098] (5)

[0099] (2.1.3) Based on the above gradient components G x and G y The gradient magnitude is calculated using equation (6). M ( x , y ) and gradient direction θ ( x , y ):

[0100] (6)

[0101] (2.2) An edge point set is generated by using a dual threshold strategy combined with connected component analysis. P .

[0102] (2.2.1) Based on gradient magnitude M ( x , y Classify each pixel:

[0103] (7)

[0104] in, T H and T L These represent the high and low thresholds, respectively, and satisfy the following conditions: 0 < T L < T H Considering both edge preservation and noise reduction effects, the calibration threshold in this embodiment is set to... T H =150、 T L =50.

[0105] (2.2.2) Weak edge pixels that satisfy the 8-neighborhood connectivity relationship with strong edge pixels are retained to obtain the edge map. E ( x , y ) is defined as:

[0106] (8)

[0107] (2.2.3) By E ( x , y Extracting edge point sets P ={( x , y ) |E ( x , y )= 1} .

[0108] (2.3) For the edge point set P Gradient-based subpixel refinement is performed on the mid-edge points.

[0109] (2.3.1) Using the current edge point ( x , y Centered on ), three points are sampled at equal intervals along the normal direction of its edge, with gradient magnitudes of respectively. M 1 , M 2 , M 3 The sub-pixel offset is calculated by fitting a quadratic parabola to the gradient profile. δ :

[0110] (9)

[0111] (2.3.2) Set the edge points ( x , y Dividing the gradient component by the gradient magnitude to perform normalization yields the unit gradient normal vector. n x , n y The calculation formula is:

[0112] (10)

[0113] (2.3.3) Generate sub-pixel level edge point set P s , where each point p s =( x s , y s )∈ P s From pixels p =( x , y )∈ PThe result is obtained by correcting formula (11):

[0114] (11)

[0115] (2.4) Separate the elliptical and linear point sets of the edge point set P.

[0116] (2.4.1) Obtain the set of parameters that make up the line using the Hough line detection algorithm. L As shown in formula (12), where, ρ The perpendicular distance from the origin to the line is... θ L For the normal to the line and x The included angle of the axis, A ( ρ,θ L ) is a vote accumulator for Hough space.

[0117] (12)

[0118] (2.4.2) The original sub-pixel set P s Divided into a set of points on a straight line and elliptic point set As shown in formula (13). Wherein, This is the distance threshold.

[0119] (13)

[0120] II. Ellipse Fitting

[0121] (1) Based on the extracted set of circular edge points Achieve fitting and optimization of segmented ellipses.

[0122] (1.1) Based on point sets The segmented ellipse is fitted using the rotating ellipse geometric model shown in formula (14).

[0123] (14)

[0124] Among them, optimization variables Ω ={ x c ,y c ,a,b,θ} represents the geometric parameters of the ellipse, ( x , y )for The pixel coordinates of the midpoint; the objective function minimizes the sum of squared residuals of all edge points relative to the standard ellipse equation after rotation transformation; constraints. a >0,b >0 is used to ensure that the length of the semi-axis is positive, satisfying the physical meaning of the geometric parameters.

[0125] (1.2) The geometric parameters of the segmented ellipse can be obtained by iteratively solving the above model using the nonlinear least squares criterion. ,in( , () represents the coordinates of the center of the dividing ellipse. , These represent the lengths of the semi-major axis and the semi-minor axis, respectively. The rotation angle of the ellipse is the difference between the local coordinate system of the ellipse and the pixel coordinate system. x The angle between axes.

[0126] (2) Based on the axial distance, the elliptical edge confusion point is divided, and the point set of the ellipticity inside and outside the calibration plate is accurately established.

[0127] (2.1) Based on the center of the ellipse ( , and rotation angle ,Will edge points in p =( x , y Transform from pixel coordinates to elliptical local coordinates, and let the transformed coordinates be represented as follows: p =( x e , y e The axial distance is defined as shown in formula (15).

[0128] (15)

[0129] (2.2) The point set for calibrating the ellipticity of the inner and outer circles of the plate is divided according to the axial distance. The outer circle point set is: ={( x e , y e )| d_af ( x e , y e )>1}, inner circle set ={( x e , y e )| d_af ( x e , y e )<1}.

[0130] (3) Perform ellipse fitting based on concentric geometric constraints on the separated inner and outer point sets.

[0131] (3.1) Set the center of the segmented ellipse and the initial value of the rotation angle. As the initial values ​​of the center and rotation angle of the inner and outer ellipses, the point sets of the outer ellipse are respectively calculated using formula (14). With the inner elliptical point set Perform ellipse fitting to obtain the initial semi-axis length of the outer ellipse. Length of the initial semi-axis of the inner ellipse .

[0132] (3.2) An ellipse fitting optimization model based on the normalized residuals of ellipse fitting is established to establish the concentric geometric constraint relationship.

[0133] (3.2.1) The normalized geometric residuals of the ellipse fitting are given in formula (16).

[0134] (16)

[0135] in, For the set of points on the edge of the ellipse, For point The axial distance to the ellipse. The residuals mentioned above are defined as the sum of squared deviations of the axial distances of all fitted points from the theoretical value of 1, used to quantify the overall fit between the point set and the elliptical model.

[0136] (3.2.2) Establish an ellipse fitting optimization model with concentric geometric constraints based on the normalized residuals of the ellipse fitting as shown in formula (17). The joint parameters to be optimized are: .

[0137] (17)

[0138] in,( x c , y c The coordinates are the shared center coordinates of the inner and outer ellipses. θ For the shared rotation angle, and These represent the semi-major axis and semi-minor axis lengths of the outer and inner ellipses, respectively.

[0139] (3.2.3) The initial parameters for model optimization are set as follows: The LM algorithm is used to perform nonlinear least squares optimization on equation (17) and iteratively update the parameters. ε J Continue until convergence, and obtain the final parameters of the inner and outer ellipses.

[0140] (18)

[0141] III. Linear Fitting

[0142] (1) Clustering of points is achieved based on the point angles in the set of points on the straight line, and the point set is divided into 6 auxiliary straight lines.

[0143] (1.1) For the sub-pixel set of straight edge The point in ( The polar angle of the point relative to the center of the ellipse can be solved using equation (19). θ pol :

[0144] (19)

[0145] in, , For point set The Middle i Each point is relative to the center of the ellipse. The coordinate difference θ pol This is the polar angle in degrees, with a value range of [-180°, 180°].

[0146] To facilitate angle matching, the polar angle is uniformly mapped to [-0°, 360°) using formula (20):

[0147] (20)

[0148] (1.2) According to The point in ( )and θ std The angular differences are assigned to the corresponding lines, points ( The corresponding line number is This allows us to obtain the set of edge points for different straight lines. ,in k =1,2,...,6, N k This represents the number of points in each group.

[0149] (2) A line fitting method based on gradient magnitude weighted SVD was designed for the line point set. Fit a straight line.

[0150] (2.1) Regarding the first k edge point set of a straight line According to formula (21), calculate the weight of each point one by one. w k,i .

[0151] (twenty one)

[0152] In the formula, M ( ) is the first k The first straight line i Points The gradient magnitude.

[0153] (2.2) No. k Weighted centroid of a set of points on a straight line C k ( x , y The definition is shown in formula (22):

[0154] (twenty two)

[0155] (2.3) Regarding the first k Construct a weighted decentralized coordinate matrix using straight lines. Q k :

[0156] (twenty three)

[0157] in, For the first k The first straight line i The weighted centroid of each point relative to this line. C k ( x , y The decentralized coordinate vector of ).

[0158] (2.4) For the matrix Q k Perform singular value decomposition (SVD) and take the right singular vector corresponding to the minimum singular value. This serves as the unit normal vector for the fitted line. Finally, the... k The straight lines are: The parameters satisfy:

[0159] (twenty four)

[0160] IV. Matching of Calibration Feature Points

[0161] (1) Extraction and hierarchical indexing of feature points in the world coordinate system

[0162] (1.1) Let the radius of the inner circle of the calibration plate in the world coordinate system be... The outer radius is The calibration plate has six evenly distributed auxiliary straight lines whose direction angles satisfy... The world coordinates of any feature point on the annulus can be expressed by formula (25):

[0163] (25)

[0164] In the formula, the radius parameter r Desirable or These correspond to the feature points on the inner and outer circles, respectively.

[0165] (1.2) Following the rule of "inner circle first, then outer circle, and within the same ring, sort all feature points and construct an ordered set of world feature points. p w This provides a unified index benchmark for subsequent feature matching.

[0166] (2) Extraction and hierarchical indexing of feature points in pixel coordinate system

[0167] (2.1) The feature points in the pixel coordinate system are obtained by solving the equations of the ellipse and the line as shown in formula (26).

[0168] (26)

[0169] (2.2) Calculate the feature points obtained from formula (26) based on formula (19). x , y Relative to the center of the ellipse ( , The polar angle is calculated, and the polar angle is mapped to the interval [-0°, 360°) using equation (20). The mapped polar angle is... θ std An auxiliary line characterizing the affiliation of the feature point.

[0170] (2.3) Transform the feature points from the pixel coordinate system to the elliptical local coordinate system, and then calculate the axial distance according to the formula (15) above. . Belongs to the outer circle, The point belongs to the inner circle.

[0171] (2.4) Referencing the sorting rules of feature points in the world coordinate system, perform hierarchical indexing of feature points in the pixel coordinate system. First, divide the inner circle point set and the outer circle point set, and then sort them in ascending order of polar angle within the same circle point set to obtain the feature point set. p pix .

[0172] V. Matching of feature points in the world coordinate system and pixel coordinate system and solving camera intrinsic parameters

[0173] (1) The sorted pixel sequence set p pix World Feature Point Set pw Map feature points one-to-one by index to construct feature point matching pairs. , ).

[0174] (2) Match feature points to {( , World feature points within )} Transformed to the camera coordinate system, after normalization and second-order radial distortion mapping, and combined with the parameter vector to be optimized. Obtain the predicted pixel coordinates corresponding to the world point. Construct the least squares objective function for reprojection error and solve for the optimal parameters:

[0175] (27)

[0176] In the formula, Let i be the coordinates of the measured pixel. For the current parameter vector The calculated predicted pixel coordinates.

[0177] (3) The LM algorithm is used to iteratively minimize the reprojection error, and the optimal parameters are output after the iteration converges. That is, the optimal camera intrinsic radial distortion parameters.

[0178] To verify the accuracy of the calibration method in this embodiment, three methods—PLF, WGILA, and IACC [7,24,25]—were selected for comparison, and comparative experiments were conducted. To suppress random errors and ensure the reliability of the experimental results, the proposed method underwent 20 repeated calibration experiments under each set of working conditions, and the statistical mean was taken. The results were then compared with the calibration results of the comparative methods to analyze the overall calibration accuracy of the different methods.

[0179] To investigate the impact of the number of calibration images on the accuracy of intrinsic parameter calibration and to compare the adaptability of different methods to different sample numbers, the number of calibration images was set to 4-20, and corresponding calibration images were generated by shooting from multiple perspectives. The CircLine-Calib method in this embodiment is compared with three other methods in terms of camera intrinsic parameters. f x , f y , c x , c y The relative root mean square error (RRMSE) results are as follows: Figure 3 As shown.

[0180] Depend on Figure 3 It can be seen that as the number of calibration images gradually increases from 4 to 20, all methods solve the problem. f x ,f y , c x , c y The RRMSE of the internal parameters generally showed a decreasing and converging trend. For the focal length parameter... f x , f y ( Figure 4 (a) Figure 4 In comparison (b), the error reduction of each method in the early stage is drastic and the curve fluctuates significantly. It usually converges when the 12th calibration image is displayed. However, the method proposed in this embodiment not only has a very small error change when the number of calibration images varies from 4 to 20, but also achieves a very small error in each number of calibration images.

[0181] At the principal point coordinates c x , c y ( Figure 3 (c) Figure 3 In terms of solving the problem (d), the method in this embodiment has more prominent advantages. The RRMSE is significantly lower than that of the three comparison methods PLF, WGILA and IACC in the entire range of image quantity. The overall curve fluctuation is very small and the stability is extremely strong.

[0182] In small-sample calibration scenarios (4-8 images for calibration), the errors of various methods show significant divergence: PLF and WGILA have extremely high initial errors and slow convergence. Even compared to the IACC scheme, the method in this embodiment still exhibits faster convergence speed across these four parameters, demonstrating stronger adaptability for small-sample calibration. When the calibration sample is expanded to 20 images, the method in this embodiment... c x , c y The RRMSE values ​​were only 0.0075% and 0.0043%, respectively, which significantly improved the accuracy of principal point solving compared with other comparative methods.

[0183] To verify the contribution of each core module to the calibration accuracy, this embodiment selects 20 calibration images with optimal algorithm performance for ablation experiments, thereby avoiding the calibration accuracy fluctuation problem caused by insufficient image quantity. The CircLine-Calib calibration method proposed in this embodiment includes three core modules: gradient magnitude weighted SVD line fitting, concentric circle constraint, and axial distance constraint. This embodiment constructs multiple sets of CircLine-Calib variant models to verify the effectiveness of each module on the final calibration result. The specific ablation settings are as follows:

[0184] (1) CircLine-Calib w / o W-SVD: Remove weighted SVD line fitting strategy;

[0185] (2) CircLine-Calib w / o Circ-C: Remove concentric circle constraints;

[0186] (3) CircLine-Calib w / o Axial-C: Cancel axial distance constraint.

[0187] Table 1 presents the quantitative results of the ablation experiments in this embodiment (20 calibration images, Gaussian noise standard deviation σ = 1 pixel, each experiment repeated 20 times). The data in the table show that the complete CircLine-Calib method has the best accuracy in calibrating intrinsic and distortion parameters. f x , f y , c x , c y The RRMSE values ​​were only 0.1056%, 0.1060%, 0.0075%, and 0.0043%. k 1 , k 2 The RRMSE values ​​were 1.7000% and 4.0829%, respectively. After removing the weighted SVD line fitting module, the algorithm's calibration accuracy significantly decreased. f x , f y The RRMSE increased by 0.0366% and 0.0445% respectively compared to the complete method, and the principal point parameters... c x , c y Error synchronization increased slightly, and overall calibration stability decreased. The algorithm performance degraded most significantly after removing the concentric circle constraint. f x , f y The RRMSE increased by 0.1527% and 0.1669% respectively compared to the complete method. c x , c y The parameter errors increased significantly, and the distortion parameter errors also deteriorated substantially, having the greatest impact on the calibration results. After removing the axial distance constraint, the calibration errors of the distortion parameters involved in the camera also increased to varying degrees. f x , f y , cx , c y All of them exhibit a certain degree of error increase, resulting in a reduction in the robustness of the algorithms.

[0188] Table 1 Ablation Experiment Results

[0189]

[0190] Based on the above ablation results, the absence of any core module will reduce the calibration accuracy and noise resistance of the method in this embodiment. This proves that the weighted SVD linear fitting strategy, concentric circle constraint, and axial distance constraint are all indispensable key modules of the CircLine-Calib method, and play an important role in improving calibration accuracy and stability.

[0191] To compare the anti-interference capabilities of different methods, robustness tests were conducted using multi-gradient Gaussian noise. This embodiment selected 20 calibration images, applying Gaussian noise with a standard deviation of 0.0–2.0 px and a step size of 0.2 px to the feature point coordinates, for a total of 11 noise levels. The intrinsic parameter calibration accuracy of this embodiment was then compared with other comparative methods to quantify the differences in noise resistance performance among the algorithms. The experimental results are as follows: Figure 4 As shown.

[0192] Depend on Figure 4 As can be seen, as the standard deviation of Gaussian noise at feature points gradually increases, the RRMSE of the intrinsic parameters of all comparison methods generally shows an upward trend.

[0193] For focal length parameter f x , f y ( Figure 4 (a) Figure 4 (b) In the PLF and WGILA methods, the error increases almost sharply with increasing noise, and the error increases significantly in the high-noise section. The error of the IACC algorithm fluctuates significantly, with the overall fluctuation range being moderate. In the low-noise range (0~0.4px), the error of the method in this embodiment is slightly higher than that of PLF and WGILA. However, after the noise intensity exceeds 0.4px, the error of the method in this embodiment remains the lowest globally, and the overall curve trend is the smoothest. Even when facing high-intensity noise disturbances, it can still output stably, and its noise robustness is better than the other comparative methods.

[0194] In principal point parameters c x , c y The solution results ( Figure 4 (c) Figure 4In the method described in this embodiment, the error remains consistently low across the entire noise gradient range, with minimal curve fluctuations. When the standard deviation of the Gaussian noise reaches 2.0px, c x The relative error is only 0.026%. c y The relative error is only 0.023%. In contrast, the principal point error of the three comparison methods, PLF, WGILA and IACC, continues to rise with the increase of noise. In particular, the IACC method fluctuates violently and the local peak is obvious. This further verifies that the method of this embodiment has strong noise resistance and robustness at the principal point parameter level.

[0195] The above embodiments are merely exemplary embodiments of the present invention and are not intended to limit the present invention. Those skilled in the art can make various modifications or equivalent substitutions to the present invention within its scope and spirit, and such modifications or equivalent substitutions should also be considered to fall within the scope of protection of the present invention.

Claims

1. A camera calibration method suitable for palletizing robots, characterized in that, Includes the following steps: Step 1: The calibration plate consists of two concentric circles and six auxiliary straight lines passing through the center of the circles, evenly distributed along the circumference of the concentric circles. Step 2: Obtain the calibration board image and preprocess the image to extract the sub-pixel sets of ellipses and lines; Step 3: Combine the axial distance separation strategy with concentric constraints to complete the concentric ellipse fitting; Step 4: Solve for the auxiliary line using angle clustering and gradient magnitude-weighted SVD; Step 5: Extract the intersection points of the auxiliary line and the concentric ellipse as feature points, and establish the mapping relationship between feature points in the pixel coordinate system and points in the world coordinate system; Step 6: Based on the feature point matching pairs in the pixel coordinate system and the world coordinate system, complete the camera calibration using the perspective projection model; Step 3 specifically involves: Step 3.1: Based on the elliptical sub-pixel point set obtained in Step 2, the segmented ellipse is fitted using a rotated ellipse geometric model. The rotated ellipse geometric model is iteratively solved using a nonlinear least squares criterion to obtain the geometric parameters of the segmented ellipse. ,in( , () represents the coordinates of the center of the dividing ellipse. , These represent the lengths of the semi-major axis and the semi-minor axis, respectively. The rotation angle of the ellipse is the difference between the local coordinate system of the ellipse and the pixel coordinate system. x The angle between the axes; Step 3.2, based on the coordinates of the center of the segmented ellipse and ellipse rotation angle The edge points in the sub-pixel point set of the ellipse are transformed from pixel coordinates to the local coordinate system of the ellipse, the axial distance is calculated, and the point set of the ellipse inside and outside the calibration plate is divided according to the axial distance. Step 3.3: Set the initial values ​​of the center and rotation angle of the segmented ellipse. Using the initial values ​​of the center and rotation angle of the inner and outer ellipses, the point sets of the inner and outer ellipses obtained in step 3.2 are fitted with ellipses using the geometric model of the rotated ellipse to obtain the initial semi-axis length of the outer ellipse. Length of the initial semi-axis of the inner ellipse An ellipse fitting optimization model based on the normalized geometric residuals of ellipse fitting is established to establish concentric geometric constraint relationships. The ellipse fitting optimization model based on the normalized geometric residuals of ellipse fitting in step 3.3 specifically establishes the concentric geometric constraint relationship: The normalized geometric residual of the ellipse fitting is: ； in, For an elliptical sub-pixel set, For point The axial distance to the ellipse; Ellipse fitting optimization model that establishes concentric geometric constraints using normalized residuals. The joint parameters to be optimized are ,in, The center coordinates are shared by the inner and outer ellipses. θ For the shared rotation angle, and These represent the lengths of the semi-major axis and semi-minor axis of the outer ellipse and the inner ellipse, respectively. The initial parameters for model optimization are set to The LM algorithm is used to perform nonlinear least squares optimization on the ellipse fitting optimization model, and the parameters are updated iteratively. ε J Until convergence, the final parameters of the inner and outer ellipses are obtained: 。 2. The camera calibration method for palletizing robots according to claim 1, characterized in that, In step 1, the auxiliary lines are evenly distributed at angles of 0°, 30°, 60°, 90°, 120°, and 150° respectively; all lines pass through the center of the circle, and each line intersects the inner and outer circles at two points, forming 24 feature points with the concentric circles; the origin of the world coordinate system... o w Located at the center of concentric circles, x w The axis to the left is positive. y w The axis pointing downwards is positive. z w The axis perpendicular to the plane of the calibration plate and pointing outwards toward the camera is positive.

3. The camera calibration method for palletizing robots according to claim 1, characterized in that, Step 2 specifically involves: Step 2.1: After converting the calibration board image into a grayscale image, Gaussian filtering and edge detection are performed to obtain the edge map. Step 2.2: Extract the edge point set from the edge map obtained in Step 2.1; Step 2.3: Perform gradient-based sub-pixel refinement on the edge points in the edge point set obtained in Step 2.2; Step 2.4: Separate the elliptical and linear point sets from the middle edge points of the sub-pixel refined edge point set obtained in Step 2.

3.

4. The camera calibration method for palletizing robots according to claim 1, characterized in that, The specific formula for the geometric model of the rotating ellipse is as follows: ； Among them, optimization variables For the geometric parameters of the ellipse, ( x , y ) represents the pixel coordinates of the sub-pixel point concentration point of the ellipse.

5. A camera calibration method for palletizing robots according to claim 1, characterized in that, Step 4 specifically involves: Step 4.1, for the points in the sub-pixel set of the straight line obtained in Step 2 ( Find the polar angle of this point relative to the center of the ellipse. : ； in, , For the sub-pixel set of straight lines i Each point is relative to the center of the ellipse. The coordinate difference The polar angle is measured in degrees, with a value ranging from -180° to 180°. Mapping to [-0°, 360°) yields Based on the concentration point of sub-pixel points of the straight line and The angular differences are assigned to the corresponding straight lines and points. The corresponding line number is This yields the edge point sets of different straight lines. ,in k =1,2,...,6, N k The number of points in each group; Step 4.2, regarding the first k edge point set of a straight line Calculate the weight of each point one by one. w k,i : ； in, For the first k The first straight line i Points The gradient magnitude; No. k Weighted centroid of a set of points on a straight line C k ( x , y )for: ； Regarding the first k Construct a weighted decentralized coordinate matrix using straight lines. Q k : ； in, For the first k The first straight line i The weighted centroid of each point relative to this line. C k ( x , y The decentralized coordinate vector of ). For matrix Q k Perform singular value decomposition and take the right singular vector corresponding to the minimum singular value. As the unit normal vector of the fitted line, the final... k The straight lines are: The parameters satisfy: .

6. A camera calibration method for palletizing robots according to claim 5, characterized in that, Step 5 specifically involves: Step 5.1, let the radius of the inner circle of the calibration plate in the world coordinate system be... The outer radius is The target is marked with six evenly distributed auxiliary straight lines, whose orientation angles satisfy... θ k =30•( k -1),( k =1,2,...,6), the world coordinates of any feature point on the annulus are represented as: ； In the formula, the radius parameter r Desirable or These correspond to feature points on the inner and outer circles, respectively. Following the rule of sorting all feature points from the inner circle to the outer circle, and from smallest to largest direction angle within the same annulus, an ordered world feature point set p is constructed. w This provides a unified index benchmark for subsequent feature matching; Step 5.2: The feature points in the pixel coordinate system are obtained by solving the equations of the ellipse and the line simultaneously. The polar angle of the feature points in the pixel coordinate system relative to the center of the ellipse is calculated, and the polar angle is mapped to the interval [-0°, 360°). The feature points are then transformed from the pixel coordinate system to the local elliptical coordinate system, and the axial distance is calculated. Belongs to the outer circle, The points belong to the inner circle; following the rule of first the inner circle, then the outer circle, and then sorting the points within the same annulus by direction angle from smallest to largest, the feature points in the pixel coordinate system are indexed hierarchically. First, the inner circle point set and the outer circle point set are divided, and then the points within the same annulus point set are sorted by polar angle from smallest to largest to obtain the feature point set. p pix ; Step 5.3: Sort the set of pixel sequences. p pix World Feature Point Set p w Map feature points one-to-one by index to construct feature point matching pairs. , ).

7. A camera calibration method for palletizing robots according to claim 1, characterized in that, Step 6 specifically involves: Based on the extrinsic parameters, the feature points are matched to the {( , World feature points within )} Transformed to the camera coordinate system, after normalization and second-order radial distortion mapping, and combined with the parameter vector to be optimized. Obtain the predicted pixel coordinates corresponding to the world feature points. Construct the least squares objective function for reprojection error and solve for the optimal parameters: ； in, Let i be the coordinates of the measured pixel. For the current parameter vector The calculated predicted pixel coordinates; The LM algorithm is used to iteratively minimize the reprojection error, and the optimal parameters are output after the iteration converges. That is, the optimal camera intrinsic radial distortion parameters.

Images