Highway traffic structural body three-dimensional geometric quantity multi-layer calibration method and system based on target reference

Abstract

The invention belongs to the field of highway traffic geometric quantity measurement, and discloses a highway traffic structural body three-dimensional geometric quantity multilayer calibration method and system based on target reference, and the method comprises the steps: collecting point cloud data of a highway infrastructure structural body through a three-dimensional laser scanner, and carrying out the denoising preprocessing; on the basis of the point cloud data after denoising preprocessing, plane errors of the three-dimensional laser scanning equipment are calibrated by using a standard plane target; on the basis of the denoised point cloud data, a standard spherical target is used to calibrate the stereo error of the three-dimensional laser scanning device; based on the calibration results of the plane error and the stereo error of the three-dimensional laser scanning equipment, the coordinates of different stations are converted through a four-parameter method, and the same-name area registration is carried out by using the maximum geometric error. The reliability, the accuracy and the adaptability of the highway infrastructure scanning equipment are improved.

Description

Technical Field

[0001] The present invention belongs to the technical field of highway traffic geometry measurement, and in particular relates to a target reference-based multi-layer calibration method and system for three-dimensional geometry of highway traffic structures. Background Art

[0002] In the field of geometric measurement, existing technologies often use target-based error calibration methods, but the following key issues still exist:

[0003] 1. Severe noise interference:

[0004] Traditional denoising methods (such as mean filtering) have limited ability to remove global random noise and cannot deeply integrate with the spatial characteristics of 3D point clouds, resulting in a high rate of false deletion of valid data.

[0005] 2. Insufficient plane calibration accuracy:

[0006] Traditional plane fitting algorithms, such as the least squares method, are sensitive to the density distribution of point clouds. In non-uniform sampling scenarios, the geometric center identification error cannot meet the requirements of high-precision target calibration.

[0007] 3. Limitations of Stereo Parameter Optimization:

[0008] Existing spherical calibration methods (such as the least squares method) only optimize the single-direction residual and do not integrate multi-dimensional gradient information, resulting in significant measurement deviations between the sphere center and radius (typical error ≥ 3mm).

[0009] 4. Algorithm efficiency bottleneck:

[0010] Traditional iterative algorithms require multiple random samplings, have high computational complexity, and are difficult to support real-time processing of large-scale point clouds.

[0011] These issues hinder the development of high-precision geometric measurement technology, necessitating a calibration system that balances theoretical rigor with engineering efficiency. Therefore, research on multi-layer calibration technology with targets for 3D geometric measurement of highway structures is crucial for improving the reliability, accuracy, and adaptability of highway infrastructure scanning equipment. Summary of the Invention

[0012] To solve the problems existing in the prior art, the present invention provides a multi-layer calibration method and system for the three-dimensional geometric quantities of highway traffic structures based on target reference, so as to solve the problem of the lack of calibration technology that takes into account both theoretical rigor and engineering efficiency in the industry, and improve the reliability, accuracy and adaptability of highway infrastructure scanning equipment.

[0013] To achieve the above object, the present invention provides the following solutions:

[0014] A target-referenced multi-layer calibration method for three-dimensional geometric quantities of a highway traffic structure, the method comprising:

[0015] Use 3D laser scanners to collect point cloud data of highway infrastructure structures and perform denoising preprocessing;

[0016] Based on the denoised point cloud data, the plane error of the 3D laser scanning equipment is calibrated using a standard plane target.

[0017] Based on the denoised point cloud data, a standard spherical target is used to calibrate the stereo error of the 3D laser scanning equipment.

[0018] Based on the calibration results of the planar error and stereo error of the 3D laser scanning equipment, the coordinates of different sites are transformed using the four-parameter method, and the maximum geometric error is used to perform homonymous region registration.

[0019] Preferably, the point cloud data of the highway infrastructure structure is collected using a 3D laser scanner and subjected to denoising preprocessing, including:

[0020] Convert the three-dimensional point cloud coordinates to the Euclidean distance of the origin to construct a one-dimensional data set;

[0021] Calculate the mean of a one-dimensional data set, Z-score or Z i ;

[0022] According to the calculation results of Z-score, set the dynamic threshold;

[0023] Eliminate the ones that satisfy |Z i ∣>Threshold point.

[0024] Preferably, based on the point cloud data after denoising preprocessing, the plane error of the three-dimensional laser scanning device is calibrated using a standard plane target, including:

[0025] Get the area of ​​the area to be measured;

[0026] Calculating the point cloud density distribution, and calculating the neighborhood radius according to the area of ​​the area to be measured and the point cloud density distribution;

[0027] Calculate the initial value of the geometric center according to the neighborhood radius;

[0028] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the calibration of the plane error of the 3D laser scanning equipment;

[0029] The area of ​​the area to be measured is obtained, including:

[0030]

[0031] Calculate point cloud density distribution, including:

[0032]

[0033] Where N r is the actual number of sampling points in the area to be measured, which is determined according to the on-site measurement requirements; N tot is the theoretical number of sampling points; λ is the empirical parameter of point cloud density; r is the neighborhood radius; D is the point cloud density distribution function;

[0034] Calculating a neighborhood radius according to the area of ​​the region to be measured and the point cloud density distribution includes:

[0035]

[0036] Calculate the initial value of the geometric center based on the neighborhood radius, including:

[0037] Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0):

[0038]

[0039] Where N i is the number of sampling points in the ith sub-region within the neighborhood, and n is the number of sub-regions;

[0040] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the plane error calibration of the point cloud data, including:

[0041] Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i :

[0042]

[0043] Where, is the temporary geometric center coordinate of the kth iteration;

[0044] The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm:

[0045]

[0046] The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than the threshold ∈.

[0047] Preferably, based on the point cloud data after denoising preprocessing, the stereo error of the three-dimensional laser scanning device is calibrated using a standard spherical target, including:

[0048] Build a spherical model;

[0049] Based on the spherical model, construct a loss function;

[0050] Based on the loss function, calculate the gradient and Hessian matrix;

[0051] Based on the gradient and the Hessian matrix, iteratively updating parameters is performed to complete the calibration of the stereo error of the three-dimensional laser scanning device;

[0052] Among them, establishing a spherical model includes:

[0053] The default spherical equation is:

[0054] (xa) 2 +(yb) 2 +(zc) 2 =r 2 ;

[0055] Where (x, y, z) are the coordinates of the point on the preset sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius;

[0056] Based on the spherical model, a loss function is constructed, including:

[0057] Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is:

[0058]

[0059] Based on the loss function, the gradient and Hessian matrix are calculated, including:

[0060] Gradient vector:

[0061]

[0062] Hessian matrix:

[0063]

[0064] Based on the gradient and Hessian matrix, the parameters are iteratively updated, including:

[0065] parameter vector The update formula is:

[0066]

[0067] Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector;

[0068] The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

[0069] Preferably, based on the calibration results of the plane error and stereo error of the three-dimensional laser scanning device, the coordinates of different sites are converted by a four-parameter method, and the maximum geometric error is used to perform homonymous region registration, including:

[0070] Calibrate the four-parameter coordinates of point cloud data;

[0071] Use NDT method to optimize calibration results;

[0072] Based on the optimized structure, the regions with the same name in the point cloud data are fixed;

[0073] Among them, the four-parameter coordinates of the point cloud data are calibrated, including:

[0074] The preset point pairs of the same name between the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i ); plane coordinate calibration is performed using the four-parameter method, and the calculation formula is:

[0075]

[0076] Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount;

[0077] Optimize calibration results using NDT methods, including:

[0078] Divide the target point cloud into a fixed-size cubic grid;

[0079] Calculate the mean μ and covariance matrix Γ of the point cloud in each grid and build a probability density model;

[0080]

[0081] Based on the probability density model, the objective is optimized by maximizing the overall likelihood function;

[0082]

[0083] Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold:

[0084]

[0085] Where M k+1 and M krepresents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function, λ is the damping factor of the algorithm, which is used to balance the gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to the gradient descent, and when λ is small, the algorithm is close to the Gauss-Newton method. I is the bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems.

[0086] Based on the optimized structure, the regions with the same name in the point cloud data are fixed, including:

[0087] Calculate the maximum Euclidean distance between the source point cloud and the target point cloud:

[0088]

[0089] Threshold setting: filter to meet d i ≤ε, where the threshold

[0090] Use the KDTree structure of the Open3D library to match point clouds.

[0091] The present invention also provides a target-referenced, multi-layer calibration system for three-dimensional geometric quantities of highway traffic structures. The system is used to implement the aforementioned method. The system includes: a preprocessing module, a plane error calibration module, a stereo error calibration module, and a region matching module.

[0092] The pre-processing module is used to collect point cloud data of highway infrastructure structures using a three-dimensional laser scanner and perform denoising pre-processing;

[0093] The plane error calibration module is used to calibrate the plane error of the three-dimensional laser scanning device using a standard plane target based on the point cloud data after denoising preprocessing;

[0094] The stereo error calibration module is used to calibrate the stereo error of the three-dimensional laser scanning device using a standard spherical target based on the point cloud data after denoising preprocessing;

[0095] The region matching module is used to convert the coordinates of different sites using a four-parameter method based on the calibration results of the plane error and stereo error of the three-dimensional laser scanning device, and to perform homonymous region registration using the maximum geometric error.

[0096] Preferably, the point cloud data of the highway infrastructure structure is collected using a 3D laser scanner and subjected to denoising preprocessing, including:

[0097] Convert the three-dimensional point cloud coordinates to the Euclidean distance of the origin to construct a one-dimensional data set;

[0098] Calculate the mean of a one-dimensional data set, Z-score or Z i ;

[0099] According to the calculation results of Z-score, set the dynamic threshold;

[0100] Eliminate the ones that satisfy |Z i ∣>Threshold point.

[0101] Preferably, based on the point cloud data after denoising preprocessing, the plane error of the three-dimensional laser scanning device is calibrated using a standard plane target, including:

[0102] Get the area of ​​the area to be measured;

[0103] Calculating the point cloud density distribution, and calculating the neighborhood radius according to the area of ​​the area to be measured and the point cloud density distribution;

[0104] Calculate the initial value of the geometric center according to the neighborhood radius;

[0105] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the calibration of the plane error of the 3D laser scanning equipment;

[0106] The area of ​​the area to be measured is obtained, including:

[0107]

[0108] Calculate point cloud density distribution, including:

[0109]

[0110] Where N r is the actual number of sampling points in the area to be measured, which is determined according to the on-site measurement requirements; N tot is the theoretical number of sampling points; λ is the empirical parameter of point cloud density; r is the neighborhood radius; D is the point cloud density distribution function;

[0111] Calculating a neighborhood radius according to the area of ​​the region to be measured and the point cloud density distribution includes:

[0112]

[0113] Calculate the initial value of the geometric center based on the neighborhood radius, including:

[0114] Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0):

[0115]

[0116] Where N iis the number of sampling points in the ith sub-region within the neighborhood, and n is the number of sub-regions;

[0117] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the plane error calibration of the point cloud data, including:

[0118] Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i :

[0119]

[0120] Where, is the temporary geometric center coordinate of the kth iteration;

[0121] The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm:

[0122]

[0123] The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than the threshold ∈.

[0124] Preferably, based on the point cloud data after denoising preprocessing, the stereo error of the three-dimensional laser scanning device is calibrated using a standard spherical target, including:

[0125] Build a spherical model;

[0126] Based on the spherical model, construct a loss function;

[0127] Based on the loss function, calculate the gradient and Hessian matrix;

[0128] Based on the gradient and the Hessian matrix, iteratively updating parameters is performed to complete the calibration of the stereo error of the three-dimensional laser scanning device;

[0129] Among them, establishing a spherical model includes:

[0130] The default spherical equation is:

[0131] (xa) 2 +(yb) 2 +(zc) 2 =r 2 ;

[0132] Where (x, y, z) are the coordinates of the point on the preset sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius;

[0133] Based on the spherical model, a loss function is constructed, including:

[0134] Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is:

[0135]

[0136] Based on the loss function, the gradient and Hessian matrix are calculated, including:

[0137] Gradient vector:

[0138]

[0139] Hessian matrix:

[0140]

[0141] Based on the gradient and Hessian matrix, the parameters are iteratively updated, including:

[0142] parameter vector The update formula is:

[0143]

[0144] Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector;

[0145] The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

[0146] Preferably, based on the calibration results of the plane error and stereo error of the three-dimensional laser scanning device, the coordinates of different sites are converted by a four-parameter method, and the maximum geometric error is used to perform homonymous region registration, including:

[0147] Calibrate the four-parameter coordinates of point cloud data;

[0148] Use NDT method to optimize calibration results;

[0149] Based on the optimized structure, the regions with the same name in the point cloud data are fixed;

[0150] Among them, the four-parameter coordinates of the point cloud data are calibrated, including:

[0151] The preset point pairs of the same name between the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i); plane coordinate calibration is performed using the four-parameter method, and the calculation formula is:

[0152]

[0153] Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount;

[0154] Optimize calibration results using NDT methods, including:

[0155] Divide the target point cloud into a fixed-size cubic grid;

[0156] Calculate the mean μ and covariance matrix Γ of the point cloud in each grid and build a probability density model;

[0157]

[0158] Based on the probability density model, the objective is optimized by maximizing the overall likelihood function;

[0159]

[0160] Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold:

[0161]

[0162] Where M k+1 and M k represents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function, λ is the damping factor of the algorithm, which is used to balance the gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to the gradient descent, and when λ is small, the algorithm is close to the Gauss-Newton method. I is the bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems.

[0163] Based on the optimized structure, the regions with the same name in the point cloud data are fixed, including:

[0164] Calculate the maximum Euclidean distance between the source point cloud and the target point cloud:

[0165]

[0166] Threshold setting: filter to meet d i ≤ε, where the threshold

[0167] Use the KDTree structure of the Open3D library to match point clouds.

[0168] Compared with the prior art, the present invention has the following beneficial effects:

[0169] 1. Efficient global noise removal:

[0170] The dynamic threshold method based on Z-score, combined with the statistical characteristics of 3D Euclidean distance, improves the accuracy of outlier removal to 95% while preserving the integrity of the main distribution point cloud;

[0171] The time complexity has been effectively optimized, supporting the processing of millions of point clouds per second.

[0172] 2. High-precision plane calibration:

[0173] The regional adaptive weighted mean algorithm dynamically adjusts the neighborhood radius according to the target area, and the plane geometric center error is stabilized within 0.1mm;

[0174] The Newton iteration method is integrated with the Hessian matrix to optimize the spherical parameters. The center error is ≤0.5mm and the radius error is ≤0.05mm, which is 80% higher than the traditional method.

[0175] 3. Innovation in Multidimensional Space Theory:

[0176] A multi-level calibration framework based on the four-parameter method and NDT is proposed to achieve joint optimization of coordinate system transformation and point cloud registration, reducing the global geometric error by 73.6%.

[0177] Introducing dynamic weight function (ω=1 / (1+d 2 )), suppress the interference of outliers on the calculation of the geometric center.

[0178] 4. Standardization and universality:

[0179] The target design complies with the standards of the highway traffic industry and is suitable for various geometric models such as spheres and planes.

[0180] The algorithm is compatible with geometric measurement of infrastructure such as bridges and tunnels and has strong scalability. BRIEF DESCRIPTION OF THE DRAWINGS

[0181] In order to more clearly illustrate the technical solution of the present invention, the following briefly introduces the drawings required for use in the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0182] Figure 1 This is a flow chart of a multi-layer calibration method for three-dimensional geometric quantities of a highway traffic structure based on target reference according to an embodiment of the present invention;

[0183] Figure 2 Schematic diagram of the structure of a multi-layer calibration system for three-dimensional geometric quantities of highway traffic structures based on target reference according to an embodiment of the present invention. DETAILED DESCRIPTION

[0184] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0185] In order to make the above-mentioned objects, features and advantages of the present invention more obvious and easy to understand, the present invention is further described in detail below with reference to the accompanying drawings and specific embodiments.

[0186] Example 1

[0187] like Figure 1 As shown, this embodiment provides a target-based multi-layer calibration method for three-dimensional geometric measurement of highway traffic structures, including the following steps:

[0188] Step 1: Point cloud data preprocessing. A denoising method for point cloud data collected by a 3D laser scanner is proposed, and the Z-score method is used to improve data quality.

[0189] Step 2: Point cloud data plane error calibration. Based on the denoised point cloud data, the plane error of the 3D laser scanning equipment is calibrated using a standard plane target. The weighted mean algorithm is improved through regional adaptive optimization to improve the accuracy of geometric center recognition and reduce the plane error of the point cloud data.

[0190] Step 3: Stereoscopic error calibration of point cloud data. To address the issue of insufficient 3D geometric measurement accuracy of the scanning device, a standard spherical target is used to calibrate the stereoscopic error of the 3D laser scanning device. The Newton iteration algorithm is used to iterate the parameter vector multiple times to improve the data accuracy of the 3D point cloud of the device.

[0191] Step 4: Calibrate point cloud data from multi-site acquisition equipment. To address data mismatches and low accuracy from multi-site scanning equipment, the coordinates of different sites are converted using the four-parameter method. The maximum geometric error is used for homonymous region registration. NDT optimization is used to calibrate the geometric errors of the point cloud data, and algorithmic compensation is used to improve point cloud data accuracy.

[0192] In this embodiment, point cloud data preprocessing: global point cloud random noise preprocessing based on Z-score: a denoising method for point cloud data collected by a three-dimensional laser scanner is proposed, and abnormal points and noise points are eliminated through the Z-score method to improve data quality.

[0193] 1. 3D point cloud coordinate conversion

[0194] Convert the three-dimensional point cloud coordinates to the Euclidean distance to the origin and construct a one-dimensional data set:

[0195]

[0196] 2 Statistical parameter calculation

[0197] Mean:

[0198] Standard Deviation:

[0199] Z-score calculation:

[0200] 3 Dynamic threshold setting

[0201] Initial threshold: According to the 3σ principle of normal distribution, the initial threshold is set to 3.

[0202] Dynamic adjustment: Gradually lower the threshold according to the actual elimination effect (taking into account the distribution characteristics of the data set and the actual application status of the road traffic industry, the value is generally 2.5, and it is dynamically adjusted downward according to the elimination effect of data anomalies) until the proportion of anomalies stabilizes.

[0203] 4. Outlier removal

[0204] Eliminate the ones that satisfy |Z i ∣>Threshold point, retain the main distribution data, where Z i The Z-score value of each point.

[0205] 5. Embodiment

[0206] Step 1: Input simulated point cloud data with noise

[0207] The coordinate values ​​of the simulated 6 points are as follows:

[0208] {(1,2,3),(2,3,3),(3,4,3),(4,5,3),(5,6,3),(15,15,15)}

[0209] Step 2: Calculate the Euclidean distance from each point to the origin

[0210] d i ={3.74,4.69,5.83,7.07,8.37,25.98}

[0211] Step 3: Calculate the mean and standard deviation

[0212] μ=9.28,σ=7.61

[0213] Step 4: Calculate the Z-score and remove outliers

[0214] Z i ={-0.73,-0.60,-0.45,-0.29,-0.12,+2.19}

[0215] The threshold is set to 2.0, and point 6 (Z6=2.14) is determined to be an abnormal point and is removed.

[0216] 6 Technical Effects

[0217] Global denoising: Applicable to random noise (such as sensor speckle) and retains the main distribution point cloud.

[0218] Dynamic adaptability: The threshold can be dynamically optimized based on data distribution to avoid excessive deletion of valid data.

[0219] Computational efficiency: low time complexity, supporting real-time processing of large-scale point clouds.

[0220] In this embodiment, point cloud data plane error calibration: plane error calibration method based on regional adaptive improved weighted mean: based on the point cloud data after denoising preprocessing, the plane error of the three-dimensional laser scanning equipment is calibrated using a standard plane target, and the regional adaptive improved weighted mean algorithm is used to improve the accuracy of geometric center recognition and reduce the plane error of the point cloud data.

[0221] 1. Regional classification

[0222] According to the actual application of common structures to be tested in the highway industry, refer to the "JGJ8 Building Deformation Measurement Specification",

[0223] "DB42 / T 1496-2019 Technical Specification for Highway Slope Monitoring", "JT / T 1037-2022 Technical Specification for Highway Bridge Structure Monitoring", "JTG F80 Highway Engineering Quality Inspection and Assessment Standard", "GB / T 39410-2020 General Specification for Low-Earth Orbit Satellite-Borne GNSS Measurement Receivers", etc. The geometric dimensions of the structure monitoring refer to JTGD60, JTG3362, JTG / TD65-05, JTG / T3365-01, JTG / TD65-06, JTG / T3360-01, JTG / TD65-05, JTG5120-2021, JTG / TH21 and other specifications. The area of ​​the area to be measured S u It is divided into four types: small, medium, large and extra large. The specific classification standards are as follows:

[0224]

[0225] Region classification provides a basis for the subsequent adaptive adjustment of neighborhood radius and point cloud density.

[0226] 2. Point cloud density distribution and weighted neighborhood radius calculation

[0227] Define the point cloud density distribution function D as the basis for determining the neighborhood radius r, then:

[0228]

[0229] Where N r The actual number of sampling points in the area to be measured is determined based on on-site measurement requirements;

[0230] N tot is the theoretical number of sampling points;

[0231] λ is an empirical parameter for point cloud density. For 3D scanning of highway infrastructure structures, the recommended value range is 0.79-0.83.

[0232] Combined with regional classification, the calculation formula of neighborhood radius is:

[0233]

[0234] 3 Calculation of initial value of geometric center

[0235] Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0):

[0236]

[0237] Where N i is the number of sampling points in the ith sub-region in the neighborhood, and n is the number of sub-regions.

[0238] 4 Weight calculation and iterative optimization

[0239] Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i :

[0240]

[0241] Where, is the temporary geometric center coordinate of the kth iteration.

[0242] The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm:

[0243]

[0244] The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than a threshold ∈ (recommended ∈ = 0.1 mm).

[0245] 5 Examples

[0246] Step 1: Design of standard target

[0247] The design standard area is 300 square centimeters, and the circular plane target with standard geometric center coordinates of (0mm, 0mm) must meet the requirements of the "JT / T 1037-2022" specification.

[0248] When only the traditional scanning method is used without plane error calibration, the target geometric center positioning error is measured as:

[0249]

[0250] The coordinate positioning error is:

[0251]

[0252] Step 2: Neighborhood radius calculation

[0253] Area Type: Medium (S u =300cm 2 )

[0254] Parameter selection: λ = 0.81

[0255] Point cloud density distribution function:

[0256] Optimal neighborhood radius:

[0257] Step 3: Calculate the initial value of the geometric center

[0258]

[0259] Step 4: Iterative optimization calibration of plane geometric center error

[0260] The weighted mean algorithm is used to calibrate the plane geometric center error. The number of iterations is limited to 20 and the iteration threshold is 0.1 mm. The core code snippet is as follows:

[0261]

[0262]

[0263] In this embodiment, the final number of iterations is 5, and the detailed iteration process is as follows:

[0264] Iteration 1: X = 0.2512 mm, Y = 0.2476 mm,

[0265] Iteration 2: X = 0.1589 mm, Y = 0.1534 mm

[0266] Iteration 3: X = 0.1135 mm, Y = 0.1128 mm

[0267] Iteration 4: X = 0.1041 mm, Y = 0.1023 mm,

[0268] Iteration 5: X = 0.1023 mm, Y = 0.0987 mm

[0269] Final coordinates: (0.1023mm, 0.0987mm)

[0270] The final coordinate error is 0.1421mm, achieving sub-millimeter precision calibration and meeting the accuracy requirements of the highway transportation industry.

[0271] 6 Technical Effects

[0272] Regional adaptation: Dynamically adjust the neighborhood radius to adapt to the point cloud distribution characteristics of different areas;

[0273] Standardized implementation: Design targets in accordance with industry specifications to ensure repeatability and engineering applicability of calibration results;

[0274] Efficiency optimization: The iterative convergence speed is fast and the accuracy can reach millimeter level.

[0275] In this embodiment, the stereo error calibration of point cloud data is performed using a Newton iteration algorithm for laser scanning equipment. This method addresses the issue of insufficient 3D feature detection accuracy in laser scanning equipment by using a standard spherical target to calibrate the device's stereo error. The Newton iteration algorithm is then used to iterate the parameter vector multiple times to improve the accuracy of the device's 3D point cloud data. Before actual deployment, the spherical target is used to calibrate the 3D scanning device's stereo error to improve its accuracy.

[0276] 1. Spherical model establishment

[0277] Define the spherical equation as:

[0278] (xa) 2 +(yb) 2 +(zc) 2 =r 2

[0279] Where (x, y, z) are the coordinates of the point on the defined sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius.

[0280] 2 Loss Function Construction

[0281] Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is defined as:

[0282]

[0283] 3 Gradient and Hessian matrix calculation

[0284] Gradient vector:

[0285]

[0286] Hessian matrix:

[0287]

[0288] by For example

[0289]

[0290] The remaining elements are calculated according to the second-order partial derivative rule to construct a 4×4 symmetric matrix.

[0291] 4 Parameter iterative update

[0292] parameter vector The update formula is:

[0293]

[0294] Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector.

[0295] The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

[0296] 5 Examples

[0297] Step 1: Experimental Setup

[0298] Target specifications: Standard spherical target, theoretical radius r = 75mm.

[0299] Data acquisition: Use a 3D laser scanner to obtain target point cloud data.

[0300] Initial parameters: The initial value is calculated by taking the mean of the point cloud coordinates, and the radius is selected as the median from the point to the center of mass to improve robustness.

[0301] Step 2: Algorithm Implementation

[0302] The Python core code is as follows:

[0303]

[0304]

[0305] Step 3: Optimize parameters and calculate errors

[0306] After iterative calculation and optimization,

[0307] Geometric center and error:

[0308] P0(x0,y0,z0)=(0.4134,0.1856,0.4432)mm=φ(P0)

[0309] Radius and error:

[0310] r=74.9953mm;φ(r)=0.0047mm

[0311] The traditional error is as follows:

[0312] P0(x0,y0,z0)=(6.6385,3.3251,2.2295)mm=φ(P0)

[0313] r=78.5553mm; φ(r)=3.5553mm

[0314] 6 Technical Effects

[0315] High-precision error compensation: The geometric center and radius of the spherical target are optimized using the Newton iteration method. The final geometric center error is reduced to (0.4134, 0.1856, 0.4432), and the radius error is 0.0047mm. Compared with the traditional least squares method (error 3.5553mm), the accuracy is significantly improved.

[0316] Fast convergence: The Newton method uses second-order derivative information (Hessian matrix) to accelerate convergence, reaching the error threshold within 100 iterations on average, significantly improving calibration efficiency.

[0317] Enhanced anti-interference capability: By integrating the spatial distribution characteristics of three-dimensional point clouds, the algorithm is significantly more robust to environmental and edge point interference than traditional methods, ensuring the reliability of three-dimensional feature measurements such as road crack depth and unevenness.

[0318] Practical application advantages: After optimization, the accuracy of device point cloud data is improved, which can directly reduce the calculation deviation of three-dimensional geometric quantities and provide highly reliable calibration data support for geometric quantity measurement in the highway industry.

[0319] In this embodiment, point cloud data calibration of multi-site acquisition equipment: To address the problems of mismatched point cloud data and low point cloud accuracy of multi-site 3D scanning equipment, the coordinates of different sites are converted using the four-parameter method, the maximum geometric error is used to align the same-name areas, and the NDT optimization method is used to calibrate the geometric error of the point cloud data. The accuracy of the point cloud data is improved through algorithm compensation.

[0320] 1 Four-parameter coordinate calibration of point cloud data

[0321] In multi-site scanning, point clouds in different coordinate systems need to be unified into the global coordinate system. Assume that the same-name point pairs of the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i ). Plane coordinate calibration is performed using the four-parameter method, and the formula is as follows:

[0322]

[0323] Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount.

[0324] Four parameter calibration steps:

[0325] ① Centroid decentralization: Calculate the centroid of the source point cloud and the target point cloud And decentralized processing.

[0326] ②Calculate intermediate variables

[0327] C1=∑(x t x s +y t y s ); C2=∑(y t x s -x t y s )

[0328] C3=∑x s ; C4=∑y s ; C5=∑x t ; C6=∑y t

[0329] ③ Solve the rotation angle:

[0330] θ=arctan2(C2,C1)

[0331] ④Calculate the translation amount:

[0332] t x =(C5-C3cosθ+C4sinθ) / N

[0333] t y =(C6-C3sinθ-C4cosθ) / N

[0334] ⑤Vertical axis calibration:

[0335] pass Compensate for elevation errors.

[0336] 2. NDT calibration optimization of point cloud data

[0337] In order to improve the calibration accuracy and efficiency, the NDT (Normal Distributions Transform) method is used to optimize the initial calibration results:

[0338] Meshing: Divide the target point cloud into a fixed-size cubic grid

[0339] Probability density modeling: Calculate the mean μ and covariance matrix Γ of the point cloud within each grid:

[0340]

[0341] Optimization goal: maximize the overall likelihood function

[0342]

[0343] Iterative optimization: Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold

[0344]

[0345] Where M k+1 and M k represents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function. λ is the damping factor of the algorithm, which is used to balance gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to gradient descent. When λ is small, the algorithm is close to Gauss-Newton method. I is a bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems.

[0346] 3. Fixed the same-named areas of point cloud data

[0347] In order to reduce manual intervention, an automatic calibration strategy of regions with the same name based on maximum geometric error is proposed:

[0348] Distance calculation: Calculate the maximum Euclidean distance between the source point cloud and the target point cloud

[0349]

[0350] Threshold setting: filter to meet d i ≤ε, where the threshold

[0351] KDTree acceleration: Use the KDTree structure of the Open3D library to quickly match point clouds and improve computing efficiency.

[0352] 4. On-site application of project implementation

[0353] A standard road section at a testing site in Beijing was selected for multi-site scanning testing:

[0354] Equipment configuration: FARO Focus X Series medium- and long-range laser scanner, laser wavelength 1550nm, ambient temperature on that day was 16.8°C.

[0355] Data acquisition: Four scanning stations were deployed to cover a 10m × 10m × 3m test space to acquire raw point cloud data.

[0356] Preprocessing: Remove noise points using the Z-score method.

[0357] Plane calibration: Combined with a standard plane target, the device's planar error is calibrated using a regional adaptive weighted average algorithm.

[0358] Stereo calibration: Combined with a standard spherical target, the Newton iteration method is used to calibrate the stereo error of the device.

[0359] Multi-site equipment calibration: The coordinate system is calibrated using the four-parameter method to perform coarse registration of multi-site equipment data; and point cloud alignment is optimized through NDT to perform fine registration of multi-site equipment data.

[0360] 5 Examples

[0361] Step 1: Experimental Setup

[0362] Test object: Use laser scanning equipment to scan a crack disease on the test road section

[0363] Parameter comparison: The geometric parameters of the crack were measured using traditional methods before calibration. A high-level standard instrument (Renishaw XL-80 laser interferometer, certified by the National Institute of Metrology, CDlx2024-01568) was used to calculate the measurement error of the traditional method. The area calculation error was 0.0181m 2 , the circumference calculation error is 18.1805mm, and the depth calculation error is 5.3902mm.

[0364] NDT parameter settings: voxel size is set to 0.01, maximum number of iterations is set to 100, and Open3D parallel acceleration is enabled.

[0365] Step 2: Algorithm Implementation

[0366] The Python core code is as follows:

[0367]

[0368] Step 3: Comparison of measured parameters

[0369] The pit-compensated morphology is obtained after data preprocessing, plane error calibration, stereo error calibration and multi-site error calibration.

[0370] The area calculation error is 0.0031m 2 (82.91% optimized compared to before calibration), the circumference calculation error is 1.2326mm (93.22% optimized compared to before calibration), and the depth calculation error is 0.7725mm (accuracy is 85.67% optimized compared to before calibration).

[0371] 6 Technical Effects

[0372] Accuracy verification: The calibrated data meets the detection accuracy requirements of geometric quantities in the highway transportation industry, and the error of feature geometric quantity detection is reduced.

[0373] Improved efficiency: Multi-site registration time (approximately 8 minutes) is shortened to 60% of traditional methods.

[0374] In summary, 1. Multi-level error calibration system based on target:

[0375] Planar calibration uses a regional adaptive improved weighted mean algorithm; stereo calibration uses the Newton iteration method; multi-site calibration integrates the four-parameter method (coarse alignment) and NDT (fine alignment) to achieve global point cloud alignment.

[0376] 2. Automated data processing process:

[0377] Noise removal uses the Z-score method; the regions with the same name are filtered by dynamic thresholding. Automatic filtering.

[0378] 3. Standardization project implementation plan:

[0379] The target design complies with highway traffic industry standards and is suitable for geometric calibration of various types of infrastructure structures and various geometric models such as spheres and planes.

[0380] Example 2

[0381] like Figure 2 As shown, the present invention also provides a target-based multi-layer calibration system for three-dimensional geometric measurement of highway traffic structures. The system is used to implement the aforementioned method. The system includes: a preprocessing module, a plane error calibration module, a stereo error calibration module, and a region matching module;

[0382] A pre-processing module is used to collect point cloud data of highway infrastructure structures using a 3D laser scanner and perform denoising pre-processing;

[0383] The plane error calibration module is used to calibrate the plane error of the 3D laser scanning equipment using a standard plane target based on the point cloud data after denoising preprocessing;

[0384] The stereo error calibration module is used to calibrate the stereo error of the 3D laser scanning device using a standard spherical target based on the point cloud data after denoising preprocessing;

[0385] The region matching module is used to convert the coordinates of different sites using the four-parameter method based on the calibration results of the plane error and stereo error of the 3D laser scanning equipment, and to perform homonymous region registration using the maximum geometric error.

[0386] In this embodiment, a 3D laser scanner is used to collect point cloud data of a highway infrastructure structure and perform denoising preprocessing, including:

[0387] Convert the three-dimensional point cloud coordinates to the Euclidean distance of the origin to construct a one-dimensional data set;

[0388] Calculate the mean of a one-dimensional data set, Z-score or Z i ;

[0389] According to the calculation results of Z-score, set the dynamic threshold;

[0390] Eliminate the ones that satisfy |Z i ∣>Threshold point.

[0391] In this embodiment, based on the point cloud data after denoising preprocessing, the plane error of the 3D laser scanning device is calibrated using a standard plane target, including:

[0392] Get the area of ​​the area to be measured;

[0393] Calculating the point cloud density distribution, and calculating the neighborhood radius according to the area of ​​the area to be measured and the point cloud density distribution;

[0394] Calculate the initial value of the geometric center according to the neighborhood radius;

[0395] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the calibration of the plane error of the 3D laser scanning equipment;

[0396] The area of ​​the area to be measured is obtained, including:

[0397]

[0398] Calculate point cloud density distribution, including:

[0399]

[0400] Where N r is the actual number of sampling points in the area to be measured, which is determined according to the on-site measurement requirements; N tot is the theoretical number of sampling points; λ is the empirical parameter of point cloud density; r is the neighborhood radius; D is the point cloud density distribution function;

[0401] Calculating a neighborhood radius according to the area of ​​the region to be measured and the point cloud density distribution includes:

[0402]

[0403] Calculate the initial value of the geometric center based on the neighborhood radius, including:

[0404] Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0):

[0405]

[0406] Where N i is the number of sampling points in the ith sub-region within the neighborhood, and n is the number of sub-regions;

[0407] Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the plane error calibration of the point cloud data, including:

[0408] Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i :

[0409]

[0410] Where, is the temporary geometric center coordinate of the kth iteration;

[0411] The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm:

[0412]

[0413] The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than the threshold ∈.

[0414] In this embodiment, based on the point cloud data after denoising preprocessing, the stereo error of the 3D laser scanning device is calibrated using a standard spherical target, including:

[0415] Build a spherical model;

[0416] Based on the spherical model, construct a loss function;

[0417] Based on the loss function, calculate the gradient and Hessian matrix;

[0418] Based on the gradient and the Hessian matrix, iteratively updating parameters is performed to complete the calibration of the stereo error of the three-dimensional laser scanning device;

[0419] Among them, establishing a spherical model includes:

[0420] The default spherical equation is:

[0421] (xa) 2 +(yb) 2 +(zc) 2 =r 2 ;

[0422] Where (x, y, z) are the coordinates of the point on the preset sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius;

[0423] Based on the spherical model, a loss function is constructed, including:

[0424] Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is:

[0425]

[0426] Based on the loss function, the gradient and Hessian matrix are calculated, including:

[0427] Gradient vector:

[0428]

[0429] Hessian matrix:

[0430]

[0431] Based on the gradient and Hessian matrix, the parameters are iteratively updated, including:

[0432] parameter vector The update formula is:

[0433]

[0434] Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector.

[0435] The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

[0436] In this embodiment, based on the calibration results of the planar error and stereo error of the 3D laser scanning device, the coordinates of different sites are converted using the four-parameter method, and the maximum geometric error is used to perform homonymous region registration, including:

[0437] Calibrate the four-parameter coordinates of point cloud data;

[0438] Use NDT method to optimize calibration results;

[0439] Based on the optimized structure, the regions with the same name in the point cloud data are fixed;

[0440] Among them, the four-parameter coordinates of the point cloud data are calibrated, including:

[0441] The preset point pairs of the same name between the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i ); plane coordinate calibration is performed using the four-parameter method, and the calculation formula is:

[0442]

[0443] Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount;

[0444] Optimize calibration results using NDT methods, including:

[0445] Divide the target point cloud into a fixed-size cubic grid;

[0446] Calculate the mean μ and covariance matrix Γ of the point cloud in each grid and build a probability density model;

[0447]

[0448] Based on the probability density model, the objective is optimized by maximizing the overall likelihood function;

[0449]

[0450] Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold:

[0451]

[0452] Where M k+1 and M krepresents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function. λ is the damping factor of the algorithm, which is used to balance gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to gradient descent. When λ is small, the algorithm is close to Gauss-Newton method. I is a bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems.

[0453] Based on the optimized structure, the regions with the same name in the point cloud data are fixed, including:

[0454] Calculate the maximum Euclidean distance between the source point cloud and the target point cloud:

[0455]

[0456] Threshold setting: filter to meet d i ≤ε, where the threshold

[0457] Use the KDTree structure of the Open3D library to match point clouds.

[0458] The embodiments described above are merely descriptions of preferred embodiments of the present invention and are not intended to limit the scope of the present invention. Without departing from the spirit of the present invention, various modifications and improvements made to the technical solutions of the present invention by persons skilled in the art should fall within the scope of protection defined by the claims of the present invention.

Claims

1. A multi-layer calibration method for three-dimensional geometric quantities of highway traffic structures based on target reference, characterized in that: The method comprises: Use 3D laser scanners to collect point cloud data of highway infrastructure structures and perform denoising preprocessing; Based on the denoised point cloud data, the plane error of the 3D laser scanning equipment is calibrated using a standard plane target. Based on the denoised point cloud data, a standard spherical target is used to calibrate the stereo error of the 3D laser scanning equipment. Based on the calibration results of the planar error and stereo error of the 3D laser scanning equipment, the coordinates of different sites are transformed using the four-parameter method, and the maximum geometric error is used to perform homonymous region registration.

2. The method according to claim 1, characterized in that Use 3D laser scanners to collect point cloud data of highway infrastructure structures and perform denoising preprocessing, including: Convert the three-dimensional point cloud coordinates to the Euclidean distance of the origin to construct a one-dimensional data set; Calculate the mean of a one-dimensional data set, Z-score or Z i ; According to the calculation results of Z-score, set the dynamic threshold; Eliminate the ones that satisfy |Z i ∣>Threshold point.

3. The method according to claim 1, characterized in that Based on the denoised point cloud data, the plane error of the 3D laser scanning equipment is calibrated using a standard plane target, including: Get the area of ​​the area to be measured; Calculating the point cloud density distribution, and calculating the neighborhood radius according to the area of ​​the area to be measured and the point cloud density distribution; Calculate the initial value of the geometric center according to the neighborhood radius; Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the calibration of the plane error of the 3D laser scanning equipment; The area of ​​the area to be measured is obtained, including: Calculate point cloud density distribution, including: Where N r is the actual number of sampling points in the area to be measured, which is determined according to the on-site measurement requirements; N tot is the theoretical number of sampling points; λ is the empirical parameter of point cloud density; r is the neighborhood radius; D is the point cloud density distribution function; Calculating a neighborhood radius according to the area of ​​the region to be measured and the point cloud density distribution includes: Calculate the initial value of the geometric center based on the neighborhood radius, including: Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0): Where N i is the number of sampling points in the ith sub-region within the neighborhood, and n is the number of sub-regions; Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the plane error calibration of the point cloud data, including: Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i : Where, is the temporary geometric center coordinate of the kth iteration; The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm: The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than the threshold ∈.

4. The method according to claim 3, characterized in that Based on the denoised point cloud data, a standard spherical target is used to calibrate the stereo error of the 3D laser scanning device, including: Build a spherical model; Based on the spherical model, construct a loss function; Based on the loss function, calculate the gradient and Hessian matrix; Based on the gradient and the Hessian matrix, iteratively updating parameters is performed to complete the calibration of the stereo error of the three-dimensional laser scanning device; Among them, establishing a spherical model includes: The default spherical equation is: (x-a) 2 +(y-b) 2 +(z-c) 2 =r 2 ; Where (x, y, z) are the coordinates of the point on the preset sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius; Based on the spherical model, a loss function is constructed, including: Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is: Based on the loss function, the gradient and Hessian matrix are calculated, including: Gradient vector: Hessian matrix: Based on the gradient and Hessian matrix, the parameters are iteratively updated, including: parameter vector The update formula is: Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector; The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

5. The method according to claim 4, characterized in that Based on the calibration results of the plane error and stereo error of the 3D laser scanning equipment, the coordinates of different sites are converted using the four-parameter method, and the maximum geometric error is used to perform homonymous region registration, including: Calibrate the four-parameter coordinates of point cloud data; Use NDT method to optimize calibration results; Based on the optimized structure, the regions with the same name in the point cloud data are fixed; Among them, the four-parameter coordinates of the point cloud data are calibrated, including: The preset point pairs of the same name between the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i ); plane coordinate calibration is performed using the four-parameter method, and the calculation formula is: Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount; Optimize calibration results using NDT methods, including: Divide the target point cloud into a fixed-size cubic grid; Calculate the mean μ and covariance matrix Γ of the point cloud in each grid and build a probability density model; Based on the probability density model, the objective is optimized by maximizing the overall likelihood function; Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold: Where M k+1 and M k represents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function, λ is the damping factor of the algorithm, which is used to balance the gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to the gradient descent, and when λ is small, the algorithm is close to the Gauss-Newton method. I is the bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems. Based on the optimized structure, the regions with the same name in the point cloud data are fixed, including: Calculate the maximum Euclidean distance between the source point cloud and the target point cloud: Threshold setting: filter to meet d i ≤ε, where the threshold Use the KDTree structure of the Open3D library to match point clouds.

6. A target-referenced, multi-layer calibration system for three-dimensional geometric quantities of highway traffic structures, the system being used to implement the method of any one of claims 1 to 5, characterized in that: The system includes: a pre-processing module, a plane error calibration module, a stereo error calibration module, and a region matching module; The pre-processing module is used to collect point cloud data of highway infrastructure structures using a three-dimensional laser scanner and perform denoising pre-processing; The plane error calibration module is used to calibrate the plane error of the three-dimensional laser scanning device using a standard plane target based on the point cloud data after denoising preprocessing; The stereo error calibration module is used to calibrate the stereo error of the three-dimensional laser scanning device using a standard spherical target based on the point cloud data after denoising preprocessing; The region matching module is used to convert the coordinates of different sites using a four-parameter method based on the calibration results of the plane error and stereo error of the three-dimensional laser scanning device, and to perform homonymous region registration using the maximum geometric error.

7. The system according to claim 6, characterized in that Use 3D laser scanners to collect point cloud data of highway infrastructure structures and perform denoising preprocessing, including: Convert the three-dimensional point cloud coordinates to the Euclidean distance of the origin to construct a one-dimensional data set; Calculate the mean of a one-dimensional data set, Z-score or Z i ; According to the calculation results of Z-score, set the dynamic threshold; Eliminate the ones that satisfy |Z i ∣>Threshold point.

8. The system according to claim 6, wherein: Based on the denoised point cloud data, the plane error of the 3D laser scanning equipment is calibrated using a standard plane target, including: Get the area of ​​the area to be measured; Calculating the point cloud density distribution, and calculating the neighborhood radius according to the area of ​​the area to be measured and the point cloud density distribution; Calculate the initial value of the geometric center according to the neighborhood radius; Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the calibration of the plane error of the 3D laser scanning equipment; The area of ​​the area to be measured is obtained, including: Calculate point cloud density distribution, including: Where N r is the actual number of sampling points in the area to be measured, which is determined according to the on-site measurement requirements; N tot is the theoretical number of sampling points; λ is the empirical parameter of point cloud density; r is the neighborhood radius; D is the point cloud density distribution function; Calculating a neighborhood radius according to the area of ​​the region to be measured and the point cloud density distribution includes: Calculate the initial value of the geometric center based on the neighborhood radius, including: Construct a neighborhood with a neighborhood radius r, including N point cloud data, and calculate the initial geometric center coordinates C0 = (x0, y0): Where N i is the number of sampling points in the ith sub-region within the neighborhood, and n is the number of sub-regions; Based on the initial value of the geometric center, the geometric center coordinates are iteratively optimized through the weighted mean algorithm to complete the plane error calibration of the point cloud data, including: Based on the distance d from the point cloud data to the temporary geometric center i , define the weight function ω i : Where, is the temporary geometric center coordinate of the kth iteration; The geometric center coordinates C = (x, y) are iteratively optimized using the weighted mean algorithm: The iteration termination condition is that the absolute value of the difference between two adjacent coordinates is less than the threshold ∈.

9. The system according to claim 8, characterized in that Based on the denoised point cloud data, a standard spherical target is used to calibrate the stereo error of the 3D laser scanning device, including: Build a spherical model; Based on the spherical model, construct a loss function; Based on the loss function, calculate the gradient and Hessian matrix; Based on the gradient and the Hessian matrix, iteratively updating parameters is performed to complete the calibration of the stereo error of the three-dimensional laser scanning device; Among them, establishing a spherical model includes: The default spherical equation is: (x-a) 2 +(y-b) 2 +(z-c) 2 =r 2 ; Where (x, y, z) are the coordinates of the point on the preset sphere, (a, b, c) are the three-dimensional coordinates of the center of the sphere, and r is the radius; Based on the spherical model, a loss function is constructed, including: Minimize the sum of the squares of the distances from all point clouds to the fitting sphere. The loss function is: Based on the loss function, the gradient and Hessian matrix are calculated, including: Gradient vector: Hessian matrix: Based on the gradient and Hessian matrix, the parameters are iteratively updated, including: parameter vector The update formula is: Where, is the updated parameter vector, is the parameter vector before updating, H is the Hessian matrix, is the gradient vector; The iteration termination condition is: radius error Δr ≤ 0.003m or the maximum number of iterations is reached.

10. The system according to claim 9, characterized in that Based on the calibration results of the plane error and stereo error of the 3D laser scanning equipment, the coordinates of different sites are converted using the four-parameter method, and the maximum geometric error is used to perform homonymous region registration, including: Calibrate the four-parameter coordinates of point cloud data; Use NDT method to optimize calibration results; Based on the optimized structure, the regions with the same name in the point cloud data are fixed; Among them, the four-parameter coordinates of the point cloud data are calibrated, including: The preset point pairs of the same name between the source point cloud and the target point cloud are: Source point cloud coordinates: P t,i =(x i ,y i ), target point cloud coordinates: P s,i =(x s,i ,y s,i ); plane coordinate calibration is performed using the four-parameter method, and the calculation formula is: Among them, a is the scale transformation parameter, θ is the rotation angle, t x and t y is the translation amount; Optimize calibration results using NDT methods, including: Divide the target point cloud into a fixed-size cubic grid; Calculate the mean μ and covariance matrix Γ of the point cloud in each grid and build a probability density model; Based on the probability density model, the objective is optimized by maximizing the overall likelihood function; Update the transformation matrix using the Levenberg-Marquardt algorithm until the difference between two adjacent transformations is less than the threshold: Where M k+1 and M k represents the transformation parameters for the kth and k+1th iterations, represents the sum of the log-likelihood function gradients of all grid cells or data points, is the sum of the second-order derivatives of the log-likelihood function, λ is the damping factor of the algorithm, which is used to balance the gradient descent and Gauss-Newton method. When λ is large, the algorithm is close to the gradient descent, and when λ is small, the algorithm is close to the Gauss-Newton method. I is the bit matrix, which is combined with λ to ensure the matrix is ​​reversible and avoid singularity problems. Based on the optimized structure, the regions with the same name in the point cloud data are fixed, including: Calculate the maximum Euclidean distance between the source point cloud and the target point cloud: Threshold setting: filter to meet d i ≤ε, where the threshold Use the KDTree structure of the Open3D library to match point clouds.

Images