Urban aerial image and LiDAR point cloud registration method using multivariate features

Abstract

The invention discloses an urban aerial image and LiDAR point cloud registration method using multivariate features, and belongs to the field of photogrammetry and remote sensing, and the method comprises the steps: carrying out the preprocessing of an image, generating a sparse feature point cloud and a three-dimensional feature line, extracting the contour of an image building roof through a deep learning method, and carrying out the recognition of the image building roof; extracting feature lines and feature surfaces of the LiDAR point cloud by using an RANSAC (Random Sample Consensus) algorithm and a PCA (Principal Component Analysis) statistical analysis method; global coarse registration is carried out on sparse feature point clouds and LiDAR point clouds generated by images through an ICL algorithm based on line and surface features; image block adjustment calculation is carried out in combination with point, line and surface multivariate control information constraints of the LiDAR point cloud, and local distortion of the image block is inhibited; and the global coarse registration and the local distortion optimization are repeatedly and iteratively carried out until the variable quantity of the error in the image orientation meets a preset requirement. According to the method, point, line and surface features are combined to form multivariate constraints, and a self-adaptive weighting strategy is adopted, so that the problem of high-precision robust registration of images and LiDAR point clouds in an urban scene can be well solved.

Description

Technical Field

[0001] The present invention belongs to the field of photogrammetry and computer vision, and specifically relates to a method for aligning urban aerial images and LiDAR point clouds using multivariate features. The key technology of the present invention is the high-precision automatic registration of aerial images and LiDAR point clouds that combines multivariate feature control information of urban scene points, lines, and surfaces. Background Art

[0002] With the accelerating pace of urbanization, the demand for detailed three-dimensional spatial information in cities is rapidly increasing. Aerial photogrammetry, a traditional method for acquiring 3D spatial information, has long relied heavily on high-precision field control points to achieve precise geometric orientation. However, the high cost, long acquisition cycle, and low degree of automation of field control points have become major bottlenecks hindering the efficient and intelligent development of aerial image processing.

[0003] The rapid development and widespread application of LiDAR technology has provided new opportunities for high-precision urban geospatial information. LiDAR can directly acquire high-precision, high-density three-dimensional point cloud data. The spatial positioning accuracy of the point cloud is usually within 10 centimeters, which can effectively replace traditional field control points as a source of control information. However, due to the differences in sensor modality and acquisition time, there are significant modal differences and spatial deviations between aerial imagery and LiDAR point clouds. In addition, complex urban scenes are rich in geometric features. Although traditional single-type feature registration methods (such as ICP and ICL) can be used in the coarse registration stage, they lack stability and robustness in complex scenes, especially in the effective use of surface features.

[0004] Existing "cloud-controlled" aerial photogrammetry methods for LiDAR point clouds mostly focus on the registration of point features or preliminary line features, neglecting the more stable surface features in urban environments, particularly the matching of building rooftop features. Furthermore, existing technologies underutilize artificial intelligence, particularly deep learning methods, for feature extraction and registration, resulting in a low degree of automation and inefficiency in the feature extraction process. Therefore, how to integrate advanced artificial intelligence methods to target the complex characteristics of urban scenes, mine and effectively utilize line and surface features, and achieve precise, efficient, and robust automatic registration of aerial imagery with LiDAR point clouds has become a critical issue that urgently needs to be addressed. Summary of the Invention

[0005] In order to solve the problems of insufficient geometric registration accuracy and poor robustness of single features between urban aerial images and LiDAR point clouds in the existing technology, the purpose of the present invention is to provide a high-precision urban aerial image and LiDAR point cloud registration method that integrates multiple features and artificial intelligence technology, especially for complex urban scenes. The method focuses on the fine matching process of line and surface features. First, deep learning and statistical analysis methods are used to extract rich point, line and surface features from images and LiDAR point clouds respectively. Then, the iterative nearest neighbor line (ICL) algorithm based on line and surface constraints is used to achieve global coarse registration of image point clouds and LiDAR point clouds. Then, the multiple features of urban scenes are further combined, and local distortion is eliminated through regional block adjustment, thereby significantly improving the accuracy and robustness of registration.

[0006] The present invention is mainly achieved through the following technical solutions:

[0007] A method for registering urban aerial images and LiDAR point clouds using multivariate features includes the following steps:

[0008] Step 1: Preprocess the image to generate sparse feature point clouds and image 3D feature lines, and use an improved deep learning model to extract the building roof outline information of the image;

[0009] Step 2: Using the statistical feature analysis method of random sampling consensus algorithm and principal component analysis, line feature and surface feature information are extracted from the LiDAR point cloud;

[0010] Step 3: Using the improved ICL algorithm, we use the correspondence between the image 3D feature lines and the LiDAR 3D feature lines, and add the projection constraints of the image roof outline and the LiDAR building plane to perform global coarse registration;

[0011] Step 4: Based on step 3, the multivariate features of points, lines, and surfaces extracted from the urban scene are used as joint control information to construct an extended energy equation for image block adjustment optimization to correct local deformation caused by camera distortion and uneven image connection strength.

[0012] Step 5: Iterate steps 3 and 4 multiple times until the change in the error in image positioning meets the preset requirements.

[0013] Optionally, in step 1, the image is preprocessed to generate a sparse feature point cloud and image three-dimensional feature lines, including: using the positioning and posture information obtained by the POS system to perform preliminary calculations of internal and external orientation elements of the image, and combining a feature matching algorithm to generate a sparse feature point cloud and image three-dimensional feature lines, wherein the sparse feature point cloud is generated by matching the same-name feature points between images, and the three-dimensional feature lines are obtained by matching the same-name line segments between images and performing three-dimensional reconstruction;

[0014] The improved deep learning model is used to extract building roof outline information from images, including: using a deep convolutional neural network model improved based on U-Net, in which the encoder uses ResNet-50 to replace the traditional convolution module to enhance feature expression capabilities, and the decoder introduces an attention mechanism to focus on building edge details and reduce background interference;

[0015] The attention mechanism calculation process is expressed as:

[0016]

[0017] in, represents the attention weight, is the high-level feature passed by the encoder, is the low-level feature passed to the decoder, 、 is the learnable convolution kernel weight, is the bias, is the Sigmoid activation function;

[0018] The binary cross entropy loss function is used during training:

[0019]

[0020] in, is the true label, is the predicted label output by the network, is the total number of samples. By minimizing the loss function, the network model can effectively learn and stably output high-precision building outlines.

[0021] Optionally, in step 2, extracting line feature and surface feature information from the LiDAR point cloud using a statistical feature analysis method based on a random sampling consensus algorithm and principal component analysis includes:

[0022] The RANSAC algorithm is used to perform a preliminary fit on the linear structure in the point cloud data. The specific implementation includes: randomly selecting a minimum number of point sets from the LiDAR point cloud, fitting the initial line segment equation, and then defining a distance threshold. All points with a distance less than the threshold are considered as supporting inliers of the model. This process is iterated until the number of inliers reaches the maximum, thereby determining the optimal linear model. The specific three-dimensional straight line equation is expressed as:

[0023]

[0024] Where, is the coordinate of a point on the line, and the direction vector is ;

[0025] After the RANSAC algorithm optimizes and obtains the optimal parameters of the line model, the PCA method is used to further analyze the spatial distribution characteristics of the inner point set and determine the most stable three-dimensional feature line:

[0026]

[0027] Where, represents the coordinates of the inner points of the line feature, represents the centroid of the interior points, is the number of interior points; the matrix The principal eigenvector corresponds to the maximum eigenvalue, and its direction is the most stable direction of the line segment;

[0028] Three non-collinear points are randomly selected to construct the initial plane equation, and then the plane equation is optimized iteratively to make more points meet the condition that the distance from the plane is less than a given threshold. The expression of the plane equation is:

[0029]

[0030] in, is the normal vector of the plane, is a constant term;

[0031] Each iteration calculates the perpendicular distance from the point to the plane:

[0032]

[0033] in, is the coordinate of any point in the LiDAR point cloud. The goal of iterative optimization is to minimize the negative value of the number of points that meet the threshold distance to obtain the optimal plane model with the most inliers.

[0034] Optionally, in step 3, the improved ICL algorithm is used to perform a global coarse registration using the correspondence between the image 3D feature lines and the LiDAR 3D feature lines, and adding the projection constraint between the image roof outline and the LiDAR building plane, including:

[0035] The 3D feature lines of the image extracted in step 1 are matched with the 3D feature lines of the LiDAR extracted in step 2 to establish an initial correspondence between the two sets of data. The rigid transformation parameters are then iteratively solved using the ICL algorithm. The core objective function of the improved ICL algorithm is to minimize the matching error between the 3D feature lines.

[0036] The specific formula is described as follows:

[0037] Define the image feature line set as , the LiDAR feature line set is { }, there is a preliminary correspondence between the two sets of characteristic lines, and the optimization objective function of the ICL algorithm is expressed as:

[0038]

[0039] Where, and Represent the three-dimensional feature lines corresponding to the image and LiDAR, represents the rotation matrix, represents the translation vector, Indicates the number of corresponding characteristic line pairs;

[0040] By iteratively solving the above equations, we can obtain the optimal and , complete the rigid registration of feature lines between the image and LiDAR point cloud;

[0041] Image contour features , the building plane features corresponding to LiDAR The currently obtained transformation matrix is ​​projected onto the spatial plane where the image contour is located to form a projection contour line , thus forming the surface contour projection error constraint function:

[0042]

[0043] in, Indicates the number of matching pairs between image contours and LiDAR surface features.

[0044] Furthermore, by combining the two constraints of line and surface, the line and surface feature constraints are integrated into a unified optimization framework, and the final overall objective function is expressed as:

[0045]

[0046] In the formula, the weight factor Used to adjust the relative strength of surface feature constraints and line feature constraints.

[0047] Optionally, in step 4, the multivariate features of points, lines, and surfaces extracted from the urban scene are used as joint control information to construct an extended energy equation for image block adjustment optimization, wherein the specific implementation of the point feature constraint includes:

[0048] The principle of the shortest distance between points and surfaces is used to select the "control points" with the same name as each sparse feature point in the LiDAR point cloud. For any sparse feature point , search for the LiDAR point cloud within the spherical neighborhood with a radius of r, and record the nearest neighbor point cloud set as , For the The number of LiDAR point clouds within the neighborhood of sparse feature points; according to the point cloud collection All points in the space are fitted using the least squares method to obtain the spatial plane , with sparse feature points To plane The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud; the geometric control constraint of the point feature It can be expressed as follows

[0049]

[0050] Where, is the prior weight, is the function for calculating the distance between spatial points and surfaces, is the average ground resolution of the image, used to convert Normalize to image unit.

[0051] Furthermore, the specific implementation of the line feature constraint in step 4 includes:

[0052] The principle of the shortest distance between points and lines is used to select the "control points" with the same name as each sparse feature point in the LiDAR point cloud; let the three-dimensional feature lines of the image and LiDAR point cloud extracted in step 1 and step 2 be and , and projected onto the horizontal plane as and , The point on (t=1,2,…,n), where n is The number of points on the image feature line; the plane projection points on the image feature line To the straight line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud;

[0053] Geometric control constraints for line features It is expressed as follows:

[0054]

[0055] Where, is the plane projection point on the image feature line To LiDAR projection line The shortest Euclidean distance, is the average ground resolution of the image, It is The weight corresponding to the line feature.

[0056] Furthermore, the specific implementation of the face feature constraint in step 4 includes:

[0057] Let the building roofs in the image and LiDAR point cloud extracted in step 1 and step 2 be and ,Will Projected onto the horizontal plane, it is recorded as For plane The edge of the contour is marked with a positioning point every certain number of pixels. (s=1,2,…,u), where u is the number of positioning points on the edge of the contour; To plane The nearest edge contour line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud;

[0058] Geometric control constraints for face features Expressed as:

[0059]

[0060] Where, Indicates that in the horizontal plane, the image contour points To LiDAR projected contour edge segment The shortest Euclidean distance, is the average ground resolution of the image, It is The weight of each face feature.

[0061] Furthermore, in step 4, the energy equation of the classic block adjustment is expanded to construct the energy equation shown in Equation (13):

[0062]

[0063] Where, is the image point residual term defined in the classic block adjustment, For the The orientation parameters of the image, For its Observation of image points on an image, Projection functions defined for collinearity condition equations; It is a geometric constraint formed by combining point, line and surface feature information.

[0064] Furthermore, appropriate prior weights are assigned to the geometric constraints formed by the joint point, line, and surface feature information to maximize the registration accuracy. The specific implementation includes:

[0065] For point and surface feature constraints, use the fitted plane and The local coplanarity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, and , taking the spatial plane fitting error as an indicator to measure coplanarity;

[0066] For line feature constraints, use the fitted line The local collinearity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, , taking the straight line fitting error as an indicator of collinearity;

[0067] When the fitting error of the spatial plane and line is small, the corresponding LiDAR control information is given a larger prior weight; otherwise, a smaller prior weight is given.

[0068] The present invention has the following advantages: (1) By combining point, line and surface features to form multivariate constraints, the problem of high-precision robust registration of images and LiDAR point clouds is solved, thereby achieving high-precision geometric orientation of aerial images of urban scenes without field control points; (2) The use of an adaptive weighting strategy can solve the problem of inconsistent accuracy and reliability of multivariate control information, thereby maximizing the registration accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0069] Figure 1 is a flow chart of the method of the present invention;

[0070] Figure 2 This is a schematic diagram of point feature constraints of the present invention;

[0071] Figure 3 Schematic diagram of line feature constraints of the present invention;

[0072] Figure 4 Schematic diagram of surface feature constraints of the present invention;

[0073] Figure 5 Schematic diagram of sparse feature points of different types of images of the present invention;

[0074] Figure 6 This is a schematic diagram of the before and after registration of aerial images and LiDAR point clouds according to the present invention. In the figure:

[0075] (a) Typical region 1 before registration, (b) Typical region 1 after registration;

[0076] (c) Typical region 2 before registration, (d) Typical region 2 after registration;

[0077] (e) Typical region 3 before registration, (f) Typical region 3 after registration. DETAILED DESCRIPTION

[0078] The technical solution of the present invention will be further specifically described below through embodiments and in conjunction with the accompanying drawings.

[0079] The present invention provides a registration method for urban aerial images and LiDAR point clouds using multivariate features. First, artificial intelligence and statistical analysis methods are used to extract rich geometric features of points, lines, and surfaces in urban scenes from aerial images and LiDAR point clouds respectively. Secondly, an ICL algorithm based on line and surface features is used to achieve efficient coarse registration between aerial image feature point clouds and LiDAR point clouds. Then, the extended regional block adjustment energy equation is constructed based on the extracted multivariate features. Local distortion is eliminated by combining point, line, and surface control information, ultimately achieving high-precision and robust registration of images and point clouds. The overall technical route is shown in the attached figure. Figure 1 As shown, the core process includes the following steps:

[0080] Step 1: Preprocess the image to generate sparse feature point clouds and image 3D feature lines, and use an improved deep learning model to extract building roof outline information in the image;

[0081] In step 1, the present invention first uses POS-assisted aerotriangulation technology to perform geometric preprocessing on the original aerial image. This process uses the positioning and attitude information obtained by the POS system to perform preliminary calculations of the image's internal and external orientation elements. This process, combined with a feature matching algorithm, generates a sparse feature point cloud and three-dimensional feature lines. The sparse feature point cloud is generated by matching the same-name feature points between images, while the three-dimensional feature lines are obtained by matching the same-name line segments between images and performing three-dimensional reconstruction. This provides the initial geometric basis for the subsequent registration process.

[0082] Specifically, for any two-dimensional feature point on an image , its corresponding three-dimensional coordinates The classical collinear equations of photogrammetry are solved to establish the geometric correspondence between the two-dimensional points in the image and the three-dimensional points in the object space, so as to obtain the preliminary three-dimensional geometric information of the image.

[0083] Secondly, to accurately extract the unique rooftop features of buildings in urban scenes, this paper uses deep learning methods for image semantic segmentation. Specifically, it employs a deep convolutional neural network model based on an improved U-Net. The encoder in the network architecture replaces the traditional convolutional module with a ResNet-50 to enhance feature representation. The decoder incorporates an attention gate to focus on building edge details and reduce background interference, effectively improving the accuracy of automatic extraction of urban building outlines.

[0084] During the forward propagation of the network, the feature map The calculation process of the spatial attention mechanism can be expressed as:

[0085]

[0086] in, represents the attention weight, is the high-level feature passed by the encoder, is the low-level feature passed to the decoder, 、 is the learnable convolution kernel weight, is the bias, is the Sigmoid activation function. The attention mechanism significantly enhances the feature response of the target area (building roof outline) through this formula, weakens the features of non-target areas, and significantly improves the network's contour extraction performance.

[0087] The data required for network training comes from the ISPRS public standard dataset and the annotated dataset built by our research team. The binary cross-entropy loss function (BCE) is used in the training process:

[0088]

[0089] in, is the true label, is the predicted label output by the network, is the total number of samples. By minimizing the loss function, the network model can effectively learn and stably output high-precision building outlines.

[0090] Step 2: Use the statistical feature analysis method of random sampling consensus algorithm (RANSAC) and principal component analysis (PCA) to automatically extract line feature and surface feature information from the LiDAR point cloud;

[0091] In step 2, in order to achieve high-precision automatic registration of images and LiDAR point clouds, the present invention uses the RANSAC algorithm and PCA analysis method to automatically and efficiently extract line features and surface features (especially building plane features) from LiDAR point cloud data.

[0092] First, during the extraction of 3D feature lines, the RANSAC algorithm is used to perform a preliminary fit to the linear structure in the point cloud data. Specifically, this is accomplished by randomly selecting a minimum set of points from the LiDAR point cloud and fitting the initial line segment equation. A distance threshold is then defined (for example, if the image ground resolution is 10 cm, the threshold is 15 cm). All points whose distance to the line segment is less than the threshold are considered inliers for the model. This process is iterated until the number of inliers reaches a maximum, thereby determining the optimal linear model. The specific 3D line equation is expressed as:

[0093]

[0094] Where, is the coordinate of a point on the line, and the direction vector is After the RANSAC algorithm optimizes and obtains the optimal parameters of the line model, the PCA method is used to further analyze the spatial distribution characteristics of the inner point set and determine the most stable three-dimensional feature line. Specifically, PCA is used to solve the eigenvalues ​​and eigenvectors of the covariance matrix of the inner point set:

[0095]

[0096] Where, represents the coordinates of the inner points of the line feature, represents the centroid of the interior points, is the number of interior points. The principal eigenvector corresponds to the maximum eigenvalue, and its direction is the most stable direction of the segment.

[0097] Secondly, this paper uses a similar RANSAC strategy to extract planar features from LiDAR point clouds, especially for prominent planar areas such as urban building roofs. The specific method is: randomly select three non-collinear points to construct an initial plane equation, and then iteratively optimize the plane equation to make more points meet the condition that the distance from the plane is less than a given threshold. The expression of the plane equation is:

[0098]

[0099] in, is the normal vector of the plane, Is a constant term. Each round of iteration calculates the vertical distance from the point to the plane:

[0100]

[0101] in, are the coordinates of any point in the LiDAR point cloud. The iterative optimization goal is to minimize the number of negative points within the threshold distance to obtain the optimal plane model with the most inliers. With the plane parameters thus determined, a large number of stable planar features, such as building roof planes, can be automatically and efficiently extracted from the LiDAR point cloud.

[0102] Step 3: Using the improved ICL algorithm, we exploit the correspondence between the 3D feature lines of the image and the 3D feature lines of the LiDAR, and add the projection constraints of the image roof outline and the LiDAR building plane to achieve fast and robust global coarse registration.

[0103] In step 3, the present invention uses an improved ICL algorithm, that is, by matching the 3D feature lines of the point cloud with the 3D feature lines of the image, the optimal rigid transformation matrix between the two sets of data is determined. Specifically, the 3D feature lines of the image extracted in step 1 are matched with the 3D feature lines of the LiDAR extracted in step 2 to establish an initial correspondence between the two sets of data. Then, the ICL algorithm is used to iteratively solve the rigid transformation parameters (rotation matrix and translation vectors The core objective function of the improved ICL algorithm is to minimize the matching error between three-dimensional feature lines. The specific formula is as follows:

[0104] Define the image feature line set as , the LiDAR feature line set is { }, there is a preliminary correspondence between the two sets of characteristic lines, and the optimization objective function of the ICL algorithm is expressed as:

[0105]

[0106] Where, and Represent the three-dimensional feature lines corresponding to the image and LiDAR, represents the rotation matrix, represents the translation vector, Represents the number of corresponding line feature pairs. By iteratively solving the above equations, the optimal and , complete the rigid registration of feature lines between the image and the LiDAR point cloud.

[0107] Secondly, the present invention further introduces surface feature constraints to improve the stability of the registration results. Specifically, the image building contour features extracted based on deep learning in step 1 are combined with the LiDAR building plane features extracted in step 2 to perform contour projection constraint matching. The specific implementation method is:

[0108] Image contour features , the building plane features corresponding to LiDAR The currently obtained transformation matrix is ​​projected onto the spatial plane where the image contour is located to form a projection contour line . This forms the surface contour projection error constraint function:

[0109]

[0110] in, Represents the number of matching pairs between image contours and LiDAR surface features. This formula further eliminates registration deviations caused by modal differences and systematic errors, improving the algorithm's adaptability and stability in complex urban environments.

[0111] Combining the above two constraints of line and surface, the present invention integrates the line and surface feature constraints into a unified optimization framework, and the final overall objective function is expressed as:

[0112]

[0113] In the formula, the weight factor It is used to adjust the relative strength of surface feature constraints and feature line constraints to ensure that the coarse registration process can obtain stable and robust results in different urban scenes.

[0114] Step 4: Using the multivariate features of points, lines, and surfaces extracted from urban scenes as joint control information, an extended energy equation is constructed to perform image block adjustment optimization to correct for camera distortion and local deformation caused by uneven image connection strength.

[0115] In step 4, the image network is adjusted by combining the multivariate control information constraints of the point, line, and surface features of the urban scene to solve the camera distortion parameters, image internal orientation elements, and image external orientation elements. This suppresses the local deformation and distortion of the image network caused by factors such as camera distortion, inconsistent image connection strength, and the resulting error accumulation, and obtains a more reliable registration result.

[0116] The process of point feature constraint is as follows: Figure 2 As shown in the figure, the point-to-plane distance principle is used to select the "control point" with the same name as each sparse feature point in the LiDAR point cloud. , search for the LiDAR point cloud within the spherical neighborhood with a radius of r, and record the nearest neighbor point cloud set as , For the The number of LiDAR point clouds within the neighborhood of sparse feature points. All points in the space are fitted using the least squares method to obtain the spatial plane , with sparse feature points To plane The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud.

[0117] Unlike field control points, there is no strict physical correspondence between image sparse feature points and LiDAR control information. Therefore, they can only be used as geometric constraints and cannot be added as control points to the regional block adjustment solution. Geometric control constraints of point features It can be expressed as follows:

[0118]

[0119] Where, is the prior weight, is the function for calculating the distance between spatial points and surfaces, is the average ground resolution of the image, used to convert Normalize to image unit.

[0120] The process of line feature constraint is as follows: Figure 3 As shown in the figure, the point-to-line distance principle is used to select the "control point" with the same name as each sparse feature point in the LiDAR point cloud. Let the three-dimensional feature lines of the image and LiDAR point cloud extracted in step 1 and step 2 be and , projected onto the horizontal plane (for plane registration) and recorded as and , The point on (t=1,2,…,n), where n is The number of points on the image feature line. To the straight line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud.

[0121] Similarly, there is no strict physical correspondence between the 3D feature lines of the image and the feature lines of the LiDAR point cloud, and the line feature control information can only be used as a geometric constraint. It can be expressed as follows:

[0122]

[0123] In formula (11), is the plane projection point on the image feature line To LiDAR projection line The shortest Euclidean distance, is the average ground resolution of the image (used for distance normalization), and The first The weight corresponding to the line feature.

[0124] The process of surface feature constraint is as follows: Figure 4 As shown, let the building roofs in the images and LiDAR point clouds extracted in step 1 and step 2 be and ,Will Projected onto the horizontal plane, it is recorded as For plane The edge of the contour is marked with a positioning point every certain number of pixels. (s=1,2,…,u), where u is the number of anchor points on the edge of the contour (including corner points). To plane The nearest edge contour line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud.

[0125] Similarly, there is no strict physical correspondence between the characteristic surfaces of image buildings and the characteristic surfaces of LiDAR point clouds, and the surface feature control information can only be used as geometric constraints. It can be expressed as follows:

[0126]

[0127] In formula (12), Indicates that in the horizontal plane, the image contour points To LiDAR projected contour edge segment The shortest Euclidean distance, is the average ground resolution of the image, The first The weight of each face feature.

[0128] The energy equation of the classic block adjustment is expanded by adding control constraints combining point, line, and surface features, and constructing the energy equation shown in Equation (13).

[0129]

[0130] Where, is the image point residual term defined in the classic block adjustment, For the The orientation parameters of the image, For its Observation of image points on an image, Projection functions defined for collinearity condition equations; It is a geometric constraint formed by combining point, line and surface feature information.

[0131] Sparse image feature points originate from diverse ground objects, such as buildings, ground, and vegetation. These corresponding real-world scenes are complex, and their neighboring LiDAR point sets possess varying local characteristics. This leads to varying errors when performing line and plane fitting on LiDAR point clouds, resulting in varying accuracy in LiDAR control information. Therefore, it is necessary to assign appropriate prior weights to LiDAR control information to maximize registration accuracy.

[0132] For point and surface feature constraints, use the fitted plane and The local coplanarity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, and , the spatial plane fitting error is used as an indicator to measure the coplanarity. For line feature constraints, the fitting line is used The local collinearity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, , the straight line fitting error is used as an indicator to measure collinearity. Figure 5 As shown in the figure, when the fitting error of spatial planes and lines is small, the corresponding LiDAR control information is given a larger prior weight; otherwise, a smaller prior weight is given. During the adjustment calculation, based on the theory of posterior variance component estimation, the posterior weights of various control information are automatically estimated. At the same time, control observations with large differences from the prior weights are eliminated to ensure the accuracy and consistency of the regional block adjustment results.

[0133] Step 5: Repeat steps 3 and 4 iteratively until the change in image positioning error meets the set requirements (less than 0.001 pixels), ultimately achieving high-precision and robust registration.

[0134] In step 5, the registration parameters are gradually optimized by iteratively implementing steps 3 and 4 multiple times until the registration accuracy stabilizes within the set threshold range, that is, the change in image positioning error is less than 0.001 pixels, ultimately ensuring high accuracy and robustness of the registration.

[0135] Example 1

[0136] A certain area was selected as an example, and an aerial orthophoto with a resolution of 0.05 meters and 27.065 LiDAR point cloud data per square meter were collected. The final effect of the registration method of the present invention is as shown in the attached figure. Figure 6 shown.

[0137] Especially for complex urban scenes, the present invention can focus on the fine matching process of line and surface features, and achieve high-precision, high-automation and strong robust registration of urban images and LiDAR data without the need for field control points.

[0138] The specific embodiments described herein are merely illustrative of the spirit of the present invention. Persons skilled in the art may make various modifications, additions, or substitutions to the described specific embodiments without departing from the spirit of the present invention or exceeding the scope of the appended claims.

Claims

1. A method for registering urban aerial images and LiDAR point clouds using multivariate features, characterized in that: The following steps are involved: Step 1: Preprocess the image to generate sparse feature point clouds and image 3D feature lines, and use an improved deep learning model to extract the building roof outline information of the image; Step 2: Using the statistical feature analysis method of random sampling consensus algorithm and principal component analysis, line feature and surface feature information are extracted from the LiDAR point cloud; Step 3: Using the improved ICL algorithm, we use the correspondence between the image 3D feature lines and the LiDAR 3D feature lines, and add the projection constraints of the image roof outline and the LiDAR building plane to perform global coarse registration; Step 4: Based on step 3, the multivariate features of points, lines, and surfaces extracted from the urban scene are used as joint control information to construct an extended energy equation for image block adjustment optimization to correct local deformation caused by camera distortion and uneven image connection strength. Step 5: Iterate steps 3 and 4 multiple times until the change in the error in image positioning meets the preset requirements.

2. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 1, wherein: In step 1, the image is preprocessed to generate a sparse feature point cloud and image 3D feature lines, including: using the positioning and attitude information obtained by the POS system to perform preliminary calculations of the internal and external orientation elements of the image, and combining a feature matching algorithm to generate a sparse feature point cloud and image 3D feature lines, wherein the sparse feature point cloud is generated by matching the same-name feature points between images, and the 3D feature lines are obtained by matching the same-name line segments between images and performing 3D reconstruction; The improved deep learning model is used to extract building roof outline information from images, including: using a deep convolutional neural network model improved based on U-Net, in which the encoder uses ResNet-50 to replace the traditional convolution module to enhance feature expression capabilities, and the decoder introduces an attention mechanism to focus on building edge details and reduce background interference; The attention mechanism calculation process is expressed as: ； in, represents the attention weight, is the high-level feature passed by the encoder, is the low-level feature passed to the decoder, 、 is the learnable convolution kernel weight, is the bias, is the Sigmoid activation function; The binary cross entropy loss function is used during training: ； in, is the true label, is the predicted label output by the network, is the total number of samples. By minimizing the loss function, the network model can effectively learn and stably output high-precision building outlines.

3. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 1, wherein: In step 2, the statistical feature analysis method using the random sampling consensus algorithm and principal component analysis is used to extract line feature and surface feature information from the LiDAR point cloud, including: The RANSAC algorithm is used to perform a preliminary fit on the linear structure in the point cloud data. The specific implementation includes: randomly selecting a minimum number of point sets from the LiDAR point cloud, fitting the initial line segment equation, and then defining a distance threshold. All points with a distance less than the threshold are considered as supporting inliers of the model. This process is iterated until the number of inliers reaches the maximum, thereby determining the optimal linear model. The specific three-dimensional straight line equation is expressed as: ； Where, is the coordinate of a point on the line, and the direction vector is ; After the RANSAC algorithm optimizes and obtains the optimal parameters of the line model, the PCA method is used to further analyze the spatial distribution characteristics of the inner point set and determine the most stable three-dimensional feature line: ； Where, represents the coordinates of the inner points of the line feature, represents the centroid of the interior points, is the number of interior points; the matrix The principal eigenvector corresponds to the maximum eigenvalue, and its direction is the most stable direction of the line segment; Three non-collinear points are randomly selected to construct the initial plane equation, and then the plane equation is optimized iteratively to make more points meet the condition that the distance from the plane is less than a given threshold. The expression of the plane equation is: ； in, is the normal vector of the plane, is a constant term; Each iteration calculates the perpendicular distance from the point to the plane: ； in, is the coordinate of any point in the LiDAR point cloud. The goal of iterative optimization is to minimize the negative value of the number of points that meet the threshold distance to obtain the optimal plane model with the most inliers.

4. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 1, wherein: In step 3, the improved ICL algorithm is used to perform a global coarse registration using the correspondence between the image 3D feature lines and the LiDAR 3D feature lines, and adding the projection constraints of the image roof outline and the LiDAR building plane, including: The 3D feature lines of the image extracted in step 1 are matched with the 3D feature lines of the LiDAR extracted in step 2 to establish an initial correspondence between the two sets of data. The rigid transformation parameters are then iteratively solved using the ICL algorithm. The core objective function of the improved ICL algorithm is to minimize the matching error between the 3D feature lines. The specific formula is described as follows: Define the image feature line set as , the LiDAR feature line set is { }, there is a preliminary correspondence between the two sets of characteristic lines, and the optimization objective function of the ICL algorithm is expressed as: ； Where, and Represent the three-dimensional feature lines corresponding to the image and LiDAR, represents the rotation matrix, represents the translation vector, represents the number of corresponding line feature pairs; By iteratively solving the above equations, we can obtain the optimal and , complete the rigid registration of line features between the image and the LiDAR point cloud; Image contour features , the building plane features corresponding to LiDAR The currently obtained transformation matrix is ​​projected onto the spatial plane where the image contour is located to form a projection contour line , thus forming the surface contour projection error constraint function: ； in, Indicates the number of matching pairs between image contours and LiDAR surface features.

5. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 4, wherein: Combining the two constraints of line and surface, the line and surface feature constraints are integrated into a unified optimization framework, and the final overall objective function is expressed as: ； In the formula, the weight factor Used to adjust the relative strength of surface feature constraints and line feature constraints.

6. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 1, wherein: In step 4, the multivariate features of points, lines, and surfaces extracted from the urban scene are used as joint control information to construct an extended energy equation for image block adjustment optimization, wherein the specific implementation of the point feature constraint includes: The principle of the shortest distance between points and surfaces is used to select the "control points" with the same name as each sparse feature point in the LiDAR point cloud. For any sparse feature point , search for the LiDAR point cloud within the spherical neighborhood with a radius of r, and record the nearest neighbor point cloud set as , For the The number of LiDAR point clouds within the neighborhood of sparse feature points; according to the point cloud collection All points in the space are fitted using the least squares method to obtain the spatial plane , with sparse feature points To plane The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud; the geometric control constraint of the point feature It can be expressed as follows ； Where, is the prior weight, is the function for calculating the distance between spatial points and surfaces, is the average ground resolution of the image, used to convert Normalize to image unit.

7. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 6, wherein: The specific implementation of the line feature constraint in step 4 includes: The principle of the shortest distance between points and lines is used to select the "control points" with the same name as each sparse feature point in the LiDAR point cloud; let the three-dimensional feature lines of the image and LiDAR point cloud extracted in step 1 and step 2 be and , and projected onto the horizontal plane as and , The point on (t=1,2,…,n), where n is The number of points on the image feature line; the plane projection points on the image feature line To the straight line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud; Geometric control constraints for line features It is expressed as follows: ； Where, is the plane projection point on the image feature line To LiDAR projection line The shortest Euclidean distance, is the average ground resolution of the image, It is The weight corresponding to the line feature.

8. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 7, wherein: The specific implementation of the face feature constraint in step 4 includes: Let the building roofs in the image and LiDAR point cloud extracted in step 1 and step 2 be and ,Will Projected onto the horizontal plane, it is recorded as For plane The edge of the contour is marked with a positioning point every certain number of pixels. (s=1,2,…,u), where u is the number of positioning points on the edge of the contour; To plane The nearest edge contour line The foot point of the perpendicular is used as the "control point" of the same name in the LiDAR point cloud; Geometric control constraints for face features Expressed as: ； Where, Indicates that in the horizontal plane, the image contour points To LiDAR projected contour edge segment The shortest Euclidean distance, is the average ground resolution of the image, It is The weight of each face feature.

9. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 8, wherein: In step 4, the energy equation of the classic block adjustment is expanded to construct the energy equation shown in Equation (13): ； Where, is the image point residual term defined in the classic block adjustment, For the The orientation parameters of the image, For its Observation of image points on an image, Projection functions defined for collinearity condition equations; It is a geometric constraint formed by combining point, line and surface feature information.

10. The method for registering urban aerial images and LiDAR point clouds using multivariate features according to claim 9, wherein: Assign appropriate prior weights to the geometric constraints formed by the joint point, line, and surface feature information to maximize the registration accuracy. The specific implementation includes: For point and surface feature constraints, use the fitted plane and The local coplanarity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, and , taking the spatial plane fitting error as an indicator to measure coplanarity; For line feature constraints, use the fitted line The local collinearity of the LiDAR point set at the time of is used to adaptively determine the prior weight of the LiDAR control information, that is, , taking the straight line fitting error as an indicator of collinearity; When the fitting error of the spatial plane and line is small, the corresponding LiDAR control information is given a larger prior weight; otherwise, a smaller prior weight is given.

Images