Comprehensive Detection of Gas Plumes from Multibeam Water Column Images with Minimisation of Noise Interferences

Multibeam echosounder systems (MBES) can record backscatter strengths of gas plumes in the water column (WC) images that may be an indicator of possible occurrence of gas at certain depths. Manual or automatic detection is generally adopted in finding gas plumes, but frequently results in low efficiency and high false detection rates because of WC images that are polluted by noise. To improve the efficiency and reliability of the detection, a comprehensive detection method is proposed in this paper. In the proposed method, the characteristics of WC background noise are first analyzed and given. Then, the mean standard deviation threshold segmentations are respectively used for the denoising of time-angle and depth-angle images, an intersection operation is performed for the two segmented images to further weaken noise in the WC data, and the gas plumes in the WC data are detected from the intersection image by the morphological constraint. The proposed method was tested by conducting shallow-water and deepwater experiments. In these experiments, the detections were conducted automatically and higher correct detection rates than the traditional methods were achieved. The performance of the proposed method is analyzed and discussed.


Introduction
New generation multibeam echosounder systems (MBES), such as the Kongsberg EM710/EM302 [1][2][3] and Teledyne Reson 7125/8125 [4,5] can collect water column (WC) data. WC data include the full acoustic information from the MBES transducer to the seafloor [1] and provide an effective way to find underwater objects such as fish schools [4,[6][7][8], wrecks [9][10][11], eelgrass [12], gas plumes [13,14], and internal ocean waves [15,16]. Gas plumes may be an indicator of the presence of gas hydrate. Methane hydrate is a clear energy source with potential importance as an energy source with environmental benefits [17]. The search for methane hydrate by detecting gas plumes from WC data is a new and low-cost alternative method to the traditional seismic exploration [2,18].
Gas plumes are formed by gas leakage [19][20][21]. Gas leakage may be modulated by temperature and pressure. The driving force of leakage is the production of methane at depth and resulting overpressure, possibly linked to fluid migration, and the buoyancy of the gas. Gas hydrate dissolution is another driver that is highly susceptible to temperature increases. Methane gas bubbles usually dissolve completely during their rise. Whether methane gas bubbles reach the sea surface depends on the water depth and initial bubble size. Many bubbles dissolve after rising a few dozens of meters [22].
combining BSs of all sectors in the ping. The above process can be conducted by commercial software for MBES data processing such as CARIS HIPS or self-developed software.
Typically, WC data are shown as the along-or across-track images. Figure 1D depicts the along-track image of a beam; the along-track image is used to rapidly find targets when the beam goes through the targets. Furthermore, accurately obtaining the shapes and sizes of targets by using the along-track image is infeasible because of their only being scanned by a beam [1].
Time-angle (T-A) and depth-across track (D-T) images are two classic across-track images [1] and are often used for detecting targets and determining the shapes of the targets. In a T-A image ( Figure 1A), the abscissa and the ordinate denote beam serial and sampling point numbers. In a D-T image ( Figure 1B), the two axes indicate swath width and vertical distance from the transducer to its nadir. The T-A and D-T images permit that the location and shape of a target can be determined accurately by the D-T image but not the T-A image. A T-A image has smaller size than the corresponding D-T image. Therefore, the T-A images are typically used for rapidly detecting possible targets, and the D-T images are adopted in recognizing targets and determining their shapes [1]. The depth-angle (D-A) image ( Figure 1C), which defines the abscissa and the ordinate as a beam serial number (corresponding to beam angle) and depth, is also adopt to reduce the size of the D-T image and improve the efficiency of the threshold segmentation. In Figure 1, minimum slant range (MSR) is defined as the minimum radial distance between MBES transducer and seafloor. Because of the interaction of side lobe and seabed echoes, the water column data beyond MSR are seriously polluted. Therefore, in this paper, the detection of gas plume is performed within the MSR. Typically, WC data are shown as the along-or across-track images. Figure 1D depicts the alongtrack image of a beam; the along-track image is used to rapidly find targets when the beam goes through the targets. Furthermore, accurately obtaining the shapes and sizes of targets by using the along-track image is infeasible because of their only being scanned by a beam [1].
Time-angle (T-A) and depth-across track (D-T) images are two classic across-track images [1] and are often used for detecting targets and determining the shapes of the targets. In a T-A image ( Figure 1A), the abscissa and the ordinate denote beam serial and sampling point numbers. In a D-T image ( Figure 1B), the two axes indicate swath width and vertical distance from the transducer to its nadir. The T-A and D-T images permit that the location and shape of a target can be determined accurately by the D-T image but not the T-A image. A T-A image has smaller size than the corresponding D-T image. Therefore, the T-A images are typically used for rapidly detecting possible targets, and the D-T images are adopted in recognizing targets and determining their shapes [1]. The depth-angle (D-A) image ( Figure 1C), which defines the abscissa and the ordinate as a beam serial number (corresponding to beam angle) and depth, is also adopt to reduce the size of the D-T image and improve the efficiency of the threshold segmentation. In Figure 1, minimum slant range (MSR) is defined as the minimum radial distance between MBES transducer and seafloor. Because of the interaction of side lobe and seabed echoes, the water column data beyond MSR are seriously polluted. Therefore, in this paper, the detection of gas plume is performed within the MSR.

Characteristics of Noise in the WC Image
A MBES transducer recorded echoes including gas plumes, noise caused by the WC environment, interferences from beam side lobes, third-party sonars, and vessels, as well as other targets [1]. Therefore, the interpretation of the image in term of gas plumes is prone to errors. The distributions of noise in the WC image are depicted in Figure 2. The distributions are as follows: 1. Sea surface noise caused by vessel and MBES transducer affects the BSs of the top of a ping. 2. Layering noise is a phenomenon in the WC image and is formed by the echoes from the non-gas plume targets, such as marine organisms, inanimate matter, and the inhomogeneous sea structure [37]. Generally, the noise distributes horizontally in different water layers and exhibits a strong BS in deep layers [26]. 3. A side lobe effect is produced by the interference of the side lobes and target echoes. The beam pattern of a conventional Mills Cross MBES is defined by a main lobe and side lobes. The scattering signal of targets can be extended given the presence of the side lobes, and seabed echoes contaminate the WC data at any slant range beyond the closest distance of approach to the seabed [1]. The side lobe effect near a seabed is called the minimum slant range (MSR) effect. The MSR effect implies that the gas plume can be detected efficiently only within the MSR.

Characteristics of Noise in the WC Image
A MBES transducer recorded echoes including gas plumes, noise caused by the WC environment, interferences from beam side lobes, third-party sonars, and vessels, as well as other targets [1]. Therefore, the interpretation of the image in term of gas plumes is prone to errors. The distributions of noise in the WC image are depicted in Figure 2. The distributions are as follows:

1.
Sea surface noise caused by vessel and MBES transducer affects the BSs of the top of a ping.

2.
Layering noise is a phenomenon in the WC image and is formed by the echoes from the non-gas plume targets, such as marine organisms, inanimate matter, and the inhomogeneous sea structure [37]. Generally, the noise distributes horizontally in different water layers and exhibits a strong BS in deep layers [26].

3.
A side lobe effect is produced by the interference of the side lobes and target echoes. The beam pattern of a conventional Mills Cross MBES is defined by a main lobe and side lobes. The scattering signal of targets can be extended given the presence of the side lobes, and seabed echoes contaminate the WC data at any slant range beyond the closest distance of approach to Sector configuration leads to the difference of noise distributions in a WC image. Because of low emission frequency and long measurement distance, noise in the outside sectors is higher than that in the middle ones.

5.
Background noise, which is formed by the echoes from turbulence, particles scattering, non-buoyant microbubbles and so on in water mass [38,39], has the characteristic of discrete distribution in the WC image.
These effects appear with the gas plumes in the WC images and challenge the traditional threshold segmentation methods, such as the 2D Otsu and iterative threshold segmentation methods mentioned in Section 1 in the gas plume detection. Furthermore, the distribution differences of noise in the different types of WC images demonstrated in Figure 1 imply the necessity of removing noise by combining the T-A, D-A, and D-T images. 4. Sector configuration leads to the difference of noise distributions in a WC image. Because of low emission frequency and long measurement distance, noise in the outside sectors is higher than that in the middle ones. 5. Background noise, which is formed by the echoes from turbulence, particles scattering, nonbuoyant microbubbles and so on in water mass [38,39], has the characteristic of discrete distribution in the WC image.
These effects appear with the gas plumes in the WC images and challenge the traditional threshold segmentation methods, such as the 2D Otsu and iterative threshold segmentation methods mentioned in Section 1 in the gas plume detection. Furthermore, the distribution differences of noise in the different types of WC images demonstrated in Figure 1 imply the necessity of removing noise by combining the T-A, D-A, and D-T images.

Detection of Gas Plumes
Gas plume detection consists of separating the gas plume from the noise in a WC image. Figure 2 illustrates that noise mostly distribute along the depth and horizontal directions and have the BS differences in various sectors in a WC image. Therefore, a detection method based on the characteristics of gas plume and noise is presented in this section. For a sector, three main steps are involved in the method: Step 1: Threshold segmentation for the T-A and D-A images The average beam intensity curves of the T-A and D-A images are first calculated, and then the threshold segmentations are conducted for the two images through the mean standard deviation threshold segmentation method to weaken the noise.
Step 2: Intersection calculation of the two segmented images The segmented D-A image is converted into the T-A image. Then an intersection operation is performed for the two segmented T-A images to maintain the common targets, remove the noncommon residual noise, and obtain an intersection image.
Step 3: Morphological constraint of the gas plume The T-A intersection image is converted into the D-T image to obtain nearly real shapes and positions of targets in the WC image, and then the morphological characteristics of the gas plumes are used as the thresholds to detect the gas plumes in the WC images.
The above process is depicted in Figure 3. The same process is used for the other sectors, and the target detection in a ping WC image can be performed.

Detection of Gas Plumes
Gas plume detection consists of separating the gas plume from the noise in a WC image. Figure 2 illustrates that noise mostly distribute along the depth and horizontal directions and have the BS differences in various sectors in a WC image. Therefore, a detection method based on the characteristics of gas plume and noise is presented in this section. For a sector, three main steps are involved in the method: Step 1: Threshold segmentation for the T-A and D-A images The average beam intensity curves of the T-A and D-A images are first calculated, and then the threshold segmentations are conducted for the two images through the mean standard deviation threshold segmentation method to weaken the noise.
Step 2: Intersection calculation of the two segmented images The segmented D-A image is converted into the T-A image. Then an intersection operation is performed for the two segmented T-A images to maintain the common targets, remove the non-common residual noise, and obtain an intersection image.
Step 3: Morphological constraint of the gas plume The T-A intersection image is converted into the D-T image to obtain nearly real shapes and positions of targets in the WC image, and then the morphological characteristics of the gas plumes are used as the thresholds to detect the gas plumes in the WC images.
The above process is depicted in Figure 3. The same process is used for the other sectors, and the target detection in a ping WC image can be performed.

Mean Standard Deviation Threshold Segmentation Algorithm
Threshold segmentation has a central position in image segmentations given the simplicity and efficiency of this process [40], in which the key is selecting a threshold. The mean standard deviation threshold segmentation algorithm is adopted in this research considering the characteristics of the noise and gas plumes in the WC images. The threshold, μ + kσ, is determined in a combined manner by the mean μ and standard deviation σ of the BSs of a sector. The k is a constant within 1-3. The determination of k is discussed in Section 5.1. The identification of the interesting targets from the noise and interferences by using μ + kσ follows the principle below: where f(x, y) is the WC data or image, and g(x, y) is the segmented binary image.

Threshold Segmentation for T-A Images
The threshold segmentation must be handled independently in different sectors given the differences in transmission frequencies differences. By using a sector as an example, if the matrix of the BSs of sampling points is M (a × b), then the segmentation process is depicted as follows: 1. The mean is calculated by averaging the BSs of each row and a column vector with a rows is obtained. 2. The mean BS is subtracted from the raw BSs of each row, and a new matrix M1 (a × b) is acquired. Figure 4 implies that the raw BSs do not obey the normal distribution, whereas M1 does. Therefore, the process is called normalization of the T-A image. 3. The threshold μ+kσ is obtained by calculating the mean and the standard deviation σ of M1. 4. M1 is smoothened to highlight the targets in the T-A image, and M2 (a × b) is obtained. In the smoothing process, the convolution smoothing is adopted for each column of M1. If a column vector is d, and the weighted window is w, then we obtain the following equation: where d is the column vector of M1, N is the length of d, and w is the weighted window function with window size M. The head and end of the sequence after the convolution should be removed to ensure N invariability before and after convolution. A Hanning function is used to eliminate

Mean Standard Deviation Threshold Segmentation Algorithm
Threshold segmentation has a central position in image segmentations given the simplicity and efficiency of this process [40], in which the key is selecting a threshold. The mean standard deviation threshold segmentation algorithm is adopted in this research considering the characteristics of the noise and gas plumes in the WC images. The threshold, µ + kσ, is determined in a combined manner by the mean µ and standard deviation σ of the BSs of a sector. The k is a constant within 1-3. The determination of k is discussed in Section 5.1. The identification of the interesting targets from the noise and interferences by using µ + kσ follows the principle below: where f (x, y) is the WC data or image, and g(x, y) is the segmented binary image.

Threshold Segmentation for T-A Images
The threshold segmentation must be handled independently in different sectors given the differences in transmission frequencies differences. By using a sector as an example, if the matrix of the BSs of sampling points is M (a × b), then the segmentation process is depicted as follows: The mean is calculated by averaging the BSs of each row and a column vector with a rows is obtained.

2.
The mean BS is subtracted from the raw BSs of each row, and a new matrix M 1 (a × b) is acquired. Figure 4 implies that the raw BSs do not obey the normal distribution, whereas M 1 does. Therefore, the process is called normalization of the T-A image. 3.
The threshold µ + kσ is obtained by calculating the mean and the standard deviation σ of M 1 .

4.
M 1 is smoothened to highlight the targets in the T-A image, and M 2 (a × b) is obtained. In the smoothing process, the convolution smoothing is adopted for each column of M 1. If a column vector is d, and the weighted window is w, then we obtain the following equation: where d is the column vector of M 1 , N is the length of d, and w is the weighted window function with window size M. The head and end of the sequence after the convolution should be removed to ensure N invariability before and after convolution. A Hanning function is used to eliminate the high-frequency interference from a non-periodic continuous signal and is thus adopted as the weighted window w:

5.
Target or noise is diagnosed by comparing the value of each sampling point in M 2 with µ + kσ.
If the value is more than µ + kσ, then the sampling point is determined as a target point and retained. Otherwise, the sampling point is determined as noise and removed.
All sectors are processed by using these five steps until the segmentation is completed.
Sensors 2017, 17, 2755 6 of 17 the high-frequency interference from a non-periodic continuous signal and is thus adopted as the weighted window w: 5. Target or noise is diagnosed by comparing the value of each sampling point in M2 with μ + kσ. If the value is more than μ + kσ, then the sampling point is determined as a target point and retained. Otherwise, the sampling point is determined as noise and removed.
All sectors are processed by using these five steps until the segmentation is completed.  Figure 5 depicts the process. The figure illustrates the T-A image formed by the raw WC data within the MSR. Moreover, the high noise appear in the two outside sectors, the arc-shaped layering noise distributes in the different water layers, and the side lobe effects accompany the gas plume and distribute horizontally in the different water layers. In Figure 5B, the mean BS curves of the three sectors reflect the noise. The layering noise is reflected in the BS variations of the different sampling numbers, whereas the side lobe effects are shown in the black rectangles and red circle. Figure 5C indicates that the side lobe effects and background noise are mostly removed after the threshold segmentation, whereas the layering noise is weakened.

Threshold Segmentation for D-A Image
The D-A image differs from the corresponding T-A image in the displayed form of a ping fan. In the D-A image, the layering noise and side lobe effects are presented as horizontal and arc-shaped  Figure 5 depicts the process. The figure illustrates the T-A image formed by the raw WC data within the MSR. Moreover, the high noise appear in the two outside sectors, the arc-shaped layering noise distributes in the different water layers, and the side lobe effects accompany the gas plume and distribute horizontally in the different water layers. In Figure 5B, the mean BS curves of the three sectors reflect the noise. The layering noise is reflected in the BS variations of the different sampling numbers, whereas the side lobe effects are shown in the black rectangles and red circle. Figure 5C indicates that the side lobe effects and background noise are mostly removed after the threshold segmentation, whereas the layering noise is weakened. the high-frequency interference from a non-periodic continuous signal and is thus adopted as the weighted window w: 5. Target or noise is diagnosed by comparing the value of each sampling point in M2 with μ + kσ. If the value is more than μ + kσ, then the sampling point is determined as a target point and retained. Otherwise, the sampling point is determined as noise and removed.
All sectors are processed by using these five steps until the segmentation is completed.  Figure 5 depicts the process. The figure illustrates the T-A image formed by the raw WC data within the MSR. Moreover, the high noise appear in the two outside sectors, the arc-shaped layering noise distributes in the different water layers, and the side lobe effects accompany the gas plume and distribute horizontally in the different water layers. In Figure 5B, the mean BS curves of the three sectors reflect the noise. The layering noise is reflected in the BS variations of the different sampling numbers, whereas the side lobe effects are shown in the black rectangles and red circle. Figure 5C indicates that the side lobe effects and background noise are mostly removed after the threshold segmentation, whereas the layering noise is weakened.

Threshold Segmentation for D-A Image
The D-A image differs from the corresponding T-A image in the displayed form of a ping fan. In the D-A image, the layering noise and side lobe effects are presented as horizontal and arc-shaped

Threshold Segmentation for D-A Image
The D-A image differs from the corresponding T-A image in the displayed form of a ping fan. In the D-A image, the layering noise and side lobe effects are presented as horizontal and arc-shaped distributions in different depths, respectively. Therefore, the layering noise that remains partly in the segmented T-A image can be removed efficiently from the D-A image by using the threshold segmentation. If the raw BS matrix of a sector in the D-A image is N (a 1 × b 1 ), then the threshold segmentation follows similar steps as described in Section 3.2.1. All sectors are segmented by the five steps, and the segmentation is performed in the D-A image. Figure 6 demonstrates the threshold segmentation process of the D-A image. The raw ping WC data used in Figure 5 is extracted to form the D-A image ( Figure 6A). The D-A image also indicates that high noise is found in the two outside sectors. The layering noise distributes horizontally, and the side lobe effects display as the arc-shaped distributions. The distributions of the two types of noise are reversed compared with the distributions depicted in Figure 5A. The noise can also be effectively reflected in the mean BS curves ( Figure 6B). The layering noise and the background noise are mostly removed after the threshold segmentation ( Figure 6C), whereas the side lobe effects are weakened. distributions in different depths, respectively. Therefore, the layering noise that remains partly in the segmented T-A image can be removed efficiently from the D-A image by using the threshold segmentation. If the raw BS matrix of a sector in the D-A image is N (a1 × b1), then the threshold segmentation follows similar steps as described in Section 3.2.1. All sectors are segmented by the five steps, and the segmentation is performed in the D-A image. Figure 6 demonstrates the threshold segmentation process of the D-A image. The raw ping WC data used in Figure 5 is extracted to form the D-A image ( Figure 6A). The D-A image also indicates that high noise is found in the two outside sectors. The layering noise distributes horizontally, and the side lobe effects display as the arc-shaped distributions. The distributions of the two types of noise are reversed compared with the distributions depicted in Figure 5A. The noise can also be effectively reflected in the mean BS curves ( Figure 6B). The layering noise and the background noise are mostly removed after the threshold segmentation ( Figure 6C), whereas the side lobe effects are weakened.

Intersection Operation for the Two Segmented Images
The background noise can be mostly removed via the above threshold segmentation for the T-A and D-A images. The segmentation for the T-A image mainly removes the side lobe effects but remains the residual layering noise in the segmented T-A image, whereas the segmentation for the D-A image eliminates the layering noise but keeps the residual side lobe effects in the segmented D-A image. The two segmentations indicate complementarity. Therefore, an intersection operation for the segmented T-A and D-A images may be performed to eliminate most residual noise from the two segmented images. Because the segmented binary D-A image differs from the segmented binary T-A image, the segmented binary D-A image need to be converted to a new binary T-A image for the intersection operation. The intersection operation of the segmented and the new T-A images can be expressed as: where G (x, y) is the intersection image of the two T-A images, g1 (x, y) is the segmented T-A image, and g2 (x, y) is the new T-A image. 0 means noise and 1 means target. The two images exhibit the same size m × n, (i, j)(m, n).

Intersection Operation for the Two Segmented Images
The background noise can be mostly removed via the above threshold segmentation for the T-A and D-A images. The segmentation for the T-A image mainly removes the side lobe effects but remains the residual layering noise in the segmented T-A image, whereas the segmentation for the D-A image eliminates the layering noise but keeps the residual side lobe effects in the segmented D-A image. The two segmentations indicate complementarity. Therefore, an intersection operation for the segmented T-A and D-A images may be performed to eliminate most residual noise from the two segmented images. Because the segmented binary D-A image differs from the segmented binary T-A image, the segmented binary D-A image need to be converted to a new binary T-A image for the intersection operation. The intersection operation of the segmented and the new T-A images can be expressed as: where G (x, y) is the intersection image of the two T-A images, g 1 (x, y) is the segmented T-A image, and g 2 (x, y) is the new T-A image. 0 means noise and 1 means target. The two images exhibit the same size m × n, (i, ⊂j)(m, n). In the intersection operation, if a target appears at the same position of the two segmented T-A images, then this target will remain in the final intersection image; otherwise, it will be removed as a noise. Therefore, the intersection operation retains the common targets and removes the residual noise in the two images.

Gas Plume Detection Based on Morphology Characteristics
At this point, most of the noise in the WC data are eliminated, and the targets remain in a ping image after the abovementioned process. However, some noise whose BSs are high remain in the intersection image because the mean standard deviation threshold segmentation method is based on BS feature of gas plumes. Therefore, we need to use other characteristics such as positional and morphological features of gas plumes to detect and recognize gas plumes from the D-T image of intersection image: (1) Position features Gas plumes that rise from the seafloor usually dissolve during the rising process, and cannot reach the sea surface. We can set a depth threshold to reduce the detection range. Setting the threshold needs to consider depth variations of gas plumes in a detected water area and completely remove the effect of residual water surface noise. Generally, the threshold is set as more than 5 m underneath the MBES transducer.
(2) Morphological feature The flare-shape of gas plumes is different from noise [41]. We can use morphological features such as height, area, and width to distinguish gas plumes from noise. The characteristics of gas plumes mentioned in Section 1 show that the heights of gas plumes are various, and higher than those of noise. The flare-shape gas plumes often present nearly linear distribution. We can define a linear structuring element B = strel('line', len, deg), then conduct first closing and second opening operations [40] to image A, namely the D-T image of intersection image in MATLAB: where • and • are opening and closing operations and ⊕ and denote the dilation and erosion operations, respectively.
The strel class in MATLAB is used to define morphological structuring element, strel ('line' len, deg) creates a linear structuring element, deg specifies the angle (in degrees) of the line, and len is approximately the distance between the centers of the structuring element members at opposite ends of the line. The deg is set to 90, and len is set by the minimum detection height that refers the detection requirement and resolution of the WC data in the study.
The area of a flare-shaped gas plume is generally greater than that of noise. Area threshold, namely the minimum detection area (height × width), are set by referring to the detection requirement and resolution of the WC data. The minimum detection height and width can be calculated quantitatively following the description in Section 5.4.

Experiment I: Detections of Gas Plumes from Shallow Water EM710 WC Data
A shallow-water MBES measurement was performed in the water area with 0.5 km × 1.5 km, roughly 160 m of water depth, and dense gas plumes in June 2012, in the South China Sea to verify the proposed method. In the MBES measurement, a Kongsberg EM 710 device was adopted. A total of 120 • of the swath coverage, three sectors, 128 beams in a ping fan, 6944.44 Hz of sampling frequency were used in a beam measurement. The dual-swath mode was set in EM 710 WC data acquisition. The dual-swath mode sets different transmission frequencies for different sectors and adjacent pings. For an odd ping, the frequencies of the middle and two outside sectors are 81 and 73 kHz, correspondingly, whereas are 97 and 89 kHz, respectively, for the even ping. A 1.6 km-length surveying line with roughly 500 m swath width, which covers the water area, was designed and 1444 pings of WC data were recorded by the Kongsberg Seafloor Information System (SIS) data acquisition software developed by Kongsberg Maritime Company at about 2.2 m of ping interval. By decoding the binary WC data files with the suffix *.wcd and performing data processing by the process mentioned in Section 2.1, WC images of recorded pings are formed by self-developed software. Based on these parameters of EM 710, the size of the detected minimum target is defined as 0.2 m in height, 2.6 m in width, and 2 m 2 in area, correspondingly by referring to Section 5.4 for later detections of gas plumes.
The 112th ping with a gas plume is extracted for the detection via the proposed method, and the k in the threshold µ + kσ is set as 1.5. The process is depicted in Figure 7. The results are presented in Figure 7F. In this figure, slight amounts of noise remain. This phenomenon is caused by the threshold segmentation algorithm. Most noise is removed because of its BS that is less than the given threshold, whereas the retained noise demonstrate higher BSs than the threshold. The morphological constraint (detailed in Section 3.4) effectively compensates the limitations of the algorithm and detects the gas plumes, as depicted in Figure 7H. The 112th ping with a gas plume is extracted for the detection via the proposed method, and the k in the threshold μ + kσ is set as 1.5. The process is depicted in Figure 7. The results are presented in Figure 7F. In this figure, slight amounts of noise remain. This phenomenon is caused by the threshold segmentation algorithm. Most noise is removed because of its BS that is less than the given threshold, whereas the retained noise demonstrate higher BSs than the threshold. The morphological constraint (detailed in Section 3.4) effectively compensates the limitations of the algorithm and detects the gas plumes, as depicted in Figure 7H. All the 1444 pings in the surveying line are processed by using the same steps as mentioned above. The total correct detection rate γ and the gas plumes correct detection rate γgas_plumes are calculated to quantitatively assess the performance of the proposed method and defined as [42]: All the 1444 pings in the surveying line are processed by using the same steps as mentioned above. The total correct detection rate γ and the gas plumes correct detection rate γ gas_plumes are calculated to quantitatively assess the performance of the proposed method and defined as [42]: γ gas_plumes = (N gas_plumes /N manual_gas_plumes ) × 100% (7) where N all is the ping number of the entire survey line, N manual_gas_plumes is the number of pings that contain gas plumes and is obtained by manual inspection, N gas_plumes is the number of pings with gas plumes and is detected by the proposed method, and N c is the correct detection number of pings, including gas plumes and null pings by the proposed method. In the experiment, the proposed method fulfills the automatic detection and achieves γ ≈ 86% and γ gas_plumes ≈ 86%, thereby indicating that the proposed method is efficient in detecting gas plumes from the seriously polluted shallow-water WC images. Based on the detection, the distributions of the gas plumes in the surveying area are illustrated in Figure 8. The 3D shape of each gas plume is also drawn by the 3D coordinates of all the echoes from the detected gas plumes. The gas plumes distribute densely in the area and have heights that vary from several meters to approximately 50 m. where Nall is the ping number of the entire survey line, Nmanual_gas_plumes is the number of pings that contain gas plumes and is obtained by manual inspection, Ngas_plumes is the number of pings with gas plumes and is detected by the proposed method, and Nc is the correct detection number of pings, including gas plumes and null pings by the proposed method. In the experiment, the proposed method fulfills the automatic detection and achieves γ ≈ 86% and γgas_plumes ≈ 86%, thereby indicating that the proposed method is efficient in detecting gas plumes from the seriously polluted shallowwater WC images. Based on the detection, the distributions of the gas plumes in the surveying area are illustrated in Figure 8. The 3D shape of each gas plume is also drawn by the 3D coordinates of all the echoes from the detected gas plumes. The gas plumes distribute densely in the area and have heights that vary from several meters to approximately 50 m.

Experiment II: Detection of Gas Plume from Deepwater EM122 WC Data
A deep-water MBES measurement was performed in the water area with 3.0 km × 40 km, roughly 1500 m of water depth and sparse gas plumes in March 2016, in the South China Sea to further verify the proposed method. In the measurement, a Kongsberg EM 122 unit was adopted. The swath coverage of 90°, six sectors, 288 beams in a ping fan, 67.34 Hz of sampling frequency in a beam measurement, dual-swath mode are set in EM 122 WC data acquisition. In the dual-swath mode, the frequencies of sectors 1-6 are 11, 12.8, 12.2, 12.95, 12.35, and 11.75 kHz correspondingly for an odd ping, and are 11.3, 13.1, 12.5, 13.25, 12.65, and 12.05 kHz for the even ping. A 40 km-length surveying line with roughly 3000 m swath width, covering the water area, was designed and 5960 pings of WC data were recorded by Kongsberg SIS at about 13 m of ping interval. By a similar process as used in the shallow-water MBES WC data, the BSs, gray levels and positions of echoes in a ping are obtained and WC images of recorded pings are formed. The height, width, and area of the detected minimum target are defined as 11.3 m, 10.5 m, and 230 m 2 , respectively, by referring these performance parameters of EM 122 for later detections of gas plumes based on morphological constraints.
The 2513th ping with a gas plume is extracted for the detection through the proposed method and the k in the threshold μ + kσ is set as 1.2. A process similar to the one shown in Figure 7 is depicted in Figure 9. The layering phenomenon of noise becomes more significant, and the intensities of noise are higher than those in shallow-water images. The same process is adopted in the detection; the results ( Figure 9H) display the detected gas plume clearly and accurately.

Experiment II: Detection of Gas Plume from Deepwater EM122 WC Data
A deep-water MBES measurement was performed in the water area with 3.0 km × 40 km, roughly 1500 m of water depth and sparse gas plumes in March 2016, in the South China Sea to further verify the proposed method. In the measurement, a Kongsberg EM 122 unit was adopted. The swath coverage of 90 • , six sectors, 288 beams in a ping fan, 67.34 Hz of sampling frequency in a beam measurement, dual-swath mode are set in EM 122 WC data acquisition. In the dual-swath mode, the frequencies of sectors 1-6 are 11, 12.8, 12.2, 12.95, 12.35, and 11.75 kHz correspondingly for an odd ping, and are 11.3, 13.1, 12.5, 13.25, 12.65, and 12.05 kHz for the even ping. A 40 km-length surveying line with roughly 3000 m swath width, covering the water area, was designed and 5960 pings of WC data were recorded by Kongsberg SIS at about 13 m of ping interval. By a similar process as used in the shallow-water MBES WC data, the BSs, gray levels and positions of echoes in a ping are obtained and WC images of recorded pings are formed. The height, width, and area of the detected minimum target are defined as 11.3 m, 10.5 m, and 230 m 2 , respectively, by referring these performance parameters of EM 122 for later detections of gas plumes based on morphological constraints.
The 2513th ping with a gas plume is extracted for the detection through the proposed method and the k in the threshold µ + kσ is set as 1.2. A process similar to the one shown in Figure 7 is depicted in Figure 9. The layering phenomenon of noise becomes more significant, and the intensities of noise are higher than those in shallow-water images. The same process is adopted in the detection; the results ( Figure 9H) display the detected gas plume clearly and accurately. We achieved γ ≈ 99%, γgas_plumes ≈ 80% when detect all 5960 pings ( Figure 10). γ is higher than γgas_plumes because of only 30 pings including gas plumes, and γgas_plumes is lower than that in Experiment I because of more complicated noise and interferences. The threshold segment method and detection process are similar to those used in Experiment I. The result shows that the proposed method is efficient in the shallow-water and deepwater detections of gas plumes.

Determination of k
k has a significant effect on the threshold segmentation. The confidence coefficient of a detected gas plume is high when the k is higher, but missing small gas plumes is possible. Conversely, some noise may be diagnosed mistakenly as gas plumes when k is too small. k is a constant that varies from 1 to 3. The noise is mainly produced by the marine organism, suspended solids, and so on. The distributions of marine organisms vary with marine temperature profile, demonstrate dense in shallow water and sparse in deep water ( Figure 11). The suspended solids, influenced by buoyance, mainly distribute in shallow water. The above results in that the noise in deepwater WC images differ from that in shallow water WC images. This phenomenon has been proven by the above experiments We achieved γ ≈ 99%, γ gas_plumes ≈ 80% when detect all 5960 pings ( Figure 10). γ is higher than γ gas_plumes because of only 30 pings including gas plumes, and γ gas_plumes is lower than that in Experiment I because of more complicated noise and interferences. The threshold segment method and detection process are similar to those used in Experiment I. The result shows that the proposed method is efficient in the shallow-water and deepwater detections of gas plumes. We achieved γ ≈ 99%, γgas_plumes ≈ 80% when detect all 5960 pings ( Figure 10). γ is higher than γgas_plumes because of only 30 pings including gas plumes, and γgas_plumes is lower than that in Experiment I because of more complicated noise and interferences. The threshold segment method and detection process are similar to those used in Experiment I. The result shows that the proposed method is efficient in the shallow-water and deepwater detections of gas plumes.

Determination of k
k has a significant effect on the threshold segmentation. The confidence coefficient of a detected gas plume is high when the k is higher, but missing small gas plumes is possible. Conversely, some noise may be diagnosed mistakenly as gas plumes when k is too small. k is a constant that varies from 1 to 3. The noise is mainly produced by the marine organism, suspended solids, and so on. The distributions of marine organisms vary with marine temperature profile, demonstrate dense in shallow water and sparse in deep water ( Figure 11). The suspended solids, influenced by buoyance, mainly distribute in shallow water. The above results in that the noise in deepwater WC images differ from that in shallow water WC images. This phenomenon has been proven by the above experiments

Determination of k
k has a significant effect on the threshold segmentation. The confidence coefficient of a detected gas plume is high when the k is higher, but missing small gas plumes is possible. Conversely, some noise may be diagnosed mistakenly as gas plumes when k is too small. k is a constant that varies from 1 to 3. The noise is mainly produced by the marine organism, suspended solids, and so on. The distributions of marine organisms vary with marine temperature profile, demonstrate dense in shallow water and sparse in deep water ( Figure 11). The suspended solids, influenced by buoyance, mainly distribute in shallow water. The above results in that the noise in deepwater WC images differ from that in shallow water WC images. This phenomenon has been proven by the above experiments and implies that k should vary with water depths. Many segmentation experiments for detecting gas plumes in shallow-water and deepwater WC images are conducted under different k to obtain an appropriate k. Figure 11 shows the three classic results at k as 1, 1.5, and 2 in shallow-water experiments and shows that k = 1.5 is the best. The result at k as 1 includes considerable noise because of the calculated threshold being too small, whereas the excessively large threshold (k as 2) leads to the segmented gas plume scattering. Therefore, k is set as 1.5 in the shallow-water detections. Similar experiments are performed for the deepwater WC images, and 1.2 is selected as the optimum k. and implies that k should vary with water depths. Many segmentation experiments for detecting gas plumes in shallow-water and deepwater WC images are conducted under different k to obtain an appropriate k. Figure 11 shows the three classic results at k as 1, 1.5, and 2 in shallow-water experiments and shows that k = 1.5 is the best. The result at k as 1 includes considerable noise because of the calculated threshold being too small, whereas the excessively large threshold (k as 2) leads to the segmented gas plume scattering. Therefore, k is set as 1.5 in the shallow-water detections. Similar experiments are performed for the deepwater WC images, and 1.2 is selected as the optimum k. Figure 11. Selection of k in the threshold segmentation.

Threshold Segmentation Methods
Four representative ping T-A images in Experiment I, A-D, are extracted and segmented via the proposed, whole-fan μ + kσ threshold segmentation, 2D Otsu, and iterative methods. The results are shown in Figure 12. The four images are seriously polluted by noise. The highlighted linear-shaped noise caused by the side lobe effects and marked with the rectangles appears clearly in A and B but unclearly in C and D. The arc-shaped layering noise marked with arc lines and induced by the inhomogeneous sea structure and distributions of marine organisms are displayed in each ping image. The segmentation results show that the 2D Otsu and iterative methods nearly cannot distinguish the gas plumes from their background. The whole-fan μ+kσ method has better performance than the two methods, but the results remain unsatisfactory. The proposed method achieved the best segmentations because the 2D Otsu and iterative methods ignore the characteristics of the noise and handle these T-A images as normal optical image. Thus, the background noise is only weakened, whereas considerable noise remains in the segmented images. The discrete arcshaped layering noises are still in the segmented T-A image, although the whole-fan μ+kσ segmentation method segments the images along the beam directions and weakens the linear-shaped side lobe effects. The reason is that the method disregards the frequency and noise differences among adjacent sectors, and the method performs the threshold segmentation in the whole fan with an invariable threshold. However, the proposed method nearly eliminates all the noise and clearly presents the gas plumes in the last binary T-A images.
The gas plumes in Experiments I and II are detected by using the four methods and combining the intersection operation and morphological constraint. The results are listed in Table 1. The proposed method performs better than the other methods do.

Threshold Segmentation Methods
Four representative ping T-A images in Experiment I, A-D, are extracted and segmented via the proposed, whole-fan µ + kσ threshold segmentation, 2D Otsu, and iterative methods. The results are shown in Figure 12. The four images are seriously polluted by noise. The highlighted linear-shaped noise caused by the side lobe effects and marked with the rectangles appears clearly in A and B but unclearly in C and D. The arc-shaped layering noise marked with arc lines and induced by the inhomogeneous sea structure and distributions of marine organisms are displayed in each ping image. The segmentation results show that the 2D Otsu and iterative methods nearly cannot distinguish the gas plumes from their background. The whole-fan µ + kσ method has better performance than the two methods, but the results remain unsatisfactory. The proposed method achieved the best segmentations because the 2D Otsu and iterative methods ignore the characteristics of the noise and handle these T-A images as normal optical image. Thus, the background noise is only weakened, whereas considerable noise remains in the segmented images. The discrete arc-shaped layering noises are still in the segmented T-A image, although the whole-fan µ + kσ segmentation method segments the images along the beam directions and weakens the linear-shaped side lobe effects. The reason is that the method disregards the frequency and noise differences among adjacent sectors, and the method performs the threshold segmentation in the whole fan with an invariable threshold. However, the proposed method nearly eliminates all the noise and clearly presents the gas plumes in the last binary T-A images.
The gas plumes in Experiments I and II are detected by using the four methods and combining the intersection operation and morphological constraint. The results are listed in Table 1. The proposed method performs better than the other methods do.

Necessity of Segmenting T-A and D-A Images
The 39th ping data with obvious noise are extracted for the three experiments below to explain the necessity of segmenting T-A and D-A images. In Experiment 1, we perform the threshold segmentation for the raw T-A image ( Figure 13B) and obtain the segmented T-A image ( Figure 13B1), transform it into the D-T image ( Figure 13B2), and detect a target by using the morphological constraint ( Figures 13B3). In Experiment 2, the same process is conducted for the raw D-A image ( Figure 13C), and the corresponding results are shown in C1-C3. In Experiment 3, we perform an intersection operation for Figures 13(B1,C1) and obtain the intersection image (Figure 13(D1)). Then, we adopt the same transformation and constraint as those in Experiment 1 and do not detect any target.
In the raw D-T image ( Figure 13A), we cannot find any targets, such as gas plumes. The side lobe effects with horizontal distribution and layering noise with arc-shaped distribution can be seen in the T-A image, respectively; the side lobe effects with arc-shaped distribution and layering noise with horizontal distribution can also be observed in the D-A image. Based on different description forms of side lobe effects and layering noise in T-A and D-A images, the threshold segmentation for T-A image completely removes the side lobe effect and attenuates the layering noise, whereas the

Necessity of Segmenting T-A and D-A Images
The 39th ping data with obvious noise are extracted for the three experiments below to explain the necessity of segmenting T-A and D-A images. In Experiment 1, we perform the threshold segmentation for the raw T-A image ( Figure 13B) and obtain the segmented T-A image ( Figure 13B1), transform it into the D-T image ( Figure 13B2), and detect a target by using the morphological constraint ( Figure 13B3). In Experiment 2, the same process is conducted for the raw D-A image ( Figure 13C), and the corresponding results are shown in C1-C3. In Experiment 3, we perform an intersection operation for Figure 13(B1,C1) and obtain the intersection image (Figure 13(D1)). Then, we adopt the same transformation and constraint as those in Experiment 1 and do not detect any target.
In the raw D-T image ( Figure 13A), we cannot find any targets, such as gas plumes. The side lobe effects with horizontal distribution and layering noise with arc-shaped distribution can be seen in the T-A image, respectively; the side lobe effects with arc-shaped distribution and layering noise with horizontal distribution can also be observed in the D-A image. Based on different description forms of side lobe effects and layering noise in T-A and D-A images, the threshold segmentation for T-A image completely removes the side lobe effect and attenuates the layering noise, whereas the layering noise is eliminated completely, and the side lobe effects are weakened in the D-A image segmentation. The segmentation results are shown in Figure 13(B1,C1). The detection results (Figure 13(B3,C3)) from the T-A or D-A images are incorrect, thereby indicating that a single T-A or D-A image is insufficient in detecting gas plumes. In the intersection image, the side lobe effect and layering noise are eliminated completely. Only certain scattering noise remains in the image. We do not detect any target in Figure 13(D3) after the morphological constraint. The result is consistent with the manual detection from the raw D-T image and further shows that the T-A and D-A images must be combined to detect the gas plume from a ping fan.
layering noise is eliminated completely, and the side lobe effects are weakened in the D-A image segmentation. The segmentation results are shown in Figures 13(B1,C1). The detection results (Figures 13(B3,C3)) from the T-A or D-A images are incorrect, thereby indicating that a single T-A or D-A image is insufficient in detecting gas plumes. In the intersection image, the side lobe effect and layering noise are eliminated completely. Only certain scattering noise remains in the image. We do not detect any target in Figure 13(D3) after the morphological constraint. The result is consistent with the manual detection from the raw D-T image and further shows that the T-A and D-A images must be combined to detect the gas plume from a ping fan.

Limitations of the Proposed Method
The proposed method achieves enhanced detection results in the above experiments but may be affected by the following factors: (1) Shape and size of the gas plume The detection capability of the proposed method depends on the MBES performance parameters, water depth, and beam angle. In a beam ray, a small target can be detected by at least two echoes. We use three echoes to detect the small target and ensure the reliability of the detection. The detection resolution, namely, the height and width of a small target can be determined by: where l is the distance of the three echoes in the beam ray, Li is the length of ith beam ray, and i and θi+1 are the ith and I + 1th beam angles, correspondingly. Li can be calculated by the water

Limitations of the Proposed Method
The proposed method achieves enhanced detection results in the above experiments but may be affected by the following factors: (1) Shape and size of the gas plume The detection capability of the proposed method depends on the MBES performance parameters, water depth, and beam angle. In a beam ray, a small target can be detected by at least two echoes. We use three echoes to detect the small target and ensure the reliability of the detection. The detection resolution, namely, the height and width of a small target can be determined by: where l is the distance of the three echoes in the beam ray, L i is the length of ith beam ray, and θ i and θ i+1 are the ith and I + 1th beam angles, correspondingly. L i can be calculated by the water depth D divided by cosθ i . Therefore, the height, width, and area of the detected minimum target are set as 0.2 m, 2.6 m, and 2 m 2 , respectively, in Experiment I and 11.3 m, 10.5 m, and 230 m 2 , correspondingly, in Experiment II. The proposed method is inefficient for the detections of gas plumes less than the minimum target. The problem can be solved by the increase of sampling rate.
(2) MBES WC data The proposed method is verified by EM 710 and EM122 WC data and needs to be tested further by other MBES WC data with different performances. In the experiments of this research, the WC data measured by the multi-sector, and dual-swath MBES are adopted. The proposed method is also suitable for detecting targets from the same MBES WC and single-frequency and single-sector WC data.
(3) Similar targets as gas plumes The similar targets as gas plumes, such as fishes or fish schools and suspended solids may affect detection of gas plumes. Literature [41,43,44] provided fish-gas discriminate protocols, quantification of seep-related methane gas emissions, and a means to quantify and assess fishes or fish schools in the presence of gas bubbles which may be useful for distinguish gas plume from fishes. In this paper, the BS and morphological features of different targets are used as thresholds to distinguish gas plumes from others. In the extreme case that the threshold parameters of gas plumes and similar targets are close, the diagnosis may be fulfilled by further mining the characteristics of gas plumes and similar targets.
(4) Targets outside of the MSR In the proposed method, the WC data outside of the MSR is eliminated because of the serious MSR effect, and only the WC data in the MSR are used for the detection. The WC data in the MSR are also adopted in the studies depicted in the literature [1,26,27]. Existing studies and the proposed method may become ineffective in detecting the targets outside of the MSR.

Conclusions and Suggestions
The detection method proposed in this paper detects gas plumes from the polluted WC images via three main steps including the mean standard deviation threshold segmentation algorithm for two small-size T-A and D-A images, intersection operation to the two segmented T-A images, and morphological constraint to the intersection image. The proposed method fulfills the automatic detection and achieves higher correct detection rate than the traditional threshold segmentation methods, even for the WC data seriously polluted.
In the applications of the proposed method, some key points should be noted. (1) The parameter k has a significant effect on the threshold segmentation and needs to be determined by experiments in the detected water area; (2) Our method is also appropriate for the MBES containing only one sector; (3) Only the gas plume within the MSR is considered for now, and the WC data outside the MSR is ignored; (4) The morphologic constraints (area, height, width) need to be chosen carefully, which can be calculated by Section 5.4(1) or properly specified.