A method for observing ocean waves based on binocular cameras

Abstract

This invention proposes a method for observing ocean waves based on a binocular camera, belonging to the field of image processing technology. The method includes: S1: The binocular camera captures ocean wave images, generating ocean wave point cloud data at the current time t; S2: The ocean wave point cloud data in S1 is converted to ocean wave point cloud data in a geodetic coordinate system; S3: The ocean wave data in S2 undergoes quality screening and optimization, including: selecting the region of interest, quality screening, converting noisy data points to null values, null value filling, and Gaussian filtering; the quality screening classifies the ocean wave data into three levels: excellent, good, and poor; S4: Based on the excellent and good quality ocean wave point cloud data, the wave height, wavelength, and period at that time are calculated; S5: S1-S4 are repeated to obtain multiple frames of ocean wave point cloud data at different times, and the effective wave height, wavelength, period, and wave spectrum are calculated. This method can fill in null values, making the overall data more reasonable and the obtained elevation data more accurate.

Description

Technical Field

[0001] This invention belongs to the field of image processing technology, specifically relating to a method for observing ocean waves based on a binocular camera. Background Technology

[0002] Ocean waves are a common marine phenomenon, and the observation of ocean wave data is of great significance for marine disaster prevention and mitigation.

[0003] Wave-measuring buoys are small buoy systems used to measure ocean waves, employing in-situ measurement technology. They can remotely measure parameters such as wave period, direction, height, directional spectrum, and power spectrum. Common methods include GPS wave measurement and gravity-based wave measurement. However, this method has limitations due to power supply constraints and the difficulty of deployment in distant waters.

[0004] Radar remote sensing technology transmits radio waves to the sea surface using radar, receives the echoes, and uses digital images to simulate sea surface conditions and calculate wave parameters. This method has limitations: insufficient image data simulation, a wide sampling range, and limited resolution. It is generally used for medium- to high-scale waves and has limited applicability to small-scale waves.

[0005] Satellite wave measurement technology, similar to radar, remotely monitors ocean waves using satellite equipment. Two existing satellite remote sensing devices are used to observe sea surface wind and wave information: radar altimeters measure the significant wave height; and synthetic aperture radar measures and retrieves the wave direction spectrum and significant wave height. The limitations of this method include the impact of long-distance equipment on accuracy, significant weather influence, non-continuous timing sampling, and resolution limitations.

[0006] In summary, while there are many current methods for measuring ocean waves, each has limited applicability. Deploying buoys in the open ocean is difficult, and radar and satellites are not suitable for small-scale wave measurements. Utilizing binocular vision technology to observe ocean waves can effectively solve this problem; however, existing binocular vision technologies are still in the research stage and are inefficient, unable to perform real-time observation and parameter inversion. To address this issue, there is an urgent need to improve existing technologies and overcome the shortcomings of traditional observation methods and current binocular vision observation techniques. Summary of the Invention

[0007] To address the aforementioned technical problems, this invention proposes a wave observation method based on a binocular camera.

[0008] To achieve the above objectives, the present invention employs the following technical content: a wave observation method based on a binocular camera, comprising the following steps:

[0009] S1: The binocular camera captures images of ocean waves and generates wave point cloud data at the current time t [X]. ct ,Y ct Z ct ]; Wave point cloud data [X ct ,Yct Z ct Established using the left camera coordinate system;

[0010] S2: Based on the IMU inertial measurement unit, the wave point cloud data in the left camera coordinate system [X ct ,Y ct Z ct Convert ocean wave point cloud data to geodetic coordinate system [X] w ,Y w Z w ];

[0011] S3: Wave data in the geodetic coordinate system of S2 [X] w ,Y w Z w Quality screening and optimization are performed, including: selecting the region of interest, quality screening, converting noisy data points to null values, null padding, and Gaussian filtering.

[0012] Quality screening will use wave data [X] w ,Y w Z w The quality of wave data is categorized as: Excellent, Good, and Poor; poor quality wave data is discarded. w ,Y w Z w Only high-quality and good-quality wave data are retained. w ,Y w Z w Proceed to the next step;

[0013] S4: Based on high-quality and good-quality wave point cloud data [X] w ,Y w Z w ] Calculate the wave height, wave wavelength, and wave period at that moment;

[0014] S5: Repeat S1-S4, with the time progressing in each cycle, to obtain multiple frames of wave point cloud data at different times, and calculate the effective wave height, effective wave wavelength, effective wave period, and wave spectrum based on the multiple frames of wave point cloud data at different times.

[0015] Furthermore, step S3, which involves selecting the region of interest and performing quality screening, includes the following steps:

[0016] S3.1: Select the region of interest, using the following formula:

[0017] X1 = int(x1*Height)

[0018] X2 = int(x2 * Height)

[0019] Y1 = int(y1*Width)

[0020] Y2 = int(y2*Width)

[0021] Where x1, x2, y1, y2 are the row and column ratios of the region of interest, Height and Width are the resolutions of the images captured by the camera, i.e., the number of pixels in rows and columns, and X1, X2, Y1, and Y2 are the index values ​​of the selected wave point cloud matrix data, i.e., the wave point cloud data of interest is extracted according to [X1:X2,Y1:Y2].

[0022] S3.2: Count the number of noise and missing values ​​in the region of interest; set quality screening rules and perform screening; based on the screening results, determine the quality of the wave point cloud data [X]. w ,Y w Z w To assess the quality of the wave point cloud data, discard the poor-quality data and retain the non-poor-quality data for subsequent steps.

[0023] S3.2 specifically includes the following steps:

[0024] S3.2.1: Obtain wave elevation data for each data point in the region of interest;

[0025] Suppose there are S data points in the region of interest. These S data points include S1 data points that contain elevation data and S2 data points that do not contain elevation data; S = S1 + S2; the data points that do not contain elevation data are called "null values".

[0026] The wave elevation data of the i-th data point is denoted as Z. i Let i ∈ [1, S1]. Calculate the noise value of the i-th data point using the following formula:

[0027]

[0028] In the formula, μ and σ are the mean and standard deviation of the wave elevation data of S1 "data points", respectively; A i Let A be the noise value of the i-th data point; if |A i If |>3, then the i-th data point is denoted as a "noise data point"; traverse the S data points to obtain the number of "noise data points" among the S1 "data points", denoted as M, where M∈[1,S1];

[0029] While traversing the S data points, we also obtain the number of null values ​​S2 mentioned above;

[0030] S3.2.2: Calculate the ratio of the total number of blank values ​​and noise points to the total number of data points in the region of interest. Let this ratio be P, and the formula is:

[0031] P = (M + S²) / S

[0032] S3.2.3: Set the quality screening rule Judge to judge the wave point cloud data [X] w ,Y w Z w The quality of ];

[0033] The formula for judgment is:

[0034]

[0035] When Judge is -1, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w The quality is poor; this frame of wave point cloud data should be removed.

[0036] When Judge is 1, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w Excellent quality;

[0037] When Judge is 0, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w [Good quality; retain the high-quality and good-quality wave point cloud data for subsequent calculations.]

[0038] Furthermore, step S3 also includes the following steps:

[0039] S3.3: If the wave point cloud data [X] w ,Y w Z w The quality assessment is rated as excellent or good, and the wave point cloud data is optimized. w ,Y w Z w Noisy data points and null values;

[0040] S3.3 specifically includes the following steps:

[0041] S3.3.1: Use row and split dilution factors to dilute the wave point cloud data [X] w ,Y w Z w [Dilution, the diluted wave point cloud data [X] w ,Y w Z w The "noise data points" in the image are converted to null values.

[0042] Set the diluted [X] w ,Y wZ w There are r1 null values ​​and r2 non-null values ​​in the data; r1 includes the null values ​​in the original wave point cloud data and the total number of null values ​​converted from "noise data points";

[0043] Let the coordinates of the b-th null value be (x... ib ,y jb );

[0044] S3.3.2: Using null coordinates (x ib ,y jb Centered on a point (R), the search range is defined by its side length R. If a valid value exists within the search range, the missing value needs to be filled in; otherwise, the missing value is not filled in. The valid value is the data point that "contains elevation data" mentioned above.

[0045] The filling process in step S3.3.2 includes the following steps:

[0046] S3.3.2.1: Calculate the Euclidean distance from the null value to each non-null value;

[0047] The formula is:

[0048]

[0049] In the formula, (x f ,y f () represents the coordinates of the f-th non-empty value; f∈[1,r2];

[0050] S3.3.2.2: Calculate the weight of the null value;

[0051] The formula is:

[0052]

[0053] In the formula, p is the set weight index;

[0054] S3.3.2.3: Filler value at empty spaces, the formula is:

[0055]

[0056] In the formula, Z valid_values This is the selected valid non-null Z-coordinate data. When calculating the X and Y axis data, simply replace it with the corresponding axis data.

[0057] S3.3.2.4: Gaussian filtering;

[0058] Gaussian function:

[0059]

[0060] Construct a kernel matrix of size (2k+1)×(2k+1) and calculate the weight of each matrix position:

[0061]

[0062] Where i and j are the coordinates of each matrix data relative to the kernel center, i = -k, -k+1, ..., k-1, k, j = -k, -k+1, ..., k-1, k;

[0063] To preserve information from subsequent elevation data, the Gaussian kernel is normalized to ensure that the sum of the filter kernels is 1.

[0064]

[0065] The original wave elevation matrix Z w Expand the matrix by k pixels symmetrically around the perimeter to generate an expanded matrix Z. wk Original wave elevation matrix Z w That is, the wave elevation data after the gaps are filled;

[0066] For the wave elevation data matrix Z w Each raw data point Z w Using (x, y) as the kernel center, calculate the smoothed value:

[0067]

[0068] The Gaussian kernel is applied to the corresponding positions of the extended matrix, and the results are multiplied point by point and then summed. Finally, the boundary portion of the extension is removed, and the smoothed result with the same size as the original data is retained for subsequent operations.

[0069] Gaussian filtering is also missing here.

[0070] Furthermore, calculate the wave point cloud data [X] at the current time t. w ,Y w Z w The wave height, wavelength, and period of ocean waves include the following steps:

[0071] S4.1: Transfer wave point cloud data [X] w ,Y w Z w The ocean waves in the image are considered as a matrix (M) with a uniform X-axis direction. t N t ,3); then matrix (M) t N t Further compression yields (M) t The vector G of ,1) t The zero-point crossing identification method is used to extract vector G. t The maximum and minimum elevation values ​​are denoted as H.max and H min；

[0072] S4.1: Calculate the wave height at the current time t;

[0073] The formula is:

[0074] H t =H max -H min ;

[0075] In the formula, H t Let t be the wave height at the current moment;

[0076] S4.2: Calculate the wavelength of the ocean wave at the current time t;

[0077] Get H max and H min Then, based on vector G t The index in the file allows you to query the corresponding H. max and H min The data of the ocean waves on the X and Y axes; for the uncompressed data matrix, there are M data points along the Y direction. t There are N data points, with corresponding values ​​along the X direction. t Quantitative data, assuming H max The corresponding data on the X and Y axes are [X max_1 X max_2 …、X max_Nt ]、[Y max_1 Y max_2 …、Y max_Nt ], H min The corresponding data on the X and Y axes are [X min_1 X min_2 …、X min_Nt ]、[Y min_1 Y min_2 …、Y min_Nt Therefore, the wavelength of the ocean wave at the current time t can be calculated; the formula is:

[0078]

[0079] Calculate ocean wave point cloud data [X] w ,Y w Z w The average wavelength of ocean waves is considered as the wave wavelength at the current time t; the formula is:

[0080]

[0081] In the formula, L t Let t be the wavelength of the ocean wave at the current time.

[0082] S4.3: Calculate the period of the ocean wave at the current time t;

[0083] There are two scenarios:

[0084] Case 1: For deep water waves with a depth greater than or equal to half the wavelength;

[0085] The formula for calculating the wave period is:

[0086]

[0087] In the formula, c represents the wave speed; g is the acceleration due to gravity; L t T is the wavelength of the ocean wave at the current time t; t It is the wave period at the current time t;

[0088] Scenario 2: For shallow water waves where the water depth is less than half the wavelength;

[0089] The formula for calculating the wave period is:

[0090]

[0091] In the formula, c represents the wave speed; g is the acceleration due to gravity; L t T is the wavelength of the ocean wave at the current time t; h is the water depth; T t It is the wave period at the current time t.

[0092] Furthermore, calculating the significant wave height, significant wavelength, significant period, and wave spectrum includes the following steps:

[0093] S5.1: Continue to cycle S1-S4 multiple times. The time will advance with each cycle. The binocular camera continuously captures images. The computer connected to the binocular camera continuously generates wave point cloud data at different times. The wave height, wave wavelength and wave period data of the wave point cloud data at different times are continuously obtained through calculation.

[0094] Suppose a total of n sets of wave data are generated. Sort the wave heights of the n sets of data in descending order, and then add the sorted data into set A in the same order. Set A is represented as:

[0095] A={{H (1) L (1) T (1)},{H (2) L (2) T (2)},…,{H (n) L (n) T (n)}}, where H (1) ≥H (2) ≥…≥H (n) ;

[0096] set up Select {{H (1) L (1) T (1)},{H (2) L (2) T (2)},…,{H (K) L (K) T (K)}};

[0097] S5.2: Calculate the significant wave height, significant wavelength, and significant period of ocean waves;

[0098] Effective wave height H s The formula is:

[0099]

[0100] Effective wavelength L s The formula is:

[0101]

[0102] Effective period T s The formula is:

[0103]

[0104] Further, in step S5, calculating the wave spectrum includes the following steps:

[0105] S5.3.1: Set the loop S1-S4 multiple times to generate u-frame wave point cloud data;

[0106] S5.3.2: In the u-frame wave point cloud data, select a fixed point whose position remains unchanged in each frame of wave point cloud data. Obtain the elevation of the same location in each frame of wave point cloud data, denoted as H. x1 H x2 ... H xu The generation time for each frame of wave point cloud data mentioned above is set to: t b1 t b2 ... t bu ;but:

[0107] t b1 to t bu The formula for the wave spectrum of a time series is:

[0108]

[0109] In the formula, ω is the angular frequency of the ocean wave. P(ω) is the power of the ocean wave; τ is the time offset.

[0110] The formula for τ is:

[0111]

[0112] By adopting the above technical solution, the present invention achieves the following beneficial effects: The present invention uses a binocular camera to collect wave elevation data, which is a continuous, non-contact, and high-density observation method. It can continuously observe waves without affecting their movement and capture more detailed wave information. Simultaneously, the use of an IMU (Inertial Measurement Unit) ensures that the data captured remains reliable even when the camera moves or slightly shifts, making it possible to mount the binocular camera observation platform on a moving platform. Furthermore, addressing issues such as missing data and noisy data acquired by the binocular camera, quality screening and optimization algorithms are used to ensure and improve data quality, resulting in higher accuracy in subsequent calculations of wave spectra and statistical parameters. This also enables the overall algorithm to achieve fully automatic and highly efficient data observation and parameter inversion. Attached Figure Description

[0113] Figure 1 This is a flowchart of the method;

[0114] Figure 2 It is a compressed wave matrix and index diagram. Detailed Implementation

[0115] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, 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.

[0116] like Figure 1 As shown, a method for observing ocean waves based on a binocular camera specifically includes the following steps:

[0117] S1: The binocular camera captures images of ocean waves and generates wave point cloud data;

[0118] Use a binocular camera aimed at the wave area, avoiding direct sunlight, and conduct observations in well-lit environments.

[0119] After the binocular camera is turned on, the left and right cameras of the binocular camera simultaneously capture images of the waves. The computer generates wave point cloud data through existing methods such as image acquisition, image filtering, automatic correction, stereo matching, and 3D reconstruction.

[0120] S2: Based on the IMU inertial measurement unit, the generated wave point cloud data in the left camera coordinate system is converted to wave point cloud data in the geodetic coordinate system.

[0121] The wave point cloud data generated in S1 is wave point cloud data in the coordinate system of a stereo camera. Different stereo cameras generate different wave point cloud data for the same wave because their coordinate systems are different. Therefore, it is necessary to convert the wave point cloud data in the stereo camera coordinate system into wave point cloud data in a unified geodetic coordinate system to facilitate measurement and evaluation under the same standard.

[0122] S2 specifically includes the following steps:

[0123] S2.1: Determine the rotation matrix R t .

[0124] In a binocular camera, the coordinate system of the left camera is based on its optical center as the origin, with the optical axis as the Z-axis, the direction pointing to the optical center of the right camera as the X-axis, and the Y-axis determined according to the right-hand rule.

[0125] The IMU (Inertial Measurement Unit) built into the left camera initializes a geodetic coordinate system with the left camera's optical center as the origin and the Z-axis parallel to and opposite to gravity. The IMU calculates the deflection angles of the left camera coordinate system to the geodetic coordinate system along the Z, X, and Y axes, denoted as θ, ... ω. The IMU (Inertial Measurement Unit) can also calculate the translation vector. The left camera rotates around the Z, X, and Y axes by angles θ. ω will yield the corresponding rotation matrices R1, R2, and R3; the calculation formula is:

[0126]

[0127] Thus, the rotation matrix R of the left camera is obtained.

[0128] The formula is:

[0129] R = R1R2R3

[0130] S2.2: The wave point cloud data generated in S1 is converted to the geodetic coordinate system.

[0131] Set the current time t, at which a frame of wave point cloud data in the left camera coordinate system is generated [X]. ct ,Y ct Z ct ], convert it to point cloud data in the geodetic coordinate system [X w ,Y w Z w ].

[0132] The formula is:

[0133]

[0134] In equation (1), [X ct ,Y ct Z ct [X] represents the wave point cloud data in the left camera coordinate system; w ,Y w Z w [This represents the wave point cloud data for the corresponding frame in the geodetic coordinate system, where Z...] w This indicates the actual elevation data of the ocean waves.

[0135] S3: The wave point cloud data at time t in the geodetic coordinate system of S2 [X] w ,Y w Z w The quality is then screened and optimized. This process includes: selecting the region of interest, quality screening, missing value imputation, and Gaussian filtering.

[0136] S3 specifically includes the following steps:

[0137] S3.1: Obtain the wave point cloud data [X] for the current frame at time t. w ,Y w Z w Select the region of interest.

[0138] X1 = int(x1*Height)

[0139] X2 = int(x2 * Height)

[0140] Y1 = int(y1*Width)

[0141] Y2 = int(y2*Width)

[0142] Where x1, x2, y1, y2 (<1) are the row and column ratios of the region of interest, Height and Width are the resolutions of the images captured by the camera, i.e., the number of pixels in rows and columns, and X1, X2, Y1, and Y2 are the index values ​​of the selected wave point cloud matrix data, i.e., the wave point cloud data of the region of interest that is cropped into a rectangle according to [X1:X2,Y1:Y2].

[0143] Selecting the region of interest can reduce the amount of data required for subsequent operations and facilitate quick identification of the wave point cloud data [X]. w ,Y w Z w The quality of the region of interest is determined by human selection.

[0144] S3.2: Count the number of noise and missing values ​​in the region of interest; set quality filtering rules and perform filtering; based on the filtering results, determine the quality of the wave point cloud data in this frame [X]. w ,Y w Z w The quality of ].

[0145] S3.2 specifically includes the following steps:

[0146] S3.2.1: Obtain wave elevation data for each “data point” in the region of interest.

[0147] Suppose there are S data points in the region of interest. These S data points include S1 data points that contain elevation data and S2 data points that do not contain elevation data. S = S1 + S2. The data points that do not contain elevation data are considered "null values".

[0148] The wave elevation data of the i-th data point (i∈[1,S1]) is denoted as Z. i The noise value of the i-th data point is calculated using the following formula:

[0149]

[0150] In equation (2), μ and σ are the mean and standard deviation of wave elevation data for S1 “data points”, respectively, and A i Let |A| be the noise value of the i-th data point. i If |>3, then the i-th data point is denoted as a "noise data point"; by comparing the thresholds as described above, we traverse the S data points to obtain the number of noise data points in the S1 "data points", which is set to M, where M∈[1,S1].

[0151] While traversing the S data points, we also obtain the number of null values ​​S2 mentioned above.

[0152] S3.2.2: Calculate the ratio of the total number of blank values ​​and noise points to the total number of data points in the region of interest. This ratio is denoted as P, and the formula is:

[0153] P = (M + S²) / S

[0154] S3.2.3: Set the quality screening rule Judge to judge the wave point cloud data [X] w ,Y w Z w The quality of ].

[0155] The formula for the judgment rule is:

[0156]

[0157] When Judge is -1, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w The quality is poor, and this frame of wave point cloud data should be discarded; that is, if the wave point cloud data of this frame at the current time t is [X], then it should be discarded. w ,Y w Z w If the quality is poor, it will not be used for subsequent calculations of effective wave height, effective wavelength, and effective period;

[0158] When Judge is 1, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w Excellent quality;

[0159] When Judge is 0, it indicates that the wave point cloud data of this frame [X] is... w ,Y w Z w Good quality.

[0160] This means that if the wave point cloud data of the current frame at time t [X] w ,Y w Z w The values ​​of "excellent" and "good" will be used in subsequent calculations of effective wave height, effective wavelength, and effective period.

[0161] If the wave point cloud data [X] w ,Y w Z w The quality assessment is rated as excellent or good, and the wave point cloud data is optimized. w ,Y w Z w Noise data points and null values ​​in ] .

[0162] S3.3 specifically includes the following steps:

[0163] S3.3.1: Use row and column dilution factors to dilute the wave point cloud data [X] w ,Y w Z w The diluted wave point cloud data [X] is then used to dilute the data. w ,Y w Z w The "noise data points" in the image are converted to null values.

[0164] Set the diluted wave point cloud data [X] w ,Y w Z wThere are r1 null values ​​and r2 non-null values. Here, r1 includes the null values ​​in the original wave point cloud data and the total number of null values ​​converted from "noise data points".

[0165] Let the coordinates of the b-th null value be (x... ib ,y jb b∈[1,r1).

[0166] Row and column dilution factors are existing techniques in image processing. For example, if both row and column dilution factors are 0.5, it means that the wave point cloud data [X] will be diluted to a minimum. w ,Y w Z w Discarding half of the data points in both rows and columns can reduce the overall amount of data.

[0167] S3.3.2: Using null coordinates (x ib ,y jb Centered on a point (R), the search range is defined by its side length R. If a valid value exists within the search range, the missing value needs to be filled in; otherwise, the missing value is not filled in. The valid value is the data point that "contains elevation data" mentioned above.

[0168] The process of filling in null values ​​in S3.3.2 includes the following steps:

[0169] S3.3.2.1: Calculate the Euclidean distance from the null value to each non-null value.

[0170] The formula is:

[0171]

[0172] In the formula, (x f ,y f ) represents the coordinates of the f-th non-empty value; f∈[1,r2].

[0173] S3.3.2.2: Calculate the weight of the null value to each non-null value.

[0174]

[0175] In the formula, p is the set weight index.

[0176] S3.3.2.3: Filler value at empty spaces, the formula is:

[0177]

[0178] In the formula, Z valid_values This is the selected valid non-null Z-coordinate data. When calculating the X and Y axis data, simply replace it with the corresponding axis data.

[0179] S3.4: Use a Gaussian filter to filter the wave point cloud data in S5 that has been filled with missing values.

[0180] The wave point cloud data in S3.4 is diluted by rows and columns as described in the aforementioned S3.3.

[0181] The wave elevation data after gap filling is smoothed using Gaussian smoothing filtering; the Gaussian function is defined in two-dimensional space as:

[0182]

[0183] Where σ is the standard deviation of the Gaussian distribution; x and y are the offsets of the filter center.

[0184] Therefore,

[0185] Construct a kernel matrix of size (2k+1)×(2k+1) and calculate the weight of each matrix position:

[0186]

[0187] Where i and j are the coordinates of each matrix data relative to the kernel center, i = -k, -k+1, ..., k-1, k, j = -k, -k+1, ..., k-1, k.

[0188] To preserve information from subsequent elevation data, the Gaussian kernel is normalized to ensure that the sum of the filter kernels is 1.

[0189]

[0190] To prevent loss of edge data, the original wave elevation matrix Z... w Expand the matrix by k pixels symmetrically around the perimeter to generate an expanded matrix Z. wk Original wave elevation matrix Z w That is, the wave elevation data after the gaps are filled.

[0191] For the wave elevation data matrix Z w Each raw data point Z w Using (x, y) as the kernel center, calculate the smoothed value:

[0192]

[0193] The Gaussian kernel is applied to the corresponding positions of the extended matrix, and the results are multiplied point by point and then summed. Finally, the boundary portions of the extension are removed, leaving a smoothed result with the same size as the original data for subsequent operations.

[0194] S4: Due to the noise and gaps in the data acquired by the binocular camera, the feasibility and accuracy of subsequent parameter calculations will be affected. Therefore, the wave point cloud data [X] obtained from S3 after quality screening and optimization is used as the basis.w ,Y w Z w Calculate the wave height H at the current time t. t Wave wavelength L t And wave cycle T t。

[0195] S4 specifically includes the following steps:

[0196] S4.1: Since ocean waves propagate unidirectionally along their direction of propagation, the characteristics of ocean waves in a certain direction do not change with the position in that direction. That is to say, in ocean wave point cloud data, ocean waves propagate along the Y-axis, their elevation changes along the Z-axis, and ocean waves are considered uniform along the X-axis. Based on this, we can assume that the ocean wave point cloud data [X] in S3, after quality screening and optimization, is... w ,Y w Z w The wave in the image is a (M) t N t The matrix M is 3). t N represents the length of the ocean wave along the Y direction at the current time t. t This represents the length of the ocean wave along the X-axis at the current time t, such as... Figure 2 As shown.

[0197] Then the matrix (M) t N t Further compression yields (M) t The vector G of ,1) t G t This refers to the curve that changes along the Y-axis.

[0198] S4.2: The zero-point crossing method is used to extract vector G. t The maximum and minimum elevation values ​​are denoted as H. max and H min Therefore, it is possible to calculate the wave point cloud data [X] at the current time t. w ,Y w Z w The wave height of the waves in ] is H t =H max -H min .

[0199] S4.3: Calculate the wavelength of the ocean waves at the current time t.

[0200] Get H max and H min Then, based on vector G t The index in the file allows you to query the corresponding H. max and H minThe data for ocean waves on the X and Y axes. For example... Figure 2 As shown, for an uncompressed data matrix, there are M data points along the Y direction. t There are N data points, with corresponding values ​​along the X direction. t Quantitative data, assuming H max The corresponding data on the X and Y axes are [X max_1 X max_2 …、X max_Nt ]、[Y max_1 Y max_2 …、Y max_Nt ], H min The corresponding data on the X and Y axes are [X min_1 X min_2 …、X min_Nt ]、[Y min_1 Y min_2 …、Y min_Nt This allows us to calculate the wavelength of the ocean wave at the current time t. The formula is:

[0201]

[0202] Calculate ocean wave point cloud data [X] w ,Y w Z w The average wavelength of ocean waves is considered as the wave wavelength at the current time t; the formula is:

[0203]

[0204] S4.4: Calculate the period of the ocean waves at the current time t.

[0205] We will discuss two scenarios:

[0206] Case 1: For deep water waves, i.e., when the water depth is greater than or equal to half the wavelength.

[0207] The formula for calculating the wave period is:

[0208]

[0209] In the formula, c represents the wave speed; g is the acceleration due to gravity (9.81 m / s²). 2 );L t T is the wavelength of the ocean wave at the current time t; t It is the wave period at the current time t.

[0210] Case 2: For shallow water waves, that is, when the water depth is less than half the wavelength.

[0211] The formula for calculating the wave period is:

[0212]

[0213] In the formula, c represents the wave speed; g is the acceleration due to gravity (9.81 m / s²). 2 );L t H is the wave wavelength at the current time t; H is the water depth.

[0214] S5: Calculate the effective wave height, effective wavelength, and effective period; calculate the wave spectrum.

[0215] In the aforementioned S1-S4, the wave height H in the point cloud data at the current time t is calculated. t Wave wavelength L t And wave cycle T t If only the above data H is used... t L t and T t Representing the effective wave height, effective wavelength, and effective period can lead to significant errors; therefore, the traditional "one-third wave method" is needed to calculate the effective wave height, effective wavelength, and effective period.

[0216] Specifically, the following steps are included:

[0217] S5.1: Continue to cycle through S1-S4 multiple times. Each cycle will advance the time. The binocular camera will continuously capture images. The computer connected to the binocular camera will continuously generate wave point cloud data at different times. The wave height, wave wavelength and wave period of the wave point cloud data at different times will be continuously obtained through calculation.

[0218] For example: loop through operations S1-S4 once to obtain the wave height H at the next time step t+1. t+1 Wave wavelength L t+1 And wave cycle T t+1 ;

[0219] Repeat the S1-S4 operation once more to obtain the wave height H at time t+2. t+2 Wave wavelength L t+2 And wave cycle T t+2 ;

[0220] And so on...

[0221] Assume that the operations S1-S4 were looped n-1 times; combined with the original wave data at time t, there are n sets of data; each set of data includes wave height, wave wavelength and wave period.

[0222] Sort the wave heights of the n sets of data from largest to smallest, and then add the sorted data into set A in the same order.

[0223] The set A is represented as:

[0224] A={{H (1) L (1) T (1)},{H (2) L (2) T (2)},…,{H (n) L (n) T (n)}}, where H (1) ≥H (2) ≥…≥H (n)

[0225] Select the first third (rounded down) of data from set A; assume the number of selected data is K. In the formula, This indicates rounding down to the nearest integer. That is, selecting {{H} (1) L (1) T (1)},{H (2) L (2) T (2)},…,{H (K) L (K) T (K)}}

[0226] S5.2: Calculate the significant wave height, significant wavelength, and significant period of ocean waves.

[0227] Effective wave height H s The formula is:

[0228]

[0229] Effective wavelength L s The formula is:

[0230]

[0231] Effective period T s The formula is:

[0232]

[0233] Effective wave height H s Effective wavelength L s and effective period T s This provides the required wave data, thus completing the purpose of the plan.

[0234] S5.3: Calculate the wave spectrum.

[0235] include:

[0236] S5.3.1: Assume that after looping S1-S4 u times, u frames of ocean wave point cloud data are generated, and the generation times of the u frames of ocean wave point cloud data are t. b1 t b2 ... t bu The u-frame wave point cloud data in step S5.3 can be the same as or different from the n-frame wave point cloud data in step S5.2.

[0237] S5.3.2: In the u-frame wave point cloud data, select a fixed point whose position remains unchanged in each frame of wave point cloud data. Obtain the elevation of the same position (fixed point) in each frame of wave point cloud data and set it as: H x1 H x2 ... H xu .

[0238] At this point, t can be calculated. b1 to t bu The formula for the wave spectrum of a time series is:

[0239]

[0240] In the formula, ω is the angular frequency of the ocean wave. P(ω) is the power of the ocean wave; τ is the time offset.

[0241] The formula for τ is:

[0242]

[0243] The wave spectrum is also wave data that needs to be observed, thus completing the purpose of the plan.

[0244] Based on the above-described preferred embodiments of the present invention, and through the foregoing description, those skilled in the art can make various changes and modifications without departing from the inventive concept. The technical scope of this invention is not limited to the contents of the specification, but must be determined according to the scope of the claims.

Claims

1. A method for observing ocean waves based on a binocular camera, characterized in that, Includes the following steps: S1: The binocular camera captures images of ocean waves and generates wave point cloud data at the current time t. Wave point cloud data Established using the left camera coordinate system; S2: Based on the IMU (Inertial Measurement Unit), the wave point cloud data in the left camera coordinate system is transferred. Ocean wave point cloud data converted to geodetic coordinate system ; S3: Wave data in geodetic coordinate system in S2 Perform quality screening and optimization; Quality screening and optimization include: selecting regions of interest, quality screening, converting noisy data points to null values, null padding, and Gaussian filtering; Quality screening will include wave data The quality is categorized as: Excellent, Good, and Poor; poor-quality wave data is then discarded. Only high-quality and good wave data are retained. Proceed to the next step; S4: Based on high-quality and good wave point cloud data Calculate the wave height, wave wavelength, and wave period at that moment; S5: Repeat S1-S4, with each iteration progressing the time, to obtain multiple frames of wave point cloud data at different times. Based on these multiple frames, calculate the effective wave height, effective wave wavelength, effective wave period, and wave spectrum. Step S3, which involves selecting the region of interest and quality filtering, includes the following steps: S3.1: Select the region of interest, using the following formula: ； in, The row and column ratio of the region of interest. These represent the resolution of the image captured by the camera, i.e., the number of pixels in rows and columns. The index value of the selected wave point cloud matrix data, i.e., according to Extract the wave point cloud data of interest; S3.2: Count the number of noise and missing values ​​in the region of interest; set quality screening rules and perform screening; based on the screening results, evaluate the wave point cloud data. To improve the quality of wave point cloud data, discard low-quality data and retain high-quality data for subsequent steps. S3.2 specifically includes the following steps: S3.2.1: Obtain wave elevation data for each data point in the region of interest; assuming there are a total of S One "data point" S The "data points" include S One data point "containing elevation data" and S Two data points that "do not contain elevation data"; S = S 1+S2; Data points that "do not contain elevation data" are "null values"; No. i The wave elevation data for each data point is denoted as Z. i ，i ∈[1, S 1], calculate the first i The noise value for each data point is calculated using the following formula: ； In the formula, μ and σ They are respectively S The mean and standard deviation of wave elevation data for a single "data point"; A i For the first i The noise value of each data point; if Then the first i Each data point is denoted as a "noise data point"; traverse... S 1 data point, obtain S The number of "noise data points" in a single "data point" is denoted as . M indivual, M ∈[1, S 1]; Traversal S At the same time as obtaining the data points, the number of the aforementioned null values ​​is also obtained. S 2; S3.2.2: Calculate the ratio of the sum of blank values ​​and noise points to the total number of data points in the region of interest, and denote this ratio as... P The formula is: ； S3.2.3: Set quality screening rules Judge Determine the wave point cloud data The quality; The formula for judgment is: ； when Judge A value of -1 indicates that the wave point cloud data in this frame is... The quality is poor, and this frame of wave point cloud data should be removed. when Judge A value of 1 indicates that the wave point cloud data in this frame is... Excellent quality; when Judge A value of 0 indicates that the wave point cloud data for that frame is... Good quality; retain the high-quality and good-quality wave point cloud data for subsequent calculations.

2. The wave observation method based on a binocular camera according to claim 1, characterized in that, Step S3 also includes the following steps: S3.3: If the wave point cloud data The quality assessment is rated as excellent or good, and the wave point cloud data is optimized. Noisy data points and null values; S3.3 specifically includes the following steps: S3.3.1: Use row and column dilution factors to dilute the wave point cloud data Dilution, the diluted wave point cloud data The "noise data points" in the data should be converted to null values. Set dilution The CCP r One null value, r Two non-null values; r 1. This includes the total number of null values ​​in the original wave point cloud data, as well as the total number of null values ​​converted from "noise data points"; Let the first b The coordinates of the null values ​​are ; S3.3.2: Using null coordinates Centered on the data point, with side length R representing the search range; if a valid value exists within the search range, the missing value needs to be filled; if no valid value exists within the search range, the missing value is not filled; the valid value is the data point "containing elevation data" mentioned above. The filling process in step S3.3.2 includes the following steps: S3.3.2.1: Calculate the Euclidean distance from the null value to each non-null value; the formula is: ； In the formula, Indicates the first f Coordinates of non-empty values; f ∈[1, r 2]; S3.3.2.2: Calculate the weight of the null value to each non-null value; The formula is: ； In the formula, p The set weight index; S3.3.2.3: Filler value at empty spaces, the formula is: ； In the formula, This is the selected valid non-null Z-coordinate data. When calculating the X and Y axis data, simply replace it with the corresponding axis data. S3.3.2.4: Gaussian filtering; Gaussian function: ； Construct a Given a kernel matrix of size [size missing], calculate the weight of each matrix position: ； in These are the coordinates of each matrix data point relative to the kernel center. , ; To preserve information from subsequent elevation data, the Gaussian kernel is normalized to ensure that the sum of the filter kernels is 1. ； Original wave elevation matrix Symmetrical expansion around the four sides k 1 pixel, generate an extended matrix Original wave elevation matrix That is, the wave elevation data after the gaps are filled; Wave elevation data matrix Each raw data point Using the kernel center, calculate the smoothed value: ； The Gaussian kernel is applied to the corresponding positions of the extended matrix, and the results are multiplied point by point and then summed. Finally, the boundary portion of the extension is removed, and the smoothed result with the same size as the original data is retained for subsequent operations.

3. The wave observation method based on a binocular camera according to claim 2, characterized in that, In step S4, the wave point cloud data at the current time t is calculated. The wave height, wavelength, and period of ocean waves include the following steps: S4.1: Transfer wave point cloud data The ocean waves can be viewed as a matrix with a uniform X-axis. Then the matrix Further compression yields vector G t The zero-point crossing identification method is used to extract vectors. G t The maximum and minimum elevation values ​​are denoted as follows: H max and H min； S4.1: Calculate the wave height at the current time t; The formula is: ； In the formula, H t Let t be the wave height at the current moment; S4.2: Calculate the wavelength of the ocean wave at the current time t; get H max and H min Then, based on the vector G t The index in the file allows you to query the corresponding... H max and H min The data of the sea waves on the X and Y axes; for the uncompressed data matrix, there are a total of [data points] along the Y direction. M t There are data points, with corresponding values ​​along the X direction. N t Quantitative data, assuming H max The corresponding data on the X and Y axes are as follows: , H min The corresponding data on the X and Y axes are Therefore, the current time can be calculated. t The wavelength of ocean waves; the formula is: ； Calculate ocean wave point cloud data The average wavelength of ocean waves is considered as the wave wavelength at the current time t; the formula is: ； S4.3: Calculate the period of the ocean wave at the current time t; There are two scenarios: Case 1: For deep water waves with a depth greater than or equal to half the wavelength; The formula for calculating the wave period is: ； In the formula, c The wave speed is represented by g; g is the acceleration due to gravity. L t It is the wavelength of the ocean wave at the current time t; T t It is the current moment. t The wave cycle at that time; Scenario 2: For shallow water waves where the water depth is less than half the wavelength; The formula for calculating the wave period is: ； In the formula, c The wave speed is represented by g; g is the acceleration due to gravity. L t It is the wavelength of the ocean wave at the current time t; h For water depth; T t It is the current moment. t The wave cycle at that time.

4. The wave observation method based on a binocular camera according to claim 3, characterized in that, Step S5, calculating the significant wave height, significant wavelength, significant period, and wave spectrum includes the following steps: S5.1: Continue to cycle S1-S4 multiple times. The time will advance with each cycle. The binocular camera continuously captures images. The computer connected to the binocular camera continuously generates wave point cloud data at different times. The wave height, wave wavelength and wave period data of the wave point cloud data at different times are continuously obtained through calculation. Suppose a total of n sets of wave data are generated. Sort the wave heights of the n sets of data in descending order, and then add the sorted data into a set in the same order. A In the middle; the set A Represented as: ,in ; set up Select ; S5.2: Calculate the significant wave height, significant wavelength, and significant period of ocean waves; Significant wave height H s The formula is: ； Effective wavelength L s The formula is: ； Effective period T s The formula is: 。 5. The wave observation method based on a binocular camera according to claim 4, characterized in that, Step S5, calculating the wave spectrum includes the following steps: S5.3.1: Set the loop S1-S4 to run multiple times, generating... u Frame of ocean wave point cloud data; S5.3.2: In the u-frame wave point cloud data, select a fixed point whose position remains unchanged in each frame of wave point cloud data. Obtain the elevation of this same location in each frame of wave point cloud data, and denot it as . The generation time for each frame of wave point cloud data mentioned above is set as follows: but: t b1 arrive t bu The formula for the wave spectrum of a time series is: ； In the formula, The angular frequency of the ocean waves; The power of the ocean waves; It's a time offset; The formula is: 。

Images