A method for constructing a building story-level entity by fusing laser point clouds and joint surveying

Abstract

The application provides a kind of laser point cloud fusion joint surveying and mapping method for constructing building storey and household level entity, comprising the following steps: graphic contour conversion into point set;Generating graphic center point set, combining the attributes of completion acceptance measurement results and real estate measurement results, based on spatial relationship for semantic fusion, and generating center point by calculating the surface graphic element of real estate measurement;Multi-modal point cloud data fusion, make up and correct building contour point;Building contour graphic segmentation extraction, the point cloud semantic information converted by joint surveying and mapping results is used as the label of point cloud segmentation extraction, the fused point cloud data is segmented, and the building edge contour and internal partition household line are extracted, and the building storey, layer and household entity are reconstructed and generated;By converting into point cloud, multi-modal point cloud data is formed by using joint surveying and mapping data, airborne laser radar, tilt photography and other data are fused, and the precise construction of building storey, layer and even household entity model is realized.

Description

Technical field

[0001] The present invention relates to the technical field of point cloud maps, and in particular to a method for constructing a building entity at the floor and household level by fusing laser point clouds and jointly mapping. [Background Technology]

[0002] Urban modeling and refined urban management have become research hotspots in recent years. Semantic extraction systems for urban facilities can provide important support for urban modeling and management. Realistic 3D is a digital virtual space that provides a realistic, three-dimensional, and temporal representation of human production, living, and ecological spaces. City-level real-world 3D is not only a crucial component of the comprehensive construction of a real-world 3D China, but also serves as a unified, authoritative 3D spatiotemporal information foundation for the development of urban digital twins, refined management, emergency response, and 3D management of natural resources. With the continuous development and maturity of oblique photography and laser point cloud 3D modeling technologies, they have become the two mainstream technologies for real-world 3D city construction in China.

[0003] While oblique photography 3D models excel at high precision, realism, omnidirectional perception, and a high degree of process automation, they suffer from significant limitations in modeling weak, fine, and distant surface textures. Laser point clouds, on the other hand, excel at high penetration, precision, efficiency, and diverse acquisition platforms and methods. However, they also suffer from a lack of texture features, a heavy workload for in-house mapping, and an incomplete understanding of fine features. Constructing building entities through the fusion of oblique photography with airborne and vehicle-borne point cloud data can lead to missing data in complex urban environments with diverse building shapes, styles, and heights, due to occlusion between buildings and other entities, impacting the quality of the 3D model. Furthermore, without internal spatial segmentation data for each floor or household, the constructed buildings can only be granular enough to capture the entire building, failing to meet the needs of applications such as refined urban management, emergency disaster mitigation, and 3D natural resource management. If you want to achieve a three-dimensional model of a building with granularity at the level of each floor and each household, detailed modeling through indoor laser scanning is not only costly and time-consuming, but also faces the problem of data loss due to occlusion.

[0004] During urban construction, various competent departments have accumulated a wealth of specialized building-related data through the management of construction projects. This includes joint surveying and mapping (real estate surveying, completion acceptance surveying, and civil air defense surveying), and other graphic data. These graphic data not only depict buildings with clear outlines and boundaries, but also boast high spatial accuracy and rich semantic attributes. They contain measurement information such as building height and elevation that is suitable for constructing 3D models. Although this data, drawn due to multi-faceted business management rules, does not fully coincide with the real-world 3D building models in terms of spatial overlap, it can be used to correct and compensate for data loss caused by occlusion in real-world 3D data collection, resulting in reduced accuracy and incomplete spatial shape representation.

[0005] At present, these thematic business data have not been applied to the construction of real-life three-dimensional fine models. How to quickly integrate the existing thematic business data with the collected tilt models and laser point cloud data to construct low-cost, fine-grained, and high-spatial-precision building entities of buildings, floors, and households has become a common problem in the industry. It is urgent to propose a new technical method to solve the above problems. [Summary of the invention]

[0006] The present invention provides a method for constructing entity models of building buildings, floors and even households by fusing laser point cloud with joint surveying and mapping. The method converts the joint surveying and mapping data into point clouds to form multimodal point cloud data, which is then integrated with airborne laser radar, oblique photography and other data to achieve the precise construction of entity models of building buildings, floors and even households. The method can not only achieve the construction of finer-grained and more accurate entities, solve the problem of automatic adaptation of geometric semantics, attribute semantics and expression granularity when constructing entity models of fine urban buildings, but also make full use of the rich thematic data of natural resources to achieve the entity construction of real-life three-dimensional fine buildings.

[0007] In order to achieve the above-mentioned purpose, the technical solutions adopted are as follows:

[0008] A method for constructing building-level entities by laser point cloud fusion and joint surveying and mapping includes the following steps:

[0009] Step S1. Convert the graphic outline into a point set

[0010] Convert the spatial surface primitives of real estate and completion acceptance measurement results into a set of points representing their contour features;

[0011] Step S2. Generate primitive center point set

[0012] The attributes of the completion acceptance measurement results and the real estate measurement results are semantically integrated based on spatial relationships, and the center points of the surface elements of the real estate measurement are calculated;

[0013] Assign the fused semantic attribute information to the center point. The center point attributes will contain the floor, room number, floor height and floor elevation value attribute information.

[0014] Step S3. Multimodal point cloud data fusion, compensation and correction of building outline points

[0015] Using a multimodal point cloud matching algorithm and a boundary contour constraint algorithm model, the point cloud generated above is fused with the point cloud generated by radar field scanning and the point cloud generated by the tilt model, and noise points are removed to correct the building outline points generated by the airborne radar point cloud.

[0016] Step S4. Building outline segmentation and extraction

[0017] The point cloud semantic information converted from the joint surveying and mapping results is used as the label for point cloud segmentation and extraction. The fused point cloud data is segmented, and the building edge contours and internal stratification and household lines are extracted to reconstruct the building, floor, and household entities.

[0018] Furthermore, the step S1 specifically further includes the following steps:

[0019] Step S1.1: First, in order to convert the joint surveying and mapping metadata into a 3D point cloud, it is necessary to calculate the polygonal approximate area edge for the metadata;

[0020] Step S1.2, then, convert the polygon base points into a contour point set, and calculate the center point for each polygon to construct a center point set to store the entity semantic attributes of the current polygon;

[0021] Among them, the Halcon operator is used to calculate the polygon matrix R, the rows and columns of each matrix represent the coordinates of the key points of the variable sequence, and each polygon represents a connected area in the image metadata;

[0022] Although the measurement error can be ignored due to the high accuracy of the surveying and mapping metadata, in order to adapt to a wider range of applications, the Halcon operator introduces an error parameter, the Hausdorff distance h, to calculate the maximum distance between the polygon area and the closed area of ​​the primitive, so that the two areas can intersect as much as possible, as shown in formula (1):

[0023]

[0024] In formula (1), U represents the area of ​​the survey and mapping metadata, and d(·) calculates the Euclidean distance between the two sets. In this way, the non-empty subset of the metric space itself is transformed into the metric space to facilitate subsequent calculations.

[0025] More importantly, the smaller the calculated error parameter is, the more contour key points are generated in the vector graphics. Therefore, setting the error parameter as small as possible generates more key points and generates a denser point cloud.

[0026] Furthermore, the step S2 specifically further includes the following steps:

[0027] Step S2.1: In order to convert the polygons of the plane into spatial points, the spatial information they possess must be considered;

[0028] The joint mapping metadata attribute information contains accurate three-dimensional information. The center point of each polygon is calculated, and the elevation information in the corresponding mapping metadata attribute is assigned to the metadata, which is then converted into the center point of the corresponding key point polygon, and the elevation information is used as the attribute value of the key point.

[0029] Step S2.2: Divide the two-dimensional plane into a grid and fill the key points in the two-dimensional grid. However, since it is impossible to ensure that every point is on a regular grid point, it is necessary to resample the plane points and generate a completely regular dense grid point pattern.

[0030] The resampling problem is defined as an assignment problem. Each grid point is assigned a key point for filling, and the total movement cost is required to be minimized. The specific processing method is shown in formula (2):

[0031]

[0032] The purpose of formula (2) is to optimize a set of permutation matrices W g and W i l , to uniquely transform the global feature F G and local features Assigned to a predefined 2D grid of points where W g Relative to F G The mapping, Redundant grids represent grid points with only some of the grid points Form a corresponding relationship, N G is the dimension of the permutation matrix and the grid; Formula (2) ensures that in the global and local neighborhoods, only a subset of grid points form a one-to-one context point, and all missing matching grid points are filled by the nearest neighbors. The matched points are mapped to the three-dimensional space according to their elevation values ​​to form a dense point cloud.

[0033] Furthermore, the step S3 specifically further includes the following steps:

[0034] The vertex set and edge set of the building outline point cloud are used to combine with the floor point cloud of the jointly surveyed building to perform building outline correction.

[0035] Step S3.1, first, use the point cloud P scanned by the airborne laser radar s and point cloud P generated by joint mapping elements t Construct graph model G respectively s and graph model G t ;

[0036] The vertices in the two graphs are matched using the Hungarian algorithm. During the matching, the number of vertices in the point cloud contour graph must be greater than the number of vertices in the joint mapping primitives. The predicted points are linked to the primitive points using the Euclidean distance.

[0037] Step S3.2, use graph matching algorithm to find point cloud P s and point cloud P t matching relationship;

[0038] Use multi-layer perceptron to extract point cloud P s Point features and encode them into a one-dimensional point sequence;

[0039] For the point cloud P t , since it is a diffuse grid point formed by filling neighboring points, the structured feature sequence is extracted using circular convolution, and the first-order graph matching similarity matrix is ​​directly obtained by multiplying the first-order feature points of formula (3):

[0040] M p =U 1 U 2T ......(3)

[0041] In order to obtain the second-order similarity, we first connect the two nodes corresponding to each edge to construct the second-order feature matrix of the two graphs, as shown in formula (4):

[0042]

[0043] Among them, formula (3) is for the first-order similarity, U 1 and U 2 are the feature vectors of each node, respectively, and the superscripts indicate different data sources;

[0044] Formula (4) For the second-order similarity, X and Y are edge feature matrices, which are constructed by setting the descriptor of the cth edge starting from the i-th node and ending at the j-th node to the concatenation matrix of the descriptors of the two nodes, W 2is a learnable parameter, so the second-order similarity has a learnable matching function. Then, the graph matching problem is approximated as the problem of finding the eigenvector corresponding to the maximum eigenvalue of the similarity matrix, which is calculated using power iteration.

[0045] The matching results of the graph are represented by a double random matrix, which intuitively reflects the possibility of establishing a matching relationship between any pair of nodes.

[0046] First, the matrix is ​​normalized by column and then by row. The double random matrix is ​​calculated through multiple iterations. However, since the values ​​of elements in the same row and column are not much different, in order to increase the obvious difference and facilitate the subsequent loss function calculation, the matrix is ​​multiplied by a constant and then processed using Softmax to obtain the probability matrix of each candidate node match. This matrix still satisfies the double random property. Finally, the loss function is calculated using the offset vector-based method to calculate the offset between the matching points. Then, the weighted sum is taken according to the probability matrix to predict the offset vector and minimize the difference between the predicted vector and the true offset. The loss function is shown in formula (5):

[0047]

[0048] Among them, d i and Denote the predicted and actual offsets from the far point to the correct matching point, respectively. Φ(·) is a penalty term. In the final loss function calculation process, the discarded points will not affect the calculation of the loss function. Through the above method, the airborne lidar point cloud and the joint mapping building layer point cloud are matched and fused. The fused point cloud is projected onto a two-dimensional plane to extract the layered information.

[0049] The boundary outline of the building layer should satisfy the manifold assumption, where each vertex can only be connected to two adjacent edges;

[0050] Given the edge set generated in the previous step, select a subset of edges to calculate by using the constrained integer programming formula; i ∈{0, 1} is defined as a binary variable to describe the edge e i Whether ∈E is activated as part of the boundary state; the set of activated edges constitutes the polygonal boundary of the layer entity;

[0051] In addition, an energy function is used to constrain the generation quality of the boundary contour, as shown in formula (6):

[0052] U(x)=(1-λ)U fidelity (x)+λU complexity (x).......(6)

[0053] In formula (6), U fidelity(x) describes the consistency of state x with the input data, U complexity (x) measures the complexity of the output boundary shape, both of which are in the interval [0, 1], and λ∈(0, 1); U fidelity (x) Evaluate the coherence between the activation edge and the boundary shape. Since the boundary shape divides the entire space into inner and outer domains, most points are located within the boundary shape and are constrained by fidelity, as shown in formula (7):

[0054] U fidelity (x) = βU points (x)+(1-β)U walls (x)......(7)

[0055] The first term measures the percentage of input points surrounded by the boundary to retain different proportions of internal cells occupied by the input points, where β∈(0, 1); while the second term is added to constrain the edge by measuring the overlap between the edge and the actual wall;

[0056] Furthermore, in order to control the complexity of the final boundary shape, as shown in formula (8):

[0057]

[0058] In formula (8), |V| is the set of vertices to be calculated. This term returns a polygon with a small number of corner vertices. Γ is a discriminant used to determine the vertex v i Whether ∈V is a corner vertex; however, since the building facade may contain slopes, it is necessary to calculate the average height of each point and then compare it with the single-layer height value in the entity attributes of the joint surveying and mapping element layer;

[0059] Then, RANSAC is used to detect planes with inconsistent z-axis from the point cloud, the average height of the support point set corresponding to each plane is calculated, the points with average height lower than the height of a single layer are merged, and the face polygons are merged.

[0060] Furthermore, the step S4 specifically further includes the following steps:

[0061] In step S4.1, the joint mapping primitives are used to give instance labels for the layer units. Therefore, the layer segmentation problem is transformed into a Markov random field problem through the standard form of the energy problem. The energy function is shown in formula (9):

[0062]

[0063] In formula (9), γ∈(0,1) is a balance term, E is all adjacent internal faces, j and k represent specific faces, and m is the total number of faces; D(l k) assigns a label to the closed surface that is consistent with the layer number instance label; the paired item V(l j , l k ) first returns a flat cell with low complexity, and secondly, layers it by the shape within the boundary;

[0064] The final layer segmentation solves the above-mentioned Markov problem through graph cut optimization technology combined with surveying and mapping metadata. First, adjacent faces are merged into large polygonal faces, each representing a layer. Then, the floor plan is simplified by deleting vertices not connected to collinear edges to achieve accurate entity segmentation of building facade layers.

[0065] Advantages of the present invention:

[0066] This invention addresses the main technical bottlenecks in constructing detailed real-life 3D models in complex urban scenes, which are low data semantic accuracy, difficulty in automated modeling, and high point cloud annotation costs due to incomplete single data angles and insufficient data attributes.

[0067] Data from joint surveying and mapping (real estate surveying, completion acceptance surveying, and civil air defense surveying) is converted into point clouds, creating multimodal point cloud data that is then integrated with airborne LiDAR and oblique photography data to achieve the precise construction of physical models of buildings, floors, and even households. The real estate and completion acceptance survey results from the joint surveying and mapping are converted into 3D point cloud data, which is then optimized and constrained with the measured airborne LiDAR point cloud. The optimized data is then fitted into vector graphics and associated with semantic attributes to segment entities at different levels. This method supplements the missing airborne LiDAR point cloud data, also correcting the spatial shape and accuracy of buildings and enriching their semantic attributes. This method integrates thematic data with LiDAR point clouds to construct detailed physical models of buildings, floors, and households at the city level, fully leveraging the enormous value of natural resource inventory data. This method not only enables the construction of finer-grained and more accurate entities, resolving the challenges of automatically adapting geometric semantics, attribute semantics, and expression granularity when constructing detailed urban building physical models, but also fully leverages the rich thematic data already available on natural resources to achieve the physical construction of detailed 3D buildings in real scenes.

[0068]

Attached drawings

[0069] Figure 1 It is an implementation flow chart of the present invention. [Specific implementation method]

[0070] The present invention will be further described below with reference to specific examples.

[0071] A method for constructing building floor and household entity by integrating laser point cloud and joint surveying and mapping, such as Figure 1 As shown, the following steps are included:

[0072] Step S1. Convert the graphic outline into a point set

[0073] Convert the spatial surface primitives of real estate and completion acceptance measurement results into a set of points representing their contour features;

[0074] Step S1.1: First, in order to convert the joint surveying and mapping metadata into a 3D point cloud, it is necessary to calculate the polygonal approximate area edge for the metadata;

[0075] Step S1.2, then, convert the polygon base points into a contour point set, and calculate the center point for each polygon to construct a center point set to store the entity semantic attributes of the current polygon;

[0076] Among them, the Halcon operator is used to calculate the polygon matrix R, the rows and columns of each matrix represent the coordinates of the key points of the variable sequence, and each polygon represents a connected area in the image metadata;

[0077] Although the measurement error can be ignored due to the high accuracy of the surveying and mapping metadata, in order to adapt to a wider range of applications, the Halcon operator introduces an error parameter, the Hausdorff distance h, to calculate the maximum distance between the polygon area and the closed area of ​​the primitive, so that the two areas can intersect as much as possible, as shown in formula (1):

[0078]

[0079] In formula (1), U represents the area of ​​the survey and mapping metadata, and d(·) calculates the Euclidean distance between the two sets. In this way, the non-empty subset of the metric space itself is transformed into the metric space to facilitate subsequent calculations.

[0080] More importantly, the smaller the calculated error parameter is, the more contour key points are generated in the vector graphics. Therefore, setting the error parameter as small as possible generates more key points and generates a denser point cloud.

[0081] Step S2. Generate primitive center point set

[0082] The attributes of the completion acceptance measurement results and the real estate measurement results are semantically integrated based on spatial relationships, and the center points of the surface elements of the real estate measurement are calculated;

[0083] Assign the fused semantic attribute information to the center point. The center point attributes will contain the floor, room number, floor height and floor elevation value attribute information.

[0084] Step S2.1: In order to convert the polygons of the plane into spatial points, the spatial information they possess must be considered;

[0085] The joint mapping metadata attribute information contains accurate three-dimensional information. The center point of each polygon is calculated, and the elevation information in the corresponding mapping metadata attribute is assigned to the metadata, which is then converted into the center point of the corresponding key point polygon, and the elevation information is used as the attribute value of the key point.

[0086] Step S2.2: Divide the two-dimensional plane into a grid and fill the key points in the two-dimensional grid. However, since it is impossible to ensure that every point is on a regular grid point, it is necessary to resample the plane points and generate a completely regular dense grid point pattern.

[0087] The resampling problem is defined as an assignment problem. Each grid point is assigned a key point for filling, and the total movement cost is required to be minimized. The specific processing method is shown in formula (2):

[0088]

[0089] The purpose of formula (2) is to optimize a set of permutation matrices W g and W i l , to uniquely transform the global feature F G and local features Assigned to a predefined 2D grid of points where W g Relative to F G The mapping, Redundant grids represent grid points with only some of the grid points Form a corresponding relationship, N G is the dimension of the permutation matrix and the grid; Formula (2) ensures that in the global and local neighborhoods, only a subset of grid points form a one-to-one context point, and all missing matching grid points are filled in by the nearest neighbors. The matched points are mapped to the three-dimensional space according to their elevation values ​​to form a dense point cloud;

[0090] Step S3. Multimodal point cloud data fusion, compensation and correction of building outline points

[0091] Using a multimodal point cloud matching algorithm and a boundary contour constraint algorithm model, the point cloud generated above is fused with the point cloud generated by radar field scanning and the point cloud generated by the tilt model, and noise points are removed to correct the building outline points generated by the airborne radar point cloud.

[0092] The vertex set and edge set of the building's outer contour point cloud are used to combine with the floor point cloud of the jointly surveyed building to perform building outer contour correction;

[0093] Step S3.1, first, use the point cloud P scanned by the airborne laser radar sand point cloud P generated by joint mapping elements t Construct graph model G respectively s and graph model G t ;

[0094] The vertices in the two graphs are matched using the Hungarian algorithm. During the matching, the number of vertices in the point cloud contour graph must be greater than the number of vertices in the joint mapping primitives. The predicted points are linked to the primitive points using the Euclidean distance.

[0095] Step S3.2, use graph matching algorithm to find point cloud P s and point cloud P t matching relationship;

[0096] Use multi-layer perceptron to extract point cloud P s Point features and encode them into a one-dimensional point sequence;

[0097] For the point cloud P t , since it is a diffuse grid point formed by filling neighboring points, the structured feature sequence is extracted using circular convolution, and the first-order graph matching similarity matrix is ​​directly obtained by multiplying the first-order feature points of formula (3):

[0098] M p =U 1 U 2T ......(3)

[0099] In order to obtain the second-order similarity, we first connect the two nodes corresponding to each edge to construct the second-order feature matrix of the two graphs, as shown in formula (4):

[0100]

[0101] Among them, formula (3) is for the first-order similarity, U 1 and U 2 are the feature vectors of each node, respectively, and the superscripts indicate different data sources;

[0102] Formula (4) For the second-order similarity, X and Y are edge feature matrices, which are constructed by setting the descriptor of the cth edge starting from the i-th node and ending at the j-th node to the concatenation matrix of the descriptors of the two nodes, W 2 is a learnable parameter, so the second-order similarity has a learnable matching function; then, the graph matching problem is approximated as the problem of finding the eigenvector corresponding to the maximum eigenvalue of the similarity matrix, where power iteration is used for calculation; then, the graph matching problem is approximated as the problem of finding the eigenvector corresponding to the maximum eigenvalue of the similarity matrix, where power iteration is used for calculation;

[0103] The matching results of the graph are represented by a double random matrix, which intuitively reflects the possibility of establishing a matching relationship between any pair of nodes;

[0104] First, the matrix is ​​normalized by column and then by row. The double random matrix is ​​calculated through multiple iterations. However, since the values ​​of elements in the same row and column are not much different, in order to increase the obvious difference and facilitate the subsequent loss function calculation, the matrix is ​​multiplied by a constant and then processed using softmax to obtain the probability matrix of each candidate node match. This matrix still satisfies the double random property. Finally, the loss function is calculated using the offset vector-based method to calculate the offset between the matching points. Then, the weighted sum is taken according to the probability matrix to predict the offset vector and minimize the difference between the predicted vector and the true offset, as shown in formula (5):

[0105]

[0106] Among them, d i and Denote the predicted and actual offsets from the far point to the correct matching point, respectively. Φ(·) is a penalty term. In the final loss function calculation process, the discarded points will not affect the calculation of the loss function. Through the above method, the airborne lidar point cloud and the joint mapping building layer point cloud are matched and fused. The fused point cloud is projected onto a two-dimensional plane to extract the layered information.

[0107] The boundary outline of the building layer should satisfy the manifold assumption, where each vertex can only be connected to two adjacent edges;

[0108] Given the edge set generated in the previous step, select a subset of edges to calculate by using the constrained integer programming formula; i ∈{0, 1} is defined as a binary variable to describe the edge e i Whether ∈E is activated as part of the boundary state; the set of activated edges constitutes the polygonal boundary of the layer entity;

[0109] In addition, an energy function is used to constrain the generation quality of the boundary contour, as shown in formula (6):

[0110] U(x)=(1-λ)U fidelity (x)+λU complexity (x).......(6)

[0111] In formula (6), U fidelity (x) describes the consistency of state x with the input data, U complexity (x) measures the complexity of the output boundary shape, both of which are in the interval [0, 1], and λ∈(0, 1); U fidelity(x) Evaluate the coherence between the activation edge and the boundary shape. Since the boundary shape divides the entire space into inner and outer domains, most points are located within the boundary shape and are constrained by fidelity, as shown in formula (7):

[0112] U fidelity (x) = βU points (x)+(1-β)U walls (x)......(7)

[0113] The first term measures the percentage of input points surrounded by the boundary to retain different proportions of internal cells occupied by the input points, where β∈(0, 1); while the second term is added to constrain the edge by measuring the overlap between the edge and the actual wall;

[0114] Furthermore, in order to control the complexity of the final boundary shape, as shown in formula (8):

[0115]

[0116] In formula (8), |V| is the set of vertices to be calculated. This term returns a polygon with a small number of corner vertices. Γ is a discriminant used to determine the vertex v i Whether ∈V is a corner vertex; however, since the building facade may contain slopes, it is necessary to calculate the average height of each point and then compare it with the single-layer height value in the entity attributes of the joint surveying and mapping element layer;

[0117] Then, RANSAC is used to detect planes with inconsistent z-axis from the point cloud, and the average height of the support point set corresponding to each plane is calculated. The points with an average height lower than the height of a single layer are merged, and the face polygons are merged.

[0118] Step S4. Building outline segmentation and extraction

[0119] The semantic information of the point cloud converted from the joint surveying and mapping results is used as the label for point cloud segmentation and extraction. The fused point cloud data is segmented, and the edge contours of the buildings and the internal layering and household division lines are extracted to reconstruct the building entities, floors, and households.

[0120] In step S4.1, the joint mapping primitives give instance labels for the layer units. Therefore, the layer segmentation problem is transformed into a Markov random field problem through the standard form of the energy problem, as shown in formula (9):

[0121]

[0122] In formula (9), γ∈(0,1) is a balance term, E is all adjacent internal faces, j and k represent specific faces, and m is the total number of faces; D(l k) assigns a label to the closed surface that is consistent with the layer number instance label; the paired item V(l j , l k ) first returns a flat cell with low complexity;

[0123] Second, layering by shape within the boundary;

[0124] The final layer segmentation solves the above-mentioned Markov problem through graph cut optimization technology combined with surveying and mapping metadata. First, adjacent faces are merged into large polygonal faces, each representing a layer. Then, the floor plan is simplified by deleting vertices not connected to collinear edges to achieve accurate entity segmentation of building facade layers.

[0125] In this embodiment, data such as joint surveying and mapping (real estate surveying and mapping, completion acceptance measurement, and civil air defense measurement) are converted into point clouds to form multimodal point cloud data, which is then integrated with airborne lidar, oblique photography, and other data to achieve accurate construction of building models of buildings, floors, and even households; the real estate and completion acceptance measurement results in the joint surveying and mapping results are converted into three-dimensional point cloud data, and mutually optimized and constrained with the measured airborne lidar point clouds. The optimized data is fitted into vector graphics, and after being associated with semantic attributes, entities of different levels are segmented; this is used to supplement the missing airborne laser point cloud data, and also realizes the method of correcting the spatial shape and accuracy of buildings and enriching the semantic attributes of buildings, realizing the fusion of thematic data and laser point clouds to construct city-level fine building models of buildings, floors, and households, and giving full play to the huge value of natural resource stock data.

[0126] At the same time, it can not only achieve the construction of finer-grained and more precise entities, solve the problem of automatic adaptation of geometric semantics, attribute semantics, and expression granularity when constructing urban fine building entity models, but also make full use of the rich thematic data of natural resources to realize the entity construction of real-life three-dimensional fine buildings.

[0127] The above-described embodiments are merely preferred embodiments of the present invention and are not intended to limit the scope of implementation of the present invention. Except for the cases listed in the specific embodiments, all equivalent changes made according to the methods and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A method for constructing building floor and household entity by laser point cloud fusion and joint surveying and mapping, characterized in that: The following steps are involved: Step S1. Convert the graphic outline into a point set Convert the spatial surface primitives of real estate and completion acceptance measurement results into a set of points representing their contour features; Step S2. Generate primitive center point set The attributes of the completion acceptance measurement results and the real estate measurement results are semantically integrated based on spatial relationships, and the center points of the surface elements of the real estate measurement are calculated; Assign the fused semantic attribute information to the center point. The center point attributes will contain the floor, room number, floor height and floor elevation value attribute information. Step S3. Multimodal point cloud data fusion, compensation and correction of building outline points Using a multimodal point cloud matching algorithm and a boundary contour constraint algorithm model, the point cloud generated above is fused with the point cloud generated by radar field scanning and the point cloud generated by the tilt model, and noise points are removed to correct the building outline points generated by the airborne radar point cloud. Step S4. Building outline segmentation and extraction The point cloud semantic information converted from the joint surveying and mapping results is used as the label for point cloud segmentation and extraction. The fused point cloud data is segmented, and the building edge contours and internal stratification and household lines are extracted to reconstruct the building, floor, and household entities.

2. The method for constructing building-floor-unit-level entities by laser point cloud fusion and joint mapping according to claim 1 is characterized in that: The step S1 specifically further includes the following steps: Step S1.1: First, in order to convert the joint surveying and mapping metadata into a 3D point cloud, it is necessary to calculate the polygonal approximate area edge for the metadata; Step S1.2, then, convert the polygon base points into a contour point set, and calculate the center point for each polygon to construct a center point set to store the entity semantic attributes of the current polygon; Among them, the Halcon operator is used to calculate the polygon matrix R, the rows and columns of each matrix represent the coordinates of the key points of the variable sequence, and each polygon represents a connected area in the image metadata; Although the measurement error can be ignored due to the high accuracy of the surveying and mapping metadata, in order to adapt to a wider range of applications, the Halcon operator introduces an error parameter, the Hausdorff distance h, to calculate the maximum distance between the polygon area and the closed area of ​​the primitive, so that the two areas intersect as much as possible, as shown in formula (1): In formula (1), U represents the area of ​​the survey and mapping metadata, and d(·) calculates the Euclidean distance between the two sets. In this way, the non-empty subset of the metric space itself is transformed into the metric space to facilitate subsequent calculations. More importantly, the smaller the calculated error parameter is, the more contour key points are generated in the vector graphics. Therefore, setting the error parameter as small as possible generates more key points and generates a denser point cloud.

3. The method for constructing building-floor-unit-level entities by laser point cloud fusion and joint mapping according to claim 1 is characterized in that: The step S2 specifically further includes the following steps: Step S2.1: In order to convert the polygons of the plane into spatial points, the spatial information they possess must be considered; The joint mapping metadata attribute information contains accurate three-dimensional information. The center point of each polygon is calculated, and the elevation information in the corresponding mapping metadata attribute is assigned to the metadata, which is then converted into the center point of the corresponding key point polygon, and the elevation information is used as the attribute value of the key point. Step S2.2: Divide the two-dimensional plane into a grid and fill the key points in the two-dimensional grid. However, since it is impossible to ensure that every point is on a regular grid point, it is necessary to resample the plane points and generate a completely regular dense grid point pattern. The resampling problem is defined as an assignment problem. Each grid point is assigned a key point for filling, and the total movement cost is required to be minimized. The specific processing method is shown in formula (2): The purpose of formula (2) is to optimize a set of permutation matrices W g and W i l , to uniquely transform the global feature F G and local features Assigned to a predefined 2D grid of points where W g Relative to F G The mapping, Redundant grids represent grid points with only some of the grid points Form a corresponding relationship, N G is the dimension of the permutation matrix and the grid; Formula (2) ensures that in the global and local neighborhoods, only a subset of grid points form a one-to-one context point, and all missing matching grid points are filled by the nearest neighbors. The matched points are mapped to the three-dimensional space according to their elevation values ​​to form a dense point cloud.

4. The method for constructing building-floor-unit-level entities by laser point cloud fusion and joint mapping according to claim 1 is characterized in that: The step S3 specifically further includes the following steps: The vertex set and edge set of the building's outer contour point cloud are used to combine with the floor point cloud of the jointly surveyed building to perform building outer contour correction; Step S3.1, first, use the point cloud P scanned by the airborne laser radar s and point cloud P generated by joint mapping elements t Construct graph model G respectively s and graph model G t ; The vertices in the two graphs are matched using the Hungarian algorithm. During the matching, the number of vertices in the point cloud contour graph must be greater than the number of vertices in the joint mapping primitives. The predicted points are linked to the primitive points using the Euclidean distance. Step S3.2, use graph matching algorithm to find point cloud P s and point cloud P t matching relationship; Use multi-layer perceptron to extract point cloud P s Point features and encode them into a one-dimensional point sequence; For the point cloud P t , since it is a diffuse grid point formed by filling neighboring points, the structured feature sequence is extracted using circular convolution, and the first-order graph matching similarity matrix is ​​directly obtained by multiplying the first-order feature points of formula (3): M p =U 1 U 2T ...... (3) In order to obtain the second-order similarity, we first connect the two nodes corresponding to each edge to construct the second-order feature matrix of the two graphs, as shown in formula (4): Among them, formula (3) is for the first-order similarity, U 1 and U 2 are the feature vectors of each node, respectively, and the superscripts indicate different data sources; Formula (4) For the second-order similarity, X and Y are edge feature matrices, which are constructed by setting the descriptor of the cth edge starting from the i-th node and ending at the j-th node to the concatenation matrix of the descriptors of the two nodes, W 2 is a learnable parameter, so the second-order similarity has a learnable matching function. Then, the graph matching problem is approximated as the problem of finding the eigenvector corresponding to the maximum eigenvalue of the similarity matrix, which is calculated using power iteration. The matching results of the graph are represented by a double random matrix, which intuitively reflects the possibility of establishing a matching relationship between any pair of nodes; First, the matrix is ​​normalized by column and then by row. The double random matrix is ​​calculated through multiple iterations. However, since the values ​​of elements in the same row and column are not much different, in order to increase the obvious difference and facilitate the subsequent loss function calculation, the matrix is ​​multiplied by a constant and then processed using softmax to obtain the probability matrix of each candidate node match. This matrix still satisfies the double random property. Finally, the loss function is calculated using the offset vector-based method to calculate the offset between the matching points. Then, the weighted sum is taken according to the probability matrix to predict the offset vector and minimize the difference between the predicted vector and the true offset. The loss function is shown in formula (5): Among them, d i and Denote the predicted and actual offsets from the far point to the correct matching point, respectively. Φ(·) is a penalty term. In the final loss function calculation process, the discarded points will not affect the calculation of the loss function. Through the above method, the airborne lidar point cloud and the joint mapping building layer point cloud are matched and fused. The fused point cloud is projected onto a two-dimensional plane to extract the layered information. The boundary outline of the building layer should satisfy the manifold assumption, where each vertex can only be connected to two adjacent edges; Given the edge set generated in the previous step, select a subset of edges to calculate by using the constrained integer programming formula; i ∈{0, 1} is defined as a binary variable to describe the edge e i Whether ∈E is activated as part of the boundary state; the set of activated edges constitutes the polygonal boundary of the layer entity; In addition, an energy function is used to constrain the generation quality of the boundary contour, as shown in formula (6): U(x)=(1-λ)U fidelity (x)+λU complexity (x)...... (6) In formula (6), U fidelity (x) describes the consistency of state x with the input data, U complexity (x) measures the complexity of the output boundary shape, both of which are in the interval [0, 1], and λ∈(0, 1); U fidelity (x) Evaluate the coherence between the activation edge and the boundary shape. Since the boundary shape divides the entire space into inner and outer domains, most points are located within the boundary shape and are constrained by fidelity, as shown in formula (7): U fidelity (x)=βU points (x)+(1-β)U walls (x)...... (7) The first term measures the percentage of input points surrounded by the boundary to retain different proportions of internal cells occupied by the input points, where β∈(0, 1); while the second term is added to constrain the edge by measuring the overlap between the edge and the actual wall; Furthermore, in order to control the complexity of the final boundary shape, as shown in formula (8): In formula (8), |V| is the set of vertices to be calculated. This term returns a polygon with a small number of corner vertices. Γ is a discriminant used to determine the vertex v i Whether ∈V is a corner vertex; however, since the building facade may contain slopes, it is necessary to calculate the average height of each point and then compare it with the single-layer height value in the entity attributes of the joint surveying and mapping element layer; Then, RANSAC is used to detect planes with inconsistent z-axis from the point cloud, the average height of the support point set corresponding to each plane is calculated, the points with average height lower than the height of a single layer are merged, and the face polygons are merged.

5. The method for constructing building-floor-unit-level entities by laser point cloud fusion and joint surveying and mapping according to claim 1 is characterized in that: The step S4 specifically further includes the following steps: In step S4.1, the joint mapping primitives are used to give instance labels for the layer units. Therefore, the layer segmentation problem is transformed into a Markov random field problem through the standard form of the energy problem. The energy function is shown in formula (9): In formula (9), γ∈(0,1) is a balance term, E is all adjacent internal faces, j and k represent specific faces, and m is the total number of faces; D(l k ) assigns a label to the closed surface that is consistent with the layer number instance label; the paired item V(l j , l k ) first returns a flat cell with low complexity, and secondly, layers it by the shape within the boundary; The final layer segmentation solves the above-mentioned Markov problem through graph cut optimization technology combined with surveying and mapping metadata. First, adjacent faces are merged into large polygonal faces, each representing a layer. Then, the floor plan is simplified by deleting vertices not connected to collinear edges to achieve accurate entity segmentation of building facade layers.

Images