Tunnel point cloud adaptive elliptic cylinder filtering method based on affine constraint

Abstract

The invention provides a tunnel point cloud adaptive elliptic cylinder filtering method based on affine constraint, and relates to the technical field of point cloud data processing, and the method comprises the steps: extracting a tunnel center line, and setting a cross section; initializing parameters of the elliptic cylinder model by using principal component analysis; setting an initial screening threshold value, and dynamically estimating an optimal threshold value by adopting an RANSAC (Random Sample Consensus) algorithm; constructing a multi-objective optimization function comprising a data fitting error term, a parameter smoothing constraint term and an affine transformation constraint term, and optimizing model parameters; and filtering the point cloud data based on the optimized model and the dynamic threshold. Affine transformation constraint and self-adaptive parameter constraint mechanisms are innovatively introduced, the constraint threshold is dynamically adjusted according to the local curvature and the deformation degree, a special recognition and processing mechanism is developed for a special area, the geometric continuity of the model and the adaptability to the complex environment are remarkably improved, and the method is suitable for large-scale popularization and application. The separation effect of noise points and structural points is improved, and a reliable data basis is provided for tunnel three-dimensional reconstruction and deformation analysis.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of point cloud data processing, and particularly relates to a tunnel point cloud adaptive elliptic cylinder filtering method based on affine constraint. BACKGROUND

[0002] In the field of tunnel point cloud data processing, accurate filtering is a key prerequisite for realizing tunnel three-dimensional reconstruction and deformation monitoring. With the wide application of three-dimensional laser scanning technology in tunnel engineering, how to extract effective information from high-noise and high-density tunnel point cloud data and accurately identify tunnel structural surfaces has become a research hotspot. However, the existing technology still has many limitations in processing complex tunnel environments.

[0003] The existing tunnel point cloud processing methods mainly include filtering methods based on statistical analysis, geometric feature extraction methods and model fitting methods. The filtering methods based on statistical analysis, such as the statistical outlier removal (SOR) model, usually use fixed parameters to identify and remove noise points. However, this kind of method is sensitive to parameter selection, and it is difficult to adapt to the dynamic changes of tunnel cross-sectional shape and size, which may lead to under-filtering or over-filtering problems. Especially when dealing with local deformation or special structures (such as shelters, cross passages, etc.) of tunnels, the method with fixed parameters often cannot accurately identify the structural characteristics, resulting in loss or misjudgment of useful information.

[0004] The geometric feature extraction method, including the method based on normal vector analysis and local curvature estimation, can better capture the geometric properties of the tunnel surface, but the calculation efficiency is low when dealing with high-noise and high-density point cloud data. In addition, this kind of method has high requirements for the uniformity of point cloud density. In the actual collected tunnel point cloud data, due to the factors such as occlusion and scanning angle, the problem of uneven density often affects the accuracy of feature extraction.

[0005] The model fitting method, such as the random sample consensus (RANSAC) algorithm and the Hough transform, identifies the tunnel structure by fitting a predefined geometric model (such as a cylinder or an elliptic cylinder). This kind of method works well in ideal conditions, but when facing complex tunnel structures (such as double-line tunnels, variable cross-section tunnels, etc.), a single model is difficult to accurately describe the geometric shape of the entire tunnel. At the same time, these methods usually use a global optimization strategy, which has high computational complexity and is difficult to meet the real-time processing requirements of large-scale point cloud data.

[0006] The prior art, such as the Chinese patent application with the publication number CN119090728A, discloses a self-adaptive elliptical cylinder filtering method based on tunnel point cloud data. The method combines the RANSAC algorithm with dynamic threshold estimation, and uses PCA initialization and iterative least squares method to optimize the parameters of the elliptical cylinder model. However, the method still has the following shortcomings: on the one hand, the method lacks an effective affine transformation constraint mechanism, making it difficult to ensure the rationality and continuity of the parameter changes between adjacent sections. Traditional methods often use fixed thresholds or simple smoothing processing for constraints, which cannot adapt to the characteristics of curvature changes in different regions, especially in curved and deformed sections of the tunnel, and fixed threshold constraints are prone to over-constraint or under-constraint problems.

[0007] On the other hand, the optimization of model parameters in the prior art often uses a single error metric (such as the distance from the point to the model), lacking a multi-objective optimization strategy that considers global geometric continuity. This leads to a lack of reasonable variation in model parameters in the longitudinal direction when processing long tunnels or tunnels under complex geological conditions, affecting the accuracy and reliability of the overall model.

[0008] In addition, the prior art lacks a special recognition and processing mechanism for special areas such as portal sections, making it difficult to cope with changes in the geometric characteristics of these areas. At the same time, the existing methods have a single means of parameter smoothing transition processing, which cannot effectively handle sudden changes in section parameters, affecting the stability and reliability of the model.

[0009] In summary, the prior art lacks a comprehensive solution that considers model fitting accuracy, geometric continuity, local characteristic adaptability, and computational efficiency when processing tunnel point cloud data, and there is an urgent need for a new method that can balance these aspects to improve the overall performance of tunnel point cloud data processing and provide a more reliable data foundation for tunnel three-dimensional reconstruction and deformation monitoring. SUMMARY

[0010] The purpose of the present application is to provide an affine constraint-based tunnel point cloud adaptive elliptical cylinder filtering method, which solves the problem of the prior art lacking comprehensive consideration of model fitting accuracy, geometric continuity, and local characteristic adaptability when processing tunnel point cloud data, and provides a reliable data foundation for tunnel three-dimensional reconstruction and deformation monitoring.

[0011] To achieve the above-mentioned purpose, the present application realizes the following technical solutions: An affine constraint-based tunnel point cloud adaptive elliptical cylinder filtering method, comprising the following steps: S1: Extracting a tunnel center line based on tunnel point cloud data, and setting sections at fixed intervals along the center line axis; S2: Initializing the parameters of the elliptical cylinder model on each section using the principal component analysis method to form an elliptical cylinder model; S3: setting an initial screening threshold for preliminarily screening abnormal points of distances from all points in the cross section to the elliptic cylindrical surface, and dynamically estimating an optimal screening threshold by using a RANSAC algorithm; S4: constructing a multi-objective optimization function containing a data fitting error term, a parameter smoothing constraint term and an affine transformation constraint term, and optimizing the elliptic cylindrical model parameters by using an optimization algorithm until the elliptic cylindrical model parameters converge; S5: filtering the tunnel point cloud data based on the optimized elliptic cylindrical model and the dynamic threshold, and outputting the filtered point cloud data and the elliptic cylindrical model parameters.

[0012] Further, in step S1, the least square method or the principal curve algorithm is used to extract the tunnel center line from the preliminarily denoised tunnel point cloud data, and the fixed interval range is 0.5-3 meters.

[0013] Further, in step S2, the principal component analysis method is used to initialize the elliptic cylindrical model parameters, specifically including: calculating the covariance matrix of the cross section point cloud, performing eigenvalue decomposition on the covariance matrix, determining the directions of the major axis and the minor axis of the elliptic cylindrical model by using the eigenvectors corresponding to the maximum eigenvalue and the second largest eigenvalue, and estimating the lengths of the major axis and the minor axis according to the proportion of the eigenvalues.

[0014] Further, in step S3, the initial screening threshold is 20 cm, and the optimal screening threshold is the 95th percentile of the distances of the inner points.

[0015] Further, in step S4, the multi-objective optimization function is wherein is the sum of squares of distances of the point cloud to the elliptic cylindrical surface, is the sum of squares of differences of parameters of adjacent cross sections, is the sum of squares of differences of affine transformation matrices between adjacent cross sections and the expected transformation matrix, and λ1 and λ2 are weight coefficients.

[0016] Further, the Levenberg-Marquardt algorithm is used to minimize the multi-objective optimization function, the initial damping factor of the algorithm is 0.01, and the damping factor adjustment rate is 10.

[0017] Further, the weight coefficient λ1 is in the range of 0.5-2.0, and the weight coefficient λ2 is in the range of 0.3-1.5.

[0018] Further, in step S4, it further includes the step of dynamically adjusting the constraint threshold according to the local curvature and the deformation degree, specifically including: dynamically calculating the direction angle change threshold according to the local curvature, the basic threshold is 2°, and the maximum threshold is 8°; and dynamically calculating the axis length change threshold according to the local deformation degree, the basic threshold is 2%, and the maximum threshold is 6%.

[0019] Further, in step S4, a step of smoothing the affine transformation parameters is further included, specifically comprising: calculating the affine transformation matrix and decomposing it into rotation, scaling and translation components; smoothing the transformation parameters using Gaussian weights; and applying cubic spline interpolation to realize smooth transition of the parameters.

[0020] Further, the curvature change rate threshold is 0.2, the optimization window for local re-optimization is 5 sections, the confidence level for abnormal change identification is 97%, and the weight decay rate of the local weighted average method is 0.7.

[0021] Further, in step S5, points within the optimal screening threshold range are retained as tunnel structure points, and points exceeding the optimal screening threshold are removed as non-structure points or noise.

[0022] Further, the expected transformation matrix of the affine transformation constraint term is constructed according to the tunnel segment type, including: the expected transformation matrix of a straight segment , wherein is a unit matrix, is a perturbation coefficient, and the value is 0.01-0.05, is a perturbation matrix; the expected transformation matrix of a curved segment contains a rotation component , which is constructed according to the local curvature ; the expected transformation matrix of a deformation segment contains a scaling component , which is constructed based on the deformation rate of adjacent sections; and the expected transformation matrix of a comprehensive segment .

[0023] Further, a multi-level verification and correction mechanism is further included, specifically comprising: a physical feasibility verification checking whether the ellipse after transformation satisfies the axis length change limit and the axis ratio limit ; a conflict constraint coordination mechanism when there is a conflict between data fitting and affine transformation constraint, the conflict determination standard is , and the data fitting accuracy is prioritized; an abnormality detection and processing mechanism marking and locally re-optimizing the case of abnormal change of the model parameters after constraint.

[0024] Further, for a tunnel straight segment with a length exceeding 50 meters, a main direction stability constraint is further included, which is calculated by accumulating the offset and the constraint condition , wherein is the length of the straight segment, and a direction cumulative penalty term is added to the multi-objective optimization function , wherein is a weight coefficient, and the value is 0.2-0.5.

[0025] Compared with the prior art, the present application has the following beneficial effects: I. The method of the present application improves the geometric continuity of the model by introducing affine transformation constraints and parameter smoothing transition processing. Affine transformation constraints ensure the physical reasonableness of parameter changes between adjacent sections; adaptive parameter change thresholds dynamically adjust according to local curvature and deformation degree, ensuring deformation detection sensitivity while avoiding over-constraint; Gaussian weight smoothing and cubic spline interpolation technology ensures the smooth transition of parameters in the longitudinal direction, making the final model more realistically reflect the actual geometric shape of the tunnel.

[0026] II. For regions with significant changes in geometric characteristics such as portal sections, the special recognition and processing mechanism developed by the present application significantly improves the adaptability of the model. Special areas are identified through curvature change rate analysis; local re-optimization processing of mutation points; parameter change statistical model identifies abnormal changes; local weighted average method ensures smooth transition. These mechanisms work together to enable the model to accurately adapt to various special structures in the tunnel, improving the overall accuracy of point cloud filtering.

[0027] III. The present application constructs a comprehensive optimization objective function including data fitting error, parameter smoothing constraint and affine transformation constraint, which ensures the geometric continuity while ensuring the fitting accuracy. Compared with traditional methods, the present application considers data fitting accuracy, parameter smoothness and geometric continuity through a multi-objective optimization function, avoiding the local optimal problem caused by only focusing on the point-to-model distance, and ensuring the global optimality of the model parameters.

[0028] IV. The method of the present application dynamically adjusts the constraint threshold according to the local curvature and deformation degree, so that the direction angle and axis length change threshold can be adaptively adjusted according to the local characteristics of the tunnel, solving the problem of poor adaptability of fixed threshold constraint in different areas. This dynamic parameter constraint mechanism can more accurately capture the local characteristics of the tunnel, ensuring the fitting accuracy while avoiding overfitting or underfitting.

[0029] V. The present application significantly improves the separation effect of noise points and structure points through dynamic threshold and multi-stage optimization. The initial threshold setting provides a reasonable starting separation standard; the dynamic adjustment of the threshold in the RANSAC process adapts to the characteristics of different areas; the multi-objective optimization improves the model fitting accuracy; the final filtering stage uses the optimized model and dynamic threshold for accurate separation, greatly improving the accuracy of structure point recognition, providing a high-quality data basis for subsequent tunnel three-dimensional reconstruction and deformation analysis. BRIEF DESCRIPTION OF DRAWINGS

[0030] Figure 1 The flowchart of an embodiment of the tunnel point cloud adaptive elliptical cylinder filtering method based on affine constraints of the present application. DETAILED DESCRIPTION

[0031] The technical solutions of the present application will be described clearly and completely below with reference to the drawings. Obviously, the described embodiments are only some of the embodiments of the present application, but not all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative work fall within the scope of the present application.

[0032] In the description of the present application, it should be noted that the terms "center", "upper", "lower", "left", "right", "vertical", "horizontal", "inner", "outer" and the like indicate the orientation or positional relationship shown in the drawings, and are only for the convenience of describing the present application and simplifying the description, and do not indicate or imply that the devices or elements referred to must have a particular orientation, be constructed and operated in a particular orientation, and therefore cannot be understood as limiting the present application. In addition, the terms "first", "second", "third" are only for descriptive purposes and cannot be understood as indicating or implying relative importance.

[0033] The present application provides a tunnel point cloud adaptive elliptic cylinder filtering method based on affine constraint, which realizes high-precision point cloud filtering for complex tunnel environment through dynamic parameter adjustment, multi-objective optimization and special position processing. The embodiments of the present application are described in detail as follows: A tunnel point cloud adaptive elliptic cylinder filtering method based on affine constraint includes the following steps: S1: Extract the tunnel center line from the preliminary denoising point cloud using the least square method or the principal curve algorithm, and set the cross section at a fixed interval along the center line axis.

[0034] In practical applications, the threshold value of preliminary denoising can be set to 2-5 times the standard deviation, preferably 3 times, to ensure that the obvious abnormal points are preliminarily removed. With the extracted center line as the reference axis, set the cross section at a fixed interval along the axis, and the interval range is 0.5-3m, preferably 1-2m. In the specific implementation process, the interval size can be adjusted according to the specific characteristics of the tunnel, such as using a larger interval (about 1.5-2m) for straight sections of the tunnel; in the curved section or deformation area, the interval is appropriately reduced (about 1-1.5m) to improve the accuracy. It should be noted that the normal direction of the cross section is perpendicular to the tunnel trend as much as possible to ensure the accuracy and continuity of the model parameters.

[0035] S2: On each cross section, use the principal component analysis method to initialize the elliptic cylinder model parameters to form an elliptic cylinder model.

[0036] The principal component analysis can identify the main distribution direction of the point cloud, so as to estimate the long-short axis ratio and direction of the section. The specific steps include calculating the covariance matrix of the section point cloud, performing eigenvalue decomposition on the covariance matrix, determining the long axis and short axis directions of the ellipse by using the eigenvectors corresponding to the maximum and second maximum eigenvalues, and estimating the long-short axis length according to the proportion of the eigenvalues. In actual application, the eigenvalue filtering threshold can be set to 5% of the maximum eigenvalue, the principal axis direction tolerance is 2°, and the outlier removal threshold is 3σ. For the case that the point cloud density is uneven, a weight compensation mechanism can be used for processing to ensure that the initialization parameters are reasonable and reliable. The weight compensation mechanism can be calculated according to the local density of the point cloud, and the points in the low-density area are given a higher weight, and the points in the high-density area are given a lower weight, so as to balance the contribution of different areas to the model fitting.

[0037] S3: Set an initial filtering threshold for preliminarily screening abnormal points of distances from all points in the section to the elliptical cylindrical surface, and dynamically estimate an optimal filtering threshold by using the RANSAC algorithm.

[0038] Specifically, after initialization, the distances from all points in the section to the elliptical cylindrical surface are calculated. An initial filtering threshold d=20cm is set for preliminarily screening abnormal points. The RANSAC algorithm is used to dynamically estimate an optimal threshold to adapt to the point cloud characteristics in different areas. The RANSAC algorithm iterates 50-500 times, and each time a minimum sample set of 5-15 points is randomly selected to fit an elliptical model, the distances of all points to the model are calculated, and the number of inliers is counted. The model with the most inliers is selected as the best fitting result, and the threshold d is dynamically adjusted according to the best fitting result, and the 95% quantile of the inlier distances is taken. This dynamic threshold strategy significantly improves the adaptability of the method to local characteristic changes. In actual application, the number of RANSAC iterations can be adjusted according to the quality of the point cloud and the characteristics of the tunnel. The number of iterations can be increased in areas with more noise and reduced in areas with less noise to balance accuracy and efficiency.

[0039] S4: A multi-objective optimization function including a data fitting error term, a parameter smoothing constraint term and an affine transformation constraint term is constructed, and an optimization algorithm is used to optimize the parameters of the elliptical cylindrical model until the parameters of the elliptical cylindrical model converge; S401: Multi-objective optimization function construction The core of the application is to construct a multi-objective optimization function to realize comprehensive optimization of the model parameters, and the mathematical expression is:

[0040] wherein, is the sum of the squared distances from the point cloud to the elliptical cylindrical surface, which ensures the data fitting accuracy, and the calculation formula is:

[0041] wherein, denotes the distance from the j th point on the i th cross section to the corresponding elliptic cylinder model; is the sum of squares of the difference of parameters between adjacent cross sections, which ensures the smooth transition of parameters, and its calculation formula is:

[0042] wherein, wherein denote the direction angle, major axis length and minor axis length of the j th cross section ellipse, respectively; is the sum of squares of the difference between the affine transformation matrix between adjacent cross sections and the expected transformation matrix, which ensures geometric continuity, and its calculation formula is:

[0043] wherein, wherein denotes the actual affine transformation matrix from the j 1 th cross section to the j th cross section, denotes the expected affine transformation matrix, denotes the Frobenius norm.

[0044] The weight coefficient λ1 ranges from 0.5 to 2.0, and λ2 ranges from 0.3 to 1.5. In the embodiment of the present application, λ1=1.0 and λ2=0.8 can be taken.

[0045] S402: Calculation and decomposition process of affine transformation matrix constructed by multi-objective optimization function Affine transformation matrix calculation: When establishing the affine transformation relationship between adjacent cross sections, first, the elliptical parameters of the adjacent two cross sections S1 and S2 are extracted by PCA analysis. The long axis direction of the ellipse of the S1 cross section is the unit vector v1, and the short axis direction is the unit vector v2; the long axis direction of the ellipse of the S2 cross section is the unit vector w1, and the short axis direction is the unit vector w2.

[0046] Construction of linear equations: A·v1=w1 A·v2=w2 Let V=[v1,v2] and W=[w1,w2]. Since v1 and v2 are orthogonal unit vectors (principal and secondary axes of the ellipse), the matrix V is invertible, and the affine transformation matrix A=W·V is calculated. -1 .

[0047] Matrix decomposition process: ​Calculate matrix A^T·A and perform eigenvalue decomposition, solve the characteristic equation det(A^T·A-λI)=0, and obtain the two eigenvalues. λ 1 and λ 2 (all are positive real numbers, because A^T·A is a positive definite symmetric matrix), construct a diagonal matrix Simultaneously, the corresponding eigenvectors q1 and q2 are obtained, and orthogonalized to form an orthogonal matrix. .

[0048] Calculate the scaling matrix , represents the scaling of the ellipse along its major and minor axes, where the scaling factor is... and These are the scaling amounts along the major and minor axes, respectively.

[0049] Calculate the rotation matrix:

[0050] Calculate the rotation angle .

[0051] S403: Construction of the expected affine transformation matrix for multi-objective optimization functions Expected affine transformation matrix It is an ideal transformation matrix constructed based on the local geometric and physical properties of the tunnel, and is used as a reference standard for the actual affine transformation of adjacent sections.

[0052] Expected transformation of straight segments: For straight segments of a tunnel, the expected affine transformation is close to the identity matrix I, indicating that the change in cross-sectional shape is very small.

[0053] in I It is the identity matrix. ε For small perturbation coefficients (values ​​ranging from 0.01 to 0.05), D This is the perturbation matrix. D The introduction of this allows for minute shape variations, avoiding over-constraint.

[0054] Desired Transformation of Curved Segments: For curved segments, the desired transformation mainly manifests as a rotational transformation, constructed based on the local curvature κ.

[0055] in α The empirical proportionality coefficient between local curvature and rotation angle is determined experimentally to be in the range of 1.2-1.8, and is dynamically adjusted for different radii of curvature. For large radii of curvature (>500m), take α =1.2-1.4 For medium radius of curvature (200-500m), take...α = 1.4-1.6 Small curvature radius (<200m) takes α = 1.6-1.8 Deformation segment expected transformation: for deformation segment, the expected transformation mainly embodies as scaling transformation, based on the deformation rate of adjacent section to construct the expected scaling matrix:

[0056] Δa and Δb respectively represent the relative change rate of long axis and short axis, β is the smoothing factor (value 0.8-0.95), the deformation rate calculation formula is:

[0057]

[0058] Here and are the next section axis length estimated based on the current deformation trend.

[0059] Comprehensive segment expected transformation: for the comprehensive segment with bending and deformation, the complete expected affine transformation matrix is calculated as:

[0060] The weight of each component is dynamically adjusted according to the dominant deformation type of the current paragraph.

[0061] S404: adaptive setting of constraint conditions of multi-objective optimization function Local geometric property analysis: first calculate the local curvature of center line , where r(s) is the center line parameter equation, and the adjacent section ellipse parameter change rate and are calculated at the same time, set the curvature judgment threshold = 0.005 and deformation judgment threshold = 0.03.

[0062] Paragraph type determination: When κ < and δa < and δb < , it is determined as a straight line segment When κ ≥ and δa < and δb < , it is determined as a bending segment When κ and δa ≥ or δb ≥ ) is determined as a deformation segment When κ ≥ and δa ≥ or δb ≥ ) is determined as a comprehensive segment Constraint parameter adaptive adjustment: α Parameter adjustment formula: α =1.2+0.6·min(1, κ / 0.01) β Parameter adjustment formula: β =0.8+0.15·min(1,max( δa , δb ) / 0.05) Linear segment perturbation coefficient: ε =0.01+0.04·max( δa, δb ) /

[0063] Dynamic control of constraint strength: introduce regional adaptability weight Control constraint strength:

[0064] wherein represents the geometric variation degree of the current section, k is the sensitivity coefficient (value 5-10), and finally The constraint term is represented as .

[0065] S405: Multi-objective optimization function constructs multi-level verification and correction mechanism Physical feasibility verification: check whether the transformed ellipse satisfies the physical constraint conditions, including: Axis length change limit:

[0066] Axis ratio limit:

[0067] If the constraint causes a violation of physical feasibility, automatically reduce the constraint strength to .

[0068] Conflict constraint coordination mechanism: when there is obvious conflict between and constraint, priority is given to guarantee data fitting accuracy. The conflict determination standard is: ​

[0069] The conflict coordination method is to temporarily reduce the constraint weight by 50% and linearly restore it in the subsequent iterations.

[0070] Abnormality detection and processing mechanism: Mark the case of abnormal change in model parameters after constraints (relative change > 20%), and use local re-optimization strategy for abnormal points, reduce the constraint weight and increase the data fitting weight, after re-optimization, use cubic spline smoothing to realize smooth transition with the front and rear sections.

[0071] S406: Multi-objective optimization function construction optimization algorithm implementation In the actual optimization process, the Levenberg-Marquardt algorithm is used to minimize the optimization function, the initial damping factor is 0.01, and the damping factor adjustment rate is 10. The optimization process includes calculating the Jacobian matrix and error vector, using QR decomposition to improve the stability of the solution, multiple random initialization to avoid local optimum, and iterative updating of model parameters until convergence (relative change rate less than 0.1%) or reaching the preset iteration number (take 50 times).

[0072] The Jacobian matrix calculation can use numerical differentiation method, and for complex cross-sectional shape, analytical expression can also be considered to improve the calculation efficiency. QR decomposition can use improved Gram-Schmidt orthogonalization method to improve numerical stability.

[0073] S407: Multi-objective optimization function construction adaptive parameter constraint mechanism Dynamic threshold calculation: dynamically calculate the direction angle change threshold according to the local curvature:

[0074] Where is 2° (basic threshold), is 8° (maximum threshold), is the reference maximum curvature (value 0.015).

[0075] According to the local deformation degree, the axis length change threshold is calculated dynamically:

[0076] Where is 2% (basic threshold), is 6% (maximum threshold), is the deformation degree evaluation index, is 0.08.

[0077] Long straight line segment main direction stability constraint: for the tunnel straight line segment with length exceeding 50m, set the main direction stability constraint, calculate and constraints (Where L is the length of straight line segment), add a directional cumulative penalty term in multi-objective optimization function Where λ 3 is a weight coefficient (0.2-0.5), to ensure the directional stability of long straight line segment.

[0078] S408: Multi-objective optimization function to build smooth transition processing mechanism Parameter sequence preprocessing: in order to ensure the smooth transition of adjacent section parameters after optimization, first of all, parameter sequence preprocessing. Initialization parameter sequence Where represents the elliptical parameter of the first j section.

[0079] Using the median absolute deviation (MAD) method to identify abnormal parameters, by calculating And mark the satisfy Abnormal points, for these abnormal points using local median filtering for preliminary correction.

[0080] Gaussian weight smoothing: the window size w is set to 5%-10% of the total section number (minimum 5 sections), Gaussian weight smoothing is applied to each parameter sequence, the smoothing formula is:

[0081] Where the weight calculation uses , When the window exceeds the boundary, the mirror expansion method is used to supplement the boundary data, and for special areas (such as hole section) the smoothing coefficient is adjusted to 60% of the original value.

[0082] Cubic spline interpolation: based on the smoothed parameter control point to build cubic spline function S(t) On each interval , S(t) Indicates a cubic polynomial:

[0083] Satisfy the following conditions: Interpolation condition:

[0084] First derivative continuity:

[0085] Second derivative continuity:

[0086] Natural boundary condition:

[0087] The polynomial coefficients of all intervals are obtained by solving the tri-diagonal linear equations. Finally, an interpolation point is set every 10-20 cm between the original cross-section positions to generate a smooth and continuous parametric curve.

[0088] Rate of change constraint verification: Verify whether the rate of change between adjacent points of the interpolated parameter sequence meets the constraint condition. If there is an interval that violates the constraint, adjust the control points of that interval and re-interpolate. Iterate this process until all intervals meet the constraint condition or the maximum number of iterations (5 times) is reached, ensuring that the final generated parameter sequence meets both the smoothness requirement and the physical constraint condition.

[0089] S409: Multi-objective optimization function construction special position processing mechanism Hole segment identification and quantization standard: For regions with significant changes in geometric characteristics such as hole segments, a special identification and processing mechanism is developed. The hole segment identification uses multiple quantization indicators: Rate of curvature change:

[0090] Rate of change of elliptical axis ratio:

[0091] Rate of change of cross-sectional area:

[0092] The region that meets all these conditions is identified as a region with significant changes in geometric characteristics that requires special processing.

[0093] Local re-optimization process: Local re-optimization is performed on the marked mutation points, with a 5-section optimization window. The specific process is as follows: Determine the optimization window: Take the marked point as the center and take its front and back 2 sections to form a local area of 5 sections Adjust the constraint weight: Adjust the smoothness constraint weight λ1 to 0.7 times the original value and the affine constraint weight λ2 to 0.5 times the original value to increase the flexibility of the model Construct a local optimization objective function:

[0094] Optimization solution: Use the Levenberg-Marquardt algorithm to optimize the objective function Parameter update: Only update the parameters of the marked point (the middle section of the window), and keep the parameters of other sections in the window unchanged Statistical model anomaly change identification: Establish a statistical model of parameter changes to identify abnormal changes with a confidence level of 97%. The specific implementation is as follows: Calculate the rate of change of adjacent sections for each parameter Calculate the mean and standard deviation of the rate of change sequence Set the threshold based on a 97% confidence level:

[0095] The position where the label change rate exceeds the threshold value is an abnormal change point Local weighted average processing: the abnormal area parameters are processed using the local weighted average method, and the weight decay rate is 0.7. For the parameter sequence , the parameter calculated after Gaussian smoothing processing is:

[0096] where the weight , and σ is the smoothing window parameter, with a value range of 2 to 5. This method realizes the smooth transition of parameters through weighted average, while retaining the main change trend of the parameters, avoiding the loss of details caused by excessive smoothing.

[0097] S5: final filtering and result output Based on the optimized elliptical cylinder model and dynamic threshold, the tunnel point cloud data is filtered, and the filtered point cloud data and elliptical cylinder model parameters are output.

[0098] Specifically, the points within the threshold range are retained as tunnel structure points (the label confidence is 98%), and the points exceeding the threshold are removed as non-structure points or noise. The model parameters are finally smoothed (the smoothing coefficient is 0.3), ensuring the overall continuity and consistency of the model. Finally, the filtered point cloud data and elliptical cylinder model parameters are output, with a processing accuracy of millimeter level, meeting the accuracy requirements of tunnel three-dimensional reconstruction and deformation analysis.

[0099] S5: based on the optimized elliptical cylinder model and dynamic threshold, the point cloud is finally filtered, and the filtered point cloud data and elliptical cylinder model parameters are output.

[0100] Specifically, the points within the threshold range are retained as tunnel structure points (the label confidence is 98%), and the points exceeding the threshold are removed as non-structure points or noise. The model parameters are finally smoothed (the smoothing coefficient is 0.3), ensuring the overall continuity and consistency of the model. Finally, the filtered point cloud data and elliptical cylinder model parameters are output, with a processing accuracy of millimeter level, meeting the accuracy requirements of tunnel three-dimensional reconstruction and deformation analysis.

[0101] In actual application process, the parameter values can be adjusted according to the specific characteristics of different tunnels. For example, for ordinary highway tunnels, a larger cross-section interval (about 2m) and a more relaxed constraint threshold can be used; for high-speed railway tunnels, a smaller cross-section interval (about 1m) and a more stringent constraint threshold can be used to meet higher accuracy requirements. For old tunnels with severe deformation, the parameter change threshold can be appropriately relaxed to improve the adaptability of the model to abnormal deformation; for newly built tunnels, the geometric continuity constraint can be strengthened to obtain smoother model results.

[0102] The method of the present application has been verified in actual tunnel cases. Taking a certain highway tunnel as an example, the tunnel is about 1.2 kilometers long, the cross section is horseshoe-shaped, and point cloud data is obtained by using mobile three-dimensional laser scanning, and the total number of points is about 250 million. After using the method of the present application, the noise point recognition accuracy reaches 96.5%, the structure point retention rate reaches 98.7%, and the processing efficiency is improved by about 35% compared with the traditional method. At the same time, the geometric continuity of the model is significantly improved, and the change of adjacent section parameters is more smooth and reasonable, especially in special areas such as portal segments, the model can accurately adapt to the change of geometric shape, providing a reliable data basis for subsequent tunnel deformation analysis.

[0103] In a non-real-time processing environment and under the condition of sufficient computing resources, the method of the present application can also be extended and replaced. For example, more complex geometric models such as hyper-elliptical cylinders or B-spline surfaces can be considered to replace the elliptical cylinder model; M-estimate sample consistency (MSAC) or progressive consistent sampling method (PROSAC) can be used to replace RANSAC algorithm; genetic algorithm, particle swarm optimization or simulated annealing algorithm can be used to replace iterative least squares method; Bayesian optimization and other methods can also be introduced to automatically adjust various parameters in the algorithm. These alternative solutions may provide better performance under certain conditions, but generally increase the computational complexity, and the balance between accuracy and efficiency needs to be considered in practical applications.

[0104] In summary, the tunnel point cloud adaptive elliptical cylinder filtering method based on affine constraint provided by the present application significantly improves the accuracy, stability and adaptability to complex environments of tunnel point cloud filtering by introducing affine transformation constraint, multi-objective optimization and adaptive parameter constraint, providing a reliable data basis for tunnel three-dimensional reconstruction and deformation monitoring.

[0105] The above embodiments are only for illustrating the technical concept and characteristics of the present application, and the purpose is to enable those skilled in the art to understand the content of the present application and implement it, and cannot limit the protection scope of the present application. Any equivalent transformation or modification made according to the spirit and essence of the present application should be covered within the protection scope of the present application.

Claims

1. An affine constraint based tunnel point cloud adaptive elliptic cylinder filtering method, characterized in that, The method comprises the following steps: S1: extracting a tunnel center line based on tunnel point cloud data, and setting cross sections at fixed intervals along the center line axis; S2: initializing an elliptic cylinder model parameter on each cross section using a principal component analysis method to form an elliptic cylinder model; S3: setting an initial screening threshold for preliminary screening of abnormal points of distances of all points in the cross section to the elliptic cylinder surface, and dynamically iteratively estimating an optimal screening threshold using a RANSAC algorithm; S4: constructing a multi-objective optimization function comprising a data fitting error term, a parameter smoothing constraint term and an affine transformation constraint term, and optimizing the elliptic cylinder model parameter using an optimization algorithm until the elliptic cylinder model parameter converges; S5: filtering the tunnel point cloud data based on the optimized elliptic cylinder model and the dynamic threshold, and outputting the filtered point cloud data and the elliptic cylinder model parameter.

2. The method of claim 1, wherein the method is characterized by, In step S1, the least squares method or the principal curve algorithm is used to extract the tunnel center line from the preliminary denoised tunnel point cloud data, and the fixed interval range is 0.5-3 meters.

3. The method of claim 1, wherein the method further comprises: In step S2, the principal component analysis method is used to initialize the elliptic cylinder model parameter, which specifically includes: calculating the covariance matrix of the cross section point cloud, performing eigenvalue decomposition on the covariance matrix, determining the directions of the major axis and the minor axis of the elliptic cylinder model using the eigenvectors corresponding to the maximum eigenvalue and the second largest eigenvalue, and estimating the lengths of the major axis and the minor axis according to the ratio of the eigenvalues.

4. The affine constraint based tunnel point cloud adaptive elliptic cylinder filtering method according to claim 1, characterized in that, In step S3, the initial screening threshold is 20 cm, and the optimal screening threshold is the 95th percentile of the distances of the inner points.

5. The affine constraint based tunnel point cloud adaptive elliptic cylinder filtering method according to claim 1, characterized in that, In step S4, the multi-objective optimization function is wherein is the sum of squares of distances from the point cloud to the elliptic cylinder surface, is the sum of squares of differences of parameters of adjacent sections, is the sum of squares of differences between the affine transformation matrix between adjacent sections and the expected transformation matrix, and λ1 and λ2 are weight coefficients. The Levenberg-Marquardt algorithm is used to minimize the multi-objective optimization function, and the initial damping factor of the algorithm is 0.01, and the damping factor adjustment rate is 10; weighting factor In step S4, it also includes the step of dynamically adjusting the constraint threshold according to the local curvature and the deformation degree, which specifically includes: dynamically calculating the direction angle change threshold according to the local curvature, using the formula: 1 has a value range of 0.5-2.0, and the weighting factor dynamically calculating the axis length change threshold according to the local deformation degree, using the formula: 2 has a value range of 0.3-1.

5.

6. The affine constraint based tunnel point cloud adaptive elliptic cylinder filtering method according to claim 1, characterized in that, In step S4, it also includes the step of smoothing the affine transformation parameter, which specifically includes: calculating the affine transformation matrix and decomposing it into rotation, scaling and translation components; using Gaussian weight to smooth the transformation parameter; and applying cubic spline interpolation to realize parameter smooth transition. base threshold is 2°, maximum threshold is 8°; In step S4, it also includes the step of special processing for the region with significant changes in geometric characteristics, which specifically includes: calculating the local curvature change rate and marking the curvature abrupt change points; locally re-optimizing the marked points; establishing a statistical model of parameter changes to identify abnormal changes; and using a local weighted average method to process the parameters of the abnormal region. base threshold is 2%, the maximum threshold is 6%.

7. The method of claim 1, wherein the method further comprises: ​ 8. The method of claim 1, wherein the method is based on affine constraints. ​ 9. The affine constraint based tunnel point cloud adaptive elliptic cylinder filtering method according to claim 1, characterized in that, The affine transformation constraint term further comprises a constraint intensity dynamic control mechanism, and the constraint intensity is controlled by a region adaptive weight The control constraint intensity, and the calculation formula is Wherein The geometric change degree of the current section is represented, k is a sensitivity coefficient (5-10), and finally The constraint term is represented as .

Images