Method for measuring three-dimensional motion vector of glacier based on low-altitude oblique photography of unmanned aerial vehicle

Abstract

The present application relates to the field of three-dimensional motion vector measurement, and particularly relates to a glacier three-dimensional motion vector measurement method based on low-altitude oblique photography of a UAV. The content includes: collecting multi-angle image data to obtain a multi-temporal original observation set; preprocessing the multi-temporal original observation set to obtain a preprocessed observation set; determining a three-dimensional point set based on the preprocessed observation set; introducing a geometric-reflection coupling consistency matching algorithm based on the preprocessed observation set to register three-dimensional points in the three-dimensional point set and determine matching points in the next phase; and obtaining a glacier three-dimensional motion vector based on the three-dimensional points in the current phase and the matching points in the next phase. The present application solves the problems of unstable cross-temporal three-dimensional corresponding point matching, large matching error, and space reference drift between multi-temporal point clouds in a low-texture and strong-reflection environment of a glacier, which leads to inaccurate calculation of a glacier three-dimensional motion vector.

Description

Technical Field

[0001] This invention relates to the field of three-dimensional motion vector measurement, and more particularly to a method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles. Background Technology

[0002] With the intensification of global climate change, accurate monitoring of glacier movement has become a crucial foundation for glacier dynamics research and disaster risk assessment. Existing methods for glacier movement monitoring mainly include total station monitoring based on ground measurements, GNSS positioning measurements, and displacement inversion methods based on satellite remote sensing imagery. While ground measurement methods offer high accuracy, they suffer from long deployment cycles, limited spatial coverage, and significant challenges in high-altitude mountainous regions. Satellite remote sensing methods, although covering a wide area, are limited by spatial resolution, revisit frequency, and atmospheric influences, making it difficult to obtain high-precision three-dimensional motion information in small- to medium-scale glacier regions. In recent years, unmanned aerial vehicle (UAV) low-altitude photogrammetry technology has been increasingly applied to the reconstruction and dynamic monitoring of glacier surface morphology due to its advantages of high resolution, maneuverability, and lower cost.

[0003] However, the above-mentioned technologies still have problems in the low texture and strong reflection environment of glaciers, such as unstable matching of three-dimensional corresponding points across time phases, large matching errors, and spatial reference drift between point clouds of multiple time phases, which leads to inaccurate calculation of glacier three-dimensional motion vectors. Summary of the Invention

[0004] This invention provides a method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by UAVs, in order to solve the problems of unstable matching of three-dimensional corresponding points across time phases and large matching errors in the low-texture and high-reflection environment of glaciers, as well as the spatial reference drift between point clouds of multiple time phases, which leads to inaccurate calculation of the three-dimensional motion vector of glaciers.

[0005] The present invention provides a method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles, comprising the following steps:

[0006] S1. Acquire multi-angle image data to obtain a multi-temporal raw observation set; preprocess the multi-temporal raw observation set to obtain a preprocessed observation set; determine the three-dimensional point set based on the preprocessed observation set;

[0007] S2. Based on the preprocessed observation set, a geometric-reflection coupling consistency matching algorithm is introduced to register the three-dimensional points in the three-dimensional point set and determine the matching point for the next time phase; based on the three-dimensional points of the current time phase and the matching point of the next time phase, the three-dimensional motion vector of the glacier is obtained.

[0008] Preferably, S1 specifically includes:

[0009] The multi-temporal raw observation set includes image data, focal length, attitude rotation matrix, camera center coordinates, flight altitude, and exposure timestamps.

[0010] Preferably, S1 specifically includes:

[0011] The data in the original observation set of multiple time phases are preprocessed, specifically including image distortion correction, global pose optimization, image enhancement processing, spatial scale unification and data time synchronization, to obtain the preprocessed observation set; the preprocessed observation set includes geometrically distortion-free images, converted focal length, camera center coordinates in a unified spatial coordinate system, attitude rotation matrix in a unified spatial coordinate system and synchronization timestamp.

[0012] Preferably, S1 specifically includes:

[0013] Scale- and rotation-invariant features are extracted and matched from the geometrically undistorted images in the preprocessed observation set to construct a set of corresponding feature points. Based on the pixel coordinates of the corresponding feature points in the geometrically undistorted images, combined with the transformed focal length and the attitude rotation matrix in the unified spatial coordinate system, a sparse 3D point set in the unified coordinate system is obtained. Based on the sparse 3D point set in the unified coordinate system, multi-view stereo dense matching is performed to generate a dense 3D point cloud. Outlier removal and smoothing filtering are performed on the dense 3D point cloud to obtain the 3D point set.

[0014] Preferably, S2 specifically includes:

[0015] In the implementation of the geometry-reflection coupling consistency matching algorithm, for any 3D point in the 3D point set, the local normal vector is determined, and the normalized viewing direction vector from the camera to the 3D point is calculated; based on the local normal vector and the normalized viewing direction vector, the viewpoint cosine is calculated.

[0016] Preferably, S2 specifically includes:

[0017] In the implementation of the geometry-reflection coupling consistency matching algorithm, the projected pixel of any three-dimensional point in the three-dimensional point set is constructed, and the mean gray level and mean gradient magnitude of the neighborhood of the projected pixel are calculated; based on the mean gray level and mean gradient magnitude of the neighborhood of the projected pixel, the normalized radiance is generated.

[0018] Preferably, S2 specifically includes:

[0019] In the implementation of the geometry-reflection coupling consistency matching algorithm, a cross-temporal candidate search is performed on the three-dimensional points in combination with the search radius to generate a candidate set.

[0020] Preferably, S2 specifically includes:

[0021] In the implementation of the geometry-reflection coupling consistency matching algorithm, the normalized radiance of each candidate 3D point in the candidate set in the next time phase is calculated, and the geometry-reflection coupling cost function is constructed by combining the radiative consistency observation confidence weight. By minimizing the geometry-reflection coupling cost function, the matching point of the 3D point in the next time phase is obtained. The radiative consistency observation confidence weight is calculated based on the viewpoint cosine.

[0022] Preferably, S2 specifically includes:

[0023] Based on the 3D points in the current time phase and their matching points in the next time phase, stable control point pairs are selected. Based on the stable control point pairs, the rotation matrix and translation vector from the next time phase to the current time phase are calculated. Based on the rotation matrix and translation vector from the next time phase to the current time phase, the coordinates of the matching points are transformed, the 3D displacement vector is calculated, and converted into the glacier's 3D motion vector.

[0024] The beneficial effects of the technical solution of the present invention are:

[0025] 1. In the data acquisition and preprocessing stage, multi-angle image data were simultaneously acquired using a five-lens oblique photography array. Combined with intrinsic parameter calibration, Brown-Conrady distortion correction, bundle adjustment optimization, and time synchronization processing, a standardized observation dataset under a unified metric spatial coordinate system was established. This processing effectively eliminated lens distortion errors, attitude drift errors, and time synchronization errors, enabling forward intersection calculations of multi-view images within a unified geometric framework. This improved the consistency and spatial accuracy of 3D point coordinate solutions, providing stable basic data input for subsequent multi-temporal point cloud registration and motion analysis. 2. In the time-series matching stage, a geometry-reflection coupled consistency matching algorithm was introduced to jointly evaluate the spatial neighborhood geometric constraints and radiation consistency characteristics of 3D points. This method incorporates the viewpoint cosine formed by the local normal vector and the normalized line-of-sight vector into the weighting system, which automatically reduces the weight of grazing observations and reduces the impact of specular reflection on matching stability from the perspective of observation geometry. At the same time, logarithmic compression and gradient normalization are used to construct a radiation consistency feature for reflection suppression, which can maintain the comparability of radiation descriptions in both high-reflection and low-texture regions. Thus, a stable three-dimensional correspondence can still be established under conditions of low texture on ice surfaces, local snow cover, or strong reflection. Attached Figure Description

[0026] Figure 1 This is a flowchart of the method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles, as described in this invention. Detailed Implementation

[0027] To further illustrate the technical means and effects adopted by the present invention to achieve its intended purpose, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. 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.

[0028] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0029] The following description, in conjunction with the accompanying drawings, details the specific scheme of the glacier three-dimensional motion vector measurement method based on UAV low-altitude oblique photography provided by this invention.

[0030] See attached document Figure 1 The diagram illustrates a flowchart of a method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles, according to an embodiment of the present invention. The method includes the following steps:

[0031] S1. Acquire multi-angle image data to obtain a multi-temporal raw observation set; preprocess the multi-temporal raw observation set to obtain a preprocessed observation set; determine the three-dimensional point set based on the preprocessed observation set;

[0032] The drone flew along a pre-set route by technical personnel and simultaneously acquired multi-angle image data using a five-lens oblique photography array, obtaining the first... The measurement (the first) Phase The raw observation data includes image data, focal length, attitude rotation matrix, camera center coordinates, flight altitude, and exposure timestamp; based on the first... Phase The original observation data were used to construct the original observation set. The time phase Indicates the first The time state corresponding to each measurement; by using Zhang Zhengyou's calibration model to calibrate the camera intrinsic parameters of each lens image, the focal length is obtained, and radial and tangential distortion coefficients are generated simultaneously for subsequent preprocessing.

[0033] The five-lens tilting camera array refers to an aerial photography system consisting of a vertically downward lens and four side-view lenses arranged at fixed tilt angles (such as 45°) along the front, back, left, and right directions, capable of simultaneously acquiring multi-view image data of the same area.

[0034] Furthermore, the data in the aforementioned original observation set undergoes preprocessing, specifically including: image distortion correction: based on the Brown-Conrady distortion model, combined with radial and tangential distortion coefficients, radial and tangential distortion correction is performed on the original pixel coordinates in the image data, converting distorted pixel positions into distortion-free pixel positions to obtain a geometrically corrected image grayscale matrix; global pose optimization: the attitude rotation matrix and camera center coordinates obtained based on the airborne inertial measurement unit are globally optimized using bundle adjustment based on feature point matching, unifying them to the same metric geographic coordinate system; image enhancement processing: histogram equalization is performed on the geometrically corrected image grayscale matrix to obtain a geometrically distortion-free image, enhancing dark details and compressing the brightness of high-reflectivity areas, reducing the impact of strong reflections from the ice surface on subsequent matching; spatial scale unification: based on... After calculating the ground resolution based on the focal length and flight altitude, the pixel standard deviation is converted to the ground metric standard deviation using the error propagation formula, achieving a quantitative mapping between pixel units and ground spatial units. Then, the focal length is converted to metric units to obtain the converted focal length. Data time synchronization: Using the timestamps of each lens exposure time as a reference, linear interpolation is used to precisely align the attitude rotation matrix, unified to the same metric geographic coordinate system, and the camera center coordinates to the corresponding exposure time, completing the time synchronization of image data and attitude data. After the above processing, a preprocessed observation set is obtained, which is a standardized observation data set including geometrically distortion-free images, converted focal lengths, camera center coordinates in a unified spatial coordinate system, attitude rotation matrices, and synchronization timestamps, used for subsequent 3D reconstruction and glacier motion vector calculation. The preprocessing processes all employ mature techniques in this field and will not be elaborated upon here.

[0035] Furthermore, for the geometrically undistorted images in the preprocessed observation set, scale-invariant and rotation-invariant features are extracted using existing scale-invariant feature transformation methods, and matched using existing KD-Tree-based nearest neighbor matching methods to establish a set of corresponding feature points. Subsequently, using the transformed focal length and attitude rotation matrix in the unified spatial coordinate system obtained in the preprocessing stage as known conditions, the pixel coordinates of the corresponding feature points in the geometrically undistorted images are used as observations. Collinearity observation equations between multi-view images are established using existing bundle adjustment methods, and a sparse 3D point set in the unified coordinate system is obtained through forward intersection calculation under multi-view geometric constraints. After obtaining the sparse 3D point set, multi-view stereo dense matching is performed using the geometrically undistorted images as input. The matching process involves searching for pixel-by-pixel correspondences among multi-view images based on epipolar geometric constraints. The sum of absolute differences (SAD) between pixel values ​​within a pre-defined window is then calculated as the matching cost. The pixel with the minimum matching cost is identified as the matching pixel, and the pixel disparity is calculated from the horizontal coordinate difference of the matching pixels. Subsequently, the corresponding feature points of the pixel disparity are converted into corresponding 3D spatial coordinates using forward intersection in photogrammetry, thus generating a dense 3D point cloud. Finally, outliers are removed from the dense 3D point cloud using statistical outlier filtering, followed by bilateral filtering for smoothing. This yields a complete 3D point set in a unified metric coordinate system for subsequent multi-temporal point cloud registration and glacier 3D motion vector calculation.

[0036] S2. Based on the preprocessed observation set, a geometric-reflection coupling consistency matching algorithm is introduced to register the three-dimensional points in the three-dimensional point set and determine the matching point for the next time phase; based on the three-dimensional points of the current time phase and the matching point of the next time phase, the three-dimensional motion vector of the glacier is obtained.

[0037] Furthermore, in order to establish stable spatial correspondences between two consecutive time phases under conditions such as low ice surface texture, specular reflection, and local snow cover, and to provide input for subsequent displacement calculations, a geometry-reflection coupling consistency matching algorithm is introduced. The specific implementation process is as follows:

[0038] Time phase As the first Phase, Phase As the first Phase, No. The corresponding three-dimensional point set is denoted as For any three-dimensional point The local normal vector is determined using existing point cloud local plane fitting methods. Specifically, in a three-dimensional point set, using three-dimensional points... Centered on the point, K neighboring points within its spatial neighborhood are selected using the existing K-nearest neighbor search method to form a local neighborhood point set. Then, the spatial centroid of this neighborhood point set is calculated, and a covariance matrix is ​​constructed based on the distribution of the neighborhood points relative to the spatial centroid of the neighborhood point set. Next, the covariance matrix is ​​decomposed into eigenvalues, and the eigenvector corresponding to the smallest eigenvalue is selected as the normal direction of the local plane where the point is located, i.e., the normal vector.

[0039] Furthermore, for each shot Calculate the distance from the camera to the 3D point Normalized view direction vector:

[0040] ;in, It is the normalized line-of-sight vector; In the phase , No. The camera center coordinates of each lens in a unified spatial coordinate system; Indicates pointing from the center of the camera to a 3D point Spatial vector; This represents the Euclidean norm, used for normalization. Further definition is given in the time phase. , No. The cosine of the perspective of each shot , Local normal vector With the normalized line-of-sight vector The inner product, i.e. Its value range is This is used to reflect the tilt angle; when it is close to 0, reflection anomalies and matching instability are more likely to occur.

[0041] Furthermore, construct three-dimensional points Projected pixels ,in, Representing a three-dimensional point In time phase Lens The projected pixels, Indicates a projection operation; Indicates the first The focal length after lens conversion; Indicates phase , No. The attitude rotation matrix of each lens in a unified spatial coordinate system; simultaneously, a normalized radiance based on logarithmic compression is introduced to suppress specular highlights. Specifically, the neighborhood of the projected pixel is determined according to the specific application scenario; for example, in high-resolution aerial photography, a 5×5 or 7×7 pixel window is often used as the projected pixel. Within the neighborhood of the image, and based on the geometrically distortion-free image, statistically determine the 3D points using existing direct traversal methods. In time phase Lens The grayscale mean of the neighborhood of the projected pixel 3D points are determined using existing gradient operators (such as the Sobel operator or the Scharr operator). In time phase Lens The mean gradient magnitude of the neighborhood of the projected pixel The normalized radiation quantity is as follows:

[0042] ;in, It is a three-dimensional point In time phase Lens Normalized radiance of projected pixels; This reflects the average reflectance benchmark; This indicates that local texture intensity is introduced for energy normalization;

[0043] Furthermore, a cross-temporal candidate search is constructed; in the first phase... Phase In the middle, the three-dimensional point set after coarse registration is first used. In three-dimensional points Construct a candidate set within the spatial neighborhood: Among them, the three-dimensional point set after coarse registration It uses the existing rigid registration method based on 3D point cloud to perform coarse registration of the 3D point set, which is used to eliminate the overall spatial pose difference and provide initial spatial neighborhood constraints for subsequent fine matching; It is the first Three-dimensional points of time phase In the 3D point set after coarse temporal registration A candidate 3D point; search radius Instead of being manually specified parameters, they are adaptively determined by measurable quantities:

[0044] ;in, In a unified spatial coordinate system, the vertical lens in time phase The camera center coordinates; It is the focal length after the vertical lens is converted; , Represents the geometrically undistorted image plane of a vertical lens direction and Standard deviation of pixels in the direction; The scaling factor between the object space scale and the image space scale; It is the overall uncertainty;

[0045] Furthermore, a geometry-reflection coupling cost function is defined, and the optimal match is selected. For each candidate 3D point... Calculate its in the first Phase Corresponding lens projection pixels And calculate the normalized radiation. Geometric-reflection coupling cost function:

[0046] ;in, It is the cost function of geometry-reflection coupling; This represents the set of visible shots (obtained after occlusion culling), with a reference range of [missing information]. ; It is the first A visible camera for a 3D point The confidence weight of radiation consistency observations is determined by the first... The calculation is obtained by dividing the absolute value of the angle cosine of each visible lens by the sum of the absolute values ​​of the angle cosines of all visible lenses. It is a three-dimensional point In time phase Visible lens Normalized radiance of projected pixels; Candidate 3D points In time phase Visible lens Normalized radiance of projected pixels; It is a geometric term used to characterize candidate 3D points. Compared to the first Temporal three-dimensional points The degree of spatial offset; It is a radiometric uniformity term used to characterize a three-dimensional point under observation from multiple visible lenses. With candidate 3D points Whether the reflection suppression radiation characteristics remain consistent across time phases.

[0047] Finally, the first Three-dimensional points of time phase In the Timing matching point ,Right now .

[0048] Based on the Three-dimensional points of time phase and three-dimensional points In the Timing matching point In the two-phase data, point pairs located in stable regions (such as bedrock exposure areas or fixed markers) that are not affected by glacier movement are selected to form stable control point pairs. Based on the stable control point pairs, existing point cloud rigid body registration or coordinate frame unification methods are used to obtain the data from the first phase. Time phase to the first The rotation matrix and translation vector of the phase, and for all Perform the same coordinate transformation to eliminate the reference drift between the two time phases; then for each pair of corresponding points ( After transformation ) Calculate from Pointer transformation The three-dimensional displacement vector is determined by the line connecting two points, and its magnitude is the three-dimensional distance between the two points. The displacement is converted into a velocity vector at the time interval of two aerial photographs, thus obtaining the three-dimensional motion vector of the glacier.

[0049] In summary, a method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles has been completed.

[0050] The order of the embodiments is for illustrative purposes only and does not represent the superiority or inferiority of the embodiments. The processes depicted in the drawings do not necessarily require a specific or sequential order to achieve the desired result. In some embodiments, multitasking and parallel processing are possible or may be advantageous.

[0051] The various embodiments in this specification are described in a progressive manner. The same or similar parts between the various embodiments can be referred to each other. Each embodiment focuses on describing the differences from other embodiments.

[0052] The above embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit it. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention, and should all be included within the protection scope of the present invention.

Claims

1. A method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles (UAVs), characterized in that, Includes the following steps: S1. Acquire multi-angle image data to obtain a multi-temporal raw observation set; preprocess the multi-temporal raw observation set to obtain a preprocessed observation set; Based on the preprocessed observation set, determine the three-dimensional point set; S2. Based on the preprocessed observation set, a geometric-reflection coupling consistency matching algorithm is introduced to register the three-dimensional points in the three-dimensional point set and determine the matching points for the next time phase; the specific implementation process of the geometric-reflection coupling consistency matching algorithm is as follows: For any 3D point in the 3D point set, determine the local normal vector and calculate the normalized viewing direction vector from the camera to the 3D point; based on the local normal vector and the normalized viewing direction vector, calculate the view cosine. A projected pixel of any 3D point in the 3D point set is constructed, and the mean gray level and mean gradient magnitude of the neighborhood of the projected pixel are calculated. Based on the mean gray level and mean gradient magnitude of the neighborhood of the projected pixel, a normalized radiance is generated. Combined with the search radius, a cross-temporal candidate search is performed on the 3D points to generate a candidate set, and the normalized radiance of each candidate 3D point in the candidate set in the next temporal phase is calculated. Combined with the radiative consistency observation confidence weight, a geometry-reflection coupling cost function is constructed. By minimizing the geometry-reflection coupling cost function, the matching point of the 3D point in the next temporal phase is obtained. The radiative consistency observation confidence weight is calculated based on the viewpoint cosine. Based on the matching points of the three-dimensional points in the current time phase and the points in the next time phase, the three-dimensional motion vector of the glacier is obtained.

2. The method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles according to claim 1, characterized in that, S1 specifically includes: The multi-temporal raw observation set includes image data, focal length, attitude rotation matrix, camera center coordinates, flight altitude, and exposure timestamps.

3. The method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles according to claim 2, characterized in that, S1 specifically includes: The data in the original observation set of multiple time phases are preprocessed, specifically including image distortion correction, global pose optimization, image enhancement processing, spatial scale unification and data time synchronization, to obtain the preprocessed observation set; the preprocessed observation set includes geometrically distortion-free images, converted focal length, camera center coordinates in a unified spatial coordinate system, attitude rotation matrix in a unified spatial coordinate system and synchronization timestamp.

4. The method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles according to claim 3, characterized in that, S1 specifically includes: Scale- and rotation-invariant features are extracted and matched from the geometrically undistorted images in the preprocessed observation set to construct a set of corresponding feature points. Based on the pixel coordinates of the corresponding feature points in the geometrically undistorted images, combined with the transformed focal length and the attitude rotation matrix in the unified spatial coordinate system, a sparse 3D point set in the unified coordinate system is obtained. Based on the sparse 3D point set in the unified coordinate system, multi-view stereo dense matching is performed to generate a dense 3D point cloud. Outlier removal and smoothing filtering are performed on the dense 3D point cloud to obtain the 3D point set.

5. The method for measuring the three-dimensional motion vector of glaciers based on low-altitude oblique photography by unmanned aerial vehicles according to claim 1, characterized in that, S2 specifically includes: Based on the 3D points in the current time phase and their matching points in the next time phase, stable control point pairs are selected. Based on the stable control point pairs, the rotation matrix and translation vector from the next time phase to the current time phase are calculated. Based on the rotation matrix and translation vector from the next time phase to the current time phase, the coordinates of the matching points are transformed, the 3D displacement vector is calculated, and converted into the glacier's 3D motion vector.

Images