CN111325666B - Airborne laser point cloud processing method based on variable resolution voxel grid - Google Patents

Abstract

The invention relates to an airborne laser point cloud processing method and application based on a variable-resolution voxel grid. First, the variable-resolution point cloud compression algorithm of the present invention can compress the original laser point cloud into a variable-resolution point cloud, and the variable-resolution voxel grid obtained based on the regular grid of the variable-resolution point cloud can cover a relatively large area. The large plane range can also ensure that the central area has a high resolution, so that it can simultaneously meet the different requirements for the plane coverage and resolution of the voxel grid in the processing of different terrain scenes. Secondly, in order to be able to handle a larger size voxel grid and adopt a larger basic network structure, the present invention builds a three-dimensional semantic segmentation network with an encoding-decoding structure based on sub-flow sparse convolution. The network has high data processing speed and small memory usage. The variable-resolution voxel grid can flexibly respond to different terrain scenarios and effectively improve the robustness of the point cloud classification model.

Description

Airborne laser point cloud processing method based on variable resolution voxel grid

technical field

The invention relates to the field of remote sensing surveying and mapping, in particular to an airborne laser point cloud data processing method.

Background technique

As an important part of smart earth, 3D geographic data, its rapid acquisition technology and data post-processing methods have always been the focus of research in the field of photogrammetry and remote sensing. As a real-time and fast active measurement method, laser scanning technology can obtain high-precision 3D point clouds and attached physical properties of the measured target all day, all day, and obtain the occluded area through vegetation gaps. three-dimensional information. Research on point cloud data compression methods has been an important research topic in the past two decades. There are many classic and effective preprocessing algorithms used to realize point cloud data preprocessing, such as projection transformation and voxel division. According to the difference of the data dimensions after point cloud preprocessing, it can be divided into two categories: (1) The method based on voxel grid, which regularizes the 3D model into a voxel grid, for example, regularizes the 3D model into a size of 30×30×30. A voxel grid (Wu, 2015), or regularized to a 32×32×32 voxel grid (Maturana, 2015). In order to better index the physical data, octrees can be introduced to construct a voxel grid based on octrees (Riegler, 2017; Wang, 2017); (2) The image-based method mainly integrates three-dimensional Models are rendered into a series of 2D views or converted into geometric feature maps or panoramas, such as multi-angle views (Su, 2015; Qi, 2016; Shi, 2015), geometric feature maps (Sinha, 2016; Nannan Qin, 2018) ). However, based on the existing research methods, point cloud data preprocessing is still in the exploratory stage, mainly because of the diversity of point clouds (different point cloud spatial distribution, such as uneven density), the complexity of the scene (types of scanned data and objects) Many: such as houses, dense forests, artificial ground, etc.) and severe occlusion.

SUMMARY OF THE INVENTION

The purpose of the present invention is to provide a method for data processing in an airborne laser point cloud, which realizes point cloud compression based on variable resolution, reduces the geographic space occupied by data, and ensures the invariance of point cloud information description, It provides good preprocessing data for subsequent point cloud classification.

The technical solution of the present invention is that an airborne laser point cloud processing method based on a variable resolution voxel grid includes the following steps:

Step 1, point cloud input and block: read the original point cloud data, use the sliding window method to block the point cloud, the block size is w×h, the sliding step size is s, w=h;

Step 2: Compression strategy based on proportional transformation: Divide each block into a compressed area and an inner area, the center point of the divided block is P _c , and the radius of the largest circumcircle of the block is denoted as R ₁ , where the inner area refers to The circular area with P _c as the center and radius R ₃ , the rest is the compressed area; for the point cloud in the compressed area, the proportional linear compression method from the outside to the inside is used, that is, only the horizontal coordinates of the points are compressed, The elevation coordinates are kept in the original value; the internal area is not compressed, and the three-dimensional coordinates of the points in the inner area are processed in the original value;

Step 3, for each block, perform point cloud compression area data compression: according to the above compression strategy, perform data compression on the point cloud in the compression area;

R₂＝0.5×R₁；同时，令点P_i与点P_c之间的连线与水平线间的夹角为α，则根据相似三角形的几何关系得到如下公式：Any original point in the given compressed region is Pi (X _i , Y _i , Z _i ), and the compressed points are P′ _i (X′ _i , Y′ _i , Z′ _i ) point _{P i and} point P The intersection point between the line between ' _i and the inner area boundary circle is P _s (X _s , Y _s , Z _s ), the radius of the outer original boundary circle is R ₁ , and the radius of the outer compressed boundary circle is R ₂ , The radius of the inner region boundary circle is R ₃ , and the center point of the inner region boundary circle is P _c (X _c , Y _c , Z _c ), where,

R ₂ =0.5×R ₁ ; at the same time, let the angle between the connecting line between point _Pi and point P _c and the horizontal line be α, then the following formula can be obtained according to the geometric relationship of similar triangles:

The proportional relationship of linear compression is as follows:

代表点P_i和P_c之间的水平距离，其计算公式如下：in,

represents the horizontal distance between points P _i and P _c , and its calculation formula is as follows:

为点P_s和点P′_i之间的水平距离，计算方式与

的计算方式一样；

is the horizontal distance between point P _s and point P′ _i , calculated in the same way as

is calculated in the same way;

Combining the above formulas gives:

Transform it to get:

also,

X _s =R ₃ ×cosα+X _c

Combining the above formulas, we finally get:

In the same way, the following formula for calculating Y′ _i is obtained:

Since only the plane coordinates are compressed and changed, the elevation coordinates of the laser point before and after compression remain unchanged, namely:

Z′ _i =Z _i

After that, the compressed variable-resolution point cloud is converted into a regular grid to obtain its dual variable-resolution regular voxel grid.

Further, the specific implementation of step 1 is as follows:

Assuming that the geographic space occupied by the point cloud in the horizontal direction is X, Y respectively, the sliding step is s, the size of the interception window is w×h, and the minimum horizontal coordinates of the whole piece of data are x ₀ , y ₀ respectively, we get The block set is P, the block steps are as follows:

(1) In the X direction, slide horizontally according to the step size s, and take the point cloud of the fixed window, where the coordinates of the upper left corner and the lower right corner of the window are (x ₀ +n×s, y ₀ ), (x ₀ + n×s+w, y ₀ +h), n∈{0,1,2,…}, until all points in the X direction are divided into blocks, and the divided data is put into P;

(2) For the Y direction, slide longitudinally according to the step size s, and take the point cloud of the fixed window, where the coordinates of the upper left corner and the lower right corner of the window are (x ₀ , y ₀ +n×s), (x ₀ + w, y ₀ +n×s+h), n∈{0,1,2,…}, until all points in the Y direction are divided into blocks, and the divided data is put into P.

The present invention also provides an application of the airborne laser point cloud processing method based on the variable resolution voxel grid described in the above technical solution in point cloud classification, further comprising the following steps:

Step 4: Perform regular gridization on the compressed variable-resolution point cloud to obtain its dual variable-resolution regular voxel grid, and input the three-dimensional gridded point cloud into the pre-trained sub-flow pattern-based Sparse convolutional encoding-decoding structural semantic segmentation network to obtain the results of point cloud classification for each classification;

Step 5: Merge the classified segmented point clouds.

Further, the network structure of the encoding-decoding structure semantic segmentation network based on sub-stream sparse convolution in step 4 is as follows:

The encoding-decoding structure semantic segmentation network based on sub-stream sparse convolution consists of an encoding layer component and a decoding layer component; the encoding layer component is composed of 5 groups of convolutional layers, and the feature maps output by these 5 groups of convolutional layers The numbers are 32, 64, 96, 128, and 160, respectively. Each group of convolutional layers consists of two convolutional layers with the same structure and a downsampling convolutional layer. Each convolutional layer is preceded by a batch normalization. The structure of the decoding layer component is symmetrical with the structure of the encoding layer component, so the decoding layer component also includes 5 groups of convolutional layers, and the number of feature maps output by these 5 groups of convolutional layers is 160, 128, 96 respectively, 64 and 32; each group of convolutional layers consists of a deconvolutional layer for upsampling and two convolutional layers with the same structure; at the same time, the features output by the convolutional layer of the encoding layer component are connected by skip connection structure and decoding. The layer component performs channel merging corresponding to the output features of the convolutional layer; finally, the output of the decoding layer component is passed to a softmax layer for voxel-wise classification.

The invention directly utilizes the existing mature semantic segmentation network framework in image classification, and at the same time can avoid problems such as doping and overlapping of different types of points in the two-dimensional projection map. First, the variable-resolution point cloud compression algorithm of the present invention can compress the original laser point cloud into a variable-resolution point cloud, and the variable-resolution voxel grid obtained based on the regular grid of the variable-resolution point cloud can cover a relatively large area. The large plane range can also ensure that the central area has a high resolution, so that it can simultaneously meet the different requirements for the plane coverage and resolution of the voxel grid in the processing of different terrain scenes. Secondly, in order to be able to deal with a larger size voxel grid and adopt a larger basic network structure, the present invention builds a three-dimensional semantic segmentation network with an encoding-decoding structure based on sub-flow sparse convolution. The network has high data processing speed and small memory usage. The variable-resolution voxel grid can flexibly respond to different terrain scenarios and effectively improve the robustness of the point cloud classification model.

Description of drawings

FIG. 1 is a schematic diagram of transformation of point coordinates in a compressed area of the present invention.

FIG. 2 is the result of airborne laser point cloud data compression according to the present invention. (a), (b), (c) are examples of variable resolution compression

FIG. 3 is a semantic segmentation network based on the encoding-decoding structure of UNet-like sub-stream type sparse convolution according to the present invention.

FIG. 4 is the classification result of the airborne laser point cloud data of the present invention.

Detailed ways

The method provided by the present invention will be further described below with reference to the accompanying drawings.

An embodiment of the present invention provides an airborne laser point cloud data compression method, which specifically includes the following steps:

Step 1, point cloud input and block: read the original point cloud data, and use the sliding window method to block the point cloud. The purpose of segmenting the input point cloud based on the sliding window is to divide the data whose original range is too large into target processing units that cover enough spatial information.

Read in the original airborne laser point cloud data, and use the sliding window algorithm to divide the point cloud into blocks. Assuming that the geographic space occupied by the point cloud in the horizontal direction is X, Y respectively, the sliding step is s, the size of the interception window is w×h, and the minimum horizontal coordinates of the whole piece of data are x ₀ , y ₀ respectively, we get The block set is P, the block steps are as follows:

(3) In the X direction, slide horizontally according to the step size s, and take the point cloud of the fixed window, where the coordinates of the upper left corner and the lower right corner of the window are (x ₀ +n×s, y ₀ ), (x ₀ + n×s+w, y ₀ +h), n∈{0,1,2,…}, until all points in the X direction are divided into blocks, and the divided data is put into P.

(4) For the Y direction, slide longitudinally according to the step size s, and take the point cloud of the fixed window, where the coordinates of the upper left corner and the lower right corner of the window are (x ₀ , y ₀ +n×s), (x ₀ + w, y ₀ +n×s+h), n∈{0,1,2,…}, until all points in the Y direction are divided into blocks, and the divided data is put into P.

In the present invention, w=h.

Step 2: Compression strategy based on proportional transformation: Divide each block into a compression area and an inner area. For the point cloud in the compressed area, a proportional linear compression method from the outside to the inside is used, that is, only the horizontal coordinates of the points are used. Compression is performed, and the original value is maintained for the elevation coordinates; the original value is maintained for the three-dimensional coordinates of the points in the inner area.

The same piece of point cloud data generally has the same density. Point cloud compression compresses the external point cloud representing spatial context information to a certain range, and keeps the spatial content represented by this part of the point cloud unchanged by increasing the density of this part of the point cloud, but the compression strength cannot be too high, otherwise it will The distance between points is too close, resulting in the loss of presentation information, such as vegetation being compressed into a vertical straight line, buildings being compressed into a horizontal straight line, etc. Therefore, considering the following situation, it is necessary to set an appropriate compression ratio for the compression area of the input point cloud to ensure that the spatial distribution of points is compressed and the information loss is minimized.

This algorithm adopts the linear compression method, that is, the point cloud in the compressed area adopts the proportional linear compression method from the outside to the inside. Considering the particularity of the airborne laser point cloud, the difference between different objects is mainly determined by the elevation of the point. Therefore, when compressing, this algorithm only compresses the horizontal coordinate of the point, and maintains the original value for the elevation coordinate. The way. The processed point cloud can not only reduce the occupied space at the geographic level, but also ensure the consistency of the content of each feature in elevation. In order to ensure that the resolution of the point cloud in the central area remains unchanged, the three-dimensional coordinates of the points in the inner area are processed in a way of keeping the original values.

Step 3, data compression in the point cloud compression area: according to the above compression strategy, data compression is performed on the point cloud in the compression area.

R₂＝0.5×R₁同时，令点P_i与点P_c之间的连线与水平线间的夹角为α，则根据图1所示的相似三角形的几何关系可得如下公式：Any original point in the given compressed region is Pi (X _i ,Y _i ,Z _i ), the compressed point is P′ _i (X′ _i ,Y′ _i ,Z′ _{i )} , the point _Pi and the point The intersection point between the line between P′ _i and the inner uncompressed boundary circle is P _s (X _s , Y _s , Z _s ), the radius of the outer original boundary circle is R ₁ , and the radius of the outer compressed boundary circle is R ₂ , the radius of the inner region boundary circle is R ₃ , and the center point of the inner region boundary circle is P _c (X _c , Y _c , Z _c ). in,

R ₂ =0.5×R ₁ At the same time, let the angle between the line connecting point P _i and point P _c and the horizontal line be α, then according to the geometric relationship of similar triangles shown in Figure 1, the following formula can be obtained:

The proportional relationship of linear compression is as follows:

代表点P_i和P_c之间的水平距离，其计算公式如下：in,

represents the horizontal distance between points P _i and P _c , and its calculation formula is as follows:

为点P_s和点P′_i之间的水平距离，计算方式与

的计算方式一样。

is the horizontal distance between point P _s and point P′ _i , calculated in the same way as

are calculated in the same way.

Combining the above formulas gives:

Transform it to get:

In addition, according to the geometric relationship shown in Figure 1, the following formula can be obtained:

X _s =R ₃ ×cosα+X _c

Combining the above formulas, we finally get:

Similarly, the following formula for calculating Y′ _i can be obtained:

Since only the plane coordinates are compressed and changed, the elevation coordinates of the laser point before and after compression remain unchanged, namely:

Z′ _i =Z _i

Through the above compression algorithm, the points in the range of R ₁ -R ₂ are compressed into R ₂ -R ₃ .

Step 4: Perform regular gridization on the compressed variable resolution point cloud to obtain its dual variable resolution regular voxel grid. The 3D gridded point cloud is input into the encoding-decoding structure semantic segmentation network based on sub-flow sparse convolution, and the result of point cloud classification is obtained;

After obtaining the compressed point cloud, the data is fed into a pre-trained semantic segmentation network based on a UNet-like sub-stream-type sparse convolutional encoder-decoder structure. The network can perform feature learning and automatic classification on the input data to ensure the correctness of the classification results.

Among them, the UNet-like network consists of an encoding layer component and a decoding layer component. Among them, the coding layer component consists of 5 groups of convolutional layers. The number of feature maps output by these 5 groups of convolutional layers is 32, 64, 96, 128 and 160 respectively. Each group of convolutional layers consists of two volumes with the same structure. It consists of a convolutional layer and a downsampling convolutional layer. Each convolutional layer is preceded by a batch normalization and ReLU layer. The structure of the decoding layer component is symmetrical with that of the encoding layer component, so the decoding layer component also includes 5 groups of convolutional layers, and the number of feature maps output by these 5 groups of convolutional layers are 160, 128, 96, 64 and 32 respectively. Each set of convolutional layers consists of a deconvolutional layer for upsampling and two convolutional layers with the same structure. At the same time, the output features of the convolutional layer of the encoding layer component will also be channel merged with the output features of the corresponding convolutional layer of the decoding layer component through the skip connection structure. Finally, the output of the decoding layer component is passed to a softmax layer for voxel-wise classification.

Step 5, merge the classified sub-block point clouds, and each sub-point cloud only takes the middle uncompressed classification result;

The sub-block point cloud classification results obtained in step 4 are merged, and each sub-point cloud only takes the middle uncompressed classification result, and finally splices into a whole point cloud.

The specific embodiments described herein are merely illustrative of the spirit of the invention. Those skilled in the art to which the present invention pertains can make various modifications or additions to the described specific embodiments or substitute in similar manners, but will not deviate from the spirit of the present invention or go beyond the definitions of the appended claims range.

Claims (4)

1. The airborne laser point cloud processing method based on the variable resolution voxel grid is characterized by comprising the following steps:

step 1, point cloud input and blocking: reading in original point cloud data, and partitioning the point cloud by using a sliding window method, wherein the size of the partition is w multiplied by h, the sliding step length is s, and w is h;

step 2, compression strategy based on equal ratio transformation: dividing each block into a compressed area and an internal area, wherein the center point of each block is P_cAnd the radius of the maximum circumscribed circle of the segment is denoted as R₁Wherein the inner region is denoted by P_cAs a circle center and has a radius of R₃The remaining is a compressed area; adopting an equal proportion linear compression method from outside to inside for the point cloud in the compressed area, namely compressing the horizontal coordinate of the point only and keeping the original value for the elevation coordinate; the internal area is not compressed, and the three-dimensional coordinates of the points of the internal area adopt a processing mode of keeping original values;

step 3, performing point cloud compression area data compression for each block: performing data compression on the point cloud in the compressed area according to the compression strategy;

any original point in a given compression region is P_i(X_i,Y_i,Z_i) The point after compression is P'_i(X′_i,Y′_i,Z′_i) Point P_iAnd point P'_iThe intersection point of the connecting line between the two and the boundary circle of the inner region is P_s(X_s,Y_s,Z_s) The radius of the outer original boundary circle is R₁The radius of the boundary circle after external compression is R₂The radius of the boundary circle of the inner region is R₃The central point of the boundary circle of the inner region is P_c(X_c,Y_c,Z_c) Wherein

at the same time, let point P_iAnd point P_cThe included angle between the connecting line and the horizontal line is alpha, and then the following formula is obtained according to the geometric relationship of similar triangles:

the linear compression is scaled as follows:

wherein,

representative point P_iAnd P_cThe horizontal distance between the two is calculated as follows:

is a point P_sAnd point P'_iThe horizontal distance between the two, the calculation mode and

the calculation mode is the same;

the above formulas are combined to obtain:

transforming it can result in:

in addition to this, the present invention is,

X_s＝R₃×cosα+X_c

combining the formulas to obtain:

by the same token, the following Y is obtained_iThe calculation formula of `:

since only the plane coordinates are compressed, the elevation coordinates of the laser spot before and after compression are unchanged, namely:

Z'_i＝Z_i

and then regularly gridding the variable-resolution point cloud obtained by compression to obtain a variable-resolution regular voxel grid which is coupled with the variable-resolution point cloud.

2. The variable resolution voxel grid-based airborne laser point cloud processing method according to claim 1, characterized in that: the specific implementation of step 1 is as follows,

assuming that the geographic spaces occupied by the point cloud in the horizontal direction are X and Y respectively, the sliding step length is s, the size of the intercepting window is w X h, and the minimum value of the horizontal coordinate of the whole data block is X₀，y₀And if the obtained block set is P, the blocking steps are as follows:

(1) transversely sliding the X direction according to the step length s, and taking the point cloud of a fixed window, wherein the coordinates of the upper left corner and the lower right corner of the window are respectively (X)₀+n×s，y₀)，(x₀+n×s+w，y₀+ h), n belongs to {0,1,2, … }, until all points in the X direction are blocked, putting the blocked data into P;

(2) longitudinally sliding in the Y direction according to the step length s, and taking the point cloud of a fixed window, wherein the coordinates of the upper left corner and the lower right corner of the window are respectively (x)₀，y₀+n×s)，(x₀+w，y₀+ n × s + h), n ∈ {0,1,2, … }, until all points in the Y direction are blocked, and the blocked data is put into P.

3. The variable resolution voxel grid-based airborne laser point cloud processing method according to claim 1, characterized in that: the method also comprises the following steps of,

step 4, performing regular grid formation on the variable-resolution point cloud obtained by compression to obtain a dual variable-resolution regular voxel grid, and inputting the three-dimensional grid point cloud into a pre-trained encoding-decoding structure semantic segmentation network based on sub-stream sparse convolution to obtain a point cloud classification result of each classification;

and 5, merging the classified partitioned point clouds.

4. The method for processing the airborne laser point cloud based on the variable resolution voxel grid according to claim 3, wherein: the network structure of the encoding-decoding structure semantic segmentation network based on the sub-stream sparse convolution in the step 4 is as follows,

the coding-decoding structure semantic segmentation network based on the sub-stream sparse convolution is composed of a coding layer component and a decoding layer component; the coding layer component consists of 5 groups of convolutional layers, the number of characteristic graphs output by the 5 groups of convolutional layers is respectively 32, 64, 96, 128 and 160, each group of convolutional layers consists of two convolutional layers with the same structure and a down-sampling convolutional layer, and each convolutional layer is provided with a batch normalization layer and a ReLU layer in front; the structure of the decoding layer part is symmetrical to that of the coding layer part, so that the decoding layer part also comprises 5 groups of convolution layers, and the number of characteristic diagrams output by the 5 groups of convolution layers is respectively 160, 128, 96, 64 and 32; each group of convolution layers consists of a reverse convolution layer for up-sampling and two convolution layers with the same structure; meanwhile, the output characteristics of the convolution layer of the coding layer component are subjected to channel merging with the output characteristics of the convolution layer corresponding to the decoding layer component through a jump connection structure; finally, the output of the decoding layer component is passed to a softmax layer for voxel-by-voxel classification.

Images