Adjustment method based on massive point clouds

Abstract

The invention discloses an adjustment method based on massive point clouds, and the method comprises the following steps: S1, calibrating a distance measurement constant of a three-dimensional laser scanner, and measuring a center total station coordinate of a reference point of a prism target ball; s2, manually picking up the scanning point cloud of the area containing the prism target ball from the scanning point cloud of the three-dimensional laser scanner, and obtaining the center scanning coordinate of the reference point of the prism target ball by adopting a fixed ball radius fitting method; s3, converting all scanning point clouds of the three-dimensional laser scanner from a scanner coordinate system to a total station horizontal coordinate system by using seven-parameter coordinate conversion; s4, performing partition processing on the scanning point cloud converted to the horizontal coordinate system of the total station; s5, picking up partition control points; s6, performing overall adjustment on the partition control points; and S7, carrying out overall adjustment on the subarea scanning point cloud. The method has the advantages that the precision requirement of deformation monitoring can be met, and the technical barrier of massive point cloud adjustment is successfully broken through.

Description

Technical Field

[0001] The present invention relates to the technical field of point cloud processing, and in particular to an adjustment method based on massive point clouds. Background Art

[0002] As the core equipment in the field of engineering monitoring, the total station continues to play an important role in the field of precision single-point measurement with its electronic angle measurement, distance measurement and coordinate solution functions (typical short-distance measurement accuracy reaches ±0.5mm). With modern iterative upgrades, the equipment has achieved deep integration of the automatic target recognition (ATR) system and the wireless transmission module, and can perform all-weather unattended automated monitoring tasks. However, its measurement mode requires a preset reflective prism, and can only obtain discrete point information. It is difficult to construct the global deformation field of the structure. It is a sampling detection method, and it is difficult to accurately capture local subtle deformation features. In large-scale infrastructure monitoring, such as airport runway deformation monitoring or landslide displacement tracking, due to its point-by-point scanning characteristics, a single full-network measurement often takes more than 5 hours, and the timeliness is significantly constrained.

[0003] In contrast, 3D laser scanners demonstrate significant technical advantages in high-density spatial data acquisition. Capable of collecting millions of point clouds per second and combining them with ultra-large spatial coverage, these devices can instantly construct 3D reality models with millimeter-level accuracy. Their core monitoring value lies in global deformation analysis, enabled by millimeter-level point cloud density. This allows for both comprehensive assessment of surface deformation and the capture of local distortion characteristics.

[0004] As the theoretical foundation for total station monitoring data processing, least squares adjustment theory has consistently maintained significant algorithmic superiority in the field of static monitoring. However, when processing high-precision laser scanning data, which can contain billions of point clouds in a single scan, traditional adjustment models not only face high-dimensional computational difficulties but also carry the potential risk of memory overflow. Summary of the Invention

[0005] The purpose of the present invention is to provide an adjustment method based on massive point clouds in accordance with the deficiencies of the above-mentioned prior art. The method uses multiple prism target balls to construct a spatial monitoring benchmark network. First, the three-dimensional coordinate data of the benchmark points are collected with high precision by a total station. Then, the scanning data of a single or multiple survey stations are rigorously adjusted based on the least squares method. Then, by comparing the coordinates after multiple adjustments, the displacement of the target point and the deformation characteristics of the overall point cloud are accurately obtained. However, in the face of massive scanning monitoring point clouds, traditional adjustment methods have significant computational bottlenecks. Directly applying the least squares algorithm will not only lead to an exponential increase in computational time, but will also cause memory resource overload. To this end, the present method constructs a hierarchical adjustment model for point cloud regional adjustment and point cloud control within the region. Computational units are established through spatial grid segmentation, and a set of feature control points is selected in each partition. In the first stage, the joint adjustment of the benchmark points and the partition control points is completed to obtain the coordinates of the partition control points after adjustment. In the second stage, regional point cloud coordinate correction is implemented based on the partition control points. This hierarchical adjustment strategy effectively implements the "overall control-local refinement" computing architecture, reducing the scale of point cloud data calculation by two orders of magnitude, and successfully breaking through the scale and efficiency barriers of least squares adjustment of massive point clouds.

[0006] The purpose of the present invention is achieved by the following technical solutions:

[0007] A massive point cloud-based adjustment method comprising the following steps:

[0008] S1: Calibrate the distance measurement constant of the 3D laser scanner and determine the center total station coordinates of the prism target sphere reference point;

[0009] S2: manually picking up a scanning point cloud of an area containing a prism target ball from the scanning point cloud of the three-dimensional laser scanner, and obtaining the central scanning coordinates of a reference point of the prism target ball by using a fixed sphere radius fitting method;

[0010] S3: using seven-parameter coordinate transformation, transforming all scanning point clouds of the three-dimensional laser scanner from the scanner coordinate system to the total station horizontal coordinate system;

[0011] S4: partitioning the scanned point cloud converted to the horizontal coordinate system of the total station;

[0012] S5: Pick up partition control points;

[0013] S6: Perform overall adjustment on the partition control points;

[0014] S7: Perform overall adjustment on the partitioned scan point cloud.

[0015] In step S1, the scan point cloud of the area containing the prism target ball is

[0016] The fixed sphere radius fitting method comprises the following steps:

[0017] The radius R of the prism target ball is taken as a known value, and the coordinates of the best fitting sphere center (C x ,C y ,C z ), so that the sum of squares of the residuals from the point to the sphere is minimized:

[0018]

[0019] Among them, the single point residual formula is:

[0020]

[0021] The fitting adopts the nonlinear least squares method, which includes the following steps:

[0022] a. Linearized residual:

[0023] Residual function in C (k) Taylor expansion:

[0024]

[0025] Where, the Jacobian matrix is:

[0026]

[0027] Where,

[0028] b. Construct a linear system:

[0029] The residual equations of all points are combined to form a linear system of equations:

[0030] J (k) ΔC=L (k) Formula 6;

[0031] The least squares solution is:

[0032]

[0033] Where:

[0034]

[0035] c. Gross error elimination:

[0036] Residual standard deviation estimate:

[0037]

[0038] Aspherical point elimination criteria:

[0039]

[0040] In step S3, the method for converting all scan point clouds of the three-dimensional laser scanner into the horizontal coordinate system of the total station includes the following steps:

[0041] Assume that o-xyh is the horizontal coordinate system of the total station, o'-x'y'h' is the scanner coordinate system, the number of point clouds is n1, the number of control points is n, and the horizontal coordinates of the total station and the scanner coordinates of the control points are expressed as (x i y i h i ) T and (x' i y' i h' i ) T , where (i=1,2,…,n);

[0042] The spatial coordinate transformation model is expressed as:

[0043]

[0044] Where (x0 y0 h0) is the translation, k is the scale variable, R1(α), R2(β), and R3(γ) are the rotation matrices, which are:

[0045]

[0046] When the number of common points n is greater than 3, the above seven parameters are determined by least squares;

[0047] For each common point, the error equation is listed as:

[0048]

[0049] Where, V=(V x V y V h ) T is the correction number;

[0050] Take approximate values x0 for seven parameters 0 、y0 0 、h0 0 , α 0 , β 0 , γ 0 、k 0 , linearize Equation 11:

[0051]

[0052] Where, the partial derivative and constant term are:

[0053]

[0054]

[0055] The normal equation is composed of Equation 12, which is iterated until convergence to obtain the translation (x0 y0 h0), rotation R1(α), R2(β), R3(γ), and scale value k. Then, according to Equation 10, all scanned point clouds are converted from the scanner coordinate system to the total station horizontal coordinate system.

[0056] In step S4, set N desired partitions;

[0057] The partition processing method of the scan point cloud converted to the horizontal coordinate system of the total station includes the following steps:

[0058] a. Calculation formula for the number of cubes:

[0059]

[0060] In the formula, the length, width, and height of the given space are L, W, and H respectively, and the side length of the cube is a. represents the rounding-up operation to ensure that the cube completely covers the space;

[0061] b. Maximum edge length search:

[0062] When the number of target cubes is N desired When , find the maximum side length a max , so that:

[0063]

[0064] Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice;

[0065] c. Minimum edge length search:

[0066] Find the minimum side length a min , so that:

[0067]

[0068] Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice;

[0069] d. Calculation of theoretical side length:

[0070] Based on the volume equipartition assumption, the calculation formula for the theoretical side length a0 is:

[0071]

[0072] e. Calculation of optimal side length:

[0073] The set of candidate edge lengths is:

[0074] {a0,a min ,a max};

[0075] By comparing the number of cubes corresponding to each candidate side length with N desired The absolute difference of * :

[0076] |NN desired |;

[0077] If the differences are the same, choose the side length closest to a0;

[0078] f. Calculation of the actual number of cubes:

[0079] The final actual number of cubes N actual , through the optimal side length a * calculate:

[0080]

[0081] In step S5, the method for picking the partition control points includes the following steps:

[0082] Traverse all points and assign them to cubes, then traverse each cube point in turn to obtain the minimum bounding box of the cube points, and take the 8 corner points and the center point of the minimum bounding box, a total of 9 points, as partition control points.

[0083] In step S6, the overall adjustment method of the partition control points includes the following steps:

[0084] After partitioning and obtaining the partition control points, the target sphere reference point is used as the starting point, and all partition control points participate in the overall least squares adjustment together; during the adjustment operation, independent adjustment calculations are performed on the plane coordinates and elevation data respectively; the horizontal angle, vertical angle and slope distance corresponding to each partition control point can be solved using the survey station coordinates and the three-dimensional coordinates of the partition control points after conversion to the horizontal coordinate system of the total station.

[0085] In step S6, using plane adjustment, the observation error equations of horizontal distance and horizontal angle can be listed respectively;

[0086] The horizontal distance correction equation is: for the horizontal distance observation value S from point i to any point j measured at the i measuring station, ij :

[0087]

[0088]

[0089] (δE' ij δN' ij ) T is the coordinate of point j in the horizontal coordinate system of station i (E' ij N' ij ) T number of corrections;

[0090] Horizontal angle correction equation, for the direction observation value R from point i to any point j measured at station i ij :

[0091]

[0092] Where θ i is the orientation angle parameter of station i;

[0093] Write the error equation in matrix form:

[0094] V = AδX-L Equation 22;

[0095] The error equation coefficient matrix A, the unknown correction factor δX, and the error equation constant term L in Equation 22 are:

[0096]

[0097] δX=(δx1 δx2 … δx n ) T ,

[0098] L=(l1 l2 … l n ) T Formula 23;

[0099] In V T Under the condition of PV=min, solve equation 22;

[0100] Depend on have to:

[0101] NδX=Cformula 24;

[0102] Equation 24 is called the normal equation, and the normal equation coefficient matrix N = A T PA, the constant term of the normal equation C = A T PL, solve Equation 24:

[0103] δX=N -1 C type 25;

[0104] The unknown number is approximated by X 0 Corrected to X 0+δX, recompose the error equation and the normal equation, solve δX, correct the approximation, and iterate until convergence, that is, the correction number δX is small enough to obtain the solution of the unknown number.

[0105] In step S7, the overall adjustment method of the partitioned scanning point cloud includes the following steps:

[0106] Based on the coordinates of the partition control points before and after adjustment, the spatial coordinate transformation model of Equation 10 is used to solve the seven-parameter transformation parameters. After the partition scan point cloud is uniformly transformed into the adjustment coordinate system through the spatial coordinate transformation model of Equation 10, the overall adjustment result of the global scan point cloud is obtained.

[0107] The advantages of the present invention are that it can not only meet the accuracy requirements of deformation monitoring, but also successfully break through the technical barriers of massive point cloud adjustment. BRIEF DESCRIPTION OF THE DRAWINGS

[0108] Figure 1 Schematic diagram of the process of the adjustment method based on massive point clouds of the present invention;

[0109] Figure 2 Schematic diagram of the partition bounding box of the present invention;

[0110] Figure 3 This is a schematic diagram of the partition control points of the present invention;

[0111] Figure 4 This is a table diagram of the scanning configuration parameters of the present invention;

[0112] Figure 5 This is a coordinate table diagram of the target sphere total station of the present invention;

[0113] Figure 6 This is the on-site scanning point cloud interface diagram of the present invention (1);

[0114] Figure 7 This is the on-site scanning point cloud interface diagram of the present invention (2);

[0115] Figure 8 This is the point cloud selection diagram (1) of the target ball of the present invention;

[0116] Figure 9 This is the point cloud selection diagram (2) of the target ball of the present invention;

[0117] Figure 10 This is a table diagram of the target ball scanner fitting coordinates of the present invention;

[0118] Figure 11 This is a table diagram showing the conversion of target sphere scanner coordinates to total station coordinates according to the present invention;

[0119] Figure 12 This is a schematic diagram of point cloud partitioning in the present invention;

[0120] Figure 13 This is the point spacing distribution map before and after the overall adjustment of the first phase of the subarea of the present invention;

[0121] Figure 14 It is a table diagram of the single-weighted mean error of the overall adjustment of each period of the present invention;

[0122] Figure 15 Figure (1) is a comparison of the regions with the same name between the first and second phases of the present invention;

[0123] Figure 16 Figure (2) is the result of comparing the regions with the same name in the first and second phases of the present invention;

[0124] Figure 17 Figure (3) is the result of comparing the regions with the same name in the first and second phases of the present invention;

[0125] Figure 18 Figure (4) is the result of comparing the regions with the same name in the first and second phases of the present invention;

[0126] Figure 19 This is a table showing the average distance and standard deviation of the point cloud of the same-name area in each period of the present invention;

[0127] Figure 20 This is a table comparing the coordinates of the target sphere detection points after adjustment in each period of the present invention with the coordinates measured by the total station. DETAILED DESCRIPTION

[0128] The features of the present invention and other related features are further described in detail below through embodiments in conjunction with the accompanying drawings to facilitate understanding by those skilled in the art:

[0129] Example: Figure 1 As shown, this embodiment relates to an adjustment method based on massive point clouds, which mainly includes the following steps:

[0130] S1: Calibrate the ranging constant of the 3D laser scanner and determine the center total station coordinates of the prism target sphere reference point.

[0131] Specifically, the scanning and monitoring process requires the reference coordinates of the target sphere's reference point, which are measured using a total station. The scanning and monitoring point cloud is measured using a scanner, which uses a non-cooperative target point cloud method. Therefore, the scanner's ranging constant and multiplication constant must be measured. This calibration can be performed by a professional metrology organization or by independently constructing a standardized calibration field. The determination of this ranging constant is not the focus of this case study; all subsequent slant distance data is assumed to have been corrected for the ranging constant.

[0132] In a three-dimensional scanning monitoring system, the benchmark point is an indispensable core element, and its stability directly determines the reliability of the monitoring results. To ensure the accurate conversion between the scanner coordinate system and the total station's horizontal coordinate system, a special prism target ball must be used as the spatial benchmark. Before monitoring is implemented, the target device must be permanently fixed in a stable area, with no less than three benchmark points, and they must be distributed as evenly as possible in space. After the deployment is completed, the total station is used to perform multiple rounds of precise coordinate acquisition of the target ball prism center, thereby establishing a high-precision benchmark control network. The establishment of the total station coordinate system of the benchmark point must be completed before formal monitoring.

[0133] S2: Manually pick the scanning point cloud of the area containing the prism target ball from the scanning point cloud of the 3D laser scanner, and use the fixed sphere radius fitting method to obtain the central scanning coordinates of the prism target ball reference point.

[0134] Specifically, after the determination of the monitoring benchmark network is completed, a three-dimensional laser scanner needs to be regularly set up to implement periodic monitoring operations. The specific station setting method can be selected according to the site conditions, such as the fixed station method or the arbitrary station method. Before implementing the scan, it is necessary to comprehensively consider the monitoring needs and accuracy requirements, and reasonably configure the core parameters such as the scanning resolution, quality level and scanning coverage. During the data collection process, special attention must be paid to ensure that the scan data contains at least 3 or more calibrated high-precision prism target ball reference points, otherwise the adjustment solution will not be able to be performed normally. After the scan is completed, all point cloud data for that period are obtained.

[0135] After the first scan is completed and the scan point cloud is obtained, you need to manually pick the approximate range of the fiducial target sphere and use the fixed sphere radius fitting method to obtain the scanning coordinates of the center of the fiducial target sphere. Later, the target sphere can be automatically identified based on the initial target sphere range.

[0136] After a single scan is completed, the scan point cloud is obtained. After manually selecting an area containing the reference target sphere, the point cloud of the area is obtained, that is, the three-dimensional point cloud data is given. In the fitting process, the radius R of the sphere is taken as a known value to solve the coordinates of the best fitting sphere center (C x ,C y ,C z ), so that the sum of squares of the residuals from the point to the sphere is minimized:

[0137]

[0138] Among them, the single point residual formula is:

[0139]

[0140] The fitting was performed using a nonlinear least squares method (Gauss-Newton method), which included the following steps:

[0141] a. Linearized residual:

[0142] Residual function in C (k) Taylor expansion:

[0143]

[0144] Where, the Jacobian matrix is:

[0145]

[0146] Where,

[0147] b. Construct a linear system:

[0148] The residual equations of all points are combined to form a linear system of equations:

[0149] J (k) ΔC=L (k) Formula 6;

[0150] The least squares solution is:

[0151]

[0152] Where:

[0153]

[0154] c. Gross error elimination:

[0155] Residual standard deviation estimate:

[0156]

[0157] Aspherical point elimination criteria:

[0158]

[0159] S3: Using the seven-parameter coordinate transformation, all scanning point clouds of the 3D laser scanner are transformed from the scanner coordinate system to the total station horizontal coordinate system.

[0160] Specifically, assuming that o-xyh is the total station horizontal coordinate system, o'-x'y'h' is the scanner coordinate system, the number of point clouds is n1, the number of control points is n, and the total station horizontal coordinates and scanner coordinates of the control points are expressed as (x i y i h i ) T and (x' i y' i h' i ) T, where (i=1,2,…,n);

[0161] The spatial coordinate transformation model (the relationship between the total station horizontal coordinates and the scanner coordinates) is expressed as:

[0162]

[0163] Where (x0 y0 h0) is the translation, k is the scale variable, R1(α), R2(β), and R3(γ) are the rotation matrices, which are:

[0164]

[0165]

[0166] When the number of common points n is greater than 3, the above seven parameters are determined by least squares;

[0167] For each common point, the error equation is listed as:

[0168]

[0169] Where, V=(V x V y V h ) T is the correction number (conversion residual);

[0170] Take approximate values x0 for seven parameters 0 、y0 0 、h0 0 , α 0 , β 0 , γ 0 、k 0 , linearize Equation 11:

[0171] Where, the partial derivative and constant term are:

[0172]

[0173] The normal equation is composed of Equation 12, which is iterated until convergence to obtain the translation (x0 y0 h0), rotation R1(α), R2(β), R3(γ), and scale value k. Then, according to Equation 10, all scanned point clouds are converted from the scanner coordinate system to the total station horizontal coordinate system.

[0174] S4: Partition processing is performed on the scanned point cloud converted to the horizontal coordinate system of the total station.

[0175] Specifically, in order to reduce the amount of least squares adjustment calculation, the point cloud after conversion to the horizontal coordinate system must be processed in different regions. The entire scanning range is processed in different regions, and N is set. desired partitions;

[0176] The partition processing method of the scan point cloud converted to the horizontal coordinate system of the total station includes the following steps:

[0177] a. Calculation formula for the number of cubes:

[0178]

[0179] In the formula, the length, width, and height of the given space are L, W, and H respectively, and the side length of the cube is a. represents the rounding-up operation to ensure that the cube completely covers the space;

[0180] b. Maximum edge length search:

[0181] When the number of target cubes is N desired When , find the maximum side length a max , so that:

[0182]

[0183] Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice;

[0184] c. Minimum edge length search:

[0185] Find the minimum side length a min , so that:

[0186]

[0187] Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice;

[0188] d. Calculation of theoretical side length:

[0189] Based on the volume equipartition assumption, the calculation formula for the theoretical side length a0 is:

[0190]

[0191] e. Calculation of optimal side length:

[0192] The set of candidate edge lengths is:

[0193] {a0,a min ,a max};

[0194] By comparing the number of cubes corresponding to each candidate side length with N desired The absolute difference of * :

[0195] |NN desired |;

[0196] If the differences are the same, choose the side length closest to a0;

[0197] f. Calculation of the actual number of cubes:

[0198] The final actual number of cubes N actual , through the optimal side length a * calculate:

[0199]

[0200] S5: Pick the partition control point.

[0201] Specifically, if Figures 2-3 As shown, all points are traversed and assigned to each cube, and then each cube point is traversed in turn to obtain the minimum bounding box of the cube point. The 8 corner points (point "1" / point "p1" to point "8" / point "p8") and the center point (point "0") of the minimum bounding box are taken as the partition control points.

[0202] S6: Perform overall adjustment on the partition control points.

[0203] Specifically, after partitioning and obtaining the partition control points, the target sphere reference point is used as the starting point, and all partition control points participate in the overall least squares adjustment together; during the adjustment operation, independent adjustment calculations are performed on the plane coordinates and elevation data respectively; the horizontal angle, vertical angle and slope distance corresponding to each partition control point can be solved using the survey station coordinates and the three-dimensional coordinates of the partition control points after conversion to the horizontal coordinate system of the total station.

[0204] By using plane adjustment, the error equations of the horizontal distance and horizontal angle observation values can be listed respectively;

[0205] The horizontal distance correction equation is: for the horizontal distance observation value S from point i to any point j measured at the i measuring station, ij :

[0206]

[0207] (δE' ij δN' ij ) T is the coordinate of point j in the horizontal coordinate system of station i (E' ij N' ij ) T number of corrections;

[0208] Horizontal angle correction equation, for the direction observation value R from point i to any point j measured at station i ij :

[0209]

[0210] Where θ i is the orientation angle parameter of station i;

[0211] Write the error equation in matrix form:

[0212] V = AδX-L Equation 22;

[0213] The error equation coefficient matrix A, the unknown correction factor δX, and the error equation constant term L in Equation 22 are:

[0214]

[0215] δX=(δx1 δx2 … δx n ) T ,

[0216] L=(l1 l2 …l n ) T Formula 23;

[0217] In V T Under the condition of PV=min, solve equation 22;

[0218] Depend on have to:

[0219] NδX=C Equation 24;

[0220] Equation 24 is called the normal equation, and the normal equation coefficient matrix N = A T PA, the constant term of the normal equation C = A T PL, solve Equation 24:

[0221] δX=N -1 C type 25;

[0222] The unknown number is approximated by X 0 Corrected to X 0 +δX, recompose the error equation and the normal equation, solve δX, correct the approximation, and iterate until convergence, that is, the correction number δX is small enough to obtain the solution of the unknown number.

[0223] Height adjustment is the same as plane adjustment, so I will not go into details here.

[0224] S7: Perform overall adjustment on the partitioned scan point cloud.

[0225] Specifically, based on the coordinates of the partition control points before and after adjustment, the spatial coordinate transformation model of Formula 10 is used to solve the seven-parameter transformation parameters. After the partition scanning point cloud is uniformly transformed into the adjustment coordinate system through the spatial coordinate transformation model of Formula 10, the overall adjustment result of the global scanning point cloud is obtained.

[0226] In this example, a factory building was selected as the test site. The length and width of the entire factory building are 102 meters and 48 meters respectively. The 3D laser scanner used in the test is Faro Focus 3D S350. The parameters of the acquisition configuration are shown in Figure 4 .

[0227] The test scanned the factory area three times, with the scanner gantry moved for each scan. The scans were performed at 1 / 2 resolution and 2X quality, with the compass, GPS, and tilt compensation sensors turned off.

[0228] Four prism target balls are evenly distributed within the factory, three of which are used as reference points and one as a monitoring point to verify its measurement accuracy. The coordinates of the prism target balls are collected using the total station multi-round measurement method. Figure 5 .

[0229] After measurement, the distance constant difference between the scanner and the total station is +8mm. Before adjustment, all scanning point clouds were corrected for the distance constant.

[0230] After collecting the three-site cloud, first use PotreeConverter to build a multi-resolution octree for the point cloud, and use Potree to browse the relevant point cloud. The browsing interface is as follows: Figures 6-7 .

[0231] like Figures 8-9 As shown in the figure, the point cloud can be browsed in the pre-programmed web software. By selecting the cube frame, the spherical point cloud of the reference point (prism target ball) can be selected respectively, and the scanning coordinates of the center of the ball can be obtained by the fitting algorithm of the target ball. Figure 10 As shown in the figure, by fixing the radius of the target sphere to 5 cm, the coordinates of the three reference target points, the corresponding number of valid points and the fitting accuracy value were recorded for three scans.

[0232] After obtaining the target sphere scanning fitting center coordinates and the total station measurement center coordinates, a seven-parameter conversion algorithm is used to convert all point clouds including the target sphere from the scanner frame to the total station frame. Figure 11 The fitting coordinates of the target ball scanner in the three scans are converted into coordinates to obtain the coordinates under the total station frame.

[0233] Since the collected point cloud is arranged in a certain order and thinned out in a certain ratio, a sparse point cloud representing the entire area is obtained. The minimum bounding box of the entire point cloud is calculated, and the cube is partitioned using the scanning point partitioning algorithm. The control points of each area are obtained using the partition control point algorithm. In this experiment, the number of cube partitions is 9. According to the point cloud partition control point picking algorithm, the minimum bounding box of each partition point is traversed, and the relevant control points are obtained for each partition in turn. The top view of the 9 areas is shown as follows: Figure 12 Each region has 9 pick points. The red box in the figure is the middle point of the minimum bounding box of each region. There are 81 regional control points in total.

[0234] Using three reference points and 81 points in 9 areas, we calculate the horizontal angle, vertical angle and slope distance of all points and perform overall adjustment. The difference in the distance between the 81 points before and after adjustment is as follows: Figure 13 , all within 2.5mm.

[0235] The overall adjustment horizontal unit weighted mean error is and the vertical unit weighted mean error is Figure 14 .

[0236] The above process data shows that the fixed diameter method is highly accurate in fitting the target sphere center, with the accuracy of the three fittings exceeding 0.0007mm. Furthermore, the residuals after converting the scanned fitting coordinates to the total station coordinates are relatively small, with the maximum residual at one point being 1.6mm and the minimum residual at one point being 0.1mm.

[0237] To verify the internal coincidence accuracy of the scans, four identically named areas were manually selected within the scan area using data from three repeated scans. These areas were located at 50, 40, 30, and 20 meters from the survey station, respectively. CloudCompare's point cloud comparison algorithm was used to calculate the difference between any two point clouds. Figures 15-18 This is a comparison of the same-named regions in Phase 1 and Phase 2. The figure records the point spacing, color mapping distance, and distance change point histogram of the two scanned regions. Figure 19 The average value and standard deviation of the distance between point clouds of the same name in the scanning area are obtained by comparing the point clouds of the same name in each of the 1st to 3rd periods.

[0238] from Figure 19 We can see that the standard deviation of repeated scans of the same point cloud increases with distance. At 50 meters, the maximum average distance between the three point clouds is 3.3mm, with a maximum standard deviation of 3.4mm. At 10 meters, the maximum average distance between the two point clouds is 1.8mm, with a maximum standard deviation of 1.1mm. Combining the three-period standard deviations and applying the error propagation theorem, the error in this scan measurement is 1.7mm.

[0239] Using the J4 target ball prism as the detection point, the difference between the scanning sphere center converted to the total station frame coordinates according to this algorithm and the total station's measured coordinates is checked to verify the external coincidence accuracy of this scanning monitoring.

[0240] When arranging the target sphere, a target sphere is reserved as a calibration point. The coordinates of the target sphere measured by the total station are compared with the coordinates after the point cloud adjustment in the third phase to obtain Figure 20 results.

[0241] By comparing the coordinates with those measured by the total station and calculating the difference, it can be calculated that the error in the difference is 0.6mm.

[0242] This demonstrates that the hybrid adjustment model constructed in this embodiment for massive point cloud data addresses the technical bottleneck encountered when using traditional least-squares adjustment for massive point clouds. Preliminary experiments have demonstrated that this model achieves an internal accuracy of ±1.7mm within a 50-meter observation range and an external accuracy of ±0.6mm in a 30-meter short baseline scenario.

[0243] In the subsequent research phase, more application tests will be carried out to build a test environment with a range of more than 150 meters, and systematically evaluate the model's three-dimensional deformation monitoring accuracy in large-scale scenes.

[0244] The beneficial technical effect of this embodiment is that it not only meets the accuracy requirements of deformation monitoring, but also successfully breaks through the technical barriers of massive point cloud adjustment.

[0245] Although the above embodiments have described in detail the concepts and embodiments of the present invention with reference to the accompanying drawings, ordinary technicians in this field can recognize that various improvements and modifications can still be made to the present invention without departing from the scope of the claims, so they are not described in detail here.

Claims

1. An adjustment method based on massive point clouds, characterized by The method comprises the following steps: S1: Calibrate the distance measurement constant of the 3D laser scanner and determine the center total station coordinates of the prism target sphere reference point; S2: manually picking up a scanning point cloud of an area containing a prism target ball from the scanning point cloud of the three-dimensional laser scanner, and obtaining the central scanning coordinates of a reference point of the prism target ball by using a fixed sphere radius fitting method; S3: using seven-parameter coordinate transformation, transforming all scanning point clouds of the three-dimensional laser scanner from the scanner coordinate system to the total station horizontal coordinate system; S4: partitioning the scanned point cloud converted to the horizontal coordinate system of the total station; S5: Pick up partition control points; S6: Perform overall adjustment on the partition control points; S7: Perform overall adjustment on the partitioned scan point cloud.

2. a kind of adjustment method based on massive point cloud as claimed in claim 1, is characterized in that In step S1, The scanning point cloud of the area containing the prism target ball is The fixed sphere radius fitting method comprises the following steps: The radius R of the prism target ball is taken as a known value, and the coordinates of the best fitting sphere center (C x ,C y ,C z ), so that the sum of squares of the residuals from the point to the sphere is minimized: Among them, the single point residual formula is: The fitting adopts the nonlinear least squares method, which includes the following steps: a. Linearized residual: Residual function in C (k) Taylor expansion: Where, the Jacobian matrix: Where, b. Construct a linear system: The residual equations of all points are combined to form a linear system of equations: J (k) ΔC=L (k) Formula 6; The least squares solution is: ΔC=(J (k)T J (k) ) -1 J (k)T L (k) Formula 7; Where: c. Gross error elimination: Residual standard deviation estimate: Aspherical point elimination criteria:

3. a kind of adjustment method based on massive point cloud as claimed in claim 2, is characterized in that In step S3, the method for converting all scan point clouds of the three-dimensional laser scanner into the horizontal coordinate system of the total station includes the following steps: Assume that o-xyh is the horizontal coordinate system of the total station, o'-x'y'h' is the scanner coordinate system, the number of point clouds is n1, the number of control points is n, and the horizontal coordinates of the total station and the scanner coordinates of the control points are expressed as (x i y i h i ) T and (x' i y' i h' i ) T , where (i=1,2,…,n); The spatial coordinate transformation model is expressed as: Where (x0 y0 h0) is the translation, k is the scale variable, R1(α), R2(β), and R3(γ) are the rotation matrices, which are: When the number of common points n is greater than 3, the above seven parameters are determined by least squares; For each common point, the error equation is listed as: Where, V=(V x V y V h ) T is the correction number; Take approximate values x0 for seven parameters 0 、y0 0 、h0 0 , α 0 , β 0 , γ 0 、k 0 , linearize Equation 11: Where, the partial derivative and constant term are: The normal equation is composed of Equation 12, which is iterated until convergence to obtain the translation (x0 y0 h0), rotation R1(α), R2(β), R3(γ), and scale value k. Then, according to Equation 10, all scanned point clouds are converted from the scanner coordinate system to the total station horizontal coordinate system.

4. a kind of adjustment method based on massive point cloud as claimed in claim 3, is characterized in that In step S4, Set N desired partitions; The partition processing method of the scan point cloud converted to the horizontal coordinate system of the total station includes the following steps: a. Calculation formula for the number of cubes: In the formula, the length, width, and height of the given space are L, W, and H respectively, and the side length of the cube is a. represents the rounding-up operation to ensure that the cube completely covers the space; b. Maximum edge length search: When the number of target cubes is N desired When , find the maximum side length a max , so that: Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice; c. Minimum edge length search: Find the minimum side length a min , so that: Iterate the solution in the interval [1, max(L, W, H)] by binary search method until the convergence condition is met 10 -6 rice; d. Calculation of theoretical side length: Based on the volume equipartition assumption, the calculation formula for the theoretical side length a0 is: e. Calculation of optimal side length: The set of candidate edge lengths is: {a0,a min ,to max }; By comparing the number of cubes corresponding to each candidate side length with N desired The absolute difference of * : |N-N desired |; If the differences are the same, choose the side length closest to a0; f. Calculation of the actual number of cubes: The final actual number of cubes N actual , through the optimal side length a * calculate:

5. a kind of adjustment method based on massive point cloud as claimed in claim 4, is characterized in that In step S5, The method for picking partition control points includes the following steps: Traverse all points and assign them to cubes, then traverse each cube point in turn to obtain the minimum bounding box of the cube points, and take the 8 corner points and the center point of the minimum bounding box, a total of 9 points, as partition control points.

6. a kind of adjustment method based on massive point cloud as claimed in claim 5, is characterized in that In step S6, The global adjustment method of the partition control points includes the following steps: After partitioning and obtaining the partition control points, the target sphere reference point is used as the starting point, and all partition control points participate in the overall least squares adjustment together; during the adjustment operation, independent adjustment calculations are performed on the plane coordinates and elevation data respectively; the horizontal angle, vertical angle and slope distance corresponding to each partition control point can be solved using the survey station coordinates and the three-dimensional coordinates of the partition control points after conversion to the horizontal coordinate system of the total station.

7. a kind of adjustment method based on massive point cloud as claimed in claim 6, is characterized in that In step S6, By using plane adjustment, the error equations of the horizontal distance and horizontal angle observation values can be listed respectively; The horizontal distance correction equation is: for the horizontal distance observation value S from point i to any point j measured at the i measuring station, ij : (δE' ij δN' ij ) T is the coordinate of point j in the horizontal coordinate system of station i (E' ij N' ij ) T number of corrections; Horizontal angle correction equation, for the direction observation value R from point i to any point j measured at station i ij : Where θ i is the orientation angle parameter of station i; Write the error equation in matrix form: V=AδX-L Equation 22; The error equation coefficient matrix A, the unknown correction factor δX, and the error equation constant term L in Equation 22 are: In V T Under the condition of PV=min, solve equation 22; Depend on have to: NδX=C Equation 24; Equation 24 is called the normal equation, and the normal equation coefficient matrix N = A T PA, the constant term of the normal equation C = A T PL, solve Equation 24: δX = N -1 Type C 25; The unknown number is approximated by X 0 Corrected to X 0 +δX, recompose the error equation and the normal equation, solve δX, correct the approximation, and iterate until convergence, that is, the correction number δX is small enough to obtain the solution of the unknown number.

8. a kind of adjustment method based on massive point cloud as claimed in claim 6, is characterized in that In step S7, the overall adjustment method of the partitioned scanning point cloud includes the following steps: Based on the coordinates of the partition control points before and after adjustment, the spatial coordinate transformation model of Equation 10 is used to solve the seven-parameter transformation parameters. After the partition scan point cloud is uniformly transformed into the adjustment coordinate system through the spatial coordinate transformation model of Equation 10, the overall adjustment result of the global scan point cloud is obtained.

Images