A brick-concrete building structure length calculation method based on point cloud

Abstract

The present application relates to a kind of brick-concrete building structure length calculation method based on point cloud, comprising the following steps: wall point cloud thickness calculation, building wall splitting, wall boundary point estimation, upper boundary fitting, structure length calculation.When calculating point cloud thickness, multiple neighborhoods are randomly selected, but only part of neighborhood points, part of neighborhood are retained;Based on wall point cloud thickness, building point cloud segmentation, wall point cloud redistribution is carried out;With single wall as basic unit, wall boundary point estimation is carried out;Based on upper boundary point, selected weight iterative least square method is introduced, and the line fitting is carried out;From upper boundary line, the region of interest is created on wall, and the structure length of wall is calculated in combination with the straight line where the vertical boundary point in the region is located.The present application can accurately calculate the structure length of brick-concrete building by selecting structure length as monitoring index, combining the advantages of point cloud monitoring and the geometric structure characteristics of brick-concrete building, and can better monitor the deformation of brick-concrete building based on structure length.

Description

Technical Field

[0001] This invention belongs to the field of building structure length calculation technology, and particularly relates to a method for calculating the length of brick-concrete building structures based on point clouds. Background Technology

[0002] Brick-concrete structures are commonly used in low-rise rural buildings. These structures have relatively simple geometry, and apart from doors and windows, their walls contain no decorative appendages, making them relatively economical. However, these buildings have a short lifespan and are susceptible to the effects of ground movement. Furthermore, rural areas are major mining regions, and the structural characteristics and the impact of mining frequently cause damage to brick-concrete buildings. Therefore, deformation monitoring of these buildings is of great significance in providing a basis for assessing their health status and damage level.

[0003] For deformation monitoring of brick-concrete buildings, traditional methods primarily use total stations, but this approach is time-consuming and labor-intensive, and can only acquire a few discrete points. With the development of surveying technology, new monitoring methods have gradually emerged, represented by Differential Interferometry (DINSAR), Unmanned Aerial Vehicle (UAV) Photogrammetry, and Terrestrial Laser Scanning (TLS). Compared to DINSAR and UAV Photogrammetry, TLS can represent surface information of objects with high precision and massive point cloud data, and is widely used.

[0004] In recent years, building deformation monitoring using TLS can be divided into two categories: a) monitoring combined with an established model, and b) monitoring directly based on point clouds. In category a, Juan et al. established a Building Information Model (BIM) and calculated displacement changes during two time periods by comparing the models. Raphael et al. also established a BIM and ultimately calculated displacement and stress changes under applied gravity. Lian et al. established a Finite Element Model (FEM) and calculated changes in principal stress and shear force in the wall under coal seam mining conditions. In category b, Batur et al. performed point cloud registration and filtering, and then calculated displacement using a cloud-to-cloud approach. Suchocki et al. divided the point cloud into strips, reduced the data volume using the optimal dataset method, and then detected wall cracks. Li et al. extracted discrete deformation from the point cloud and, combined with a settlement prediction model, calculated wall subsidence, horizontal movement, and tilting deformations. Yang et al. created a triangular mesh vector model based on point clouds and combined it with the alpha shape algorithm, ultimately obtaining the wall tilt and crack depth. Based on the above, it can be observed that type A monitoring requires the use of various commercial software when building the model, applying forces, and comparing models. Type B monitoring, on the other hand, is conducted at the point cloud level, allowing for more flexible selection of monitoring methods and indicators; this type of monitoring appears to be a developing trend.

[0005] In China, building monitoring indicators include subsidence, horizontal movement, tilt, curvature, and horizontal deformation. Subsidence and horizontal movement primarily affect the water level at which the building is located, without causing structural damage. In practical applications, due to the relatively small values ​​of tilt and curvature, horizontal deformation is the main factor causing structural damage. Calculating this value typically requires the placement of manual markers on the building. For damage caused by coal mining, calculations can also be performed in conjunction with mining conditions. However, both methods yield only a small amount of horizontal deformation data and require the registration of multiple observation data into the same coordinate system, increasing manpower and resources and introducing errors. When the structural damage caused by horizontal deformation reaches a certain level, the most direct phenomenon is the appearance of cracks in the walls, reflected geometrically as a change in the length of the wall structure. Therefore, this invention directly selects structural length as the monitoring indicator and, considering the advantages of point cloud-based monitoring and the geometric characteristics of brick-concrete buildings, proposes a point cloud-based method for calculating the structural length of brick-concrete buildings. Summary of the Invention

[0006] The purpose of this invention is to provide a method for calculating the length of brick-concrete building structures based on point clouds in order to solve the above-mentioned problems.

[0007] The present invention achieves the above objectives through the following technical solutions:

[0008] A method for calculating the length of a brick-concrete building structure based on point clouds includes the following steps:

[0009] S1: Obtain the point cloud of the building wall, preprocess the point cloud of the building wall, and obtain the processed point cloud of the building wall.

[0010] S2: Set multiple neighborhoods within the processed building wall point cloud, filter the point clouds within the multiple neighborhoods to obtain effective neighborhood points, and calculate the wall point cloud thickness value.

[0011] S3: Determine a first distance threshold and a second distance threshold based on the wall point cloud thickness value; perform building point cloud segmentation on the processed building wall point cloud based on the first distance threshold to obtain the component point cloud of each single wall; and redistribute the wall point cloud of each single wall based on the second distance threshold to obtain the redistributed wall point cloud.

[0012] S4: The wall boundary points are estimated by using the improved point cloud boundary estimation method based on normals.

[0013] S5: Process the boundary points of the wall to obtain the upper boundary line. Under the intersection constraint, use the weighted iterative least squares method to fit the upper boundary line to obtain the fitted upper boundary line.

[0014] S6: Based on the fitted upper boundary line, create a region of interest on the wall, and calculate the length of the wall structure according to the straight line where the vertical boundary point in the region of interest is located.

[0015] As a further optimization of the present invention, in step S1, the point cloud of the building wall is obtained, and the point cloud of the building wall is preprocessed to obtain the processed point cloud of the building wall. The specific process is as follows:

[0016] The building is scanned to obtain a point cloud of the building walls. The point cloud of the building walls is then subjected to noise reduction, thinning, and ground filtering to obtain a processed point cloud of the building walls.

[0017] As a further optimization of the present invention, in step S2, multiple neighborhoods are set within the processed building wall point cloud, and the point clouds within the multiple neighborhoods are filtered to obtain effective neighborhood points, and the wall point cloud thickness value is calculated. The specific process is as follows:

[0018] N points are randomly selected within the processed point cloud of the building wall, and each point is assigned a neighborhood with a radius of R.

[0019] Calculate the distance d from each point in the neighborhood to the plane containing that neighborhood. Points whose distance d is less than the distance error are considered valid neighbors. Calculate the point cloud thickness TH for each neighborhood. Nei As shown in the following formula:

[0020]

[0021] Where RMSE1 represents the distance mean square error of a single neighborhood; n and m are the number of original points and the number of effective neighborhood points in a single neighborhood, respectively; [d1,d2,…,d i ,…,d n ] and [d1,d2,…,d j ,…,d m ] represents the set of distances from the original point in a single neighborhood to the plane containing that neighborhood, and the set of distances from the effective neighboring points in a single neighborhood to the plane containing that neighborhood, respectively; d mean Represents the set [d1, d2, ..., d j ,…,d m The mean distance from a valid neighboring point to the plane containing the neighboring points;

[0022] The point cloud thickness TH of each neighborhood Nei The mean value is calculated to obtain the mean thickness of the neighborhood point cloud. Neighborhoods with a thickness error less than the mean value are considered valid neighborhoods. The wall point cloud thickness value TH is then calculated. Wall As shown in the following formula:

[0023]

[0024] Where RMSE2 represents the mean square error of the neighborhood point cloud thickness, and N and M are the original number of neighborhoods and the effective number of neighborhoods, respectively. and These are the point cloud thickness sets of the original neighborhood and the point cloud thickness sets of the effective neighborhood, respectively. For set The mean point cloud thickness of the effective neighborhood in the data.

[0025] As a further optimization of the present invention, in step S3, a first distance threshold and a second distance threshold are determined based on the wall point cloud thickness value. Based on the first distance threshold, the processed building wall point cloud is segmented to obtain the component point cloud of each single-sided wall. Based on the second distance threshold, the component point cloud of each single-sided wall is redistributed to obtain the redistributed wall point cloud. The specific process is as follows:

[0026] A first distance threshold and a second distance threshold are determined based on the wall point cloud thickness value, wherein the first distance threshold is greater than the wall point cloud thickness value and the second distance threshold is equal to the wall point cloud thickness value;

[0027] The processed building wall point cloud is projected onto the XOY plane. Based on the first distance threshold and the random sampling consensus algorithm, the point cloud projected onto the XOY plane is processed to obtain the component point cloud of a single wall. The component point cloud of the single wall is deleted from the processed building wall point cloud, completing one update of the data set. The random sampling consensus algorithm processing steps and the data set update steps are repeated until the component point cloud of each single wall is obtained.

[0028] A random sampling consensus algorithm is used to process the point cloud of each single-sided wall to obtain the point cloud of the single-sided wall; the difference set is calculated between the point cloud of the component of the single-sided wall and the point cloud of the single-sided wall to obtain the difference set point cloud; the distance from each point in the difference set point cloud to each single-sided wall is calculated; each point in the difference set point cloud whose distance to the single-sided wall is less than a second distance threshold is selected and assigned to the point cloud of the single-sided wall to obtain the assigned wall point cloud of each single-sided wall.

[0029] As a further optimization of the present invention, in step S4, an improved point cloud boundary estimation method based on normals is used to estimate the wall boundary points of the allocated wall point cloud, and the specific process is as follows:

[0030] The weighted iterative least squares method is used to calculate the point cloud of the allocated wall to obtain the model parameters of the plane. Based on the model parameters of the plane, the rotation angles between the allocated wall point cloud and the XOZ and XOY planes are calculated, as shown in the following formula:

[0031] θ=X A X B / |X A ||X B |;

[0032] Where θ is the rotation angle; || is the magnitude of the model parameter, where X A X B These are the pre-defined, known model parameters for plane A and plane B, respectively.

[0033] The wall point cloud after assignment is rotated to be parallel to the XOY plane based on the rotation angle, resulting in the rotated wall point cloud;

[0034] Select a point in the rotated wall point cloud, connect the point to its neighboring points, calculate the angle between adjacent lines, and arrange them according to the size of the angle.

[0035] Calculate the difference between adjacent included angles after arrangement, and take the point where the difference is greater than the preset angle threshold as the wall boundary point.

[0036] As a further optimization of the present invention, in step S5, the wall boundary points are processed to obtain the upper boundary line. Under the intersection constraint condition, the upper boundary line is fitted using the weighted iterative least squares method to obtain the fitted upper boundary line. The specific process is as follows:

[0037] The wall boundary points are divided into inner boundary points and outer boundary points, and the boundary line of the straight line where the outer boundary point is located is obtained, including the vertical boundary line and the horizontal boundary line. The horizontal boundary line adjacent to the roof is taken as the upper boundary line.

[0038] The upper boundary line is vertically projected onto the XOY plane, and the model parameters of the upper boundary line projected onto the XOY plane are solved using the weighted iterative least squares method.

[0039] Based on the model parameters of the upper boundary line on the XOY plane, the rotation angles between the wall plane containing the upper boundary line and the XOZ and YOZ planes are calculated respectively.

[0040] Project the upper boundary line horizontally onto the XOZ or YOZ plane with a smaller rotation angle, and use the weighted iterative least squares method to solve for the model parameters of the projected upper boundary line. Lines A and B of the upper boundary line on the two planes are represented by the following two line equations:

[0041] k1x+b1=y=(z-b2) / k2;

[0042] k3x+b3=y=(z-b4) / k4;

[0043] Where k1, k2, b1, b2 are the parameters of the equation expression of line A in the spatial coordinate system, and k3, k4, b3, b4 are the parameters of the equation expression of line B in the spatial coordinate system.

[0044] Solving the equations of two lines yields their intersection points. Let points A and B be the intersection points obtained from the solution. Then, the fitted upper boundary line is shown in the following equation:

[0045]

[0046] Among them, (x a ,y a ,z a ), (x b ,y b ,z b ( ) are the coordinates of the intersection points A and B of the boundary lines, respectively;

[0047] As a further optimization of the present invention, in step S6, based on the fitted upper boundary line, a region of interest is created on the wall, and the length of the wall structure is calculated according to the straight line where the vertical boundary point of the region of interest is located. The specific process is as follows:

[0048] Set the step size to H, and translate the fitted upper boundary line downwards along the wall plane in sequence. Use the fitted upper boundary line after each translation as the center line to define the region of interest with a height of 2H.

[0049] Filter out the boundary points within the region of interest that are located on the wall itself and the boundary points located on adjacent walls;

[0050] The boundary points located on the wall itself and the boundary points located on adjacent walls are merged, and the model parameters of the vertical boundary lines in the region of interest are calculated using the weighted iterative least squares method, as shown in the following formula:

[0051] k5x + b5 = k6y + b6 = z;

[0052] Where k5, k6, b5, and b6 are the equation expression parameters of the model parameters of the vertical boundary line in the spatial coordinate system;

[0053] Let the number of translations be I. The fitted upper boundary line after the I-th translation is shown in the following formula:

[0054] (xf) / g = y = (xh-HI) / l;

[0055] Based on the model parameters of the vertical boundary line and the fitted upper boundary line after the I-th translation, and assuming that points a and b are two intersection points within the region of interest, then the wall length is:

[0056]

[0057] Where L is the length of the wall structure at the current translation position; (x a ,y a ,z a ), (x b ,y b ,z b ) are the coordinates of the intersection point a and point b within the region of interest, respectively.

[0058] The beneficial effects of this invention are as follows:

[0059] This invention selects structural length as a monitoring indicator, combines the advantages of point cloud monitoring with the geometric characteristics of brick-concrete buildings, and can calculate the structural length of brick-concrete buildings very accurately. Based on the structural length, the deformation of brick-concrete buildings can be better monitored. Attached Figure Description

[0060] Figure 1 This is a flowchart of the method of the present invention;

[0061] Figure 2 This is a point cloud diagram of the building and a single wall surface according to the present invention;

[0062] Figure 3 This is a schematic diagram of the wall boundary point estimation of the present invention;

[0063] Figure 4 This is a diagram showing the solution of the boundary line model parameters of the present invention;

[0064] Figure 5 This is a schematic diagram illustrating the calculation of the wall structure length according to the present invention;

[0065] Figure 6 This is a schematic diagram of an embodiment of the present invention;

[0066] Figure 7 This is a flowchart of the deformation monitoring process of the present invention;

[0067] Figure 8 This is a schematic diagram of the wall point cloud thickness calculation results in an embodiment of the present invention;

[0068] Figure 9 This is a flowchart illustrating intermediate results of an embodiment of the present invention;

[0069] Figure 10 This is a diagram showing the structural length deformation in an embodiment of the present invention;

[0070] Figure 11 These are actual photos taken at the measurement site in an embodiment of the present invention;

[0071] Figure 12 This is a schematic diagram of the boundary estimation results in an embodiment of the present invention;

[0072] Figure 13 This is a schematic diagram illustrating the structural length error caused by the boundary point estimation method of this invention;

[0073] Figure 14 This is a schematic diagram illustrating the structural length error caused by the weighting function of this invention. Detailed Implementation

[0074] The present application will now be described in further detail with reference to the accompanying drawings. It should be noted that the following specific embodiments are only used to further illustrate the present application and should not be construed as limiting the scope of protection of the present application. Those skilled in the art can make some non-essential improvements and adjustments to the present application based on the above application content.

[0075] Example 1

[0076] like Figure 1-14As shown, a method for calculating the length of a brick-concrete building structure based on point clouds includes the following steps:

[0077] The wall point cloud thickness is calculated by randomly selecting multiple neighborhoods, but only retaining a portion of the neighborhood points and a portion of the neighborhoods to achieve robust calculation.

[0078] Building wall segmentation is achieved in two steps, based on the wall point cloud thickness and after determining two distance thresholds: building point cloud segmentation and wall point cloud redistribution.

[0079] Wall boundary point estimation was performed using a single wall as the basic unit. This estimation method is derived from the point cloud boundary estimation method based on normals, and only requires calculating the normals once and transforming the coordinate system once.

[0080] Upper boundary fitting: Based on the upper boundary point, the weighted iterative least squares method is introduced, and the line is fitted under the constraint that two adjacent upper boundary lines intersect at a point.

[0081] The structural length is calculated by creating a region of interest on the wall that can be translated downwards from the upper boundary line. After each translation, the structural length of the wall is calculated by combining the vertical boundary point within the region with the straight line.

[0082] This invention employs the Random Sample Consensus Algorithm (RANSAC) and the Weighted Iterative Least Squares (WILS) method multiple times in calculating the structure length; therefore, these are described as fundamental principles.

[0083] Random sampling consensus algorithm

[0084] When the dataset can be clearly described by a certain model and there are no outliers or gross errors, the least squares (LS) method can be used to solve for the model parameters, which can be described by formula (1).

[0085] X = (B T PB) -1 B T PL(1)

[0086] Where X is the matrix of solved model parameters; B is the matrix of coefficients of the observation equation; L is the matrix of observations; and P is the matrix of observation weights (Pij = 0, i ≠ j). Since there are no outliers or glitches, P is taken as the identity matrix.

[0087] In reality, datasets inevitably contain outliers, and separating normal values ​​objectively and effectively is challenging. RANSAC can address this issue once the model describing the data is determined. The separation process is as follows:

[0088] Step 1: Randomly select n values ​​from the dataset and solve the model parameters using formula (1).

[0089] Step 2: Calculate the geometric distance between the data values ​​and the model, and count the number of values ​​that are less than the distance threshold.

[0090] Step 3: Repeat the above parameter calculation and count. After reaching the required number of iterations, the data values ​​that are less than the distance threshold corresponding to the maximum number are considered normal values.

[0091] In addition, RANSAC can also solve for model parameters. In this case, the model parameters corresponding to the maximum number should be taken as the final result.

[0092] Weighted Iterative Least Squares Method

[0093] Similarly, the dataset also contains some outliers. These outliers affect both LS and RANSAC when solving for model parameters. WILS, however, can mitigate the impact of outliers by using a weighting function to change their weights in the adjustment process.

[0094] Taking the IGG weighting function as an example, its expression is shown in formula (3). The main process of solving the model parameters using WILS is as follows:

[0095] Step 1: Based on the entire dataset, solve for the model parameters using formula (1).

[0096] Step 2: Calculate the residuals for each data value using formula (2). Then update the observation weight matrix according to formula (3).

[0097] Step 3: Repeat the above steps to solve the model and update the weight matrix. After reaching the required number of iterations, use the final solved model parameters as the result.

[0098] V = BX - L(2)

[0099] Where V is the residual matrix of each data value.

[0100]

[0101] Among them, v i σi represents the residual obtained from solving for the i-th data value; s1 and s2 are threshold coefficients to ensure robustness; σ0 is the standard deviation of the residual obtained from solving for the i-th data value. Here, s1 is set to 1.0 and s2 is set to 2.5.

[0102] Wall point cloud thickness calculation

[0103] When a point cloud can be described by a planar model, the point cloud thickness can be represented by the distance from the point to the plane. Based on the building point cloud observed for the first time, N points were randomly selected, and a neighborhood with a radius of R was defined for each point to calculate the wall point cloud thickness.

[0104] It is important to note that not all points in each neighborhood are included in the point cloud thickness calculation. For each neighborhood, the distance from a neighboring point to the plane containing the neighborhood is calculated. Considering the presence of point cloud noise, only points whose distance is less than the distance error are retained. In this embodiment, the distance d from a point within the neighborhood to the plane containing that neighborhood is calculated. Points whose distance d is less than the distance error are considered valid neighboring points, and the point cloud thickness TH of each neighborhood is calculated. Nei As shown in the following formula:

[0105]

[0106] Where RMSE1 represents the distance mean square error of a single neighborhood; n and m are the number of original points and the number of effective neighborhood points in a single neighborhood, respectively; [d1,d2,…,d i ,…,d n ] and [d1,d2,…,d j ,…,d m ] represents the set of distances from the original point in a single neighborhood to the plane containing that neighborhood, and the set of distances from the effective neighboring points in a single neighborhood to the plane containing that neighborhood, respectively; d mean Represents the set [d1, d2, ..., d j ,…,d m The mean distance from a valid neighboring point to the plane containing the neighboring points;

[0107] The point cloud thickness TH of each neighborhood Nei The mean value is calculated to obtain the mean thickness of the neighborhood point cloud. Neighborhoods with a thickness error less than the mean value are considered valid neighborhoods. The wall point cloud thickness value TH is then calculated. Wall As shown in the following formula:

[0108]

[0109] Where RMSE2 represents the mean square error of the neighborhood point cloud thickness, and N and M are the original number of neighborhoods and the effective number of neighborhoods, respectively. and These are the point cloud thickness sets of the original neighborhood and the point cloud thickness sets of the effective neighborhood, respectively. For set The mean point cloud thickness of the effective neighborhood in the data.

[0110] Building wall dismantling

[0111] Walls are the basic units that make up a building and also the basic units for calculating structural length. In addition to the point cloud of the building walls, the scanner also acquires point clouds of other objects within the scanning field. Furthermore, due to the disordered nature of point clouds, the point clouds of each wall are arranged irregularly, collectively forming the building point cloud. Based on the above, the decomposition of building walls consists of two steps: building point cloud segmentation and wall point cloud redistribution.

[0112] Building point cloud segmentation

[0113] The wall is approximately perpendicular to the XOY plane in three-dimensional space. Figure 2 Taking a building consisting of four walls as an example, in this plane, the geometric features of a single wall are represented by straight lines, and the geometric features of the building are represented by multiple closed straight lines connected together. Based on this feature, RANSAC is introduced into building point cloud segmentation.

[0114] Specifically, the first step is to project the building point cloud onto the XOY plane so that it can be segmented into wall elements. The second step is to perform RANSAC; the points corresponding to the maximum number of elements represent the component point cloud of a single wall. The third step is to remove the component point clouds of each single wall from the building point cloud, completing one update of the dataset. The fourth step is to repeat RANSAC and dataset updates until the component point cloud of each single wall is obtained.

[0115] Among them, B, X, and L are constructed as follows Figure 2 As shown. To prevent the straight line projected onto the XOY plane from being divided into multiple lines, the distance threshold should be greater than the wall point cloud thickness.

[0116] Wall point cloud redistribution

[0117] Using a distance threshold greater than the point cloud thickness can prevent the straight line projected onto the XOY plane from being segmented into multiple lines. However, this can also lead to over-segmentation of the building's point cloud.

[0118] Based on the order of segmentation, over-segmentation can be divided into three categories:

[0119] a. The point cloud of the first segmented single wall includes the point cloud of that single wall, but there are point clouds of adjacent walls on both sides.

[0120] b. The point cloud of the single wall that was later segmented only includes the point cloud of that single wall, but the point clouds on both sides are missing.

[0121] c. Some walls have point clouds that include the point cloud of that single wall, but one side has point clouds of adjacent walls while the other side is missing. Therefore, the point clouds of the segmented single walls need to be redistributed.

[0122] The first step is to use RANSAC to obtain the point cloud of each individual wall. The second step is to calculate the difference between the point clouds of the individual walls and the point clouds of the individual walls. The third step is to calculate the distance from each point in the difference set to each individual wall. The fourth step is to assign the point to the point cloud of the wall whose distance is less than a distance threshold.

[0123] Among them, B, X, and L are constructed as follows Figure 3 As shown. At this point, the distance threshold should be equal to the wall point cloud thickness, and the same point can be assigned to multiple single-sided walls.

[0124] Wall boundary point estimation

[0125] Currently, point cloud boundary estimation methods based on normals are widely used. This method first establishes a 3D coordinate system using the plane normals formed by point Pi in the wall point cloud and its neighboring points. Then, within the XOY plane of this coordinate system, boundary points and non-boundary points are determined based on the geometric relationship between point Pi and its neighboring points. In this method, the neighborhood for calculating the normals is usually larger than the neighborhood for determining boundary points; a larger neighborhood results in better boundary estimation but slower estimation speed.

[0126] In three-dimensional space, given the model parameters of two planes, it is possible to rotate one plane to make it parallel to the other. Assume the model parameters of planes A and B are X and X, respectively. A X B Then the rotation angle can be expressed as in formula (4). At this time, based on the rotation angle, plane A can be rotated to be parallel to plane B.

[0127] θ=X A X B / |X A ||X B | (4)

[0128] Where θ is the rotation angle; || is the modulus of the model parameter.

[0129] Considering that all points in the single-sided wall point cloud are approximately located on the same plane, and in order to reduce the given parameters and improve the estimation speed, this paper makes some changes to the estimation method. The main change is that the wall point cloud is rotated to be parallel to the XOY plane before boundary estimation. This rotation process can be divided into two steps: I. Rotate the wall point cloud around the Z-axis to be parallel to the XOZ plane, II. Rotate the wall point cloud around the X-axis to be parallel to the XOY plane. This process can be described by formula (5).

[0130]

[0131] Where (x,y,z) and (x',y',z') are the coordinates of the point before and after the rotation, respectively; θ1 and θ2 are the rotation angles in steps I and II, respectively, that is, the angles with the XOZ and XOY planes.

[0132] For the split single-sided wall, the boundary point estimation process of this method can be described as follows:

[0133] The first step is to solve for the planar model parameters using WILS. The second step is to calculate the rotation angle between the wall point cloud and the XOZ and XOY planes using formula (4). The third step is to rotate the wall point cloud to be parallel to the XOY plane using formula (5). The fourth step is as follows... Figure 3 As shown, in the rotated wall point cloud, point Pi is connected to its neighboring points. Fifth step: Calculate the angle between adjacent connecting lines and arrange them according to their size. Sixth step: Calculate the difference between adjacent angles, and use point Pi with a difference greater than the angle threshold as the boundary point. Seventh step: Perform an inverse transformation on formula (5) to restore the wall point cloud to its original posture before rotation.

[0134] upper boundary line fitting

[0135] Wall boundary points can be divided into two categories: internal boundary points and external boundary points. Internal boundary points are generated by the presence of doors and windows or point clouds and are not involved in subsequent calculations. Boundary lines are the straight lines containing external boundary points and can also be divided into two categories: vertical boundary lines and horizontal boundary lines. The horizontal boundary line adjacent to the roof, i.e., the upper boundary line, is the basis for calculating the structural length of the building.

[0136] like Figure 4 As shown, observing the geometric features of the upper boundary line reveals that the upper boundary lines of adjacent walls share the same intersection point. Furthermore, outlier values ​​are inevitably present in the estimated boundary points, and both LS and RANSAC methods struggle to avoid their influence. Therefore, this method incorporates WILS and performs upper boundary line fitting under the constraint of the intersection point.

[0137] The first step is to project the upper boundary line vertically onto the XOY plane and solve the model parameters of the line using WILS. The second step is to calculate the rotation angle between the wall plane containing the upper boundary line and the XOZ or YOZ plane using formula (4). The third step is to project the upper boundary line horizontally onto the XOZ or YOZ plane with the smaller angle, and then solve the model parameters using WILS. Wherein, B, X, and L are constructed as follows... Figure 4 As shown. At this time, the upper boundary lines A and B can be described by formulas (6) and (7) respectively. Given the equations of the two lines, their intersection points can be solved. Assuming that points A and B are the two intersection points of the upper boundary lines, the upper boundary lines can be described by formula (8).

[0138] k1x+b1=y=(z-b2) / k2(6)

[0139] k3x+b3=y=(z-b4) / k4(7)

[0140]

[0141] Among them, (x a ,y a ,z a ), (x b ,y b ,z b The coordinates of the intersection points A and B of the boundary lines are shown below.

[0142] g = (xa-xb) / (ya-yb) is the denominator of the expression on the left side;

[0143] g=xa+(xa-xb) / (ya-yb) is the part after the numerator symbol on the left side of the equation;

[0144] The same applies to h / l.

[0145] Structural length calculation

[0146] To calculate the length of the wall structure, this method translates the upper boundary line downwards along the wall plane in increments of H. Then, using the upper boundary line after each translation as the center line, a region of interest with a height of 2H is defined.

[0147] like Figure 5 As shown, in each region of interest, the first step is to extract the boundary points within the region. Since boundary point estimation is performed on a per-wall basis, boundary points on the same vertical boundary are divided into two parts, located within the boundary points of the wall itself and its adjacent walls, respectively. The second step involves merging these divided boundary points and using WILS to estimate the vertical boundary line model parameters within the region of interest. Here, B, X, and L are constructed as follows... Figure 4 As shown, the vertical boundary line can be described by formula (9). At this time, the upper boundary line A after the I-th translation can be described by formula (10). Assuming that points a and b are two intersection points within a certain region of interest, the length of the wall can be calculated by formula (11).

[0148] k5x + b5 = k6y + b6 = z (9)

[0149] (xf) / g=y=(xh-HI) / l (10)

[0150]

[0151] Where L is the length of the wall structure at the current translation position; (x a ,y a ,z a ), (x b ,yb ,z b (a) and (b) are the coordinates of point a and point b within the region of interest, respectively.

[0152] The case study using this method was conducted in Bozhou City, China. Two observations were performed during the monitoring period, approximately 11 months apart. Scanner parameters and measurement times are as follows: Figure 6 As shown in (a). For Figure 6 (b) The point cloud is preprocessed with noise reduction, thinning, and ground filtering, according to Figure 7 The structural length of the building was calculated, and deformation analysis was performed.

[0153] In this method, N=40 was randomly selected for the wall point cloud thickness. Regarding the value of R, following a trial-and-error approach, R was initially increased in 5cm increments within the range of 5-20cm. Figure 8 As shown in (a)-(d), the point cloud thickness values ​​of each neighborhood and wall increase with the increase of R. However, after R reaches 15cm, this change is not significant. Therefore, R was increased twice more at 10cm intervals. At this point, observations were made. Figure 8 (a)-(f) show that when R reaches 15cm, the thickness of the retained neighborhood point cloud does not change significantly, mostly not exceeding 3cm. Furthermore, the wall point cloud thickness calculated from the neighborhood point cloud thickness also does not change significantly, remaining at approximately 1.7cm. Therefore, this paper uses the wall point cloud thickness of 1.7cm calculated with N=40 and R=15cm as the benchmark for determining the distance threshold.

[0154] Analysis of structural length and deformation results

[0155] In building wall segmentation, the distance thresholds for building point cloud segmentation and wall point cloud redistribution should be greater than and equal to the wall point cloud thickness, respectively. Based on the calculation results above, the distance thresholds are determined to be 4.0cm and 1.7cm, respectively, resulting in the wall segmentation results as follows. Figure 9 As shown in figure a. Due to tree obstruction near wall D, the acquired point cloud of this wall has patchy gaps, resulting in some missing data in the point cloud of wall D after segmentation. However, adjacent walls have been completely segmented and their structures are intact. Subsequently, boundary point estimation was performed for the point cloud of each individual wall. Figure 9 In the estimation results shown in b, the building outline formed by the boundary points of each wall is clear. Similarly, due to the influence of observation conditions, a few noise points appear at the boundary points of wall D, but the inner and outer boundary points are still clearly distinguishable and will not affect subsequent steps. Then, based on the upper boundary points, an upper boundary line was fitted. To visually display the fitting results, the intersection points of the fitted lines were calculated. Observation Figure 9c. It can be observed that the intersection point of the upper boundary line and the vertical boundary point lie exactly on the same straight line: on the vertical boundary, which reflects that the fitted line matches the actual situation. When calculating the length of the structure, with H = 0.25m as the height of the region of interest, the upper boundary line was translated 20 times from top to bottom. The intersection points and structural lengths calculated in each region of interest are as follows: Figure 9 As shown in d and Table 1.

[0156] According to formula (15), by subtracting the structural lengths of structures with the same number of translations, the deformation information of the building at that translation position can be obtained. To visually display the deformation, the results are plotted as follows: Figure 10 As shown.

[0157] L e =L2-L1(15)

[0158] Among them, L e For structural length deformation, L1 is the structural length at a certain translation position during the first observation in two consecutive observations, and L2 is the structural length at the same translation position during the second observation.

[0159] Table 1 Calculation results of structural length

[0160]

[0161]

[0162] Observation Table 1 and Figure 10 It can be observed that for each wall, the curve shape of the structural length between the two observations shows a high degree of similarity. The deformation at various heights of wall A ranges from 34 to 48 mm. Observing its deformation trend reveals that, centered on the displacement number 10, i.e., 2.5 m from the roof, the deformation tends to increase with distance from that point. This is further supported by on-site photographs. Figure 11 Looking at wall A, the crack width gradually decreases from the roof to the ground, but the number of cracks increases from one to two. The change in structural length reflected by the cracks is basically consistent with the calculation results in this paper. The deformation at various heights of wall B is small, ranging from 0 to 12 mm. No cracks were found on the surface of this wall during the on-site inspection. The deformation at various heights of wall C ranges from -18 to 12 mm, slightly greater than that of wall B. Figure 11 As shown in Figure b, a crack appears in the wall from top to bottom, while no other cracks appear in the wall itself. The deformation at various heights of wall D ranges from 0 to 36 mm. (Combined with...) Figure 6 As can be observed in section b, this wall contains the main lighting structures of the building, such as doors and windows, and has relatively little wall structure. Furthermore, the cracks appear in areas with less wall structure. Therefore, the deformation fluctuations at various heights of this wall are small, and it does not exhibit a deformation fluctuation trend similar to that of wall A.

[0163] Comparison of time and results for boundary point estimation

[0164] For wall boundary point estimation, an adaptation was performed based on the normal-based method. In this paper, the proposed method and the normal-based method were used for boundary estimation, and their time and results were compared. The estimation time and results obtained based on the wall point cloud from the first observation are shown in Table 2. Figure 12 As shown.

[0165] Table 2 Boundary Point Estimation Time

[0166]

[0167]

[0168] As shown in Table 2, the proposed method requires approximately 1 / 10 the time compared to the normal-based method. This is because, during the computation, the proposed method avoids calculating the neighborhood normals for each point and also eliminates the need to re-estimate the coordinates of the point cloud in each neighborhood coordinate system. Only one normal calculation and one coordinate system transformation are required. Figure 12 It can be seen that both the proposed method and the normal-based method can successfully estimate the wall boundary points. To further compare the boundary point estimation results, the absolute value of the error in the structural length at the same location is calculated based on the estimation results of the two methods. Figure 13 The results show that the absolute value of the error does not exceed 8 mm, and most remain within 0-5 mm. This demonstrates that the boundary points estimated by the two methods have high consistency, and the time required by the proposed method is significantly less than that of the normal-based method. The influence of the weighting function on the calculation of structural length.

[0169] Besides IGG in this invention, Huber and Hampel are also commonly used weighting functions. Unlike the IGG weighting function, which divides the data into three segments, the Huber weighting function shown in formula (16) divides the data into two segments, while the Hampel weighting function shown in formula (17) divides the data into four segments. In this method, to verify the influence of different weighting functions, the structural length was calculated based on the boundary points in the first observation. The structural lengths calculated using Huber and Hampel as weighting functions were compared with those calculated using IGG as weighting function, and the absolute value of the error of the structural length at the same location was obtained.

[0170] like Figure 14The results show that, in the absolute error between Huber and IGG, wall A has one location with an error of approximately 6mm. In the absolute error between Hampel and IGG, wall B has one location with an error of approximately 7mm. The absolute errors at the remaining locations are mainly distributed between 0-5mm. This is primarily due to the different structures of the three weighting functions. The second segment of the Huber weighting function is divided into two parts in the IGG weighting function, and the third segment of the IGG weighting function is divided into two parts in the Hampel weighting function. However, overall, the three weighting functions have a relatively small impact on the results, all remaining within 5mm.

[0171]

[0172] Where s is the threshold coefficient for ensuring robustness, which is taken as 1.5 in this paper.

[0173]

[0174] Among them, s1, s2, and s3 are threshold coefficients to ensure robustness, which are set to 1.5, 2.0, and 3.0 respectively in this paper.

[0175] Evaluation of the accuracy of structural length calculation

[0176] To assess accuracy, additional observations were conducted on building wall A on the day of the second observation. Based on the wall point cloud thickness, and with H = 0.25m as the region of interest height, the structural length of wall A was calculated. After acquiring the deformation between the two observations, the ratio of the absolute value of the deformation to the structural length (i.e., the deformation rate) was calculated, and the results are shown in Table 3.

[0177] The second and additional observations were two consecutive observations taken on the same day, and their true deformation value should be 0 mm. Therefore, the deformation values ​​in the table above can be considered as error values, and the deformation rate as the error rate, used to evaluate the accuracy of the results. Table 3 shows that the deformation values ​​between these two consecutive observations are distributed between -6 and 6 mm. This is mainly due to random errors such as scanner observation errors and multi-station data splicing errors. Combined with the fact that the deformation rate at each location does not exceed 0.05%, it can be shown that the deformation values ​​obtained by the method presented in this paper have high accuracy and can achieve millimeter-level precision.

[0178] Table 3 Accuracy Evaluation Results

[0179]

[0180]

[0181] For brick-concrete buildings, this method proposes a structural length calculation method and demonstrates its application in deformation monitoring. Through wall point cloud thickness calculation, building wall segmentation, wall boundary point estimation, upper boundary fitting, and structural length calculation, combined with case studies and discussion, the main conclusions are as follows:

[0182] (1) The deformation obtained by the structural length calculated by this method is basically consistent with the deformation of the wall cracks.

[0183] (2) In the boundary point estimation part, both the method and the normal-based method can successfully estimate the wall boundary points, and the method takes only 1 / 10 of the time. When the estimation result is the only variable, the absolute value of the calculated structural length error does not exceed 8mm, and most of them remain in the range of 0-5mm.

[0184] (3) In the weighted iterative least squares method, the weighting function is the only variable. The structural length calculated by IGG is used as the standard, and the absolute value of the error between the Huber and Hampel calculation results and the standard is mainly distributed in the range of 0-5mm.

[0185] (4) Assuming that the deformation is 0 mm between two consecutive observations with a very short interval, the error of the structural length calculated by this method is distributed between -6 and 6 mm, with an error rate of no more than 0.05%, which can achieve mm-level accuracy.

[0186] The embodiments described above are merely examples of several implementations of the present invention, and while the descriptions are relatively specific and detailed, they should not be construed as limiting the scope of the present invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the scope of protection of the present invention.

Claims

1. A point cloud based length calculation method for a masonry structure, characterized in that, The method comprises the following steps: S1: obtaining a building wall point cloud, preprocessing the building wall point cloud, and obtaining a processed building wall point cloud; S2: setting a plurality of neighborhoods in the processed building wall point cloud, filtering the point cloud in the plurality of neighborhoods, obtaining effective neighborhood points, and calculating a wall point cloud thickness value; S3: determining a first distance threshold value and a second distance threshold value based on the wall point cloud thickness value, performing building point cloud segmentation on the processed building wall point cloud based on the first distance threshold value, obtaining a component point cloud of each single-wall, and performing wall point cloud re-distribution on the component point cloud of each single-wall based on the second distance threshold value to obtain a distributed wall point cloud; S4: performing wall boundary point estimation on the distributed wall point cloud by using an improved normal-based point cloud boundary estimation method to obtain wall boundary points; S5: processing the wall boundary points to obtain an upper boundary line, and fitting the upper boundary line under intersection constraint by using a selected weight iterative least square method to obtain a fitted upper boundary line; In the step S5, the wall boundary points are processed to obtain an upper boundary line, and the upper boundary line is fitted under intersection constraint by using a selected weight iterative least square method to obtain a fitted upper boundary line, and the specific process is as follows: The wall boundary points are divided into inner boundary points and outer boundary points to obtain a boundary line of a straight line where the outer boundary points are located, including a vertical boundary line and a horizontal boundary line, and the horizontal boundary line adjacent to a roof is taken as the upper boundary line; The upper boundary line is vertically projected onto an XOY plane, and a model parameter of the upper boundary line projected onto the XOY plane is solved by using a selected weight iterative least square method; The model parameter of the upper boundary line on the XOY plane is used to calculate a rotation angle between a wall plane where the upper boundary line is located and an XOZ plane or a YOZ plane; The upper boundary line is horizontally projected onto the XOZ plane or the YOZ plane with a smaller rotation angle, and a model parameter of a straight line after the projection of the upper boundary line is solved by using a selected weight iterative least square method, and the upper boundary line on the two planes is represented by the following two straight line equations: ； ； where k1, k2, b1 and b2 are equation expression parameters of line A in a space coordinate system, and k3, k4, b3 and b4 are equation expression parameters of line B in the space coordinate system; The two straight line equations are solved to obtain intersection points, and points A and B are intersection points obtained by solving, and the fitted upper boundary line is as follows: ； wherein (x a , y a , z a ), (x b , y b , z b ) are the coordinates of the intersection point A and point B of the boundary line, respectively; S6: based on the fitted upper boundary line, an interest region is created on a wall, and a wall structure length is calculated according to a straight line where a vertical boundary point in the interest region is located; In the step S6, based on the fitted upper boundary line, an interest region is created on a wall, and a wall structure length is calculated according to a straight line where a vertical boundary point in the interest region is located, and the specific process is as follows: A step length H is set, the fitted upper boundary line is sequentially translated downward along a wall plane, and a region with a height of 2H is defined as an interest region with the fitted upper boundary line after each translation as a center line; The boundary points in the interest region that are located in the wall itself and the boundary points that are located in adjacent walls are screened out. The boundary points located in the wall body itself and the boundary points located in the adjacent wall body are combined, and the model parameters of the vertical boundary line in the region of interest are calculated by using the selected weight iterative least square method, as shown in the following formula: ； Wherein, k5, k6, b5, b6 are the equation expression parameters of the model parameters of the vertical boundary line in the spatial coordinate system; The fitting upper boundary line after the Ith translation is as shown in the following formula: ； Based on the model parameters of the vertical boundary line and the fitting upper boundary line after the Ith translation, the points a and b are two intersection points in the region of interest, and the wall length is: ； Wherein, L is the wall structure length of the current translation position; (x a , y a , z a ), (x b , y b , z b ) are the coordinates of the intersection point a and point b in the region of interest, respectively.

2. The method for calculating the length of a brick-concrete building structure based on point clouds according to claim 1, characterized in that: In the step S1, the building wall point cloud is obtained, the building wall point cloud is preprocessed, and the processed building wall point cloud is obtained, and the specific process is as follows: The building is scanned to obtain the building wall point cloud, the building wall point cloud is subjected to noise reduction, thinning and ground filtering to obtain the processed building wall point cloud.

3. The method for calculating the length of a brick-concrete building structure based on point cloud according to claim 1, characterized in that: In the step S2, a plurality of neighborhoods are set in the processed building wall point cloud, and the point clouds in the plurality of neighborhoods are filtered to obtain effective neighborhood points, and the wall point cloud thickness value is calculated, and the specific process is as follows: N points are randomly selected in the processed building wall point cloud, and each point is set to define a neighborhood with R as the radius; The distance d of the points in the neighborhood to the plane where the neighborhood is located is calculated, and the points with the distance d less than the mean error of distance are regarded as valid neighborhood points, and the point cloud thickness TH of each neighborhood is calculated Nei As shown in the following formula: ； wherein RMSE1 represents the mean square error of the distance of the single neighborhood; n and m represent the number of original points in the single neighborhood and the number of effective neighborhood points, respectively; [d1, d2, …, d i , …, d n ] and [d1, d2, …, d j , …, d m ] represent the set of distances from the original points in the single neighborhood to the plane of the neighborhood and the set of distances from the effective neighborhood points in the single neighborhood to the plane of the neighborhood, respectively; and d mean represents the mean of the distances from the effective neighborhood points in the set [d1, d2, …, d j , …, d m ] to the plane of the neighborhood. The point cloud thickness TH of each neighborhood Nei The mean value of the neighborhood point cloud thickness is calculated, and the neighborhood with a thickness error less than the mean value is taken as an effective neighborhood, and the wall point cloud thickness value TH is calculated Wall As shown in the following formula: ； Wherein, RMSE2 represents the neighborhood point cloud thickness mean square error, N and M are respectively the original neighborhood quantity and the number of effective neighborhoods; , , …, , …, ] and [ , , …, ,…, ] are respectively the original neighborhood point cloud thickness set and the effective neighborhood point cloud thickness set; TH meanNei is the mean value of the effective neighborhood point cloud thickness in the set[ , , …, , …, ].

4. The method for calculating the length of a brick-concrete building structure based on point clouds according to claim 1, characterized in that: In the step S3, the first distance threshold and the second distance threshold are determined based on the wall point cloud thickness value, the processed building wall point cloud is segmented based on the first distance threshold to obtain the constituent point cloud of each single wall, and the constituent point cloud of each single wall is re-distributed based on the second distance threshold to obtain the distributed wall point cloud, and the specific process is as follows: The first distance threshold and the second distance threshold are determined based on the wall point cloud thickness value, the first distance threshold is greater than the wall point cloud thickness value, and the second distance threshold is equal to the wall point cloud thickness value; The processed building wall point cloud is projected onto the XOY plane, the point cloud projected onto the XOY plane is processed based on the first distance threshold and the random sample consensus algorithm to obtain the constituent point cloud of the single wall, the constituent point cloud of the single wall is separated from the processed building wall point cloud to complete one update of the data set, and the random sample consensus algorithm processing step and the data set updating step are repeated until the constituent point cloud of each single wall is obtained; The random sample consensus algorithm is used to process the constituent point cloud of each single wall to obtain the single wall point cloud, and the difference set point cloud is obtained by subtracting the constituent point cloud of the single wall from the single wall point cloud, the distance of each point in the difference set point cloud to each single wall is calculated, and the points in the difference set point cloud whose distance to the single wall is less than the second distance threshold are distributed to the point cloud of the single wall to obtain the distributed wall point cloud of each single wall.

5. The method for calculating the length of a brick-concrete building structure based on point clouds according to claim 1, characterized in that: In the step S4, the improved normal-based point cloud boundary estimation method is used to estimate the wall boundary points of the distributed wall point cloud to obtain the wall boundary points, and the specific process is as follows: The selected weight iterative least square method is used to calculate the wall point cloud after distribution, to obtain the model parameters of the plane; the rotation angles between the wall point cloud after distribution and XOZ plane and XOY plane are calculated based on the model parameters of the plane, as shown in the following formula: ； where θ is the rotation angle; || is the model parameter length, where X A , X B are the model parameters of the preset known plane A and plane B, respectively; The wall point cloud after distribution is rotated to be parallel to XOY plane based on the rotation angles, to obtain the wall point cloud after rotation; A point is selected in the wall point cloud after rotation, the point and its neighborhood points are connected, the included angles between the adjacent lines are calculated, and the angles are arranged according to the size; The difference values between the adjacent included angles after arrangement are calculated, and the points with the difference values greater than the preset angle threshold are taken as the wall boundary points.

Images