Original terrain air-ground integrated measurement method

Abstract

The invention discloses an air-ground integrated measurement method for an original terrain, and relates to the technical field of engineering surveying, and the method comprises the steps: S1, carrying multi-sensor aerial survey by an unmanned aerial vehicle, obtaining an image and POS positioning data, scanning a blind area of the unmanned aerial vehicle by using a three-dimensional laser scanner, obtaining laser point cloud data, and carrying out the scanning of the laser point cloud data; constructing a multi-source data set including optical images, POS positioning data and laser point cloud, and forming an air-ground complementary data acquisition mode; s2, geometric correction and aerial triangulation encryption are carried out on the unmanned aerial vehicle image, denoising and multi-station registration are carried out on the laser point cloud data, and the data quality is improved; s3, unifying a multi-source data space reference through coordinate system unification, feature fusion and decision-level fusion, and realizing feature complementation; s4, constructing a digital elevation model and a digital surface model based on the fused data, and generating a topographic mapping map; the method breaks through the limitation of a single technology through air-ground technology complementation, and ensures the global data coverage of complex terrains.

Description

Technical Field

[0001] The present invention relates to the technical field of engineering surveying, and particularly relates to an integrated measurement method for original terrain and open space. Background Art

[0002] In the field of engineering surveying, high-precision mapping of complex terrain areas (such as mountains, hills, canyons, etc.) has always been a technical problem faced by the industry. In such areas, the terrain often has large undulations, dense vegetation coverage, and poor accessibility in some areas. Traditional surveying methods are difficult to meet the requirements of high-precision and full-coverage of terrain data in fields such as modern engineering construction, geological disaster monitoring, and natural resource surveys.

[0003] Currently, traditional terrain surveying methods mainly include manual on-site surveying and satellite remote sensing surveying. Manual on-site surveying relies on surveying personnel to use equipment such as total stations and levels for point-by-point measurement. In complex terrain areas, surveying personnel not only need to face high-risk operations such as climbing steep mountains and crossing dense jungles, but also are limited by terrain conditions, resulting in extremely low measurement efficiency. For example, in steep slope areas with a slope exceeding 60°, it is difficult to set up manual surveying equipment, the distribution of measurement points is sparse, and it is difficult to obtain complete terrain data; in areas with dense vegetation, the instrument line of sight is easily blocked, resulting in large errors in measurement data, making it difficult to achieve high-precision measurement, and the operation cycle of manual surveying is long, making it difficult to meet the timeliness requirements of large-scale terrain mapping.

[0004] Although satellite remote sensing surveying can quickly obtain large-scale terrain images, it is limited by spatial resolution and cannot accurately capture minute geomorphic features and local terrain changes in complex terrains. For example, in mountain terrain mapping, satellite remote sensing is difficult to accurately reflect local terrain details such as valleys and gullies, and in bad weather conditions such as clouds, fog, and rain, the imaging quality seriously deteriorates, and effective measurement data cannot be obtained, resulting in difficulties in ensuring the timeliness and accuracy of the measurement.

[0005] With the development of surveying and mapping technologies, unmanned aerial vehicle (UAV) aerial surveying and three-dimensional laser scanning technologies have gradually been applied to the field of terrain surveying. UAV aerial surveying, with its advantages of flexibility and high operation efficiency, can quickly obtain large-area terrain image data. However, in complex terrain areas, due to limitations in flight altitude and perspective, there are many measurement blind spots. For example, in areas with large vertical drops such as cliffs and deep valleys, UAVs cannot approach for measurement due to safety distance restrictions, resulting in data gaps; at the same time, in environments with uneven lighting and unstable airflows, the data obtained by UAV aerial surveying is prone to defects such as geometric deformation and insufficient overlap, affecting the measurement accuracy. Three-dimensional laser scanning technology can accurately obtain the three-dimensional information of ground objects through point cloud data. However, in large-area terrain surveying, the three-dimensional laser scanner needs to frequently change measurement stations, the data acquisition efficiency is low, and the equipment operation is complex, with high professional requirements for operators, making it difficult to quickly complete full-coverage measurement of complex terrain areas.

[0006] Therefore, how to overcome the limitations of a single technology in complex terrain surveying and achieve high-precision, high-efficiency, and full-coverage terrain surveying has become a technical problem that urgently needs to be solved currently. Summary of the Invention

[0007] The purpose of the present invention is to provide an integrated method for surveying the original terrain and open spaces. By using unmanned aerial vehicle (UAV) aerial survey to quickly obtain high-resolution images and POS positioning data over a large area, covering open areas, and at the same time using a three-dimensional laser scanner to supplement the measurement of the blind spots of the UAV, a "complementary ground and air" data acquisition mode is formed, solving the limitations of a single technology in terrain adaptability and ensuring full-domain data coverage of complex terrains.

[0008] The present invention is realized through the following technical solutions:

[0009] An integrated method for surveying the original terrain and open spaces includes the following steps:

[0010] S1. Conduct aerial survey by carrying multiple sensors on a UAV to obtain high-resolution images and POS positioning data. At the same time, use a three-dimensional laser scanner to scan the blind spots of the UAV to obtain laser point cloud data, and construct a multi-source data set including optical images, POS positioning data, and laser point clouds.

[0011] S2. Perform geometric correction and aerial triangulation encryption processing on the UAV images, and perform denoising and multi-station registration processing on the laser point cloud data.

[0012] S3. Through coordinate system unification, feature fusion, and decision-level fusion, make the spatial reference of multi-source data unified and the features complementary.

[0013] S4. Based on the fused data, construct a digital elevation model and a digital surface model, and generate a topographic surveying map.

[0014] In this solution, the mobility and flexibility of the UAV are utilized to collect large-scale high-resolution optical images and POS data. Synchronously, a three-dimensional laser scanner is used to accurately supplement the point clouds in the blind spots of the UAV (such as steep slopes, dense vegetation areas, etc.), forming a full-domain data acquisition mode of "surface coverage + point precision", solving the problem of coverage blind spots of a single technology in terrain adaptability; and through geometric correction, aerial triangulation encryption of UAV images and denoising and registration of laser point clouds, combined with coordinate system unification, feature-level and decision-level fusion, the spectral texture advantages of images and the three-dimensional geometric elevation advantages of point clouds are complementary, significantly improving the terrain modeling accuracy, making up for the error defects of a single data source in plane positioning, elevation inversion, and ground object recognition. Then, relying on the automated data processing process and full-process quality control, the fieldwork time and in-field manual intervention are greatly shortened, forming an efficient closed-loop from data acquisition to result output, and finally achieving the measurement goal of "high-precision, high-efficiency, and full-coverage" under complex terrains.

[0015] In one embodiment, the route planning of UAV aerial survey is dynamically adjusted according to the terrain complexity, and the adjustment formula is:

[0016]

[0017] where H is the real-time flight altitude (m), H0 is the initial altitude, V is the real-time wind speed (m / s), k1 is the wind speed influence coefficient, E is the real-time exposure value, E0 is the initial exposure, L is the real-time light intensity, and L0 is the standard light intensity.

[0018] In this solution, the route planning of UAV aerial survey is dynamically adjusted according to the terrain complexity. Its function is to ensure the quality of aerial survey data and flight safety under complex terrains by responding to changes in environmental parameters in real time. Specifically, in areas with strong winds such as mountainous areas, the flight altitude is increased by increasing the k1 coefficient to avoid the instability of the UAV's attitude and image blurring caused by strong winds, and at the same time ensure the image coverage accuracy in different terrain areas; and the camera exposure value E is adaptively adjusted according to the real-time light intensity L. In an environment with uneven light (such as mountainous areas where shaded areas and strong light areas alternate, and water area reflection areas), overexposure or underexposure problems can be effectively avoided, and the consistency of image color and brightness can be guaranteed; the combination of the two enables the UAV to fly safely in complex terrains and obtain high-resolution image data with small geometric distortion.

[0019] In one embodiment, when the three-dimensional laser scanner performs supplementary survey in the blind area, the density of the acquired laser point cloud data is calculated and noise is removed through the following formula:

[0020]

[0021] where Density i is the local density of point i, n is the number of points within the neighborhood radius r, and Distance ij is the Euclidean distance between point i and neighborhood point j; when Density i <T D and Distance ij >T Dist , it is determined as a noise point, T D is the density threshold, and T Dist is the distance threshold.

[0022] In this solution, the calculation of local density can measure the density of the point cloud around point i. The fewer the number of points n within the neighborhood radius r, the lower the density, indicating that this point may be an isolated noise; the Euclidean distance quantifies the spatial separation degree between point i and neighborhood points. The larger the distance, the more likely it is that this point deviates from the true terrain features. When a certain point simultaneously satisfies that the density is lower than the threshold T D and the distance exceeds the threshold TDist When it is the case, it is determined as a noise point and removed, which can effectively filter out the outliers generated by environmental interference (such as flying birds, reflective objects) or equipment errors during the scanning process, and retain the high-density and closely related point clouds of the real terrain and features.

[0023] In one embodiment, in actual measurement, due to factors such as large measurement range and complex terrain, scanning is often required at multiple measurement stations, and the point cloud data of each measurement station is in different coordinate systems. Therefore, registration is required. When the three-dimensional laser scanner performs supplementary measurement in the blind area, the acquired laser point cloud data is registered at multiple stations through the following method, including:

[0024] Rough registration: Find the corresponding point pairs between different measurement stations through the geometric features of the laser point cloud, and then calculate the initial transformation matrix through these corresponding point pairs to roughly align the point cloud data of different measurement stations;

[0025] During the rough registration process, the optimal transformation matrix is evaluated and selected through the following objective function:

[0026]

[0027] where, R is the rotation matrix, t is the translation vector, are the point cloud coordinates of the same-name target points of the two measurement stations;

[0028] Fine registration: Continuously calculate the distance between corresponding points through iteration and adjust the transformation matrix to minimize the distance between corresponding points until the preset convergence condition is met;

[0029] During the fine registration process, the distance error between corresponding points is minimized through iteration:

[0030]

[0031] Until the root mean square error (RMSE) < 1mm, the spatial splicing of the multi-station laser point cloud is achieved, where, p i is the i-th point in the source point cloud, q i is the point in the target point cloud corresponding to the i-th point p i in the source point cloud.

[0032] In this solution, in the rough registration stage, the initial transformation matrix is calculated by finding the point pairs corresponding to the geometric features of the laser point cloud, which can roughly align the point cloud data of different survey stations. By evaluating and screening the optimal transformation matrix with the objective function, the most realistic matrix can be determined among many possible transformations, avoiding large deviations caused by random selection and providing a better starting state for subsequent fine registration. In the fine registration stage, the distance between corresponding points is continuously calculated iteratively and the transformation matrix is adjusted to minimize the distance error of the corresponding points until the root mean square error is less than 1 mm, enabling high-precision spatial stitching, accurately fusing the point cloud data of different survey stations in a unified coordinate system, eliminating the remaining errors after rough registration, and ensuring that the finally obtained point cloud data can accurately reflect the true terrain and ground object information.

[0033] In one embodiment, in the aerial triangulation encryption of UAV images, multi-scale extrema are detected by the Gaussian pyramid and the Difference of Gaussian (DOG) operator. The formula is:

[0034] DOG(x, y, σ) = (G(kσ) - G(σ)) * I(x, y)

[0035] where G(σ) is the Gaussian kernel function, k is the adjacent scale factor, and I(x, y) is the image gray value.

[0036] In this solution, the Gaussian pyramid generates a multi-resolution image sequence by performing Gaussian blurring and downsampling on the image at different scales, enabling the algorithm to capture multi-level features from the macroscopic terrain contour to the microscopic ground object details. The Difference of Gaussian operator highlights the scale-invariant extrema (such as corner points, edge points, etc.) by calculating the difference between adjacent scale Gaussian blurred images (i.e., the DOG response). These extrema have strong adaptability to illumination changes, viewing angle differences, and terrain undulations. This process can effectively extract stable feature points, providing sufficient and reliable homologous point pairs for subsequent image matching, significantly improving the accuracy of solving the exterior orientation elements of the image and calculating the encrypted point coordinates in aerial triangulation encryption. Especially in complex terrain areas such as mountainous areas and vegetation-covered areas, it can reduce the feature matching errors caused by scale changes and ensure the reliability of aerial triangulation encryption.

[0037] In one embodiment, the multi-source data coordinate system is unified through the following conversion formula:

[0038]

[0039] where R is the rotation matrix, ε x 、ε y 、ε z are the rotation parameters, m is the scale factor, and (X0, Y0, Z0) are the translation parameters.

[0040] In this solution, the rotation matrix R, the rotation parameters ε x 、εy , ε z , together with the scale factor m and the translation parameters (X0, Y0, Z0), constitute a rigid body transformation model, which can accurately transform the heterogeneous data (such as image coordinates and point cloud coordinates) obtained by different sensors (UAVs, laser scanners) from their respective independent coordinate systems to a unified local plane coordinate system. By performing least squares adjustment on the measured coordinates of ground control points to solve for the seven parameters, the conversion error can be controlled within millimeters (plane error ≤ 5 mm, elevation error ≤ 8 mm), ensuring seamless stitching of multi-source data in terms of spatial position and avoiding problems such as terrain misalignment and feature offset caused by inconsistent coordinate references.

[0041] In one embodiment, in the feature-level fusion, the spectral indices extracted from UAV images satisfy the following formula:

[0042]

[0043] where ρ NIR , ρ RED , ρ GREEN are the reflectances in the near-infrared, red, and green bands, and λ1 ≥ λ2 ≥ λ3 are the eigenvalues of the covariance matrix.

[0044] In this solution, the spectral indices extracted from UAV images are combined with the eigenvalues of the covariance matrix of the laser point cloud to achieve complementary advantages of the spectral information of the images and the geometric information of the point cloud. Specifically, the spectral indices, through the ratio of the reflectances in the near-infrared, red, and green bands, can sensitively identify targets such as vegetation and water bodies (e.g., NDVI > 0.3 distinguishes vegetation, NDWI > 0.5 locates water bodies), providing semantic information on the types of features; while the eigenvalues of the covariance matrix of the laser point cloud reflect the geometric properties of the terrain surface (describing the degree of terrain undulation). After cascading the two, a composite feature vector is formed, which not only retains the spectral sensitivity of the images to vegetation coverage and water body distribution but also integrates the geometric description ability of the point cloud for terrain slope and roughness, providing multi-source feature support for subsequent feature classification and terrain analysis.

[0045] In one embodiment, for the complex terrain features in the vegetation-covered area and the water-land interface area, to achieve intelligent collaboration and accuracy optimization of multi-source data, in the decision-level fusion, the digital elevation optimization formula for the vegetation-covered area is:

[0046] Z DEM = 0.7·Z laser + 0.3·Z img

[0047] where Z laser is the elevation of the ground points in the laser point cloud, and Z img is the elevation of the UAV point cloud;

[0048] The shoreline extraction in the land-water interface area needs to meet the following conditions simultaneously:

[0049] NDWI > 0.5 and |Z i -Z j | > 0.3m, where Z i and Z j are the elevations of adjacent points.

[0050] In this solution, by combining the water body spectral characteristics of UAV images and the elevation mutation characteristics of laser point clouds, it is ensured that the shoreline positioning simultaneously meets the dual constraints of spectral characteristics and terrain mutation, avoiding misjudgments of a single data source (such as misjudgment due to water surface reflection relying only on images, and missing judgment due to sparseness relying only on point clouds), and controlling the shoreline extraction accuracy within 10 cm.

[0051] In one embodiment, during the generation of the topographic surveying and mapping map, it is necessary to combine the high-density surface points of the laser point cloud and the texture points of the UAV image, and optimize the terrain details through the following formula:

[0052]

[0053] Z(S0) is the terrain attribute value at the position S0 to be predicted, n is the number of known data points used for interpolation, Z(S i ) is the terrain attribute value at the i-th known position S i ), and λ i is the weight coefficient of the i-th known point.

[0054] In this solution, the high-density surface points of the laser point cloud provide high-precision discrete elevation control data to ensure the elevation accuracy of key terrain positions, while the texture points of the UAV image supplement the large-scale continuous planar position and surface texture information through dense matching of point clouds. When the two are fused through the interpolation formula, the weight coefficient can be dynamically adjusted according to the point cloud density and image resolution, which can not only retain the high-precision characteristics of the laser point cloud in the blind area, but also use the image texture points to fill the low-density area. This operation can effectively repair the defects of a single data source, make the generated DEM smooth and natural in the plain area, retain the micro-topographic undulations in the mountain area, and make the ground object boundaries clear after the DSM is overlaid with the image texture, finally realizing the complete restoration of terrain details.

[0055] In one embodiment, to avoid problems such as texture misalignment and blurring, and improve the quality and accuracy of the model, when generating the digital surface model, the highest elevation within the grid is taken and the UAV image texture is overlaid, and the texture mapping error ≤ 1 pixel.

[0056] In summary, compared with the prior art, the present invention mainly has the following advantages and beneficial effects:

[0057] 1. The present invention constructs an air-ground complementary data acquisition mode of "area coverage + point precision" through the large-scale rapid coverage of UAV aerial survey and the precise supplementary survey of the blind areas of 3D laser scanning, solving the problem of incomplete coverage of traditional single technologies in complex terrains;

[0058] 2. Through coordinate system unification, feature-level fusion, and decision-level fusion, the present invention realizes the deep complementarity between the planar texture advantages of UAV images and the elevation geometric advantages of laser point clouds, solving problems such as elevation distortion caused by vegetation occlusion of UAVs and large stitching errors in laser scanning in traditional methods. BRIEF DESCRIPTION OF THE DRAWINGS

[0059] The drawings described herein are used to provide a further understanding of the embodiments of the present invention, form a part of this application, and do not constitute a limitation on the embodiments of the present invention. In the drawings:

[0060] Figure 1 is a flowchart of the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0061] To make the objectives, technical solutions, and advantages of the present invention clearer and more understandable, the present invention will be further described in detail below in conjunction with the embodiments and the drawings. The illustrative embodiments and descriptions of the present invention are only used to explain the present invention and do not constitute a limitation on the present invention.

[0062] Embodiment

[0063] Traditional topographic survey methods often struggle to meet the challenges of complex terrains, suffering from problems such as measurement blind areas, insufficient accuracy, and low efficiency. The following embodiments provide an integrated air-ground measurement method for the original terrain by combining UAV aerial survey and 3D laser scanning technology, realizing an air-ground complementary data acquisition mode and effectively solving these problems. The implementation process is as follows:

[0064] As Figure 1 shown, it includes the steps:

[0065] S1. Conduct aerial survey by mounting multiple sensors on a UAV to obtain high-resolution images and POS positioning data. At the same time, scan the blind areas of the UAV with a 3D laser scanner to obtain laser point cloud data, and construct a multi-source data set including optical images, POS positioning data, and laser point clouds;

[0066] Specifically, first, conduct aerial survey work using a drone equipped with multiple sensors. Among them, select a suitable drone according to factors such as the size of the measurement area, terrain complexity, and flight environment. For example, for a large open area, a fixed-wing drone with a long endurance can be selected; for a small-scale and complex terrain area, a multi-rotor drone has more advantages as it can hover flexibly and fly at low altitudes. For the multi-sensor configuration: select a camera with high pixels, a large aperture, and good color reproduction ability to obtain clear and detailed image data; during installation, the POS positioning system should ensure a rigid connection between the sensor and the drone platform to reduce the impact of vibration and displacement on the measurement accuracy.

[0067] Before conducting route planning, collect information such as the terrain, landform, and meteorology of the measurement area, and understand the obstacles and no-fly zones within the area. Generally speaking, according to the scope and shape of the measurement area, the heading overlap should be no less than 70%, and the side overlap should be no less than 60% to ensure the accuracy of subsequent data processing. Before aerial survey, it is also necessary to dynamically adjust the route planning according to the terrain complexity. For complex terrain areas, such as mountains and canyons, use the formula H = H0 + k1·V, to determine the real-time flight altitude H, where H is the real-time flight altitude (m), H0 is the initial altitude, V is the real-time wind speed (m / s), k1 is the wind speed influence coefficient, E is the real-time exposure value, E0 is the initial exposure, L is the real-time light intensity, and L0 is the standard light intensity. Through this method, high-resolution images and accurate POS positioning data can be obtained, efficiently covering open areas.

[0068] Then, use a three-dimensional laser scanner to scan the blind areas that are difficult for the drone to cover (such as steep slopes and dense vegetation areas) to obtain laser point cloud data. During the scanning process, calculate the density and remove noise from the obtained laser point cloud data. Calculate the local density of point i through the formula where n is the number of points within the neighborhood radius r. When Density i <T D and Distance ij >T Dist , it is determined as a noise point and removed. Among them, Distance ij is the Euclidean distance between point i and neighborhood point j, (x i , y i , z i ) is the three-dimensional coordinates of the i-th point in the current measuring station, (x j , y j , z j ) is the three-dimensional coordinates of the point corresponding to point i among the j-th points in the neighborhood, T D is the density threshold, T Distis the distance threshold. Finally, a multi-source dataset containing optical images, POS positioning data, and laser point clouds is constructed.

[0069] The construction process is as follows:

[0070] S2. Geometric correction and aerial triangulation encryption are performed on the UAV images, and denoising (already described in the previous text) and multi-station registration are performed on the laser point cloud data;

[0071] Specifically, geometric correction is performed on the images obtained by the UAV. The geometric correction can select existing polynomial models or collinearity equation models. In the aerial triangulation encryption link, Gaussian pyramids and differential Gaussian operators are used to detect multi-scale extrema in the following way:

[0072] DOG(x, y, σ) = (G(kσ) - G(σ)) * I(x, y)

[0073] where G(σ) is the Gaussian kernel function, k is the adjacent scale factor, and I(x, y) is the image gray value. Stable feature points are extracted by this method to provide sufficient and reliable homologous point pairs for subsequent image matching, and improve the accuracy of exterior orientation element calculation and encrypted point coordinate calculation.

[0074] Then, due to factors such as large measurement range and complex terrain, scanning often needs to be carried out at multiple measurement stations, and the point cloud data of each measurement station is in different coordinate systems. Therefore, registration is required. Multi-station registration is divided into two stages: rough registration and fine registration;

[0075] Rough registration stage: Corresponding point pairs between different measurement stations are found through the geometric features of the laser point cloud, and then the initial transformation matrix is calculated through these corresponding point pairs to roughly align the point cloud data of different measurement stations;

[0076] In the rough registration process, the optimal transformation matrix is evaluated and selected through the following objective function:

[0077]

[0078] where R is the rotation matrix and t is the translation vector, are the homologous target point cloud coordinates of two measurement stations. Through this calculation process, the initial coordinate systems between different measurement stations are aligned, and the overall position error of the point cloud is reduced from "meter level" to "centimeter level", providing a feasible initial value for fine registration.

[0079] Fine registration: Based on the rough registration result, dense matching is performed on all point clouds. By continuously iteratively calculating the distance between corresponding points and adjusting the transformation matrix, the distance between corresponding points is minimized until the preset convergence condition is met;

[0080] In the fine registration process, the distance error between corresponding points is minimized through iteration:

[0081] Among them, p i is the i-th point in the source point cloud, and q i is the point in the target point cloud corresponding to the i-th point p of the source point cloud i . Until the root mean square error (RMSE) < 1 mm, the local detail error remaining after coarse registration is solved, enabling millimeter-level accurate spatial stitching of multi-station lidar point clouds and ensuring seamless connection of point cloud boundaries.

[0082] Then enter the next step of multi-source data fusion:

[0083] S3. Through coordinate system unification, feature fusion, and decision-level fusion, the spatial reference of multi-source data is unified and features are complementary;

[0084] First, the multi-source data coordinate system is unified through the following conversion formula:

[0085]

[0086] Among them, R is the rotation matrix, ε x , ε y , ε z are rotation parameters, m is the scale factor, (X0, Y0, Z0) are translation parameters, and (X s , Y s , Z s ) are the coordinates in the original coordinate system. The heterogeneous data obtained by different sensors are accurately converted from their respective independent coordinate systems to a unified local plane coordinate system, eliminating coordinate system differences and laying a foundation for subsequent feature fusion.

[0087] To organically combine the features of different data sources to give full play to their respective advantages and improve the usability of data and the accuracy of analysis, the features of different data sources are subjected to feature fusion, and feature fusion includes spectral feature and geometric feature extraction;

[0088] Among them, in the spectral features, spectral indices are extracted for UAV image data, and the extraction of spectral indices satisfies the following formula:

[0089]

[0090] Among them, ρ NIR , ρ RED , ρ GREEN are the reflectances of the near-infrared, red, and green bands, and λ1 ≥ λ2 ≥ λ3 are the eigenvalues of the covariance matrix. These spectral indices can reflect information such as vegetation coverage and water body distribution of ground objects.

[0091] Among them, in geometric features, the aforementioned acquired lidar point cloud data can provide rich geometric features, such as the density, curvature, normal vector, etc. of the point cloud. By using the point cloud terrain analysis method to calculate the curvature of the local neighborhood of the point cloud, the degree of bending of the terrain surface can be quantitatively described, so as to distinguish different terrain features such as planes and curved surfaces. Calculating the curvature of the point cloud can distinguish different terrain features such as planes and curved surfaces.

[0092] Cascade the extracted spectral features and geometric features to form a comprehensive feature vector. Assume that the spectral feature vector is S = [s1, s2,..., s m , S m is the multi-band reflectance data of the UAV image, and the geometric feature vector is G = [g1, g2,..., g n , G n is the lidar point cloud data. Then the comprehensive feature vector F can be expressed as F = [s1, s2,..., s m , g1, g2,..., g n .

[0093] In actual operation, it is necessary to ensure that the scales and ranges of different features are similar to avoid some features having too much influence on the subsequent matching algorithm. The normalization method can be used to scale each feature value to the interval [0, 1]:

[0094]

[0095] where x is the original feature value, x min , x max are the minimum and maximum values of this feature respectively.

[0096] Then, use the feature matching algorithm, such as the descriptor-based matching algorithm. Commonly used feature descriptors include SIFT (Scale-Invariant Feature Transform), SURF (Speeded-Up Robust Features), etc. These descriptors have characteristics such as scale invariance and rotation invariance, and can accurately describe feature points in different images or point cloud data; use the feature matching algorithm to find homologous feature points between the feature descriptors of different data sources. These homologous feature points can be used as the basis for subsequent decision-level fusion for further analysis and processing.

[0097] Decision-level fusion is to make a comprehensive decision on multi-source data based on feature fusion to obtain more accurate and reliable results. For the complex terrain features in the vegetation-covered area and the water-land junction area, the lidar point cloud data can penetrate the vegetation to obtain the elevation information of the ground surface, but may be interfered by the vegetation; while the UAV image data can provide the spectral information of the vegetation, but it is difficult to accurately obtain the elevation of the ground surface. Therefore, the weighted average method can be used for elevation optimization:

[0098] Z DEN= 0.7·Z lascr + 0.3·Z img

[0099] where Z DEM is the elevation value of the final digital elevation model (DEM), Z laser is the elevation of the ground points of the lidar point cloud, and Z img is the elevation of the UAV point cloud;

[0100] In the water-land boundary area, it is necessary to comprehensively consider spectral features and elevation features to extract the shoreline. Certain threshold conditions can be set. When the shoreline extraction in the water-land boundary area satisfies NDWI > 0.5 and |Z i - Z j | > 0.3 m, where Z i and Z j are the elevations of adjacent points, it is determined as a shoreline point. In this way, the water-land boundary line can be accurately extracted.

[0101] After the multi-source data fusion is completed, the following steps are entered:

[0102] S4. Construct a digital elevation model and a digital surface model based on the fused data, and generate a topographic surveying and mapping map;

[0103] Specifically, to accurately extract ground points from the fused data, algorithms such as progressive morphological filtering or cloth simulation filtering can be used. Since the collected point data is discrete, the following interpolation method is needed to obtain the elevation values within the entire area. The formula: where Z(x, y) is the elevation value of the point to be interpolated (x, y), Z i is the elevation value of the i-th known point, and w i is the weight of the i-th known point.

[0104] After obtaining the elevation values of each grid point through interpolation calculation, organize them into a regular grid matrix to form an initial digital elevation model. In some embodiments, to further refine the digital elevation model, fracture line information (such as rivers, ridge lines, etc.) can be combined, and these positions are constrained during the interpolation process so that the digital elevation model can more accurately reflect the true characteristics of the terrain.

[0105] Finally, evaluate the accuracy of the constructed digital elevation model. If the accuracy does not meet the standard, it is necessary to readjust the interpolation parameters or collect more ground point data.

[0106] Meanwhile, also based on the previously fused multi-source data, the research area is divided into several grids. For each grid, the elevation value of the highest point is extracted from the point cloud data containing all ground objects as the elevation value of the grid to construct the basic framework of the digital surface model. The texture of the UAV image is overlaid on the digital surface model to enhance the visualization effect of the model. During the texture mapping process, it is necessary to ensure that the texture mapping error ≤ 1 pixel, which can be achieved by using a texture mapping algorithm based on feature matching. Then, post-processing operations such as smoothing and edge detection are performed on the digital surface model to remove the jagged edges and noise in the model and improve the quality of the model.

[0107] Generate a topographic survey map by combining the constructed digital elevation model and digital surface model, and optimize the topographic details through interpolation using the following formula:

[0108]

[0109] Z(S0) is the topographic attribute value at the position S0 to be predicted, n is the number of known data points used for interpolation, Z(S i ) is the topographic attribute value at the i-th known position S i , and λ i is the weight coefficient of the i-th known point;

[0110] For some areas with complex topographic changes, the number of known points can be increased or the weight coefficients can be adjusted to improve the expression accuracy of topographic details.

[0111] *************************************************

[0112] The specific embodiments described above further elaborate on the purpose, technical solutions, and beneficial effects of the present invention. It should be understood that the above are only specific embodiments of the present invention and are not used to limit the protection scope of the present invention. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. An integrated measurement method for original terrain and open space, characterized in that It includes the following steps: S1. Conduct aerial surveys by carrying multiple sensors on an unmanned aerial vehicle (UAV) to obtain high-resolution images and POS positioning data. At the same time, scan the blind area of the UAV with a three-dimensional laser scanner to obtain laser point cloud data, and construct a multi-source dataset containing optical images, POS positioning data, and laser point clouds; S2. Perform geometric correction and aerial triangulation encryption processing on the UAV images, and perform denoising and multi-station registration processing on the laser point cloud data; S3. Through coordinate system unification, feature fusion, and decision-level fusion, unify the spatial reference and complement the features of multi-source data; S4. Based on the fused data, construct a digital elevation model and a digital surface model, and generate a topographic survey map.

2. The integrated measurement method for original terrain and open space according to claim 1, wherein The route planning of UAV aerial surveys is dynamically adjusted according to the terrain complexity, and the adjustment formula is: where H is the real-time flight altitude (m), H0 is the initial altitude, V is the real-time wind speed (m / s), k1 is the wind speed influence coefficient, E is the real-time exposure value, E0 is the initial exposure, L is the real-time light intensity, and L0 is the standard light intensity.

3. The integrated measurement method for original terrain and open space according to claim 2, characterized in that, When the three-dimensional laser scanner performs supplementary measurement in the blind area, calculate and remove noise from the obtained laser point cloud data density through the following formula: Among them, Density i is the local density of point i, n is the number of points within the neighborhood radius r, and Distance ij is the Euclidean distance between point i and neighborhood point j; when Density i <T D and Distance ij >T Dist it is determined as a noise point, T D is the density threshold, and T Dist is the distance threshold.

4. The integrated measurement method for original terrain and open space according to claim 3, characterized in that When the three-dimensional laser scanner performs supplementary measurement in the blind area, perform multi-station registration on the obtained laser point cloud data through the following methods, including: Rough registration: Find corresponding point pairs between different measurement stations through the geometric features of the laser point cloud, and then calculate the initial transformation matrix through these corresponding point pairs to roughly align the point cloud data of different measurement stations; During the rough registration process, evaluate and select the optimal transformation matrix through the following objective function: where \(R\) is the rotation matrix and \(t\) is the translation vector, are the coordinates of the homologous target point clouds of the two survey stations; Fine registration: Continuously iterate to calculate the distance between corresponding points and adjust the transformation matrix to minimize the distance between corresponding points until the preset convergence condition is met; During the fine registration process, iteratively minimize the distance error between corresponding points: Until the root mean square error (RMSE) < 1 mm, enabling the spatial stitching of multi-station laser point clouds, where p i is the i-th point in the source point cloud, q i is the point in the target point cloud corresponding to the i-th point p i in the source point cloud.

5. The integrated measurement method for original terrain and open space according to claim 4, characterized in that In the aerial triangulation encryption of UAV images, detect multi-scale extrema through the Gaussian pyramid and the difference of Gaussian operator, and the formula is: DOG(x,y,σ)=(G(kσ)-G(σ))*I(x,y) where G(σ) is the Gaussian kernel function, k is the adjacent scale factor, and I(x,y) is the image gray value.

6. The integrated measurement method for original terrain and open space according to claim 5, characterized in that, The unification of the multi-source data coordinate system is carried out through the following conversion formula: where R is the rotation matrix, ε x , ε y , ε z are rotation parameters, m is the scale factor, and (X0, Y0, Z0) are translation parameters.

7. A method for integrated measurement of original terrain and open space according to claim 6, characterized in that, In the feature-level fusion, the spectral index extracted from the UAV image satisfies the following formula: n = eigenvector(cov(N(p)))[0], where ρ NIR , ρ RED , ρ GREEN are the reflectances in the near-infrared, red, and green bands, and λ1≥λ2≥λ3 are the eigenvalues of the covariance matrix.

8. A method for integrated measurement of original terrain and open space according to claim 7, characterized in that In the decision-level fusion, the digital elevation optimization formula for the vegetation-covered area is: Z DEM = 0.7·Z laser + 0.3·Z img Among them, Z laser is the elevation of the ground points of the laser point cloud, and Z img is the elevation of the UAV point cloud; The extraction of the shoreline in the water-land junction area needs to satisfy both: NDWI > 0.5 and |Z i -Z j | > 0.3m, Z i and Z j and Z are the elevations of adjacent points.

9. The integrated measurement method of original terrain and open space according to claim 8, characterized in that During the process of generating the topographic survey map, it is necessary to combine the high-density surface points of the laser point cloud and the texture points of the UAV image, and interpolate and optimize the terrain details through the following formula, and the formula is: Z(S0) is the terrain attribute value at the position S0 to be predicted, n is the number of known data points for interpolation, Z(S i ) is the terrain attribute value at the i-th known position S i , and λ i is the weight coefficient of the i-th known point.

10. A method for integrated measurement of original terrain and open space according to claim 9, characterized in that, When generating the digital surface model, take the elevation of the highest point within the grid and overlay the texture of the UAV image, and the texture mapping error ≤ 1 pixel.

Images