Multi-beam outlier automatic filtering method combining uncertainty and density clustering algorithm

Abstract

The application discloses a kind of multi-beam outlier automatic filtering methods combined with uncertainty and density clustering algorithm, belong to multi-beam sounding data processing technical field.The application constructs grid node absorbable water depth point selection model according to CUBE filtering algorithm, and gridding is carried out to sounding data;The one-dimensional clustering of node absorbable water depth value is carried out using DBSCAN clustering algorithm, a plurality of water depth hypotheses are constructed according to the difference between water depth values, then the water depth hypothesis and uncertainty of node are updated using Kalman filter;Finally, the water depth hypothesis with higher reliability is selected from a plurality of water depth hypotheses as the real water depth value of node.The application can provide a kind of multi-beam sounding data automatic cleaning algorithm selection for domestic multi-beam sounding data processing software.

Description

Technical Field

[0001] This invention belongs to the field of multibeam bathymetry data processing technology, and relates to an automatic method for cleaning up outliers in multibeam bathymetry data. Background Technology

[0002] Multibeam bathymetry (MBES) systems are a primary technology for seabed topographic surveying, and removing outliers from bathymetry data is a crucial step in obtaining high-precision results. In recent years, the volume of raw multibeam bathymetry data has exploded, and manual processing suffers from long processing times and low efficiency, severely limiting the efficiency of processing massive amounts of data. The CUBE (Combined Uncertainty and Bathymetry Estimator algorithm) is the core algorithm of the CARIS multibeam bathymetry data processing software. CARIS has continuously improved and refined the CUBE filtering algorithm (see CARISSIPS and HIPS Change List), and its effectiveness has been widely recognized in the industry. However, because the CUBE filtering algorithm is a trade secret, the internal model of the continuously improved CUBE filtering algorithm is not publicly available, making it impossible for researchers to investigate its specific implementation methods. Replicating the CUBE filtering algorithm based on existing literature cannot achieve the same results as the CARIS built-in algorithm. Faced with the urgent need to process massive amounts of multibeam bathymetry data, the long-standing use of foreign software has been subject to numerous limitations. The demand for mature and comprehensive domestically developed multibeam bathymetry data processing software is increasingly pressing. The lack of a reliable automatic outlier removal algorithm for multibeam bathymetry data hinders the localization of this software. Therefore, proposing a clear, parameter-simple, and reliable automatic multibeam bathymetry value removal algorithm and verifying its effectiveness is a major challenge that urgently needs to be addressed in the localization of multibeam bathymetry data processing software. Summary of the Invention

[0003] Based on in-depth research of the CUBE filtering algorithm, this invention fully draws on the CUBE filtering algorithm framework and process, introduces the DBSCAN clustering algorithm (Density-Based Spatial Clustering of Applications with Noise) to identify and cluster outliers, and improves the CUBE filtering algorithm's multi-depth hypothesis tracking model. It proposes a multi-beam automatic outlier filtering algorithm that combines uncertainty and density clustering: A grid node absorbable depth point selection model is constructed based on the CUBE filtering algorithm, and the bathymetry data is gridded; the DBSCAN clustering algorithm is used to perform one-dimensional clustering of the absorbable depth values ​​of the nodes, and multiple depth hypotheses are constructed based on the differences between the depth values; then, Kalman filtering is used to update the depth hypotheses and uncertainties of the nodes; finally, the depth hypothesis with higher confidence is selected from the multiple depth hypotheses as the true depth value of the node. Based on the implementation of the CUBE filtering algorithm, through comparison of manual interactive processing method, CUBE filtering algorithm and the present invention, the measured data and simulation test results show that the filtering results of the present invention are more robust, can effectively handle single and continuous outliers in multibeam bathymetry data, and have wider applicability, providing a new idea to help the development of domestic multibeam bathymetry data processing software.

[0004] To achieve the above objectives, the technical solution of the present invention is as follows:

[0005] An automatic filtering method for multi-beam outliers combining uncertainty and density clustering algorithms mainly includes six steps: estimation surface construction, water depth value clustering, centroid point calculation, water depth and uncertainty estimation, selection of optimal estimation, and filtering of raw bathymetry data (see...). Figure 1 The specific steps are as follows:

[0006] (1) Estimating the construction of the surface

[0007] The extrapolated surface is a gridded surface independent of the original bathymetry data (x, y, z). Within the plane coordinates (x, y) of the original bathymetry data, a grid spacing is set, and the grid intersections are called nodes (x, y). k ,y k ,z k (k is the number of nodes). Node planar coordinates (x...) k ,y k The water depth value (z) of the node is determined by the meshing process. k The depth is estimated from the water depth values ​​of its neighboring water depth points. When constructing the estimation surface, the absorbable water depth points of the nodes should be calculated for subsequent steps; the specific process is as follows:

[0008] First, calculate the node absorption domain Rcap using formula (1):

[0009]

[0010] Where μcapdiaca is the node capture distance ratio coefficient; μcapdismin is the node minimum absorption area; and Gendepth is the average water depth within a circle centered on the node with a radius of three times the grid.

[0011] Then, the capture radius Rinflu at the water depth point is calculated using formula (2):

[0012]

[0013] Where Grid is the grid spacing; TVUmax is the limit of vertical uncertainty allowed for different levels of hydrographic surveying in the IHO S-44 (2008 version) standard; TVU is the vertical uncertainty of the depth point; and μdistexp is the proportionality coefficient of uncertainty as distance increases.

[0014] According to IHO S-44 (2008 edition) standard, the allowable limit of vertical uncertainty for different grades of hydrographic surveys is expressed using the following formula:

[0015]

[0016] The values ​​of a and b are shown in Table 1, and d is the water depth value of the depth point.

[0017] Table 1. Minimum Standards for Hydrographic Surveying of Different Grades (IHO S-44)

[0018]

[0019] A water depth can only be absorbed by a node if the water depth point is located within the node's absorption domain and the node is also located within the water depth point's capture radius.

[0020] (2) Water depth value clustering

[0021] After step (1) is completed, the water depth points within the node's absorption domain that can be absorbed have been determined, and the nodes and water depth points have been associated on the two-dimensional plane. At this point, the node will typically have many possible water depth values. The DBSCANS clustering algorithm is used to perform one-dimensional clustering of the water depth values ​​of the water depth points that can be absorbed by the node. The good clustering results of the DBSCANS clustering algorithm form the basis for updating the node's water depth values ​​and uncertainty estimations, as well as selecting the optimal values.

[0022] The DBSCAN clustering algorithm has two key parameters (Eps, MinPts) to describe the density of the sample distribution in the neighborhood. The parameter Eps describes the Euclidean distance threshold of the neighborhood of a data point, used to calculate the number of depth points within a circular region centered on the core point and with Eps as the radius. Theoretically, Eps ranges from (0, +∞); however, based on numerous computational examples, a value of 0.8 to 1 for Eps is considered ideal. Considering the small differences in depth between multibeam bathymetry data points, the parameter Eps is set to a range of (0,1], and in practice, a value of 1 is used. A smaller Eps setting leads to an increase in the number of clusters, but a decrease in the number of depth points within each cluster. The parameter MinPts describes the threshold for the number of data points within a neighborhood of Eps for a core point. MinPts is fixed at 1 to ensure that nodes absorbing fewer depth points can be clustered normally without being identified as noise points. When data quality is poor and nodes absorb too many abnormal depth values, setting MinPts to 1 effectively distinguishes between different levels of abnormal depth values, ensuring that the number of clusters for "normal data points" is greater than that for "abnormal data points." By reasonably setting the two key parameters Eps and MinPts of the DBSCAN clustering algorithm, clusters for "normal data points" and "abnormal data points" are constructed, establishing multiple depth hypotheses to prepare for outlier removal.

[0023] (3) Calculation of the centroid

[0024] After step (2) is completed, the water depth values ​​of the water depth points within the node absorption domain are clustered into classes (the number of classes is ≥1; otherwise, it is an empty node and should be removed). The first water depth value entering the filter must represent the characteristics of that class of water depth values ​​and is called the water depth value of the centroid (including uncertainty). Considering the robustness of this invention in processing abnormal data, it is necessary to control the water depth value of the centroid. Therefore, this invention uses the water depth value of each class with the smallest difference from the average water depth value of that class as the water depth value of the centroid in that class; specifically as follows: First, calculate the average value of all water depth values ​​in each class, then calculate the difference between each water depth value in that class and the average water depth value, and select the water depth value corresponding to the smallest difference as the water depth value of the centroid (see Figure 2 ).

[0025] (4) Estimation of water depth and uncertainty

[0026] In step (3), the water depth value (including uncertainty) of the centroid point is determined. A Kalman filter is then used to estimate the water depth value and uncertainty. The estimated water depth value representing the possible water depth of that node is called the water depth assumption. The Kalman filter estimation formula is:

[0027] One-step prediction uncertainty:

[0028]

[0029] Water depth prediction:

[0030] z j [n|n-1]=z j [n-1|n-1](5)

[0031] Kalman filter gain K j [n]:

[0032]

[0033] The difference between the water depth value and the predicted value, v j [n]:

[0034] v j [n] = d j [n]-z j [n|n-1] (7)

[0035] Node depth value update:

[0036] z j [n|n]=z j [n|n-1]+K j [n]×v j [n] (8)

[0037] Node uncertainty update:

[0038]

[0039] In the above formula, σ j [n-1|n-1]、z j [n-1|n-1] represent the uncertainty and water depth value of the j-th node at time n-1 after Kalman filtering, respectively; σ j [n|n-1]、z j [n|n-1] represent the uncertainty and water depth value predicted by the j-th node at time n-1 to time n, respectively; d j [n] is the water depth entering the Kalman filter at time n; σ j [n] is d j Uncertainty of [n]; σ j [n|n] represents the uncertainty of the j-th node at time n after Kalman filtering.

[0040] This step ends when the absorbable water depth values ​​around node j are sequentially fed into the Kalman filter.

[0041] (5) Selection of the optimal valuation

[0042] Step (4) uses a Kalman filter to estimate the possible uncertainties and water depth values ​​representing the nodes. This invention selects the water depth assumption with the minimum uncertainty as the optimal estimate for the nodes. The number of clusters obtained by the DBSCAN clustering algorithm affects the number of water depth points within them, and the number of water depth points absorbed by the water depth assumption ultimately affects the final uncertainty of the water depth assumption. The selection method of this invention considers the influence of the number of clusters obtained by the DBSCAN clustering algorithm, as well as the influence of the number of water depth points in each cluster and the influence of uncertainty. It is an optimal estimate selection method that combines all three factors.

[0043] (6) Filtering of raw depth sounding data

[0044] After selecting the optimal estimate in step (5), the water depth value (z) of each node on the surface was estimated in step (1). k The estimation surface has been determined. The estimated surface is then overlaid onto the original bathymetry data to remove outliers. After the estimation in step (4), the uncertainty of the water depth value at the node is very small, and it can be considered as the water depth value closest to the true observation value obtained under existing measurement conditions. This invention uses an outlier judgment criterion based on approximate true values ​​for outlier removal, as shown in the following formula:

[0045]

[0046] in, The optimal valuation for the node finally selected in step (5); d i Let be the water depth value of the absorbable water depth point within the node's absorption domain, and i be the number of absorbable water depth points within the node's absorption domain (1≤i≤N). k is generally taken as 2 or 3.

[0047]

[0048] It is worth noting that a depth point may be absorbed by multiple nodes. When any node identifies the depth point as an outlier, the depth point is filtered out.

[0049] The beneficial effects of this invention are as follows: This invention draws upon the CUBE filtering algorithm's grid node absorbable depth point selection model and proposes an automatic outlier filtering algorithm for multibeam bathymetry, combining uncertainty and density clustering algorithms. The CARIS software has continuously improved the CUBE filtering algorithm, which is now a trade secret. The localization of multibeam bathymetry data processing software urgently requires a clear, simple, and reliable outlier cleaning algorithm for multibeam bathymetry data to process massive amounts of multibeam bathymetry data and support standard multibeam bathymetry data processing workflows. This invention can provide a choice of automatic multibeam bathymetry data cleaning algorithm for domestically produced multibeam bathymetry data processing software. Attached Figure Description

[0050] Figure 1 This is a flowchart of the method of the present invention.

[0051] Figure 2 This is a flowchart of step (3).

[0052] Figure 3 The seabed topography before and after edge beam removal; (a) before removal, (b) after removal.

[0053] Figure 4 The following is a comparison chart of the effects of Scheme 1: (a) shows the original data, (b) shows the data after manual removal of outliers, (c) shows the data after removal of outliers using the CUBE filtering algorithm, and (d) shows the effect of the present invention in removing outliers.

[0054] Figure 5(a) shows a single outlier point added manually; Figure 5(b) is a magnified view of the added area.

[0055] Figure 6(a) shows the manually added continuous outlier points; Figure 6(b) is a magnified view of the added area.

[0056] Figure 7 The results are shown for different processing methods: (a) CUBE filtering algorithm processing scheme 2 result, (b) processing scheme 2 result of the present invention, (c) CUBE filtering algorithm processing scheme 3 result, and (d) processing scheme 3 result of the present invention. Detailed Implementation

[0057] To make the problems solved by the model of this invention, the methods used, and the effects achieved clearer, the invention will be further described in detail below with reference to the accompanying drawings and experiments. It should be understood that the specific experiments described herein are merely for illustrative purposes and not for limiting the invention. Furthermore, it should be noted that, for ease of description, only the parts relevant to the invention are shown in the accompanying drawings, not all of them.

[0058] The experimental data comes from a performance comparison project of domestically produced multibeam echo sounding systems. The data is multibeam measurement data of seabed topography in the Jintang sea area of ​​Zhoushan, Zhejiang Province in 2020. The iBeam8120 domestic multibeam echo sounding system was used. This system operates at a frequency of 400kHz, with 256 beams per ping, an actual operating angle of 140°, and a depth resolution of 6mm. Installation deviation corrections were made during data acquisition, and corrections for sound velocity, tides, and positioning accuracy were performed later. The experimental data has been freed from the influence of systematic errors.

[0059] Data preparation consists of three steps. First, the horizontal and vertical uncertainties of the sounding points are calculated using the HGM model in CARIS software. Second, considering that edge beams contain a large number of clustered gross errors, leading to the loss of real seabed information and hindering the absorption of depth points by the CUBE surface nodes to estimate the actual seabed conditions, based on experience in multibeam bathymetry, a beam opening angle three times the measured depth value is selected as the experimental data, and edge beams outside the beam opening angle are removed. The multibeam ping opening angle is 140°, the coverage width is approximately 305 meters, and the average water depth is about 80 meters. Proportionally, approximately 200 sounding points in the center are selected as valid data for verification. The original data contains 50,944 sounding points; after removing edge beams, there are 39,599 sounding points. Third, a manual removal method is used to eliminate gross errors in the experimental data, preparing data for simulation experiments. Before manual removal of gross errors, there were 39,599 sounding points; after manual removal of gross errors, there were 39,392 sounding points.

[0060] Based on the six detailed steps of the multi-beam automatic outlier filtering method combining uncertainty and density clustering algorithms described in the technical solution, three implementation schemes are designed to demonstrate the effectiveness of the method and verify its adaptability in handling outliers of different shapes. An approximately square region is selected from the manually removed outlier multi-beam data as the simulation example data. The selected region contains 15 Ping data points, each Ping containing 18 sounding points, for a total of 270 sounding points. The region has a side length of approximately 15m and an area of ​​approximately 225㎡.

[0061] Option 1: Process the raw depth sounding data.

[0062] Option 2: Randomly generate 10 Gaussian noise points with a mean of 80m and a standard deviation of 3m, and randomly add them to the selected area. Any single noise point that is significantly raised above the surrounding seabed topography should be identified as a gross error and removed.

[0063] Option 3: Randomly generate 20 Gaussian noise points with a mean of 80m and add them to two consecutive Pings in the region, adding 10 noise points to each Ping. The added consecutive noise points that are significantly raised above the surrounding seabed topography should be identified as gross errors and removed.

[0064] The experimental results show that: Figure 4 As shown in (a), (b), (c), and (d), the filtering effects of the manual removal method and the method described in this invention are relatively ideal. However, the depth measurement data processed by the CUBE filtering algorithm still contains a few outliers. Comparison and analysis with the original measured data revealed that the remaining outliers in the CUBE filtering algorithm mainly appear in areas where consecutive outliers occur in the original depth measurement data. Figure 7In (a), (b), (c), and (d), both the CUBE filtering algorithm and the method described in this invention can identify and remove single outliers. However, when faced with consecutive outliers, the CUBE filtering algorithm fails to completely remove them, while the method described in this invention can completely remove consecutive outliers, demonstrating its excellent ability to clean up measurement outliers of different forms.

[0065] Finally, it should be noted that the above experiments are only used to illustrate the method of the present invention, and not to limit it; although the present invention has been described in detail, those skilled in the art should understand that modifications to the foregoing method, or equivalent substitutions of some or all of the method features, do not mean that the essence of the corresponding method deviates from the scope of the method of the present invention.

Claims

1. A method of multi-beam outlier automatic filtering combined with uncertainty and density clustering algorithm, characterized in that, The method comprises the following steps: (1) Estimation of surface construction The estimated surface is a gridded surface independent of the original sounding data (x, y, z); within the range of the original sounding data plane coordinates (x, y), a grid spacing is set, the grid intersection point is called a node (x k ,y k ,z k ), where k is the number of nodes; the node plane coordinates (x k ,y k ) are determined by gridding, and the water depth value z k of the node is estimated from the water depth values of its adjacent water depth points; when constructing the estimated surface, the absorbable water depth points of the node should be calculated, and the specific process is as follows: First, the node absorption domain Rcap is calculated using formula (1): Wherein, μcapdiaca is the node capture distance proportional coefficient; μcapdismin is the minimum node absorption domain; Gendepth is the average water depth in the circle with the node as the center and three times the grid as the radius; Then, the capture radius Rinflu of the water depth point is calculated using formula (2): Wherein, Grid is the grid spacing; TVUmax is the limit vertical uncertainty allowed by different levels of hydrographic survey in IHO S-44-2008 standard; TVU is the vertical uncertainty of the water depth point, and μdistexp is the proportional coefficient of the uncertainty with the increase of the distance; Only when the water depth point is located in the node absorption domain and the node is also located in the capture radius of the water depth point, the water depth point can be absorbed by the node; (2) Water depth value clustering After step (1), the water depth points available for absorption in the node absorption domain have been determined, and the node and the water depth point have completed the association in the two-dimensional plane; at this time, the node will obtain many water depth values, and the DBSCAN clustering algorithm is used to perform one-dimensional clustering on the water depth values of the water depth points available for node absorption; (3) Calculation of centroid point After step (2) is completed, the water depth values of the water depth points in the node absorption domain are clustered; if the number of classes is greater than or equal to 1, the node is empty and should be removed; the first water depth value entering the filter represents the characteristics of the water depth values of the class, which is called the water depth value of the centroid point; the water depth value of each class that has the smallest difference with the average water depth value of the class is taken as the water depth value of the centroid point in the class; (4) Estimation of water depth value and uncertainty In step (3), the water depth value of the centroid point is determined, and the Kalman filter is enabled to estimate the water depth value and the uncertainty; the estimated water depth value representing the possible water depth value of the node is called the water depth hypothesis; the Kalman filter estimation formula is: One-step prediction uncertainty: Water depth prediction: z j [n|n-1] = z j [n-1|n-1] (5) Kalman filter gain K j [n]: Difference between the water depth value and the predicted value v j [n]: v j [n] = d j [n] - z j [n | n - 1] (7) Node water depth value update: z j [n|n] = z j [n|n-1] + K j [n] x v j [n] (8) Node uncertainty update: wherein σ j [n-1|n-1] and z j [n-1|n-1] are the uncertainty and water depth value of the jth node after updating by Kalman filter at time n-1, respectively; σ j [n|n-1] and z j [n|n-1] are the uncertainty and water depth value of the jth node predicted at time n from time n-1, respectively; d j [n] is the water depth value entering the Kalman filter at time n; σ j [n] is the uncertainty of d j [n]; σ j [n|n] is the uncertainty of the jth node after updating by Kalman filter at time n. After the water depth values available for absorption around the node enter the Kalman filter in turn and end, this step ends; (5) Optimal estimate selection In step (4), the Kalman filter is used to estimate the possible uncertainty and water depth value of the node, and the water depth hypothesis with the smallest uncertainty is selected as the optimal estimate of the node; (6) Filtering of original sounding data Step (5) is the optimal value selection, step (1) in the estimated surface on each node depth value z k It has been determined that the estimated surface is covered on the original sounding data, and the outliers are removed; the outliers are removed using the approximate true value-based outlier judgment criterion, and the formula is as follows: wherein, is the optimal estimate finally selected by the node in step (5); d i is the water depth value of the water depth point that can be absorbed in the absorption domain of the node, i is the number of water depth points that can be absorbed in the absorption domain of the node, 1≤i≤N; k is 2 or 3; A water depth point can be absorbed by multiple nodes, and if any node identifies the water depth point as an abnormal value, the water depth point is filtered out.

2. The method of claim 1, wherein, In step (1), the limit vertical uncertainty allowed by different levels of hydrographic survey is calculated using formula: Wherein, a and b are constants, and d is the water depth value of the water depth point.

3. The method according to claim 1 or 2, characterized in that, In step (2), the DBSCAN clustering algorithm has two key parameters Eps and MinPts to describe the sample distribution tightness of the neighborhood; wherein, the value range of the parameter Eps is set to (0, 1]; and the parameter MinPts is fixed at 1, which ensures that nodes with a small number of absorbed water depth points can be normally clustered without being identified as noise points.

Images