Sea surface anomaly detection method based on sea surface roughness model using airborne LiDAR bathymetry
Through the airborne LiDAR deep-sea surface anomaly point detection method based on the sea surface roughness model, the coupled principal component analysis and ISODATA clustering algorithm are used to solve the problem of uncertainty error in sea surface images in traditional methods, and the submarine sounding accuracy is improved.
Patent Information
- Application Number
- CN202111365420.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2021-11-18
- Publication Date
- 2025-05-16
- Estimated Expiration
- 2041-11-18
AI Technical Summary
The traditional sea surface abnormal point detection method is based on sea surface images and is susceptible to uneven light, uneven pixel brightness and camera perspective, resulting in uncertainty errors and affecting the accuracy of onboard LiDAR depth sounding.
The airborne LiDAR deep-sea surface anomaly point detection method based on the sea surface roughness model is used to extract multiple waveforms and geometric feature parameters, and the sea surface roughness model is constructed using the coupled principal component analysis method, and the sea surface anomaly points are detected and identified through the ISODATA clustering algorithm.
It effectively avoids uncertain errors in sea surface images, reduces the false topography problems caused by sea surface anomalies, and improves the accuracy of airborne laser submarine depth sounding.
Smart Images

Figure CN114266986B_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the technical field of ocean surveying and mapping, and in particular relates to a method for detecting abnormal points on an airborne LiDAR bathymetric sea surface based on a sea surface roughness model. Background Art
[0002] The airborne LiDAR bathymetric system has the characteristics of high measurement accuracy, high measurement point density, high work efficiency, strong mobility, and measurement continuity. It is particularly suitable for rapid detection of complex terrains such as shallow water areas and areas near islands and reefs, and can achieve seamless splicing of the above-water and underwater terrain of the coastline. The airborne LiDAR bathymetric system uses a 532nm blue-green laser with strong water penetration ability as a means of seabed detection. When the blue-green laser reaches the sea surface, part of the laser beam returns along the incident path, and the other part of the laser beam penetrates the sea surface to reach the seabed. Due to the influence of external factors such as waves, sea breezes, and tides, waves will be generated on the sea surface. When the wind speed reaches a certain critical value, the waves will break, and a large number of water droplets and water foam will be generated at the crest of the wave, and a large number of bubbles will be generated inside and on the surface of the fluctuating water body. Sea surface anomalies refer to laser points on the sea surface caused by ocean whitecaps, breaking waves, water foam and bubbles during the airborne LiDAR bathymetry process. Their scattering effect will reduce the accuracy of water depth measurement to a certain extent. At the same time, sea surface anomalies will also change the distribution of water vapor density and salinity in the atmosphere-water interface layer. When the blue-green laser penetrates the atmosphere-water interface formed by the ocean anomalies, the laser incident angle and refraction angle will undergo serious path deviations, which will produce serious false seabed topography and affect the ALB bathymetry accuracy.
[0003] Most of the traditional sea surface anomaly detection methods are based on the sea surface image to identify the sea surface whitecaps. Since the image is easily affected by uneven lighting, uneven pixel brightness and camera perspective, the sea surface image itself has a certain uncertainty error, which will affect the accuracy of whitecap identification to a certain extent. In view of the uncertainty error existing in the traditional method, in order to reduce the influence of the ALB bathymetric error caused by sea surface anomalies such as ocean whitecaps, it is necessary to propose an airborne LiDAR bathymetric sea surface anomaly detection method based on the sea surface roughness model. The present invention uses the measured sea surface waveform and laser point cloud as the data basis, avoiding the uncertainty error of the sea surface image in the traditional method. Summary of the invention
[0004] In view of the above-mentioned technical problems existing in the prior art, the present invention proposes a method for detecting sea surface anomaly points in airborne LiDAR bathymetry based on a sea surface roughness model. The method has a reasonable design and solves the influence of uncertainty errors existing in the prior art. Based on the measured sea surface waveform and laser point cloud data, the method can effectively detect sea surface anomaly points, eliminate the false seabed terrain caused by sea surface anomaly points, and improve the accuracy of airborne laser seabed bathymetry.
[0005] In order to achieve the above object, the present invention adopts the following technical solution:
[0006] The method for detecting sea surface anomaly points using airborne LiDAR bathymetry based on a sea surface roughness model includes the following steps:
[0007] Step 1: Extract multivariate waveform and geometric feature parameters based on the acquired airborne LiDAR bathymetric sea surface data;
[0008] Step 2: Use the coupled principal component analysis method to linearly combine the extracted multivariate feature vectors into a sea surface roughness model that can express the comprehensive characteristic information of the sea surface;
[0009] Step 3: Based on the constructed waveform + geometric sea surface roughness model, the ISODATA clustering algorithm is used to detect and identify sea surface anomalies.
[0010] Preferably, in step 2, the sea surface roughness model is constructed as follows:
[0011] Step 2.1: Coupled principal component analysis is to extract m orthogonal directions from the input multivariate waveform and geometric feature vector by optimizing the information criterion, maximizing the variance of the projection data or minimizing the reconstruction error. Make the input feature vector have as large a variance as possible in these m directions; therefore, the input feature vector x∈R n Can be Zhang Cheng's m-dimensional subspace representation; among them, η=1,2,…,m; the inner product is used in the m-dimensional subspace
[0012] Step 2.2: Find the cell direction The projection of the input feature vector x along these directions has the maximum deviation, as shown in formula (19):
[0013]
[0014] In the formula, E PCA (w) is a semi-positive definite function;
[0015] Step 2.3: Set Then there is
[0016]
[0017] Where w = αc η(η=1,2,…,m;α∈R), c is the deviation matrix; when α=1, w is a unit vector, and only when w=c1, the Hessian matrix H(w) of w is semi-positive definite; therefore, w will eventually converge to the direction of c1, and the converged E PCA (w) is equal to the maximum value λ1;
[0018] Step 2.4: Remove E PCA (w) and then maximize E PCA (w), then E PCA (w) obtains the maximum value λ2 in the direction of w=c2;
[0019] Step 2.5: Repeat steps 2.1 to 2.4 until all m main directions are and the principal components of the input feature vector x All are derived;
[0020] Step 2.6: Based on the obtained principal components The input feature vector x is linearly combined to obtain a sea surface roughness model that can express the comprehensive characteristic information of the sea surface.
[0021] Preferably, in step 3, the specific method of detecting sea surface anomalies using the ISODATA clustering algorithm is as follows:
[0022] Step 3.1: Select initial parameters and distribute the sea surface roughness data to each cluster center according to the indicators; the specific steps include the following:
[0023] Step 3.1.1: Input sea surface waveform + geometric roughness data, and preset C initial cluster centers {z1,z2,z3,…,z C}, the initial position can be chosen arbitrarily;
[0024] Step 3.1.2: Preset parameters;
[0025] K—the expected number of clusters, C—the number of initial cluster centers; θ N —The minimum number of roughnesses allowed in each cluster, θ E —The maximum relative standard deviation of the distribution of each feature component within the class, greater than θ E The cluster of θ is split, C —The minimum distance between two cluster centers, less than θ C The two clusters of N are merged. T —The maximum number of “merges” that can be performed during each iteration, N S —maximum number of iterations;
[0026] Step 3.1.3: Assign the sea surface roughness to the nearest cluster S according to the nearest neighbor rule j, assuming D j =min{||xz j ||,j=1,2,…,C}, that is, ||xz j || is the smallest distance, then x∈S j ;
[0027] Step 3.1.4: If S j The number of roughness N in j <θ N , then cancel the cluster, and C is reduced by 1;
[0028] Step 3.2: Calculate the distance index function of the sea surface roughness in each cluster; specifically include the following steps:
[0029] Step 3.2.1: Modify each cluster center, as shown in formula (21):
[0030]
[0031] Step 3.2.2: Calculate each clustering domain S j The average distance between the medium roughness and each cluster center is shown in formula (22):
[0032]
[0033] Step 3.2.3: Calculate the total average distance between all sea surface roughness and their corresponding cluster centers, as shown in formula (23):
[0034]
[0035] Step 3.3: Determine the splitting, merging and iterative operations, which specifically include the following steps:
[0036] Step 3.3.1: If the number of iterations is greater than N S , then set θ C =0, go to step 3.5.1;
[0037] Step 3.3.2: If That is, if the number of cluster centers is less than or equal to half of the expected number of clusters, go to step 3.4.1 and perform a split operation on the existing clusters;
[0038] Step 3.3.3: If C ≥ 2K, do not perform splitting and go to step 3.5.1; otherwise, go to step 3.3.4;
[0039] Step 3.3.4: If When the number of iterations is an odd number, it switches to the split operation; when the number of iterations is an even number, it switches to the merge operation;
[0040] Step 3.4: Split operation, including the following steps:
[0041] Step 3.4.1: Calculate the standard deviation vector of the roughness distance in each cluster, as shown in formula (24):
[0042] σ j =(σ 1j ,σ 2j ,…,σ nj ) T (twenty four);
[0043] Among them, the components of the standard deviation vector are:
[0044]
[0045] In the above formula, i = 1, 2, ..., n is the dimension of the roughness vector; j = 1, 2, ..., C is the number of clusters; N j For S j The number of roughness in ;
[0046] Step 3.4.2: Calculate the standard deviation vector {σ j ,j=1,2,…,C} jmax ;
[0047] Step 3.4.3: For any maximum component σ jmax , j = 1, 2, ..., C, if σ jmax >θ E , and at the same time, one of the following two conditions is met:
[0048] ① And N j >2(θ N +1), that is, S j The total number of roughness in the cluster exceeds twice the minimum number of point clouds allowed;
[0049] ②
[0050] Then the class S j Split into two classes, the original cluster center z j Delete, and let C = C + 1; the centers of the two new classes z j + and z j - respectively at the original cluster center z j Based on adding and subtracting σ jmax The corresponding component aσ jmax , where 0 <a≤1;
[0051] Step 3.5: Merge operation, including the following steps:
[0052] Step 3.5.1: Calculate the distances of all cluster centers, as shown in formula (26):
[0053] D βγ =||z β -z γ ||, β=1,2,…,C-1,γ=2,3,…,C (26);
[0054] Step 3.5.2: Compare D βγ and θ C The size of D βγ <θ C The values of are arranged in ascending order according to the smallest distance, that is:
[0055] {D β1γ1 ,D β2γ2 ,…,D βLγL}, D β1γ1 <D β2γ2 <… <D βLγL (27);
[0056] Step 3.5.3: Set the distance D βlγl The two cluster centers z βl and z γl Merge and get the new center as formula (28):
[0057]
[0058] In the formula, the two merged cluster center vectors are weighted by the roughness in their cluster domains respectively, so that is the true mean vector;
[0059] Step 3.6: Re-iterate and calculate various indicators to converge the clustering results, identify and remove abnormal points on the sea surface;
[0060] If the input parameters need to be modified, go to step 3.1.1, otherwise go to step 3.1.2, calculate to the maximum number of iterations, detect and remove abnormal points on the sea surface through the roughness clustering results, and end the calculation.
[0061] Beneficial technical effects brought by the present invention:
[0062] The present invention proposes an airborne LiDAR bathymetric sea surface anomaly point detection method based on a sea surface roughness model. Compared with the prior art, the present invention constructs a sea surface roughness model based on the extracted multivariate waveform and geometric feature vector using a coupled principal component analysis method, and detects sea surface anomalies using an ISODATA clustering algorithm. The present invention uses the measured sea surface waveform and laser point cloud as data basis, avoids the uncertainty error of the sea surface image in the traditional method, solves the problem of false seabed terrain caused by the path deviation of the laser incident angle and refraction angle caused by the sea surface anomaly, and improves the accuracy of seabed bathymetry. BRIEF DESCRIPTION OF THE DRAWINGS
[0063] Figure 1 This is a flow chart of the method for detecting sea surface anomalies using airborne LiDAR bathymetry based on a sea surface roughness model according to the present invention.
[0064] Figure 2 Schematic diagram of coupled principal component analysis of two-dimensional input vectors in the present invention.
[0065] Figure 3 This is a flow chart of the ISODATA clustering algorithm used in the present invention to detect sea surface anomalies. DETAILED DESCRIPTION
[0066] The present invention is further described in detail below with reference to the accompanying drawings and specific embodiments:
[0067] The present invention provides a method for detecting sea surface anomalies using airborne LiDAR bathymetry based on a sea surface roughness model. Figure 1 As shown, the following steps are included:
[0068] Step 1: Extract multivariate waveform and geometric feature parameters based on airborne LiDAR bathymetric sea surface data.
[0069] The extraction of multivariate waveforms and geometric features of airborne LiDAR bathymetry has improved the data mining capabilities of airborne LiDAR bathymetry and provided a data basis for application areas such as sea surface anomaly detection.
[0070] Step 2: Use the coupled principal component analysis method to linearly combine the extracted multivariate feature vectors into a sea surface roughness model that can express the comprehensive characteristic information of the sea surface.
[0071] The extracted multivariate waveform and geometric features can reflect the sea surface information from different levels, but there is also information overlap and inconsistent weights. Therefore, the coupled principal component analysis algorithm is used to linearly combine the extracted waveform features and geometric features into a sea surface roughness model that can express the comprehensive characteristic information of the sea surface. Figure 2 .
[0072] In a further embodiment, step 2 specifically includes the following steps:
[0073] Step 2.1: Coupled principal component analysis is to extract m orthogonal directions from the input multivariate waveform and geometric feature vector by optimizing the information criterion, maximizing the variance of the projection data or minimizing the reconstruction error. Make the input feature vector have as large a variance as possible in these m directions; therefore, the input feature vector x∈R n Can be Zhang Cheng's m-dimensional subspace representation; among them, η=1,2,…,m; the inner product is used in the m-dimensional subspace
[0074] Step 2.2: Find the cell direction The projection of the input feature vector x along these directions has the maximum deviation, as shown in formula (19):
[0075]
[0076] In the formula, E PCA (w) is a semi-positive definite function;
[0077] Step 2.3: Set Then there is
[0078]
[0079] Where w = αc η (η=1,2,…,m;α∈R), c is the deviation matrix; when α=1, w is a unit vector, and only when w=c1, the Hessian matrix H(w) of w is semi-positive definite; therefore, w will eventually converge to the direction of c1, and the converged E PCA (w) is equal to the maximum value λ1;
[0080] Step 2.4: Remove E PCA (w) and then maximize E PCA (w), then E PCA (w) obtains the maximum value λ2 in the direction of w=c2;
[0081] Step 2.5: Repeat steps 2.1 to 2.4 until all m main directions are and the principal components of the input feature vector x All are derived;
[0082] Step 2.6: Based on the obtained principal components The input feature vector x is linearly combined to obtain a sea surface roughness model that can express the comprehensive characteristic information of the sea surface.
[0083] In the specific implementation, the coupled principal component analysis method constructs an efficient feature parameter selection mechanism, which linearly combines multi-index factors such as waveform and geometric feature vectors into a roughness model. This model can reflect most of the information of the original feature vector, and the information contained is non-repetitive, thus obtaining more scientific and effective data information.
[0084] Step 3: Based on the constructed waveform + geometric sea surface roughness model, the ISODATA clustering algorithm is used to detect and identify sea surface anomalies.
[0085] The ISODATA clustering algorithm can be used to set initial parameters, use the "merge" and "split" mechanisms, and iterate according to the initial cluster center and the number of categories set, until the sea surface waveform + geometric roughness data is classified into multiple clustering results, thereby detecting sea surface anomalies. The ISODATA clustering algorithm detects sea surface anomalies. Figure 3 .
[0086] In a further embodiment, step 3 specifically includes the following steps:
[0087] Step 3.1: Select initial parameters and distribute the sea surface roughness data to each cluster center according to the indicators; the specific steps include the following:
[0088] Step 3.1.1: Input sea surface waveform + geometric roughness data, and preset C initial cluster centers {z1,z2,z3,…,z C}, the initial position can be chosen arbitrarily;
[0089] Step 3.1.2: Preset parameters;
[0090] K—the expected number of clusters, C—the number of initial cluster centers; θ N —The minimum number of roughnesses allowed in each cluster, θ E —The maximum relative standard deviation of the distribution of each feature component within the class, greater than θ E The cluster of θ is split, C —The minimum distance between two cluster centers, less than θ C The two clusters of N are merged. T —The maximum number of “merges” that can be performed during each iteration, N S —maximum number of iterations;
[0091] Step 3.1.3: Assign the sea surface roughness to the nearest cluster S according to the nearest neighbor rule j , assuming D j =min{||xz j ||,j=1,2,…,C}, that is, ||xzj || is the smallest distance, then x∈S j ;
[0092] Step 3.1.4: If S j The number of roughness N in j <θ N , then cancel the cluster, and C is reduced by 1;
[0093] Step 3.2: Calculate the distance index function of the sea surface roughness in each cluster; specifically include the following steps:
[0094] Step 3.2.1: Modify each cluster center, as shown in formula (21):
[0095]
[0096] Step 3.2.2: Calculate each clustering domain S j The average distance between the medium roughness and each cluster center is shown in formula (22):
[0097]
[0098] Step 3.2.3: Calculate the total average distance between all sea surface roughness and their corresponding cluster centers, as shown in formula (23):
[0099]
[0100] Step 3.3: Determine the splitting, merging and iterative operations, which specifically include the following steps:
[0101] Step 3.3.1: If the number of iterations is greater than N S , then set θ C =0, go to step 3.5.1;
[0102] Step 3.3.2: If That is, if the number of cluster centers is less than or equal to half of the expected number of clusters, go to step 3.4.1 and perform a split operation on the existing clusters;
[0103] Step 3.3.3: If C ≥ 2K, do not perform splitting and go to step 3.5.1; otherwise, go to step 3.3.4;
[0104] Step 3.3.4: If When the number of iterations is an odd number, it switches to the split operation; when the number of iterations is an even number, it switches to the merge operation;
[0105] Step 3.4: Split operation, including the following steps:
[0106] Step 3.4.1: Calculate the standard deviation vector of the roughness distance in each cluster, as shown in formula (24):
[0107] σ j =(σ 1j ,σ 2j ,…,σ nj ) T (twenty four);
[0108] Among them, the components of the standard deviation vector are:
[0109]
[0110] In the above formula, i = 1, 2, ..., n is the dimension of the roughness vector; j = 1, 2, ..., C is the number of clusters; N j For S j The number of roughness in ;
[0111] Step 3.4.2: Calculate the standard deviation vector {σ j ,j=1,2,…,C} jmax ;
[0112] Step 3.4.3: For any maximum component σ jmax , j = 1, 2, ..., C, if σ jmax >θ E , and at the same time, one of the following two conditions is met:
[0113] ① And N j >2(θ N +1), that is, S j The total number of roughness in the cluster exceeds twice the minimum number of point clouds allowed;
[0114] ②
[0115] Then the class S j Split into two classes, the original cluster center z j Delete, and let C = C + 1; the centers of the two new classes z j + and z j - respectively at the original cluster center z j Based on adding and subtracting σ jmax The corresponding component aσ jmax , where 0 <a≤1;
[0116] Step 3.5: Merge operation, including the following steps:
[0117] Step 3.5.1: Calculate the distances of all cluster centers, as shown in formula (26):
[0118] Dβγ =||z β -z γ ||, β=1,2,…,C-1,γ=2,3,…,C (26);
[0119] Step 3.5.2: Compare D βγ and θ C The size of D βγ <θ C The values of are arranged in ascending order according to the smallest distance, that is:
[0120] {D β1γ1 ,D β2γ2 ,…,D βLγL}, D β1γ1 <D β2γ2 <… <D βLγL (27);
[0121] Step 3.5.3: Set the distance D βlγl The two cluster centers z βl and z γl Merge and get the new center as formula (28):
[0122]
[0123] In the formula, the two merged cluster center vectors are weighted by the roughness in their cluster domains respectively, so that is the true mean vector;
[0124] Step 3.6: Re-iterate and calculate various indicators to converge the clustering results, identify and remove abnormal points on the sea surface;
[0125] If the input parameters need to be modified, go to step 3.1.1, otherwise go to step 3.1.2, calculate to the maximum number of iterations, detect and remove abnormal points on the sea surface through the roughness clustering results, and end the calculation.
[0126] In the specific implementation, the "merge" and "splitting" mechanisms of the ISODATA algorithm are used to effectively detect the information of sea surface anomalies, and the sea surface after removing the sea surface anomalies is closer to the real calm sea surface. Therefore, the ISODATA clustering algorithm has good applicability and feasibility for sea surface anomaly detection and can improve the accuracy of airborne laser seabed bathymetry.
[0127] In summary, the present invention provides a method for detecting abnormal points on the sea surface using airborne LiDAR bathymetry based on a sea surface roughness model, the method comprising: extracting multivariate waveforms and geometric feature vectors based on the acquired airborne LiDAR bathymetry sea surface data, using a coupled principal component analysis method to linearly combine the extracted feature vectors into a sea surface roughness model that can express comprehensive feature information of the sea surface; based on the constructed sea surface roughness model, detecting and identifying abnormal points on the sea surface using an ISODATA clustering algorithm. The present invention uses measured sea surface waveforms and laser point clouds as data basis, avoids the uncertainty error of sea surface images in traditional methods, solves the problem of false seabed topography caused by path deviation of laser incident angles and refraction angles caused by abnormal points on the sea surface, and improves the accuracy of seabed bathymetry.
[0128] Of course, the above description is not a limitation of the present invention, and the present invention is not limited to the above examples. Changes, modifications, additions or substitutions made by technicians in this technical field within the essential scope of the present invention should also fall within the protection scope of the present invention.
Claims
1. An airborne LiDAR bathymetric sea surface anomaly detection method based on a sea surface roughness model, characterized by: The following steps are involved: Step 1: Extract multivariate waveform and geometric feature parameters based on the acquired airborne LiDAR bathymetric sea surface data; Step 2: Use the coupled principal component analysis method to linearly combine the extracted multivariate feature vectors into a sea surface roughness model that can express the comprehensive characteristic information of the sea surface; Step 3: Based on the constructed waveform + geometric sea surface roughness model, the ISODATA clustering algorithm is used to detect and identify sea surface anomalies.
2. The method for detecting sea surface anomalies using airborne LiDAR bathymetry based on a sea surface roughness model according to claim 1, characterized in that: In step 2, the sea surface roughness model is constructed as follows: Step 2.1: Coupled principal component analysis is to extract m orthogonal directions from the input multivariate waveform and geometric feature vector by optimizing the information criterion, maximizing the variance of the projection data or minimizing the reconstruction error. Make the input feature vector have as large a variance as possible in these m directions; therefore, the input feature vector x∈R n Can be Zhang Cheng's m-dimensional subspace representation; among them, Inner product Step 2.2: Find the cell direction The projection of the input feature vector x along these directions has the maximum deviation, as shown in formula (19): In the formula, E PCA (w) is a semi-positive definite function; Step 2.3: Set Then there is Where w = αc η (η=1,2,…,m;α∈R), c is the deviation matrix; when α=1, w is a unit vector, and only when w=c1, the Hessian matrix H(w) of w is semi-positive definite; therefore, w will eventually converge to the direction of c1, and the converged E PCA (w) is equal to the maximum value λ1; Step 2.4: Remove E PCA (w) and then maximize E PCA (w), then E PCA (w) obtains the maximum value λ2 in the direction of w=c2; Step 2.5: Repeat steps 2.1 to 2.4 until all m main directions are and the principal components of the input feature vector x All are derived; Step 2.6: Based on the obtained principal components The input feature vector x is linearly combined to obtain a sea surface roughness model that can express the comprehensive characteristic information of the sea surface.
3. The method for detecting sea surface anomalies using airborne LiDAR bathymetry based on a sea surface roughness model according to claim 1, characterized in that: In step 3, the specific method of detecting sea surface anomalies using the ISODATA clustering algorithm is as follows: Step 3.1: Select initial parameters and distribute the sea surface roughness data to each cluster center according to the indicators; the specific steps include the following: Step 3.1.1: Input sea surface waveform + geometric roughness data, and preset C initial cluster centers {z1,z2,z3,…,z C }, the initial position can be chosen arbitrarily; Step 3.1.2: Preset parameters; K—the expected number of clusters, C—the number of initial cluster centers; θ N —The minimum number of roughnesses allowed in each cluster, θ E —The maximum relative standard deviation of the distribution of each feature component within the class, greater than θ E The cluster of θ is split, C —The minimum distance between two cluster centers, less than θ C The two clusters of N are merged. T —The maximum number of "merges" that can be performed during each iteration, N S —maximum number of iterations; Step 3.1.3: Assign the sea surface roughness to the nearest cluster S according to the nearest neighbor rule j , assuming D j =min{||xz j ||,j=1,2,…,C}, that is, ||xz j || is the smallest distance, then x∈S j ; Step 3.1.4: If S j The number of roughness N in j <θ N , then cancel the cluster, and C is reduced by 1; Step 3.2: Calculate the distance index function of the sea surface roughness in each cluster; specifically include the following steps: Step 3.2.1: Modify each cluster center, as shown in formula (21): Step 3.2.2: Calculate each clustering domain S j The average distance between the medium roughness and each cluster center is shown in formula (22): Step 3.2.3: Calculate the total average distance between all sea surface roughness and their corresponding cluster centers, as shown in formula (23): Step 3.3: Determine the splitting, merging and iterative operations, which specifically include the following steps: Step 3.3.1: If the number of iterations is greater than N S , then set θ C =0, go to step 3.5.1; Step 3.3.2: If That is, if the number of cluster centers is less than or equal to half of the expected number of clusters, go to step 3.4.1 and perform a split operation on the existing clusters; Step 3.3.3: If C ≥ 2K, do not perform splitting and go to step 3.5.1; otherwise, go to step 3.3.4; Step 3.3.4: If When the number of iterations is an odd number, it switches to the split operation; when the number of iterations is an even number, it switches to the merge operation; Step 3.4: Split operation, including the following steps: Step 3.4.1: Calculate the standard deviation vector of the roughness distance in each cluster, as shown in formula (24): s j =(s 1j ,s 2j ,…,s nj ) T (24); Among them, the components of the standard deviation vector are: In the above formula, i = 1, 2, ..., n is the dimension of the roughness vector; j = 1, 2, ..., C is the number of clusters; N j For S j The number of roughness in ; Step 3.4.2: Calculate the standard deviation vector {σ j ,j=1,2,…,C} jmax ; Step 3.4.3: For any maximum component σ jmax , j = 1, 2, ..., C, if σ jmax >θ E , and at the same time, one of the following two conditions is met: ① And N j >2(θ N +1), that is, S j The total number of roughness in the cluster exceeds twice the minimum number of point clouds allowed; ② Then the class S j Split into two classes, the original cluster center z j Delete and set C = C + 1; the centers of the two new classes and They are respectively at the original cluster center z j Based on adding and subtracting σ jmax The corresponding component aσ jmax , where 0 <a≤1; Step 3.5: Merge operation, including the following steps: Step 3.5.1: Calculate the distances of all cluster centers, as shown in formula (26): D βγ =||z β -z γ ||,β=1,2,…,C-1,γ=2,3,…,C (26); Step 3.5.2: Compare D βγ and θ C The size of D βγ <θ C The values of are arranged in ascending order according to the smallest distance, that is: {D β1γ1 ,D β2γ2 ,…,D βLγL },D β1γ1 <D β2γ2 <…<D βLγL (27); Step 3.5.3: Set the distance D βlγl The two cluster centers z βl and z γl Merge and get the new center as formula (28): In the formula, the two merged cluster center vectors are weighted by the roughness in their cluster domains respectively, so that is the true mean vector; Step 3.6: Re-iterate and calculate various indicators to converge the clustering results, identify and remove abnormal points on the sea surface; If the input parameters need to be modified, go to step 3.1.1, otherwise go to step 3.1.2, calculate to the maximum number of iterations, detect and remove abnormal points on the sea surface through the roughness clustering results, and end the calculation.