CN106091972B - A kind of building change detecting method projecting dot density based on moving window - Google Patents

Abstract

The invention discloses a building change detection method based on the projected point density of a moving window, comprising the steps of: collecting and registering point cloud data in two phases; performing dimensionality reduction analysis on a stable wall; The projected points and establish the structural coordinate system; convert the point cloud data of the two phases into the structural coordinate system; use the K-means classification method to obtain the brick point cloud and mortar point cloud; project the mortar point cloud to the Z direction and the Y direction respectively; define Fixed window length L _fix and moving window length L _move , calculate point cloud line density changes along the Z direction and Y direction; find the horizontal and vertical dividing lines between bricks, calculate the coordinates of the four corners of each brick and the center of each brick ; Change detection using two corresponding brick centers. The invention automatically extracts the deformation information of the center of each brick on the wall, greatly improves the degree of automation, fully excavates the original data, and ensures that the deformation information of each brick of the masonry structure building can be effectively obtained under the premise of ensuring accuracy.

Description

A Building Change Detection Method Based on Projected Point Density of Moving Window

technical field

The invention designs a building change detection method, and in particular relates to a building change detection method based on a moving window.

Background technique

In cities, buildings are the main places for human activities, and their safety status is related to human daily life and economic activities. Therefore, it is of great significance to detect changes in buildings. Based on the change detection and repair of buildings, especially The deformation and damage of buildings caused by earthquakes has been the focus of research in the industry in recent years. In the development history of building structures, masonry structures have been widely used in the foundation structure of buildings since ancient times. Bricks, stones, masonry And adobe and other blocks, a combination of mortar (mortar, clay slurry, etc.) However, its masonry strength is low and its earthquake resistance is poor. Therefore, the research on change detection of masonry buildings has extremely important practical significance. The existing image-based change detection methods, according to the level of information processing , can be divided into pixel-based, feature-based, and target-based change detection, and the change detection of LiDAR (Light Detection And Ranging) LiDAR data has gradually begun to be studied in the field of photogrammetry and computer vision, and image-based change detection The accuracy of the method is not high enough, the target is not clear enough, and it is easily affected by the image quality. The obtained deformation information is not rich enough, and it is impossible to accurately obtain the deformation information of each brick in the masonry structure. In addition, the low degree of automation also seriously affects Efficiency of change detection.

Contents of the invention

Purpose of the invention: The purpose of the present invention is to provide a method for detecting building changes based on the density of projected points of the moving window to address the deficiencies of the prior art.

Technical solution: A method for detecting building changes based on the density of moving window projection points of the present invention comprises the following steps:

(1) Use the laser scanner system to scan the same building in two phases to obtain the point cloud data of the building surface, and set m targets around the changing building, where m≥3, the observation value of the laser scanner system is Three-dimensional coordinates and laser reflection intensity of building surface points;

(2) Using the target set in step 1), calculate the coordinate transformation parameter Z of the point cloud in two phases, and register the point cloud;

(3) Select the point cloud data of the fixed wall, use the principal component analysis method to perform dimension reduction analysis on it, obtain the eigenvectors (v _i , i=1,2,3), and calculate the projection of the station on the plane where the fixed wall is located , take it as the coordinate origin, establish the structural coordinate system, and convert the point cloud data from the station coordinate system to the structural coordinate system;

(4) Based on the intensity information of the point cloud data, the K-means clustering method is used to classify the point cloud of the changing wall, and the brick point cloud and the mortar point cloud are separated;

(5) Using the mortar point cloud obtained in step (4), project the point cloud coordinates to the Z and Y directions, define the fixed window length L _fix and the moving window length L _move , and calculate the distance along the Z direction and Y direction respectively by moving the window Direction point cloud line density change;

(6) According to the linear density variation that step (5) obtains, obtain horizontal and vertical dividing line between each brick respectively, calculate the coordinates of four corner points of each brick, set up the brick model;

(7) obtain each brick point cloud according to the brick model that step (6) obtains, calculate each brick center;

(8) Obtain deformation information according to the three-dimensional coordinates of the brick centers corresponding to the two phases obtained in step (7).

Preferably, during the two-phase scanning in step (1), the target is fixed.

Preferably, the three-dimensional coordinate conversion equation described in step (2) and step (3) is specifically as follows:

Let matrix A be the three-dimensional coordinates of the point cloud in the A coordinate system, and matrix B be the three-dimensional coordinates of the point cloud in the B coordinate system. The three-dimensional coordinate transformation equations of the two coordinate systems A and B are as follows:

(Δx, Δy and Δz represent the translation of the coordinate origin, k is the scale factor, k=0, R is the rotation matrix from the A coordinate system to the B coordinate system)

Preferably, in step (2), the calculation of the registration conversion parameters by the two-period data is specifically as follows:

The coordinate transformation parameter Z is further written as, Z=[Δx, Δy, Δz, ε _x , ε _y , ε _z , 1] The least square method is used to estimate the coordinate transformation parameter Z, and the estimate of the coordinate transformation parameter Z can be obtained as :

Z＝(A ^T Q ^-1 A) ^-1 A ^T QB

In the formula, Q is the covariance matrix of m target coordinate measurement errors in the B coordinate system, and the form is as follows:

Preferably, the specific method of converting from the station coordinate system to the structure coordinate system in the step (3) is:

Assuming the three-dimensional coordinates of the scanning point X {X _i =(xi _, y _i , _zi )|i=1,2,…,n}, construct the corresponding covariance matrix:

in, is the center of gravity coordinates of the point set, and the matrix C is subjected to principal component analysis, and the three eigenvalues λ ₁ , λ ₂ , λ ₃ can be obtained in descending order, and λ ₁ ≥ λ ₂ > λ ₃ > 0, corresponding to λ ₃ The eigenvector v ₃ of , and v ₃ is the normal vector, v ₃ is the unit vector of the X-axis of the structure coordinate system in the station coordinate system, and the Z-axis of the structure coordinate system is consistent with the station coordinate system, and the Y-axis is vertical Based on the determined XOZ plane, a right-handed coordinate system is formed, and the projection S'(x _s , y _s , z _s ) of the station coordinate S(0,0,0) on the plane where the fixed wall is located is calculated, and it is used as the structural coordinate system The origin of the coordinates, so the translation vector (Δx, Δy, Δz) = (-x _s , -y _s , -z _s ), after establishing the translation parameters and coordinate axis rotation parameters, the point cloud data after the two-phase registration Rotate to the structure coordinate system;

Preferably, the specific method for separating wall bricks and mortar by the K-means clustering method based on strength information proposed in step (4) is:

The clustering error square sum function E is used as the clustering criterion function, and the point intensity information is used as the classification attribute, where, x _ij is the jth sample of class i, m _i is the cluster center or centroid of class i, n _i is the number of samples of class i, and the K-means clustering algorithm finds k best clusters through repeated iterations. The cluster center, where k=2, assigns all n sample points to the cluster center closest to it, so that the sum of squared clustering errors E is the smallest, the process is as follows:

a, randomly designate k cluster centers m _i (i=1,2,...,k);

b, find the nearest cluster center for each sample _xi , and assign it to this class;

c. Recalculate the new centers of each cluster: N _i is the current sample number of the i-th cluster;

d, calculate the deviation,

e, if the E value converges, then return m _i (i=1,2,…,k), the algorithm terminates, otherwise return b;

Preferably, step (5) uses the mortar point cloud to calculate the line density change of the point cloud based on the window moving method. The specific method is as follows:

Let the average width of the mortar be the known quantity L _mortar , define the fixed window length L _fix and the moving window length L _move , and the three satisfy the following relationship:

L _mortgage ≈L _window +2L _move

Calculate the number of moving windows along the Z and Y directions respectively:

Where [] is rounding symbol,

Calculate the number of points in each window along the Z direction and the Y direction respectively:

n _zi (i=1,2,…,n _y ),n _yi (i=1,2,…,n _z )

Calculate the line density of points along the Z and Y directions separately:

Density_z＝(n _z(i-1) +n _zi +n _z(i+1) )/(3L _fix )(i＝2,3,…,(n _z -1))

Density_y＝(n _y(i-1) +n y +n _y(i+1) )/(3L _fix )(i ₌ 2,3,...,(n _y -1))

For each window in the Z direction and Y direction, calculate the line density change rate:

Grad(i,1)=Density_y(i)-Density_y(i-1)

Grad(i,2)=Density_y(i+1)-Density_y(i)

Preferably, described in step (6) utilizes linear density variation, calculates horizontal and vertical dividing line between each brick, and its concrete method is as follows:

For the analyzed Z direction or Y direction, the point cloud line density is significantly higher than the non-joint area at the mortar joint perpendicular to this direction, so the average density of the window in this direction is selected as the threshold ( n _total is the total number of point clouds), when the point cloud density in the window is greater than the threshold, and Grad(i,1)>0 and Grad(i,2)<0, the fixed window is the range of the mortar joint , Calculate the center line of the mortar joint through the point cloud within the range of the mortar joint, and then calculate the coordinates of the four corner points of the brick model through the center line of the joint around the brick to establish the brick model.

Beneficial effects: The method for detecting building changes based on the density of moving window projection points of the present invention has a high degree of automation on the one hand, and fully utilizes the original data on the other hand with high precision. The deformation information of each brick of the structural building can be effectively obtained.

Description of drawings

Fig. 1 is a flowchart of the present invention;

Fig. 2 is the schematic diagram of station coordinate system of the present invention;

Fig. 3 is a schematic diagram of a structural coordinate system of the present invention;

Fig. 4 is a schematic diagram of a brick model of the present invention;

Fig. 5 is a certain wall surface K mean value classification result that the present invention chooses;

Fig. 6 is a histogram of point line density variation of point cloud projection in the present invention;

Fig. 7 is the schematic diagram of each brick central point that the present invention extracts;

Figure 8 is the three-dimensional deformation information of bricks acquired by the present invention.

Detailed ways

The technical solutions of the present invention will be described in detail below, but the protection scope of the present invention is not limited to the embodiments.

As shown in Figure 1, a building change detection method based on the density of moving window projection points includes the following steps:

(1) Use the laser scanner system to scan the same building in two phases to obtain the point cloud data of the building surface, and set m targets in the fixed area around the changing building. Generally, m≥4. Considering the laser scanner It comes with an automatic leveling function, and the instrument is in a horizontal state for each measurement, so m≥3 is sufficient, and the observation values of the laser scanner system are the three-dimensional coordinates of the building surface points and the laser reflection intensity;

(2) Using the target set in step (1), calculate the coordinate transformation parameter Z of the point cloud in two phases, and register the point cloud;

(3) Select the point cloud data of the fixed wall, use the principal component analysis method to perform dimension reduction analysis on it, obtain the eigenvectors (v _i , i=1,2,3), and calculate the projection of the station on the plane where the fixed wall is located , take it as the coordinate origin O, establish a structural coordinate system, and convert the point cloud data from the station coordinate system to the structural coordinate system;

(4) Based on the intensity information of the point cloud data, the K-means clustering method is used to classify the point cloud of the changing wall, and the brick point cloud and the mortar point cloud are separated;

(5) Using the mortar point cloud obtained in step (4), project the point cloud coordinates to the Z and Y directions, define the fixed window length L _fix and the moving window length L _move , and calculate the distance along the Z direction and Y direction respectively by moving the window Direction point cloud line density change;

(6) According to the linear density variation that step (5) obtains, obtain horizontal and vertical dividing line between each brick respectively, calculate the coordinates of four corner points of each brick, set up the brick model;

(7) obtain each brick point cloud according to the brick model that step (6) obtains, calculate each brick center;

(8) Obtain deformation information according to the three-dimensional coordinates of the brick centers corresponding to the two phases obtained in step (7).

Taking "the change detection of a masonry structure building before and after an earthquake" as an example, the present invention is further elaborated:

{1} Use the Leica C10 laser scanner system to scan the building, set up the instrument at the position shown in Figure 5, and set up 4 targets in the stable area outside the experimental building, scan before and after the earthquake test respectively, and obtain two phases Laser point cloud data on the surface of buildings, observations include two types: three-dimensional coordinates, laser reflection intensity;

{2}As shown in Figure 2, with the station as the center and the vertical direction as the Z axis, use the target set in step {1} to calculate the coordinate conversion parameter Z of the two-phase point cloud, and coordinate the coordinates of the point cloud after the earthquake test Accurate to the coordinate system of the point cloud station before the seismic test, so far, the two point clouds are in the same coordinate system, so as to facilitate the deformation analysis in the later stage;

{3}As shown in Figure 3, the three coordinate axes are respectively parallel to the three adjacent faces of the building centered on a certain corner point. The point cloud data of the fixed wall surface is selected, and the dimensionality reduction analysis is carried out by principal component analysis to obtain the fixed The normal vector of the wall, that is, the X-axis unit vector of the structural coordinate system is [0.9348 0.3551 0.0011], the Z-axis unit vector of the structural coordinate system is [0 0 1], so the Y-axis unit vector is [0.3551 -0.9348 0], and the station coordinates The projection point from the origin of system coordinates to the fixed wall is [8.27663.1442 0.0095]. According to the above parameters and conversion relationship, the point cloud data is converted from the station coordinate system to the structure coordinate system;

{4} As shown in Figure 5, according to the intensity information of the point cloud data, the point cloud of the two phases of the wall is classified by the K-means clustering method, and the brick point cloud and the mortar point cloud are separated;

{5} After on-site investigation and measurement, the average width of mortar between bricks is 10-12mm. Therefore, according to L _mortar ≈ L _window + 2L _move , the fixed window length L _fix = 0.008m and the moving window length L _move = 0.002m , using the mortar point cloud obtained in step {4}, project the point cloud coordinates to the Z direction and Y direction, and calculate the point cloud line density changes along the Z direction and Y direction by moving the window, and the results are shown in Figure 6;

{6}According to the linear density change obtained in step {5}, respectively obtain the vertical and horizontal dividing lines between each brick, and calculate the coordinates of the four corner points of each brick, and establish the brick as shown in Figure 4 Model;

{7}As shown in Figure 7, according to the brick model obtained in step {6}, the point cloud data of each brick is obtained, and the center of each brick is calculated;

{8}As shown in Figure 8, according to the three-dimensional coordinates of the brick centers corresponding to the two phases obtained in step {7}, the deformation information is obtained.

Claims (7)

1. a kind of building change detection method based on moving window projection point density, it is characterized in that: comprise the following steps: 1) Use the laser scanner system to scan the same building in two phases to obtain the point cloud data of the building surface, and set m targets around the changing building, where m≥3, the observation value of the laser scanner system is the building The three-dimensional coordinates of the object surface point and the laser reflection intensity; 2) Using the target set in step 1), calculate the coordinate transformation parameter Z of the point cloud in two phases, and register the point cloud; 3) Select the fixed wall point cloud data, use the principal component analysis method to perform dimension reduction analysis on it, obtain the feature vector v _i , i = 1, 2, 3, calculate the projection of the station on the plane where the fixed wall is located, and calculate its As the coordinate origin, establish a structural coordinate system, and convert the point cloud data from the station coordinate system to the structural coordinate system; 4) Based on the intensity information of the point cloud data, the K-means clustering method is used to classify the point cloud of the changing wall, and the brick point cloud and the mortar point cloud are separated; 5) Using the mortar point cloud obtained in step 4), it is advisable to set the average width of the mortar as L _mortar , project the coordinates of the point cloud to the Z direction and the Y direction, and define the fixed window length L _fix and the moving window _length L move to satisfy: L _Mortar ≈L _window +2L _move ; determine the number of moving windows according to the maximum and minimum values of the point coordinate component Z or Y in the point cloud, and calculate the line density changes of the point cloud along the Z direction and the Y direction through the moving window; 6) Step 5) The linear density change of the obtained point cloud. For the analyzed Z direction or Y direction, the point cloud linear density is significantly higher than the non-joint area at the mortar joint perpendicular to this direction. Therefore, select this The average density of the window in the direction is used as the threshold n _total is the total number of point clouds. When the point cloud density in the window is greater than the threshold and Grad(i,1)>0 and Grad(i,2)<0, the fixed window is the range of the mortar joint. Calculate the center line of the mortar joint through the point cloud within the range of the mortar joint, and then calculate the coordinates of the four corner points of the brick model through the center line of the joint around the brick to establish the brick model; 7) Obtain each brick point cloud according to the brick model obtained in step 6), obtain the geometric center of the brick point cloud, and obtain each brick center; 8) Obtain deformation information according to the three-dimensional coordinates of the brick centers corresponding to the two phases obtained in step 7). 2. The building change detection method based on the projected point density of the moving window as claimed in claim 1, characterized in that: step 1) During the two-stage scanning process, the target is fixed. 3. the building change detection method based on moving window projection point density as claimed in claim 1, is characterized in that: step 2) and described in step 3) three-dimensional coordinate conversion equation is specifically as follows: Let matrix A be the three-dimensional coordinates of the point cloud in the A coordinate system, and matrix B be the three-dimensional coordinates of the point cloud in the B coordinate system. The three-dimensional coordinate transformation equations of the two coordinate systems A and B are as follows: Among them, Δx, Δy and Δz represent the translation of the coordinate origin, k is the scale factor, k=0, and R is the rotation matrix from the A coordinate system to the B coordinate system. 4. the building change detection method based on moving window projection point density as claimed in claim 1, is characterized in that: step 2) described by the calculation of two-period data registration transformation parameter is specifically as follows: The coordinate transformation parameter Z is further written as, Z=[Δx, Δy, Δz, ε _x , ε _y , ε _z , 1], the coordinate transformation parameter Z is estimated by the least square method, and the estimation of the coordinate transformation parameter Z can be obtained for: Z＝(A ^T Q ^-1 A) ^-1 A ^T QB In the formula, Q is the covariance matrix of m target coordinate measurement errors in the B coordinate system, and the form is as follows: 5. the building change detection method based on moving window projection point density as claimed in claim 1, is characterized in that: step 3) mentioned is transformed into principal component analysis dimensionality reduction and coordinate origin in the structure coordinate system by station coordinate system The specific method of determination is: Assuming the three-dimensional coordinates of the scanning point X {X _i =(x _i ,y _i , _zi )|i=1,2,…,n}, construct the corresponding covariance matrix: in, is the center of gravity coordinates of the point set, and the matrix C is subjected to principal component analysis, and the three eigenvalues λ ₁ , λ ₂ , λ ₃ can be obtained in descending order, and λ ₁ ≥ λ ₂ > λ ₃ > 0, corresponding to λ ₃ The eigenvector v ₃ of , and v ₃ is the normal vector, v ₃ is the unit vector of the X-axis of the structure coordinate system in the station coordinate system, and the Z-axis of the structure coordinate system is consistent with the station coordinate system, and the Y-axis is vertical Based on the determined XOZ plane, a right-handed coordinate system is formed, and the projection S'(x _s , y _s , z _s ) of the station coordinate S(0,0,0) on the plane where the fixed wall is located is calculated, and it is used as the structural coordinate system The origin of the coordinates, so the translation vector (Δx, Δy, Δz) = (-x _s , -y _s , -z _s ), after establishing the translation parameters and coordinate axis rotation parameters, the point cloud data after the two-phase registration Rotate to the structure coordinate system. 6. the building change detection method based on moving window projection point density as claimed in claim 1, is characterized in that: step 4) the concrete method that the proposed K-means clustering method based on intensity information separates metope brick and mortar for: The clustering error square sum function E is used as the clustering criterion function, and the point intensity information is used as the classification attribute, where, x _ij is the jth sample of the i-th class, m _i is the cluster or centroid of the i-th class, and n _i is the number of samples in the i-th class. The K-means clustering algorithm finds k best clusters through repeated iterations. The cluster center, where k=2, assigns all n sample points to the cluster center closest to it, so that the sum of squared clustering errors E is the smallest, the process is as follows: Step 1, randomly designate k cluster centers m _i (i=1,2,...,k); Step 2, find the nearest cluster center for each sample _xi , and assign it to this class; Step 3, recalculate the new centers of each cluster: N _i is the current sample number of the i-th cluster; Step 4, calculate the deviation, Step 5, if the E value converges, then return mi ( _i =1,2,...,k), and the algorithm terminates, otherwise, return to Step2. 7. the building change detection method based on moving window projection point density as claimed in claim 1, is characterized in that: step 5) mentions and utilizes mortar point cloud, calculates point cloud linear density change based on window moving method, its specific method as follows: Calculate the number of moving windows along the Z and Y directions respectively: Where [] is the rounding symbol, y _max and y _min are the maximum and minimum values of the Y coordinate component of the point cloud, respectively; z _max and z _min are the maximum and minimum values of the Z coordinate component of the point cloud, respectively; Calculate the number of points in each window along the Z direction and the Y direction respectively: n _zi (i=1,2,…,n _y ),n _yi (i=1,2,…,n _z ) Calculate the line density of points along the Z and Y directions separately: Density_z＝(n _z(i-1) +n _zi +n _z(i+1) )/(3L _fix )(i＝2,3,…,(n _z -1)) Density_y＝(n _y(i-1) +n y +n _y(i+1) )/(3L _fix )(i ₌ 2,3,...,(n _y -1)) For each window in the Z direction and Y direction, calculate the line density change rate: Grad(i,1)=Density_y(i)-Density_y(i-1) Grad(i,2)=Density_y(i+1)−Density_y(i).