Three-dimensional reconstruction method of ocean waves based on binocular image data

Abstract

This invention discloses a method for 3D reconstruction of ocean waves based on binocular image data. First, wave image data is acquired through a binocular acquisition platform. Next, a camera model is calibrated based on the data acquired by the binocular system. This involves calculating camera intrinsic parameters using the Zhang Zhengyou calibration method and optimizing camera extrinsic parameters based on image correlation and the spatial distribution of feature points. Rapid stereo matching is then performed on weakly textured wave images. Finally, 3D reconstruction is performed on the wave point cloud. This method enables high-precision self-calibration of the binocular measurement system, completing rapid 3D reconstruction based on binocular image data.

Description

Technical Field

[0001] The invention belongs to a method for three-dimensional reconstruction of ocean waves based on binocular image data, and relates to the field of three-dimensional reconstruction. Background Art

[0002] The ocean is the cradle of life on Earth, covering approximately 70% of the Earth's surface. Furthermore, it is an essential habitat for humans and other organisms. Ocean waves generally refer to wind-generated waves in the ocean, primarily including wind waves, swells, and offshore waves. Wave size varies significantly under varying wind speed, direction, and terrain conditions, typically with periods ranging from a few tenths of a second to tens of seconds, wavelengths from tens of centimeters to hundreds of meters, and heights from a few centimeters to over 20 meters. In unusual terrain, wave heights can reach over 30 meters. Ocean waves are closely intertwined with marine activities, including maritime transportation, offshore oil and other mineral exploration and development, marine fisheries, marine engineering, and maritime military operations. Therefore, understanding ocean wave information, such as wave height, period, and speed, is crucial for preventing marine disasters and for offshore engineering construction. In the study of ocean wave dynamics, efficiently acquiring ocean wave data is a hot topic for in-depth theoretical research and related applied studies. Currently commonly used wave 3D reconstruction instruments such as wave buoys and X-band radars are unable to collect wave data within a certain range for a long time, and therefore cannot meet the data needs for in-depth theoretical research on waves. Summary of the Invention

[0003] In response to the above problems, the present invention provides a method for three-dimensional reconstruction of ocean waves based on binocular image data, comprising four modules: ocean wave image data acquisition (1), visual measurement system calibration (2), weak-texture ocean wave image stereo matching optimization (3), and ocean wave point cloud three-dimensional reconstruction (4). The method can achieve high-precision self-calibration of the binocular measurement system and complete rapid three-dimensional reconstruction based on ocean wave binocular image data. The functions of each module in the technical solution adopted by the present invention to solve the above problems are as follows.

[0004] The processing steps of the visual measurement system calibration module are as follows:

[0005] S21. Obtaining the camera's intrinsic calibration and lens distortion model through Zhang's calibration method;

[0006] S22, using the lens distortion model obtained in step (S21), performing dedistortion processing on the wave image data captured by the camera;

[0007] S23, using a feature extraction algorithm to obtain a feature point set for the dedistorted binocular image processed in step (S22), and obtaining preliminary left and right view matching point pairs through a feature point matching algorithm;

[0008] S24, performing correlation evaluation on the matching point pairs obtained in step (S23), and eliminating error point pairs with low correlation scores;

[0009] S25, optimizing the spatial distribution of the matching points processed in step (S24), first uniformly dividing the left view image into M circular areas of equal size using the following formula,

[0010]

[0011] Where: r is the radius of the preset circular area; is the horizontal coordinate of the center of the circle; The coordinates of the matching points corresponding to the left view are stored in the Kd tree by constructing a Kd tree, and the correlation scores of the matching points are adjusted according to the following formula and stored in the Kd tree.

[0012]

[0013] in: is the correlation score of the matching point. Then, the nearest matching point to the center of each circle in each area is obtained through nearest neighbor search. The distance calculation adopts the formula shown below:

[0014]

[0015] in: is the coordinate of the feature point. This is used to achieve sparse adjustment of the feature points in the area circle. Usually, a pair of wave images cannot satisfy the requirement of obtaining matching points in every area circle. Therefore, by adding matching point information of other wave image pairs under the same camera perspective, the matching points in the area circle are filled;

[0016] S26, calculating the basic matrix and extrinsic parameter estimation based on the matching points processed in step (S25);

[0017] S27. Calculate the reprojection error of the left and right views using the extrinsic parameter matrix of step (S26) and evaluate it.

[0018] The processing steps of the weak texture wave image stereo matching optimization module are as follows:

[0019] S31. Optimizing the search space of the left and right views using the camera extrinsic matrix;

[0020] S32. Optimize the disparity calculation of the SGM algorithm by combining the disparity range estimation;

[0021] S33. GPU acceleration is used to optimize the calculation rate in the process of step (S32).

[0022] The steps of the wave point cloud 3D reconstruction module are as follows:

[0023] S41, calculating the three-dimensional spatial coordinate information of the waves by triangulation;

[0024] S42. Remove noise from the ocean wave point cloud obtained in step (S41) by performing point cloud filtering. Because a general dense matching framework is used, no assumptions are made about specific application scenarios. Therefore, the point cloud is filtered based on the smooth and continuous characteristics of ocean waves. Specifically, the total difference along the camera's Z axis between a point in the reconstructed point cloud and the surrounding four points is calculated using the following formula. All total differences are then sorted, and finally, the 3D coordinate points with large total differences are removed.

[0025]

[0026] Optionally, in the visual measurement system calibration module, the feature point extraction algorithm can be any of SIFT, ORB, and SURF. The feature point matching algorithm can be any of brute force matching, RANSAC matching, and FLANN matching. Error point pairs with low correlation scores are eliminated, meaning those with scores below 0.98. The radius of the preset circular area must be no less than 50 pixels, and the number of circular areas must be no less than 100. The average reprojection error threshold is 0.5 pixels.

[0027] Optionally, in the 3D reconstruction module of the ocean wave point cloud, if the parallax angle between two matching pixels and the 3D coordinate point is less than 20 degrees, it is considered a large error point and is removed. In the point cloud filtering, the coordinate points with large total errors account for 5% of the total number of points.

[0028] The beneficial effects of this application are:

[0029] The present invention provides a method for 3D reconstruction of ocean waves based on binocular image data. First, a binocular acquisition platform acquires ocean wave image data. Next, a camera model is calibrated based on the data acquired by the binocular system. This involves calculating camera intrinsic parameters using the Zhang Zhengyou calibration method and optimizing camera extrinsic parameters based on image correlation and the spatial distribution of feature points. Rapid stereo matching is then performed on weakly textured ocean wave images, and finally, 3D reconstruction is performed on the ocean wave point cloud. This method enables high-precision self-calibration of the binocular measurement system, completing rapid 3D reconstruction based on binocular ocean wave image data.

[0030] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following briefly introduces the drawings required for use in the embodiments. It should be understood that the following drawings only show certain embodiments of the present application and therefore should not be regarded as limiting the scope. For ordinary technicians in this field, other relevant drawings can be obtained based on these drawings without paying any creative work.

[0031] Figure 1 Flowchart of the Ocean Wave binocular reconstruction algorithm.

[0032] Figure 2 Improved local matching algorithm flowchart. DETAILED DESCRIPTION

[0033] The following will be combined with the drawings in the embodiments of this application to clearly and completely describe the technical solutions in the embodiments of this application. Obviously, the embodiments described are only part of the embodiments of this application, not all of the embodiments. Based on the embodiments in this application, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of this application.

[0034] In one embodiment, the present invention provides a method for 3D reconstruction of ocean waves based on binocular image data with high computational accuracy and high computational efficiency. The method comprises six modules: ocean wave image data acquisition, visual measurement system calibration, stereo matching optimization of weak texture ocean wave images, GPU-based stereo matching hardware acceleration, 3D reconstruction of ocean wave point clouds, and construction of a mean sea level coordinate system. Figure 1 As shown, this embodiment provides a method for 3D reconstruction of ocean waves based on binocular image data. The method includes four modules: ocean wave image data acquisition (1), visual measurement system calibration (2), stereo matching optimization of weak texture ocean wave images (3), and 3D reconstruction of ocean wave point clouds (4). The specific implementation steps are as follows.

[0035] The processing steps of the visual measurement system calibration module are as follows:

[0036] S21. Obtaining the camera's intrinsic calibration and lens distortion model through Zhang's calibration method;

[0037] S22, using the lens distortion model obtained in step (S21), performing dedistortion processing on the wave image data captured by the camera;

[0038] S23. For the dedistorted binocular image processed in step (S22), a feature extraction algorithm is used to obtain a feature point set, and a feature point matching algorithm is used to obtain preliminary left and right view matching point pairs. The feature point extraction algorithm may be SIFT, ORB, or SURF, and the feature point matching algorithm may be brute force matching, RANSAC matching, or FLANN matching.

[0039] S24, evaluate the correlation of the matching point pairs obtained in step (S23), and use the digital image correlation algorithm to calculate the correlation of the neighborhood of the feature points of the left and right views. First, the deformation of the neighborhood of the feature point is described by the second-order deformation function shown in the following formula:

[0040]

[0041] in: is the pixel coordinate of the feature point area;

[0042] is the deformation parameter to be determined. Then, the local grayscale correlation function is established by the zero-mean normalized minimum distance square standard, and the correlation score of the matching point pair is calculated by the following formula:

[0043]

[0044] Where: M is the size of the feature point neighborhood; is the deformation increment; is the grayscale mean of the image neighborhood; is the normalization function of the image neighborhood. Then, the matching points are sorted from high to low according to the correlation score using the sorting algorithm, and the matching points with a correlation score less than 0.98 are eliminated;

[0045] S25. Optimize the spatial distribution of the matching points processed in step (S24). First, evenly divide the left view image into M circular areas of equal size using the following formula:

[0046]

[0047] Where: r is the radius of the preset circular area; is the horizontal coordinate of the center of the circle; The coordinates of the matching points corresponding to the left view are stored in the Kd tree by constructing a Kd tree, and the correlation scores of the matching points are adjusted according to the following formula and stored in the Kd tree.

[0048]

[0049] in: is the correlation score of the matching point. Then, the nearest matching point to the center of each circle in each area is obtained through nearest neighbor search. The distance calculation adopts the formula shown below:

[0050]

[0051] in: is the coordinate of the feature point. This is used to achieve sparse adjustment of the feature points in the area circle. Usually, a pair of wave images cannot satisfy the requirement of obtaining matching points in every area circle. Therefore, by adding matching point information of other wave image pairs under the same camera perspective, the matching points in the area circle are filled;

[0052] S26. Calculate the basic matrix and extrinsic parameter estimation based on the matching points processed in step (S25). The essential matrix E of the left and right views can be obtained by using the eight-point method to calculate the following formula:

[0053]

[0054] in: is the camera intrinsic parameter matrix; is the expression of the matching points of the left and right images in normalized coordinates. By performing singular value decomposition on the formula shown below, the camera's external parameters can be obtained.

[0055]

[0056] in: Represents the antisymmetric matrix of the translation vector; R represents the rotation matrix;

[0057] S27. Calculate the reprojection error of the left and right views using the extrinsic parameter matrix from step (S26) and evaluate it. The reprojection error of the left and right views is calculated using the camera extrinsic parameters. If the reprojection error is greater than a set threshold of 0.5 pixels, remove the points with large single-point reprojection error. Then, refill the feature points in the area, re-estimate the camera extrinsic parameters, and recalculate the reprojection error. Repeat these steps until the camera extrinsic parameters with the smallest reprojection error are obtained.

[0058] The processing steps of the weak texture wave image stereo matching optimization module are as follows:

[0059] S31. Optimize the search space for the left and right views using the camera extrinsic matrix. Perform epipolar correction on the binocular image to satisfy the epipolar constraint. The epipolar constraint constrains the search space for matching points to a straight line within the image plane, eliminating many meaningless searches.

[0060] S32, combining the parallax range estimation to optimize the parallax calculation. Figure 2The improved local matching algorithm framework shown calculates the disparity maps of the left and right wave images and optimizes the disparity calculation by adding a disparity range estimation module. The specific steps are: first, using a local sliding window, the left and right images are encoded through the Census transform; then, the disparity range is estimated by combining the position information of the matching points; then, the cost cube of the left and right views is calculated by calculating the Hamming distance; then, the cost aggregation operation is performed and the disparity map obtained according to the "winner takes all" strategy is obtained. The cost aggregation formula is as follows:

[0061]

[0062] in: is the path cost of pixel p along path r; To match the cost value; The disparity discontinuity and smoothness judgment conditions are added. Finally, the left-right consistency detection and median filtering are added to optimize the disparity map.

[0063] S33. GPU acceleration is used to optimize the calculation rate of step (S32). Stereo matching efficiency of high-resolution images is optimized by taking advantage of the number of GPU cores. In terms of optimizing memory access, the concurrent thread bundles are set to multiples of 32. First, the GPU's video memory allocation includes: image pixel data, Census transformation data, cost aggregation data, and disparity data. The cost aggregation data is stored in disparity as the main order to facilitate continuous access when the cost data is aggregated, and the remaining data is set in row-major order. Then, when performing Census encoding, in order to avoid a large number of repeated accesses to global memory, 32x32 threads and shared memory for storing 32x32 integers are allocated to each thread block. Next, for the initial cost calculation, each thread block is allocated a number of threads equal to the disparity range (D) and shared memory three times the size of the disparity range. D is set as a multiple of 32. Each thread loops through D adjacent pixels in the left image. Before the initial cost calculation, the census transform data for these D adjacent pixels and all corresponding pixels in the right image are loaded into shared memory. Computing the initial cost of D pixels in the left image only requires reading the census transform data from the right image 2xD times. Secondly, for path aggregation, a thread block is allocated DxM threads, where M is the number of pixels the thread block is responsible for. Each D thread calculates the path cost for a pixel under all disparity values. The threads loop sequentially and calculate the path cost for all pixels along a path. Furthermore, to minimize the path cost for a pixel under all disparities, a reduction method with low time complexity is used. Finally, for disparity optimization, the aggregation cost accumulation step, optimal disparity value selection step, and sub-pixel optimization step for each path are implemented within the same kernel function. The accumulated path cost can be used directly to calculate the optimal disparity and sub-pixel disparity values, thus avoiding the need to re-access the accumulated cost array in separate kernel functions. To improve access efficiency, shared memory is used to store the accumulated cost values, preventing subsequent optimal disparity and sub-pixel disparity calculations from accessing global memory again. Furthermore, the accumulated cost values ​​do not need to be stored in separate global memory, reducing memory usage. The kernel function's thread allocation strategy is highly similar to the path aggregation algorithm, ensuring that the same thread bundle accesses the continuous cost space of the same pixel. To reduce the number of parallel reduction targets during optimal disparity selection, the same strategy of having one thread responsible for accumulating two disparity values ​​is still used. The difference is that each thread no longer needs to execute a loop, resulting in a higher degree of parallelization.

[0064] The steps of the wave point cloud 3D reconstruction module are as follows:

[0065] S41. Calculate the three-dimensional coordinate information of the waves through triangulation. Since the collinearity equation describes the collinearity of the object point, the optical center, and the image point, it is equivalent to the projection process of a three-dimensional point to a two-dimensional point. At this time, by adding the image of another view, two straight lines are formed in the three-dimensional space. The object point is the "intersection point" of these two straight lines. The "intersection point" here refers to the least squares intersection point, because these two straight lines usually do not intersect at one point due to deviation. The collinearity equation system is established by the following formula

[0066]

[0067] Then rewrite the above formula into a coefficient matrix and perform singular value decomposition. The eigenvector corresponding to the minimum eigenvalue is the coordinate of the 3D point to be determined. Finally, calculate the parallax angle between the two matching pixels and the 3D coordinate point. If the angle is less than 20 degrees, it is considered a large error point and is eliminated.

[0068] S42. Remove noise from the ocean wave point cloud obtained in step (S41) by filtering the point cloud. Because a general dense matching framework is used, no assumptions are made about specific application scenarios. Therefore, the point cloud is filtered based on the smooth and continuous characteristics of the ocean waves. Specifically, the total difference along the camera's Z axis between a point in the reconstructed point cloud and the surrounding four points is calculated using the following formula. All total differences are then sorted, and finally, the 5% of 3D coordinate points with the largest total differences are removed.

[0069] .

Claims

1. A method for 3D reconstruction of ocean waves based on binocular image data, comprising four modules: ocean wave image data acquisition (1), visual measurement system calibration (2), stereo matching optimization of weak-texture ocean wave images (3), and 3D reconstruction of ocean wave point clouds (4); The processing steps of the visual measurement system calibration module are as follows: S21. Obtaining the camera's intrinsic calibration and lens distortion model through Zhang's calibration method; S22, using the lens distortion model obtained in step (S21), performing dedistortion processing on the wave image data captured by the camera; S23, using a feature extraction algorithm to obtain a feature point set for the dedistorted binocular image processed in step (S22), and obtaining preliminary left and right view matching point pairs through a feature point matching algorithm; S24, performing correlation evaluation on the matching point pairs obtained in step (S23), and eliminating error point pairs with low correlation scores; S25, optimizing the spatial distribution of the matching points processed in step (S24), first uniformly dividing the left view image into M circular areas of equal size using the following formula; ； in: r is the radius of the preset circular area; is the horizontal coordinate of the center of the circle; The coordinates of the matching points corresponding to the left view are stored in the Kd tree by constructing a Kd tree, and the correlation scores of the matching points are adjusted according to the following formula and stored in the Kd tree; ； in: is the correlation score of the matching points, and then the nearest matching point to the center of each area circle is obtained through nearest neighbor search. The distance calculation uses the following formula; ； in: The coordinates of the feature points are used to achieve sparse adjustment of the feature points in the area circle. Usually, a pair of ocean wave images cannot satisfy the requirement of obtaining matching points in every area circle. Therefore, by adding the matching point information of other ocean wave image pairs under the same camera perspective, the matching points in the area circle are filled. S26, calculating the basic matrix and extrinsic parameter estimation based on the matching points processed in step (S25); S27, calculating the reprojection error of the left and right views using the extrinsic parameter matrix of step (S26), and evaluating it; The processing steps of the weak texture wave image stereo matching optimization module are as follows: S31. Optimizing the search space of the left and right views using the camera extrinsic matrix; S32. Optimize the disparity calculation of the SGM algorithm by combining the disparity range estimation; S33, using GPU acceleration to optimize the calculation rate in step (S32); The steps of the wave point cloud 3D reconstruction module are as follows: S41, calculating the three-dimensional spatial coordinate information of the waves by triangulation; S42. The ocean wave point cloud obtained in step (S41) is subjected to point cloud filtering to remove noise. Since a general dense matching framework is used and no assumptions are made about specific application scenarios, the point cloud is filtered using the smooth and continuous characteristics of the ocean wave. The specific implementation method is as follows: first, the total difference between a point in the reconstructed point cloud and the surrounding four points along the camera Z axis is calculated according to the following formula, then all the total differences are sorted, and finally the three-dimensional coordinate points with large total differences are eliminated; 。 2. The method for 3D reconstruction of ocean waves based on binocular image data according to claim 1, characterized in that: In the visual measurement system calibration module, the feature point extraction algorithm can be any one of SIFT, ORB, and SURF; the feature point matching algorithm can be any one of brute force matching, RANSAC matching, and FLANN matching; error point pairs with low correlation scores are eliminated, referring to matching points with scores below 0.98; the radius of the preset circular area is not less than 50 pixels, the number of circular areas is not less than 100, and the average reprojection error threshold is 0.5 pixels.

3. The method for 3D reconstruction of ocean waves based on binocular image data according to claim 1, characterized in that: In the wave point cloud 3D reconstruction module, if the parallax angle between two matching pixel points and the 3D coordinate point is less than 20 degrees, it is considered a large error point and is eliminated. In the point cloud filtering, the coordinate points with large total errors account for 5% of the total number of points.

Images