Method for quickly measuring structure based on oblique photography

Abstract

The application relates to the technical field of image data processing, and discloses a structure rapid actual measurement method based on oblique photography. Two-dimensional structure line features in oblique photography multi-view images are extracted; the spatial intersection relationship of the two-dimensional structure lines in different views is calculated based on a collinear condition equation and polar geometric constraints, and a three-dimensional space structure line skeleton is generated through back projection; a graph optimization energy function is constructed by taking the three-dimensional space structure line skeleton as a hard constraint, a distance penalty term is introduced into a cost aggregation process of multi-view stereo matching, a disparity iterative propagation is guided to be carried out in a tangent direction at a structure edge by a diffusion algorithm, and a dense three-dimensional point cloud which retains edge geometric features is generated; surface meshes are extracted based on the dense three-dimensional point cloud, and a surface fitting is carried out to output actual measurement results. The application overcomes the point cloud hole and edge deformation defects caused by dense matching failure in a weak texture area, guarantees the integrity and authenticity of edge geometric features, and improves actual measurement accuracy and reliability.

Description

Technical Field

[0001] This invention relates to the field of image data processing technology, and discloses a rapid method for measuring structures based on oblique photography. Background Technology

[0002] Oblique photogrammetry acquires a sequence of two-dimensional images of a target structure from multiple perspectives using multi-angle cameras mounted on a flight platform. In conventional structural measurement procedures, the data processing system extracts feature points from these two-dimensional images, searches for corresponding points between adjacent images based on the principle of grayscale consistency, and reconstructs a sparse three-dimensional point cloud by solving for the camera's intrinsic and extrinsic parameters and pose relationships. Subsequently, a conventional approach uses a multi-view stereo matching algorithm to estimate and fuse the sparse point cloud's depth map. Through planar scanning or cost calculation and aggregation in parallax space, a dense three-dimensional point cloud model is generated. Finally, the system extracts the structure's surface mesh based on this dense point cloud and outputs measured data such as distances and angles through surface fitting and geometric calculations.

[0003] In the process of generating dense point clouds through multi-view stereo matching, the cost aggregation step is crucial in determining the matching quality. Conventional cost aggregation algorithms are based on the assumption of local smoothness, assuming that adjacent pixels have similar disparity values ​​and that cost propagation and diffusion occur isotropically within the image plane. However, building structures typically have large areas of weak texture, such as painted walls or glass curtain walls, where grayscale gradient changes are subtle. Within these weak texture areas, the matching cost curve often appears flat or contains multiple similar minima, resulting in a lack of clear convergence direction for disparity propagation during cost aggregation. This easily leads to disparity overflow and matching chaos, causing the generated point cloud to exhibit flying spots, holes, or blurred edges at the edges of structures.

[0004] Existing multi-view stereo matching methods based on grayscale consistency suffer from disparity propagation going out of control in weakly textured regions of structures due to the lack of geometric constraints beyond grayscale information. When disparity propagation crosses the edges of structures, conventional isotropic diffusion mechanisms fail to distinguish depth jumps on either side of the edge, causing erroneous disparities from weakly textured regions to spread towards the edges, rendering dense matching completely ineffective in such areas. This results in sparse, porous, and geometrically deformed 3D point clouds generated at the structure edges, destroying the true topology and geometric features of the structure's edge contours. Consequently, the measured dimensions of the structure extracted from this dense point cloud deviate from the true physical dimensions, creating the core technical problem of measurement distortion. Summary of the Invention

[0005] The purpose of this invention is to provide a rapid method for measuring structures based on oblique photography, which can effectively solve the problems mentioned in the background art.

[0006] To achieve the above objectives, the technical solution adopted by the present invention is as follows: Rapid structural measurement methods based on oblique photography include: Extracting two-dimensional structural line features from multi-view images acquired by oblique photogrammetry; Based on the collinearity condition equation and epipolar geometric constraints, the spatial intersection relationship of the two-dimensional structural line in different views is calculated, and a three-dimensional spatial structural line skeleton is generated by back projection. Using the three-dimensional spatial structure line skeleton as a hard constraint, a graph optimization energy function is constructed. The distance penalty term of the structure line endpoints and edge points in the graph optimization energy function is introduced into the cost aggregation process of multi-view stereo matching, guiding the diffusion algorithm to perform disparity iterative propagation along the tangent direction at the structure edge, generating a dense three-dimensional point cloud that retains the geometric features of the structure edge. Based on the dense 3D point cloud, the surface mesh of the structure is extracted and surface fitting is performed to output the measured results of the structure size.

[0007] Preferably, the extraction of two-dimensional structural line features from multi-view images acquired by oblique photography includes: A multi-scale feature extraction network is constructed to extract feature maps of different receptive fields from the multi-view images; Edge response calculation is performed on the feature map to obtain an initial edge probability map; By combining the image gradient direction information, non-maximum suppression and hysteresis thresholding are performed on the initial edge probability map to filter out edge pixels whose gradient magnitude meets the preset conditions. An edge connection strategy based on local geometric consistency is adopted to perform topological connections on the edge pixels to form continuous two-dimensional structural line features.

[0008] Preferably, the step of calculating the spatial intersection relationship of the two-dimensional structural lines in different views based on the collinearity condition equation and epipolar geometric constraints, and generating a three-dimensional spatial structural line skeleton by back projection, includes: The two-dimensional structural line is parameterized as the intersection line of a three-dimensional spatial plane, and a three-plane intersection geometric model of the two-dimensional structural line and the spatial plane is established. Based on the collinearity condition equation, the pixel coordinates of the two-dimensional structural line are mapped to the spatial plane, and the epipolar constraint is used to eliminate mismatched line segments of the two-dimensional structural line in different views. Triangulation is performed on the retained matching line segments to solve the spatial parameters of the three-plane intersection geometric model, and back projection is used to generate a continuous three-dimensional spatial structure line skeleton.

[0009] Preferably, a graph optimization energy function is constructed using the three-dimensional spatial structure line skeleton as a hard constraint. The distance penalty terms for the endpoints and edge points of the structure line in the graph optimization energy function are introduced into the cost aggregation process of multi-view stereo matching, including: The endpoints and edge points of the three-dimensional spatial structure skeleton are used as vertices of the graph model. Data items are constructed using the geodesic distance between the vertices, and smoothing items are constructed using the angle between the normal vectors of adjacent vertices. The graph optimization energy function is generated by combining these elements. Extract the distance penalty term from the graph optimization energy function, and superimpose the distance penalty term as a regularization term onto the initial cost volume of multi-view stereo matching. Then, perform filtering and aggregation on the initial cost volume.

[0010] Preferably, the surface mesh of the structure is extracted based on the dense three-dimensional point cloud and surface fitting is performed to output the measured results of the structure's dimensions, including: Normal estimation and clustering segmentation are performed on the dense three-dimensional point cloud to extract different surface regions of the structure; Using the three-dimensional spatial structure skeleton as boundary constraints, the surface region is anisotropically meshed to generate a triangular mesh that fits the edge of the structure. The triangular mesh is fitted with a surface based on geometric feature preservation. The spatial geometric parameters of the fitted surface are calculated, and the measured distance and angle results of the structure are output based on the spatial geometric parameters.

[0011] Preferably, the construction of a multi-scale feature extraction network to extract feature maps of different receptive fields from the multi-view images includes: An encoder with parallel dilated convolution branches is used to perform feature encoding on the multi-view image to obtain multi-scale feature maps; The multi-scale feature map is input into the feature fusion layer, and the weight distribution of the multi-scale feature map in the spatial dimension is calculated through the spatial attention mechanism, and the weight distribution of the multi-scale feature map in the channel dimension is calculated through the channel attention mechanism. The multi-scale feature map is weighted and fused based on the weight distribution of the spatial dimension and the weight distribution of the channel dimension to output a multi-scale fused feature map.

[0012] Preferably, the step of using epipolar constraints to eliminate mismatched line segments of the two-dimensional structure line in different views includes: For the target two-dimensional structure line in the reference view, calculate the epipolar equation of the target two-dimensional structure line in the neighboring view; Determine the projection region of the envelope of the target two-dimensional structural line in the neighborhood view; Calculate the overlap ratio between the epipolar segment corresponding to the epipolar equation and the projected region to obtain the epipolar intersection-union ratio value; Two-dimensional structural lines in the neighborhood view whose epipolar intersection-union ratio is lower than a preset threshold are identified as mismatched line segments and removed, while matching line segment pairs that satisfy the epipolar geometric constraints are retained.

[0013] Preferably, the distance penalty term is superimposed as a regularization term onto the initial cost volume of multi-view stereo matching, and the initial cost volume is filtered and aggregated, including: The distance penalty term is mapped to the disparity space to construct an anisotropic diffusion tensor in the disparity domain. The neighboring pixel weights in the cost aggregation process are determined based on the anisotropic diffusion tensor, such that a larger diffusion weight is set in the tangential direction of the three-dimensional spatial structure skeleton, and a smaller diffusion weight is set in the normal direction of the three-dimensional spatial structure skeleton. The initial cost body is aggregated by semi-global matching based on the neighboring pixel weights, and the aggregated cost body is output.

[0014] Preferably, before extracting the two-dimensional structural line features from the multi-view images acquired by oblique photography, the method further includes: Obtain the initial intrinsic and extrinsic parameters and initial pose of the oblique photography camera; Project the two-dimensional structural line feature onto three-dimensional space and calculate the reprojection error of the two-dimensional structural line feature between adjacent views; A bundle adjustment optimization model based on line features is constructed with minimizing the reprojection error as the objective function. The initial intrinsic and extrinsic parameters and initial pose of the oblique photography camera are jointly iteratively solved and updated to obtain the optimized target camera parameters and target pose. Based on the target camera parameters and target pose, subsequent spatial intersection calculations are performed.

[0015] Preferably, after outputting the measured results of the structure dimensions, the method further includes: Acquire airborne lidar point cloud data of the target structure, and extract planar and line features from the airborne lidar point cloud data; The planar and line features in the airborne lidar point cloud data and the corresponding features in the dense 3D point cloud are coarsely registered and finely registered based on iterative nearest points to obtain the coordinate transformation matrix. Based on the coordinate transformation matrix, the dense 3D point cloud and the airborne lidar point cloud data are fused and compared to identify and repair the occlusion and void areas in the dense 3D point cloud, and the verified measured results are output.

[0016] Compared with the prior art, the beneficial effects of the present invention are as follows: 1. This invention extracts two-dimensional structural line features from oblique photogrammetric multi-view images and calculates their spatial intersection relationships in different views based on collinearity condition equations and epipolar geometric constraints, generating a three-dimensional spatial structural line skeleton through back-projection. Based on this, an optimization energy function is constructed using the three-dimensional spatial structural line skeleton as a hard constraint. Distance penalty terms for structural line endpoints and edge points are introduced into the cost aggregation process of multi-view stereo matching, guiding the diffusion algorithm to iteratively propagate disparity along the tangential direction at the structural edge. This mechanism provides geometric boundary conditions independent of grayscale information for disparity propagation in weakly textured regions. The distance penalty term forcibly blocks erroneous disparity propagation paths across structural edges, changing the blindness of conventional isotropic diffusion. Within weakly textured regions, disparity can be correctly aggregated along the structural line extension direction, while disparity propagation is strictly suppressed in the direction perpendicular to the edge normal. This mechanism overcomes the point cloud holes and edge deformation defects caused by failed dense matching in weakly textured and textureless regions, ensuring the integrity and authenticity of the geometric features of the structure's edges and improving the accuracy of measured structure dimensions.

[0017] 2. In extracting 2D structural line features, this invention employs a multi-scale feature extraction network combined with spatial and channel attention mechanisms to obtain a fused feature map containing rich contextual information. By combining image gradient direction information and a geometrically consistent edge connection strategy, the completeness of continuous structural line extraction is improved. During the backprojection generation of the 3D spatial structural line skeleton, the 2D structural lines are parameterized into a three-plane intersection geometric model. Epipolar constraints are used to eliminate mismatched line segments, improving the spatial localization accuracy of the 3D spatial structural line skeleton. In the cost aggregation stage, the distance penalty term is mapped to the disparity space to construct an anisotropic diffusion tensor. Diffusion weights are dynamically allocated according to the tangent and normal directions of the 3D spatial structural line skeleton; a large weight is assigned to the tangent direction to ensure surface smoothness, while a small weight is assigned to the normal direction to maintain sharp edges. These mechanisms further suppress mismatches across structural edges, reduce the redundancy of feature extraction and matching, and improve the convergence speed and overall experimental efficiency of the 3D reconstruction process. Attached Figure Description

[0018] Figure 1 This is a flowchart of the oblique photography multi-view image acquisition and preprocessing process of the present invention; Figure 2 This is a flowchart of the two-dimensional structural line feature extraction process of the present invention; Figure 3 This is a flowchart illustrating the generation process of the three-dimensional spatial structure line skeleton of the present invention. Figure 4 This is a flowchart of the graph optimization and multi-view stereo matching cost aggregation process of the present invention; Figure 5 This is a flowchart of the process for generating dense 3D point clouds and measuring the dimensions of structures according to the present invention. Figure 6 This is a flowchart of the multi-scale feature extraction and precise matching of three-dimensional structure lines of the present invention. Detailed Implementation

[0019] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0020] Please refer to Figure 1 This embodiment provides a rapid field measurement method for structures based on oblique photogrammetry. First, a multi-view image sequence of the target structure is acquired using an oblique photogrammetry system. This oblique photogrammetry system is mounted on a fixed-wing UAV platform and includes five industrial cameras with different perspectives: one vertically downward orthophoto camera and four oblique cameras facing forward, backward, left, and right respectively. The initial intrinsic parameters of each camera are obtained through laboratory calibration, including focal length, principal point coordinates, radial distortion coefficient, and tangential distortion coefficient. The initial extrinsic parameters of each camera are obtained through a combined navigation system of a global navigation satellite system and an inertial measurement unit, including the rotation matrix and translation vector at each exposure time.

[0021] The acquired multi-view images are preprocessed, including radiometric correction and distortion correction. Radiometric correction is used to eliminate illumination differences between different cameras and sensor response inhomogeneities, while distortion correction is used to eliminate radial and tangential distortions of the camera lens, resulting in corrected multi-view images.

[0022] refer to Figure 2 Two-dimensional structural line features are extracted from the corrected multi-view image. First, the image gradient of each pixel is calculated, with the gradient magnitude and direction obtained using the Sobel operator. For a pixel at a specified coordinate, both its horizontal and vertical gradients are calculated through convolution operations. To calculate the horizontal gradient, a three-row, three-column horizontal convolution kernel is used to convolve with the gray values ​​of the pixel and its eight neighboring pixels. The values ​​in the first row of this kernel, from left to right, are -1, 0, and +1; the values ​​in the second row are -2, 0, and +2; and the values ​​in the third row are -1, 0, and +1. The result is the horizontal gradient value of that pixel. Similarly, to calculate the vertical gradient, a three-row, three-column vertical convolution kernel is used to convolve with the corresponding pixel's gray value. The values ​​in the first row of this kernel, from left to right, are -1, -2, and -1; the values ​​in the second row are all zero; and the values ​​in the third row are +1, +2, and +1. The result is the vertical gradient value of that pixel. Convolution operation refers to the calculation process of multiplying the convolution kernel with the corresponding image region element by element and then summing the results.

[0023] For a pixel at a specified coordinate, its gradient magnitude is calculated as follows: First, obtain the horizontal and vertical gradient values ​​of the pixel. Then, square both gradient values, sum them, and finally take the square root of the sum. The resulting value is the gradient magnitude of the pixel. Its gradient direction is calculated as follows: Divide the vertical gradient value of the pixel by the horizontal gradient value to obtain the quotient. Then, calculate the arctangent of the quotient. The result is the gradient direction of the pixel.

[0024] Non-maximum suppression (NMS) is applied to the calculated gradient magnitude map. For each pixel, its gradient magnitude is compared with the gradient magnitudes of its two adjacent pixels along the gradient direction. If the gradient magnitude of the current pixel is not the maximum of these three pixels, its gradient magnitude is set to 0; otherwise, its gradient magnitude is retained. NMS can refine edges to a single pixel width.

[0025] Hysteresis thresholding is applied to the gradient magnitude map after non-maximum suppression. Two thresholds are set: a high threshold and a low threshold. and low threshold ,in Gradient magnitude greater than Pixels with gradient magnitudes between [values ​​missing] are marked as strong edge pixels. and Pixels between these points are marked as weak edge pixels, and those with gradient magnitudes less than [a certain value] are [marked as weak edge pixels]. The pixels with weak edges are set to 0. Then, all weak edge pixels are traversed; if a strong edge pixel exists in its 8-neighborhood, it is retained as an edge pixel; otherwise, it is set to 0. Hysteresis thresholding can remove false edges caused by noise while preserving true weak edges.

[0026] The selected edge pixels are topologically connected to form continuous 2D structure line features. An edge connection strategy based on local geometric consistency is employed. Starting with each unvisited edge pixel, neighboring edge pixels are searched along the direction perpendicular to its gradient direction. If the difference between the gradient direction of a neighboring pixel and the gradient direction of the current pixel is less than a preset angle threshold, the neighboring pixel is connected to the current edge segment. This process is repeated until no neighboring pixels meeting the criteria are found, forming a complete 2D structure line. This process is performed on all edge pixels to obtain all 2D structure line features.

[0027] refer to Figure 3Based on the collinearity condition equation and epipolar geometric constraints, the spatial intersection relationship of two-dimensional structure lines in different views is calculated, and a three-dimensional spatial structure line skeleton is generated by back projection. The collinearity condition equation describes the geometric relationship that the spatial point, the camera optical center, and the image point are all on the same straight line. For any spatial point, its imaging coordinates on the image plane can be calculated using camera parameters and spatial coordinates: the three-dimensional coordinates of the spatial point are transformed to the camera coordinate system using a rotation matrix and a translation vector, and the horizontal and vertical coordinate values ​​are divided by the depth coordinate value and then multiplied by the camera focal length to obtain the horizontal and vertical pixel coordinates of the image point.

[0028] Epipolar geometric constraints describe the geometric relationship between two views. For an image point in the reference view... Its corresponding image point in the neighborhood view It must be located at the polar line Above. Under epipolar geometric constraints, the corresponding point of an image point in the reference view in the neighboring view must fall on the corresponding epipolar line. This epipolar line can be calculated using the fundamental matrix and the coordinates of the reference image point. The fundamental matrix is ​​used to characterize the epipolar geometric relationship between two views. Its calculation process is as follows: first, obtain the inverse matrix of the intrinsic parameter matrix of the reference view and the transpose inverse matrix of the intrinsic parameter matrix of the neighboring view; then, combine the relative rotation matrix and the antisymmetric matrix of the translation vector between the two views, and perform matrix multiplication operations in sequence to finally obtain the fundamental matrix.

[0029] The 2D structure line is parameterized as the intersection of three-dimensional spatial planes. For each 2D structure line in each view, it corresponds to a 3D plane that passes through the camera's optical center and contains the 2D structure line. The equation of this plane can be obtained by solving for the coordinates of the camera's optical center and the coordinates of two points on the 2D structure line. For three matching 2D structure lines in three different views, the three corresponding spatial planes intersect at a single spatial line, which is the 3D spatial structure line.

[0030] The pixel coordinates of the 2D structure line are mapped to the corresponding spatial plane based on the collinearity condition equation. For each pixel on the 2D structure line, a ray originating from the camera's optical center can be obtained according to the collinearity condition equation. The intersection of this ray with the spatial plane is the spatial point corresponding to that pixel. Epipolar constraints are used to eliminate mismatched segments of the 2D structure line in different views. For a 2D structure line in the reference view, its epipolar line in the neighboring views is calculated. If the distance between a 2D structure line in the neighboring views and the epipolar line is greater than a preset distance threshold, it is identified as a mismatched segment and eliminated.

[0031] Triangulation is performed on the retained matching line segments to solve for the spatial parameters of the three-plane intersection geometric model. For three matching two-dimensional structural lines, the equations of their corresponding three spatial planes can be solved simultaneously to obtain a system of linear equations. By solving this system of linear equations, the parameters of the three-dimensional spatial structural lines can be obtained, including the direction vector and the coordinates of a point on the line. This process is repeated for all matching two-dimensional structural lines to obtain all three-dimensional spatial structural lines. These structural lines are then connected to form a continuous three-dimensional spatial structural line skeleton.

[0032] refer to Figure 4 A graph optimization energy function is constructed using a 3D spatial structural skeleton as a hard constraint. The vertices of the graph model consist of the endpoints and edge points of the 3D spatial structural skeleton, and the edges consist of the connections between adjacent vertices. The total energy of the graph optimization is composed of two parts: the energy of the data items and the energy of the smoothing items. During calculation, the energy values ​​of the data items and the smoothing items, multiplied by their corresponding weighting coefficients, are added together to obtain the total energy value. The weighting coefficients are used to adjust the proportion of influence of the data items and the smoothing items in the total energy, achieving a balance between fitting accuracy and smoothness.

[0033] The data item measures the tightness of the fit between the vertices of the graph model and the skeleton of the 3D spatial structure. Its value is the sum of the distance values ​​corresponding to all vertices. For each vertex in the graph, the straight-line distance from that vertex to the nearest point on the skeleton of the 3D spatial structure is calculated, and the distance values ​​of all vertices are added together to obtain the total energy value of the data item.

[0034] The smoothing term measures the consistency of the normal vectors of adjacent vertices in a graph model. Its value is the sum of the squared differences in the normal vectors of all adjacent vertex pairs. For each pair of adjacent vertices that are connected, the difference between the normal vectors of the two vertices is calculated, and the square of the magnitude is taken. This value is then added to all adjacent vertex pairs to obtain the total energy value of the smoothing term.

[0035] By minimizing the graph optimization energy function, the vertex positions of the 3D spatial structure skeleton are optimized, resulting in a more accurate and smoother 3D spatial structure skeleton.

[0036] The distance penalty term for the endpoints and edge points of the structure line in the graph optimization energy function is introduced into the cost aggregation process of multi-view stereo matching. First, the initial cost volume for multi-view stereo matching is calculated. The initial cost value for stereo matching is calculated using a normalized cross-correlation method. For each pixel in the reference view and each candidate disparity, the calculation follows these rules: A fixed-size window centered on the current pixel is taken. The average gray value of the reference view pixels within the window and the average gray value of the pixels within the window after offsetting the corresponding disparity in the neighboring views are calculated. The average gray value of each corresponding pixel in the two windows is subtracted from the average gray value, and the results are multiplied digit by digit to obtain the numerator. The sum of the squares of the gray value differences between the two windows is calculated and the square root is taken. The two square root results are then multiplied to obtain the denominator. The numerator is divided by the denominator to obtain the normalized cross-correlation coefficient. Finally, the coefficient is subtracted from the value 1 to obtain the initial matching cost value of the current pixel under that disparity.

[0037] Extract the distance penalty term in the energy function to optimize the image. For each pixel in the reference view... and every possible parallax According to parallax Calculate the three-dimensional spatial point corresponding to this pixel. The coordinates of the point are then calculated in three-dimensional space. Shortest distance to the nearest three-dimensional spatial structure skeleton The value of the distance penalty term is divided into two cases based on the shortest distance from the spatial point to the structural line skeleton: when the shortest distance is less than or equal to the preset distance threshold, the distance penalty term is zero and no penalty is applied; when the shortest distance is greater than the preset distance threshold, the distance penalty term is the penalty coefficient multiplied by the distance value exceeding the threshold. The greater the distance exceeding the threshold, the higher the corresponding penalty value.

[0038] The distance penalty term is added as a regularization term to the original initial cost volume of multi-view stereo matching to obtain the regularized cost volume. Specifically, the calculation method is as follows: for the initial cost value corresponding to each pixel and each candidate disparity, the value of the corresponding distance penalty term multiplied by the regularization weight coefficient is added, and the result is the regularized cost value of that pixel under that disparity.

[0039] The regularized cost volume is filtered and aggregated. A semi-global matching algorithm is used for cost aggregation, and cost propagation is performed in multiple directions. During cost propagation, a guided diffusion algorithm performs disparity iterative propagation along the tangential direction at the structural edges. For each pixel... Based on the tangent and normal directions of the nearby 3D spatial structural skeleton, the propagation weights in different directions are adjusted. A larger propagation weight is set in the tangent direction so that parallax can be correctly aggregated along the extension direction of the structural line; a smaller propagation weight is set in the normal direction to block the propagation of erroneous parallax across the edge of the structure.

[0040] The disparity map is calculated based on the aggregated cost volume. For each pixel... The disparity that minimizes the aggregation cost is selected as the optimal disparity for that pixel. Post-processing of the disparity map includes left-right consistency checking, sub-pixel interpolation, and median filtering. Left-right consistency checking removes mismatched pixels from the disparity map, sub-pixel interpolation improves the accuracy of the disparity map, and median filtering removes noise from the disparity map.

[0041] refer to Figure 5 The disparity maps of multiple views are fused to generate a dense 3D point cloud. For each view's disparity map, the disparity of each pixel is converted into the coordinates of a 3D point based on the camera's intrinsic and extrinsic parameters. The 3D points of all views are merged, and duplicate points are removed to obtain a dense 3D point cloud.

[0042] The surface mesh of the structure is extracted based on a dense 3D point cloud, and surface fitting is performed to output the measured dimensions of the structure. First, the normal vector of the dense 3D point cloud is estimated. Principal component analysis is used to calculate the normal vector of each point. For each point whose normal vector is to be calculated, several nearest neighboring points are selected, and a covariance matrix is ​​constructed according to the following rules: First, the average coordinates of all neighboring points are calculated to obtain the average neighborhood coordinates; then, the average coordinates are subtracted from the coordinates of each neighboring point to obtain the corresponding difference vector; each difference vector is multiplied by its transpose to obtain the corresponding matrix; the matrices corresponding to all neighboring points are summed, and finally, the sum is divided by the total number of neighboring points to obtain the final covariance matrix.

[0043] For covariance matrix Eigenvalue decomposition yields three eigenvalues. and the corresponding three feature vectors Minimum eigenvalue corresponding feature vector That is, the normal vector of that point.

[0044] Clustering and segmentation are performed on dense 3D point clouds to extract different surface regions of structures. A region growing algorithm is used for clustering and segmentation. First, a seed point is selected, with a high confidence level for its normal vector. Then, starting from the seed point, neighboring points that meet the following conditions are added to the same region: the spatial distance between the neighboring point and the seed point is less than a preset distance threshold, and the angle between the normal vectors of the neighboring points and the seed point is less than a preset angle threshold. This process is repeated until no neighboring points that meet the conditions can be found, forming a surface region. This process is repeated for all unassigned points to obtain all surface regions.

[0045] Anisotropic meshing is performed on each surface region using a 3D spatial structural skeleton as boundary constraints. The mesh is refined near the structural lines to ensure geometric accuracy at the edges, while it is sparsed in flat regions to reduce computation. The Delaunay triangulation algorithm is used to generate triangular meshes. The constraints of the 3D spatial structural skeleton are introduced during the triangulation process to ensure that the edges of the triangular mesh coincide with the structural lines as much as possible, generating triangular meshes that fit the edges of the structure.

[0046] A surface fitting based on geometric feature preservation is performed on the triangular mesh. B-spline surfaces are used for fitting. For each surface region, the vertices of the triangular mesh are used as control points to construct a B-spline surface. During the fitting process, constraints from a 3D spatial structural line skeleton are introduced to ensure the continuity and sharpness of the surface at the structural lines.

[0047] Calculate the spatial geometric parameters of the fitted surface. For planar surfaces, calculate their plane equation, normal vector, area, and perimeter; for curved surfaces, calculate their curvature, Gaussian curvature, and mean curvature. Based on these geometric parameters, calculate the measured distances and angles of the structure. For example, calculate the distance between two parallel planes as the width of the structure, and calculate the angle between two intersecting planes as the angle of the structure. Table 1 shows a comparison of the reprojection errors before and after bundle adjustment optimization based on line features.

[0048] Table 1 Comparison of reprojection errors before and after beam adjustment optimization based on line features Table 1 shows a comparison of the reprojection errors of each camera before and after the beam adjustment optimization based on line features. As can be seen from Table 1, the reprojection errors of each camera are significantly reduced after optimization, with the average reprojection error decreasing from 2.36 pixels to 0.43 pixels. This indicates that the camera parameters and poses have been effectively optimized, providing an accurate geometric basis for subsequent 3D reconstruction.

[0049] This embodiment realizes a complete process for rapid measurement of structures based on oblique photogrammetry. By extracting two-dimensional structural line features and back-projecting them to generate a three-dimensional spatial structural line skeleton, it is introduced as a hard constraint into the cost aggregation process of multi-view stereo matching. This guides the parallax propagation along the tangent direction of the structural line, blocks the propagation of erroneous parallax across the edge of the structure, generates a dense three-dimensional point cloud that preserves the geometric features of the structure's edge, and finally outputs accurate measured results of the structure's dimensions.

[0050] In a preferred embodiment, reference Figure 6A multi-scale feature extraction network was constructed to extract feature maps from multi-view images with different receptive fields. This network employs an encoder-decoder structure, with the encoder containing three parallel dilated convolution branches. Each branch has a different dilation rate (1, 2, and 4), corresponding to different receptive field sizes. Dilated convolution can expand the receptive field of the convolution kernel without increasing the number of parameters or computational cost, thus capturing feature information at different scales.

[0051] The first dilated convolution branch has a dilation rate of 1, corresponding to a standard convolution, with a receptive field size of 3×3, used to capture local detail features. The second dilated convolution branch has a dilation rate of 2 and a receptive field size of 5×5, used to capture mid-scale features. The third dilated convolution branch has a dilation rate of 4 and a receptive field size of 9×9, used to capture large-scale contextual features. Each dilated convolution branch contains two convolutional layers, each followed by a batch normalization layer and a ReLU activation function.

[0052] The multi-scale feature maps output from the three dilated convolution branches are input into the feature fusion layer. The feature fusion layer employs spatial attention and channel attention mechanisms to calculate the weight distribution of the multi-scale feature maps in the spatial and channel dimensions, respectively, and then performs weighted fusion of the multi-scale feature maps based on the weight distribution.

[0053] Spatial attention mechanisms are used to capture the differences in importance of feature maps across spatial dimensions. The computation process is as follows: First, the input feature maps are subjected to average pooling and max pooling respectively, resulting in two spatial feature maps. Then, these two feature maps are concatenated along the channel dimension, passed through a 7×7 convolutional layer and a sigmoid activation function, to obtain the spatial attention map. The generation process of the spatial attention map is as follows: First, the input feature maps are subjected to average pooling and max pooling respectively, resulting in two single-channel spatial feature maps; these two feature maps are concatenated along the channel dimension, passed through a 7x7 convolutional layer for feature fusion, and then numerically mapped using an activation function to finally obtain the spatial attention map. The activation function constrains the output value to the interval between zero and one, representing the feature weights at the corresponding spatial locations.

[0054] Channel attention mechanisms are used to capture the differences in importance of feature maps across different channels. The calculation process is as follows: First, the input feature map is subjected to average pooling and max pooling respectively, resulting in two channel-dimensional feature vectors. Then, these two feature vectors are input into a shared multilayer perceptron, which contains two fully connected layers. The first fully connected layer reduces the dimension of the feature vector to 1 / 16 of its original size, and the second fully connected layer restores the dimension to its original size. The outputs of the two multilayer perceptrons are summed and passed through a sigmoid activation function to obtain the channel attention map. The channel attention map is used to characterize the importance of each channel in the feature map. Its generation process is as follows: First, the input feature map is subjected to global average pooling and global max pooling respectively, resulting in two one-dimensional channel feature vectors. These two feature vectors are then input into the same shared multilayer perceptron for feature transformation, and the two transformed results are summed bit by bit. Finally, an activation function maps the values ​​to the interval between zero and one, ultimately obtaining the weights corresponding to each channel, i.e., the channel attention map.

[0055] The spatial attention map and channel attention map are multiplied by the input feature map respectively to obtain the spatially weighted feature map and the channel-weighted feature map. Then, these two weighted feature maps are added together to obtain the multi-scale fused feature map.

[0056] Edge response calculation is performed on the multi-scale fused feature map to obtain an initial edge probability map. A 1×1 convolutional layer is used to compress the number of channels of the multi-scale fused feature map to 1, and then a sigmoid activation function is applied to obtain the edge probability value of each pixel, forming the initial edge probability map.

[0057] Non-maximum suppression and hysteresis thresholding are applied to the initial edge probability map based on image gradient direction information. For each pixel in the initial edge probability map, its edge probability value is compared with the two adjacent pixels along the gradient direction. If the edge probability value of the current pixel is not the maximum of these three pixels, its edge probability value is set to 0; otherwise, its edge probability value is retained. Then, hysteresis thresholding is applied to the edge probability map after non-maximum suppression to filter out edge pixels whose edge probability values ​​meet preset conditions.

[0058] An edge connection strategy based on local geometric consistency is employed to topologically connect edge pixels, forming continuous two-dimensional structural line features. Starting from each unvisited edge pixel, neighboring edge pixels are searched along the direction perpendicular to its gradient direction. The angle between the gradient direction of the neighboring pixel and the current pixel, as well as the spatial distance, are calculated. If the angle is less than a preset angle threshold and the distance is less than a preset distance threshold, the neighboring pixel is connected to the current edge segment. Simultaneously, the geometric consistency between the neighboring pixel and the current edge segment is checked. If the position of the neighboring pixel is consistent with the extension direction of the current edge segment, it is added; otherwise, the search stops. This process is repeated until no neighboring pixels meeting the conditions can be found, forming a complete two-dimensional structural line.

[0059] The two-dimensional structural line is parameterized as the intersection line of a three-dimensional spatial plane, establishing a three-plane intersection geometric model of the two-dimensional structural line and the spatial plane. For a two-dimensional structural line in the reference view... It corresponds to an optical center in three-dimensional space that is accessed through a reference camera. And includes plane .flat The equation can be obtained by and The two points on and The solution is obtained. Similarly, for two matching two-dimensional structure lines in the neighborhood view... and They each correspond to a plane in three-dimensional space. and Three planes , and Intersect at a straight line in space This straight line is the three-dimensional spatial structure line.

[0060] The three-plane intersection geometric model solves for the common intersection line, i.e., the three-dimensional spatial structure line, by solving the simultaneous equations of the three spatial planes. The three spatial planes correspond to the projection planes formed by the two-dimensional structure lines in the three views and the optical centers of their respective cameras. Each plane satisfies the standard three-dimensional plane equation form, i.e., the linear combination of the plane parameters and the coordinates of the spatial points equals zero. The solution set of the three plane equations is the spatial line obtained by the intersection of the three planes.

[0061] The pixel coordinates of the 2D structure line are mapped to the corresponding spatial plane based on the collinearity condition equation. For each pixel on the 2D structure line, a ray originating from the camera's optical center can be obtained according to the collinearity condition equation. The intersection of this ray with the spatial plane is the spatial point corresponding to that pixel. Epipolar constraints are used to eliminate mismatched line segments of the 2D structure line in different views.

[0062] For the target 2D structure line in the reference view, calculate its epipolar equation in the neighboring view. For the two endpoints of the target 2D structure line... and Calculate their epipolar lines in the neighborhood view respectively. and The polar segment of the target 2D structure line in the neighborhood view is... and The line segments between them.

[0063] Determine the projection region of the target 2D structure line's envelope bounding box in the neighboring view. The envelope bounding box of the target 2D structure line is the smallest rectangle containing the structure line. Based on the homography matrix between the reference view and the neighboring view, project the four vertices of the envelope bounding box into the neighboring view to obtain the projection region. .

[0064] The epipolar intersection-union ratio (IUGR) is calculated by determining the overlap ratio between the epipolar segment and the projected area. The IUGR measures the degree of overlap between the epipolar segment and the projected area. It is calculated as follows: first, find the intersection of the epipolar segment and the projected area, and calculate the area of ​​the intersection; then, find the union of the epipolar segment and the projected area, and calculate the area of ​​the union; finally, divide the intersection area by the union area to obtain the IUGR value.

[0065] Set the polar line intersection ratio threshold. The polar line intersection ratio is lower than Two-dimensional structural lines in the neighborhood view are identified as mismatched line segments and removed, while matching line segment pairs that satisfy the epipolar geometric constraints are retained.

[0066] Triangulation is performed on the retained matching line segments to solve for the spatial parameters of the three-plane intersection geometric model. For the three matching two-dimensional structural lines, the equations of their corresponding three spatial planes can be solved simultaneously to obtain a system of linear equations. This system of linear equations is solved using singular value decomposition to obtain the direction vectors of the three-dimensional spatial structural lines. and a point on a straight line The coordinates of a point on a line are given. The three-dimensional spatial structure line is expressed parametrically. Given the coordinates of a reference point on the line and the direction vector of the line, the coordinates of any point on the line can be calculated as follows: multiply the direction vector by an arbitrary real-valued parameter, then add it to the coordinates of the reference point. The result is the coordinates of the point on the line corresponding to the parameter position.

[0067] Perform this operation on all matching 2D structural lines to obtain all 3D spatial structural lines. These 3D spatial structural lines are then stitched together according to their connection relationships to form a continuous 3D spatial structural line skeleton. The parameter configurations for each layer of the multi-scale feature extraction network are shown in Table 2.

[0068] Table 2 Parameter configuration of each layer of the multi-scale feature extraction network Table 2 shows the parameter configurations of each layer in the multi-scale feature extraction network. As can be seen from Table 2, the network contains three parallel dilated convolutional branches, each with a different dilation rate, capable of capturing feature information at different scales. The feature fusion layer combines spatial attention and channel attention mechanisms, adaptively adjusting the weights of features at different scales and in different channels, thus improving the expressive power of edge features.

[0069] This embodiment constructs a multi-scale feature extraction network and combines spatial and channel attention mechanisms to obtain a fused feature map containing rich contextual information, thereby improving the completeness and accuracy of 2D structural line feature extraction. By establishing a three-plane intersection geometric model and combining epipolar intersection-union ratio to eliminate mismatched line segments, the spatial positioning accuracy of the 3D spatial structural line skeleton is improved.

[0070] In a preferred embodiment, the endpoints and edge points of the three-dimensional spatial structure skeleton are used as vertices of the graph model. Data terms are constructed using the geodesic distances between vertices, and smoothing terms are constructed using the angle between the normal vectors of adjacent vertices. These are combined to generate a graph optimization energy function. The graph model's vertex set... This includes all endpoints and edge points of the 3D spatial structure line skeleton. Edge points are points sampled at equal intervals along the 3D spatial structure line skeleton; the sampling interval is determined based on the length and accuracy requirements of the structure line. (This refers to the edge set of the graph model.) It includes the connection relationships between adjacent vertices, which are vertices that are directly connected on the skeleton of the three-dimensional spatial structure.

[0071] Data Items This measure of the fit between a vertex and the skeleton of a 3D spatial structure is defined as the sum of the Euclidean distances from each vertex to the nearest point on its corresponding 3D spatial structure skeleton. For a vertex... The nearest point on its corresponding three-dimensional spatial structure skeleton The distance from a point to a line can be calculated using the formula for distance from a point to a line. For a three-dimensional line expressed in parametric form, given the coordinates of a reference point on the line and the direction vector of the line, the distance from any point in space to the line is calculated as follows: First, calculate the vector formed by the point and the reference point of the line. Then, perform a three-dimensional cross product operation on this vector and the direction vector of the line, and calculate the magnitude of the cross product vector. Finally, divide this magnitude by the magnitude of the direction vector itself. The final value is the distance from the point to the line. Here, the cross product operation refers to the standard cross product calculation of three-dimensional vectors.

[0072] Smoothing terms The normal vector consistency between adjacent vertices is measured and defined as the sum of the squares of the differences between the normal vectors of adjacent vertices. The normal vector of a vertex is calculated using the tangent directions of the 3D spatial structure skeleton near that vertex. For vertices on the 3D spatial structure skeleton... Its tangent direction is the direction vector of the structure line where the vertex is located. The normal vector is any vector perpendicular to the tangent direction. In this embodiment, the vector perpendicular to the tangent direction and with the smallest angle to the Z-axis of the world coordinate system is chosen as the normal vector.

[0073] Combining data terms and smoothing terms to generate a graph for optimizing the energy function: in, The weighting coefficients for the smoothing term are used to balance the contributions of the data terms and the smoothing term.

[0074] The vertex positions of the 3D spatial structure skeleton are optimized by minimizing the graph optimization energy function. The gradient descent method is used to solve this optimization problem, iteratively updating the vertex positions until the change in the energy function is less than a preset threshold or the maximum number of iterations is reached. The optimized 3D spatial structure skeleton is more accurate and smoother, providing reliable geometric constraints for subsequent cost aggregation.

[0075] The distance penalty term in the energy function of the image optimization is extracted and superimposed as a regularization term onto the initial cost volume of multi-view stereo matching. This initial cost volume is then filtered and aggregated. The distance penalty term measures the distance between a point in 3D space and the skeleton of the 3D spatial structure; the greater the distance, the larger the penalty. For each pixel in the reference view... and every possible parallax According to parallax Calculate the three-dimensional spatial point corresponding to this pixel. The coordinates of the point are then calculated in three-dimensional space. Shortest distance to the nearest three-dimensional spatial structure skeleton Distance penalty item and Proportional.

[0076] A distance penalty term is mapped to the disparity space to construct an anisotropic diffusion tensor within the disparity domain. This anisotropic diffusion tensor describes the diffusion intensity in different directions during cost aggregation. For each pixel, an anisotropic diffusion tensor is constructed by combining the tangent and normal directions of the nearby 3D spatial structure line skeleton. The diffusion tensor consists of two superimposed components: the tangent component is the outer product of the tangent diffusion coefficient and the tangent vector with its own transpose, and the normal component is the outer product of the normal diffusion coefficient and the normal vector with its own transpose. The tangent diffusion coefficient is greater than the normal diffusion coefficient, achieving stronger diffusion along the structure line direction and weaker diffusion perpendicular to the structure line direction.

[0077] The neighboring pixel weights in the cost aggregation process are determined based on the anisotropic diffusion tensor. During cost aggregation, for each pixel... Its neighboring pixels The weights are determined by the diffusion tensor. and pixels and The weights of neighboring pixels are determined by the spatial distance between them. During the cost aggregation process, the weights of neighboring pixels are calculated using a Gaussian function: first, the coordinate difference vector between the current pixel and its neighboring pixels is obtained; then, this difference vector is multiplied on the left by the inverse matrix of the diffusion tensor, and then multiplied on the right by the difference vector itself to obtain a quadratic result; this result is divided by the square of twice the standard deviation of the spatial domain, and the negative value is taken; then, the exponent of this negative value is calculated with the natural constant as the base, and finally, the weight value corresponding to the neighboring pixel is obtained.

[0078] because The diffusion coefficient is larger in the tangential direction, and the corresponding weight is also larger, which allows the parallax to be correctly aggregated along the direction of the structure line; the diffusion coefficient is smaller in the normal direction, and the corresponding weight is also smaller, which blocks the propagation of erroneous parallax across the edge of the structure.

[0079] The initial cost volume is aggregated using semi-global matching based on the weights of neighboring pixels, resulting in an aggregated cost volume. The semi-global matching algorithm propagates the cost in multiple directions; the cost propagation formula for each direction is as follows: in, For the direction of dissemination, For pixels In parallax Along the direction The aggregate value, The penalty coefficient is 1 for a change in parallax. This is the penalty coefficient for parallax changes greater than 1.

[0080] In this embodiment, the neighboring pixel weights are... Introducing the cost propagation formula, we adjust the penalty coefficients in different directions. For propagation along the tangential direction, we decrease the penalty coefficient. and This allows parallax to propagate more easily along the tangent direction; for propagation along the normal direction, the penalty coefficient is increased. and This makes it difficult for parallax to propagate along the normal direction.

[0081] Based on the aggregated cost volume, a disparity map is calculated to generate a dense 3D point cloud that preserves the geometric features of the structure's edges. Normal vector estimation and clustering segmentation are performed on the dense 3D point cloud to extract different surface regions of the structure. Principal component analysis is used to calculate the normal vector of each point, and a region growing algorithm is employed for clustering segmentation. During the region growing process, constraints from a 3D spatial structural line skeleton are introduced, ensuring that the segmentation boundary coincides with the structural line, thus improving the accuracy of surface region segmentation.

[0082] Using a 3D spatial structural skeleton as boundary constraints, anisotropic meshing is performed on the surface region to generate triangular meshes that fit the edges of the structure. The Delaunay triangulation algorithm is employed, inserting the 3D spatial structural skeleton as constraint edges into the triangular mesh during the meshing process, ensuring that the edges of the triangular mesh coincide with the structural lines as much as possible. In regions near the structural lines, the mesh edge length is reduced and the mesh density is increased to maintain geometric accuracy at the edges; in flat regions, the mesh edge length is increased and the mesh density is increased to reduce computational cost.

[0083] A surface fitting based on geometric feature preservation is performed on the triangular mesh. The spatial geometric parameters of the fitted surface are calculated, and the measured distance and angle results of the structure are output based on these parameters. The moving least squares method is used for surface fitting. During the fitting process, constraints from a three-dimensional spatial structural line skeleton are introduced to ensure the continuity and sharpness of the surface at the structural lines. The fitting formula for the moving least squares method is: in, As basis functions, For coefficients, The number of basis functions.

[0084] coefficient By minimizing the weighted least squares error, we obtain: in, For the weight function, For sampling points, The function value at the sampling point. This represents the number of sampling points.

[0085] In this embodiment, the weight function used for surface fitting adopts a Gaussian function. The specific calculation rule is as follows: First, calculate the magnitude of the coordinate difference vector between the point to be fitted and the sampling point. Square the magnitude as the numerator, and use the square of the kernel width as the denominator. Divide the numerator by the denominator, take the negative value of the ratio, and then raise the power of this negative value to the natural constant. The final value obtained is the weight corresponding to that sampling point. The kernel width is used to control the rate at which the weight decays with distance.

[0086] During the fitting process, for sampling points located on the skeleton of the three-dimensional spatial structure line, the value of their weight function is increased, so that the fitted surface is closer to these points, ensuring the geometric accuracy at the structure line.

[0087] Airborne LiDAR point cloud data of the target structure is acquired, and planar and line features are extracted from the data. Airborne LiDAR point cloud data possesses high accuracy and reliability, providing accurate geometric information about the structure. A random sample consensus algorithm is used to extract planar features; for each plane, its equation and normal vector are calculated. A Hough transform-based method is used to extract line features; for each line, its direction vector and the coordinates of a point on the line are calculated.

[0088] Coarse registration of planar and line features from airborne LiDAR point cloud data with corresponding features from dense 3D point clouds, and fine registration based on iterative nearest point (RIB) are performed to obtain coordinate transformation matrices. Coarse registration employs a sample consistency initial registration algorithm, which obtains the initial coordinate transformation matrix by matching feature descriptors. Fine registration uses an iterative nearest point algorithm, which iteratively finds corresponding points between two point clouds, minimizes the distance between corresponding points, and optimizes the coordinate transformation matrix.

[0089] Based on a coordinate transformation matrix, dense 3D point cloud data is fused and compared with airborne LiDAR point cloud data to identify and repair occluded and void areas in the dense 3D point cloud. Due to occlusion issues in oblique photography, some void areas may exist in the dense 3D point cloud. By comparing with the airborne LiDAR point cloud data, these void areas are identified, and then the corresponding areas in the airborne LiDAR point cloud data are used for repair. The repaired point cloud is more complete and accurate, providing a reliable data foundation for subsequent dimensional measurements.

[0090] Output the verified measured results. Based on the fused point cloud, surface mesh extraction and surface fitting are performed again to calculate the measured distances and angles of the structure. The measured results are compared with the design drawings to verify the accuracy of the measured results. The weight allocation of the anisotropic diffusion tensor in different directions is shown in Table 3.

[0091] Table 3 Weight allocation of anisotropic diffusion tensor in different directions Table 3 shows the weight distribution of the anisotropic diffusion tensor in different directions. As can be seen from Table 3, the diffusion weight is largest (1.00) when the direction is aligned with the tangent direction, and smallest (0.00) when the direction is aligned with the normal direction. The diffusion weight gradually decreases as the angle between the direction and the tangent direction increases. This weight distribution allows parallax to propagate sufficiently along the tangent direction of the structural lines, while propagating almost nothing in the normal direction, effectively blocking the propagation of erroneous parallax across the structural edges and maintaining the sharpness of the structural edges.

[0092] This embodiment constructs an anisotropic diffusion tensor and dynamically allocates diffusion weights based on the tangent and normal directions of the 3D spatial structural skeleton, further suppressing mismatch phenomena across structural edges. By fusing and comparing dense 3D point cloud data with airborne lidar point cloud data, occluded and void areas are repaired, improving the completeness and accuracy of the measured results.

Claims

1. A rapid method for measuring structures based on oblique photography, characterized in that, include: Extracting two-dimensional structural line features from multi-view images acquired by oblique photogrammetry; Based on the collinearity condition equation and epipolar geometric constraints, the spatial intersection relationship of the two-dimensional structural line in different views is calculated, and a three-dimensional spatial structural line skeleton is generated by back projection. Using the three-dimensional spatial structure line skeleton as a hard constraint, a graph optimization energy function is constructed. The distance penalty term of the structure line endpoints and edge points in the graph optimization energy function is introduced into the cost aggregation process of multi-view stereo matching, guiding the diffusion algorithm to perform disparity iterative propagation along the tangent direction at the structure edge, generating a dense three-dimensional point cloud that retains the geometric features of the structure edge. Based on the dense 3D point cloud, the surface mesh of the structure is extracted and surface fitting is performed to output the measured results of the structure size.

2. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, The extraction of two-dimensional structural line features from multi-view images acquired by oblique photography includes: A multi-scale feature extraction network is constructed to extract feature maps of different receptive fields from the multi-view images; Perform edge response calculations on the feature map to obtain an initial edge probability map; By combining the image gradient direction information, non-maximum suppression and hysteresis thresholding are performed on the initial edge probability map to filter out edge pixels whose gradient magnitude meets the preset conditions. An edge connection strategy based on local geometric consistency is adopted to perform topological connections on the edge pixels to form continuous two-dimensional structural line features.

3. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, The process of calculating the spatial intersection relationship of the two-dimensional structural lines in different views based on the collinearity condition equation and epipolar geometric constraints, and generating a three-dimensional spatial structural line skeleton through back projection, includes: The two-dimensional structural line is parameterized as the intersection line of a three-dimensional spatial plane, and a three-plane intersection geometric model of the two-dimensional structural line and the spatial plane is established. Based on the collinearity condition equation, the pixel coordinates of the two-dimensional structural line are mapped to the spatial plane, and the epipolar constraint is used to eliminate mismatched line segments of the two-dimensional structural line in different views. Triangulation is performed on the retained matching line segments to solve the spatial parameters of the three-plane intersection geometric model, and back projection is used to generate a continuous three-dimensional spatial structure line skeleton.

4. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, Using the three-dimensional spatial structure line skeleton as a hard constraint, a graph optimization energy function is constructed. The distance penalty terms for the endpoints and edge points of the structure line in the graph optimization energy function are introduced into the cost aggregation process of multi-view stereo matching, including: The endpoints and edge points of the three-dimensional spatial structure skeleton are used as vertices of the graph model. Data items are constructed using the geodesic distance between the vertices, and smoothing items are constructed using the angle between the normal vectors of adjacent vertices. The graph optimization energy function is generated by combining these elements. Extract the distance penalty term from the graph optimization energy function, and superimpose the distance penalty term as a regularization term onto the initial cost volume of multi-view stereo matching. Then, perform filtering and aggregation on the initial cost volume.

5. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, Based on the dense 3D point cloud, the surface mesh of the structure is extracted and surface fitting is performed. The measured dimensions of the structure are then output, including: Normal estimation and clustering segmentation are performed on the dense three-dimensional point cloud to extract different surface regions of the structure; Using the three-dimensional spatial structure skeleton as boundary constraints, the surface region is anisotropically meshed to generate a triangular mesh that fits the edge of the structure. The triangular mesh is fitted with a surface based on geometric feature preservation. The spatial geometric parameters of the fitted surface are calculated, and the measured distance and angle results of the structure are output based on the spatial geometric parameters.

6. The rapid field measurement method for structures based on oblique photography according to claim 2, characterized in that, The construction of a multi-scale feature extraction network to extract feature maps of different receptive fields from the multi-view images includes: An encoder with parallel dilated convolution branches is used to perform feature encoding on the multi-view image to obtain multi-scale feature maps; The multi-scale feature map is input into the feature fusion layer, and the weight distribution of the multi-scale feature map in the spatial dimension is calculated through the spatial attention mechanism, and the weight distribution of the multi-scale feature map in the channel dimension is calculated through the channel attention mechanism. The multi-scale feature map is weighted and fused based on the weight distribution of the spatial dimension and the weight distribution of the channel dimension to output a multi-scale fused feature map.

7. The rapid field measurement method for structures based on oblique photography according to claim 3, characterized in that, The method of using epipolar constraints to eliminate mismatched line segments of the two-dimensional structure line in different views includes: For the target two-dimensional structure line in the reference view, calculate the epipolar equation of the target two-dimensional structure line in the neighboring view; Determine the projection region of the envelope of the target two-dimensional structural line in the neighborhood view; Calculate the overlap ratio between the epipolar segment corresponding to the epipolar equation and the projected region to obtain the epipolar intersection-union ratio value; Two-dimensional structural lines in the neighborhood view whose epipolar intersection-union ratio is lower than a preset threshold are identified as mismatched line segments and removed, while matching line segment pairs that satisfy the epipolar geometric constraints are retained.

8. The rapid field measurement method for structures based on oblique photography according to claim 4, characterized in that, The distance penalty term is superimposed as a regularization term onto the initial cost volume of multi-view stereo matching, and the initial cost volume is filtered and aggregated, including: The distance penalty term is mapped to the disparity space to construct an anisotropic diffusion tensor in the disparity domain; The neighboring pixel weights in the cost aggregation process are determined based on the anisotropic diffusion tensor, such that a large diffusion weight is set in the tangential direction of the three-dimensional spatial structure skeleton, and a small diffusion weight is set in the normal direction of the three-dimensional spatial structure skeleton. The initial cost body is aggregated by semi-global matching based on the neighboring pixel weights, and the aggregated cost body is output.

9. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, Before extracting the two-dimensional structural line features from the multi-view images acquired by oblique photography, the method further includes: Obtain the initial intrinsic and extrinsic parameters and initial pose of the oblique photography camera; Project the two-dimensional structural line feature onto three-dimensional space and calculate the reprojection error of the two-dimensional structural line feature between adjacent views; A bundle adjustment optimization model based on line features is constructed with minimizing the reprojection error as the objective function. The initial intrinsic and extrinsic parameters and initial pose of the oblique photography camera are jointly iteratively solved and updated to obtain the optimized target camera parameters and target pose. Based on the target camera parameters and target pose, subsequent spatial intersection calculations are performed.

10. The rapid field measurement method for structures based on oblique photography according to claim 1, characterized in that, Following the output of the measured structural dimensions, the following is also included: Acquire airborne lidar point cloud data of the target structure, and extract planar and line features from the airborne lidar point cloud data; The planar and line features in the airborne lidar point cloud data and the corresponding features in the dense 3D point cloud are coarsely registered and finely registered based on iterative nearest points to obtain the coordinate transformation matrix. Based on the coordinate transformation matrix, the dense 3D point cloud and the airborne lidar point cloud data are fused and compared to identify and repair the occlusion and void areas in the dense 3D point cloud, and the verified measured results are output.

Images