Side-Scan Sonar Image Mosaic Using Couple Feature Points with Constraint of Track Line Positions

To obtain large-scale seabed surface image, this paper proposes a side-scan sonar (SSS) image mosaic method using couple feature points (CFPs) with constraint of track line positions. The SSS geocoded images are firstly used to form a coarsely mosaicked one and the overlapping areas between adjacent strip images can be determined based on geographic information. Inside the overlapping areas, the feature point (FP) detection and registration operation are adopted for both strips. According to the detected CFPs and track line positions, an adjustment model is established to accommodate complex local distortions as well as ensure the global stability. This proposed method effectively solves the problem of target ghosting or dislocation and no accumulated errors arise in the mosaicking process. Experimental results show that the finally mosaicked image correctly reflects the object distribution, which is meaningful for understanding and interpreting seabed topography.


Introduction
Acoustic waves move quite efficiently through water and can travel over great distances to conduct ocean observation and monitoring [1].Side-scan sonar (SSS) is the most widely used for remote sensing in water environment.SSS transducer built in a towfish is usually towed behind a surveying vessel by a cable, emits a wide-angle beam and receives the seabed echoes at fixed time intervals in a linear path to form the seabed image simultaneously on port and starboard sides [1][2][3].SSS images can vividly reflect the seabed targets, topography features and sediment distributions [1,2], which has important applications in underwater target detection [4,5], terrain reconstruction [6,7], underwater navigation and localization [8] and many other fields.A single strip SSS image can only reflect the local-regional seabed characteristics.To reflect the seabed overall, several stripes must be mosaicked together to obtain the large-scale seabed image.At present, the mosaic methods are mainly based on image geographic information or image common features, which are, respectively, named the geocoding method and feature-based method.
The former method is to construct image mosaics in the same geographic frame according to the image geocoding information [3,9].During the surveying using SSS, the track line positions can be accurately obtained through Short Base Line (SBL) or Ultra Short Base Line (USBL) measurement [10,11].However, these instruments are rarely installed in most towfishes because of high cost.In such a situation, the track line positions are usually reckoned through the vessel positions provided by GNSS and the cable length.Next, the reckoned towfish positions are combined with heading data and attitude data to be used to calculate every SSS image pixel coordinates with the classical flat bottom assumption [1][2][3]12].After these two calculation steps, the SSS images are geocoded and can be mosaicked in the same geographic frame.However, the wave and currents, the vessel velocity and the cable length variation may result in inaccurate towfish location [13,14], which will damage the following SSS image geocoding process.Moreover, when experiencing a complicated seafloor topography, inaccurate heading or attitude data, the calculated SSS image pixel coordinates may also be inaccurate [15,16].As a result, the same target may have different coordinates in different SSS geocoded images and there may exist target ghosting or dislocation in the mosaicked result.
The latter method is based on the image registration technique.Through detecting the common features of the overlapping area in different strip images, a relationship model can be established between the reference and sensed image.According to the relationship model, the sensed image will be transformed into the reference coordinate system and fused with the reference one to form a mosaicked image.As the common features are adopted to connect the reference and sensed image, the coordinates of feature points (FPs) in the sensed image will be transformed to be accordant with that of the reference one.Thus, the problem of target ghosting or dislocation will not happen.There have been several methods developed in this field, including using the thin plate spline (TPS) to build the relationship model between multi-temporal SSS images in the overlapping area [17].Daniel et al. used the highlight area and shadow zone of sonar images to obtain the common features; however, this only works for images with no orientation or scale change [18].Cyril et al. used the correlation ratio and mutual information (CR and MI) of adjacent image intensity to obtain the position information of common features.The performance of the method will be influenced by the complex noise from ocean environment and variational gains [19].Meanwhile, intensity based techniques are also bound to fail for homogeneous areas with low entropy content [20].Zhao et al. adopted the Speeded-Up Robust Features (SURF) algorithm to detect the couple feature points (CFPs) of adjacent strips and used a rotate translation model in the mosaic process, which takes advantage of the feature consistence in different SSS images [21].Ye et al. used nonlinear diffusion filtering to build a nonlinear scale space for SSS image feature matching, which improved noise-reduction performance [22].Tao et al. used the SURF algorithm with prior position knowledge, which could reduce the mismatching rate and the computation time [23].However, Zhao, Ye and Tao's studies are mainly focused on two adjacent images and ignore the question that the accumulated errors of the image coordinates will arise in the following mosaic process when using the former mosaicked one as the reference.As a result, the finally mosaicked image will deviate from the actual measurement range and cannot reflect the true seabed surface.
To solve the mentioned problems above, this study proposes a novel SSS image mosaic method using CFPs with constraint of track line positions.The paper is organized as follows.Section 2 describes the proposed SSS image mosaic method in detail.Section 3 designs the experiments to verify the proposed method and theoretically analyzes the results.Section 4 discusses more about the performance of the proposed method.Section 5 draws the conclusions according to the experiments and discussions.

SSS Image Mosaic Method
The proposed SSS image mosaic method is mainly divided into two steps, which are, respectively, named coarsely mosaicking and refined adjustment.The first step is to construct the image mosaic in the same geographic frame according to the SSS geocoded image pixel coordinates.The second step is to establish an adjustment model to revise the coarsely mosaicked result using the common features of the overlapping areas between different SSS images and the corresponding track line positions.These two steps are described in detail below.
the data decoding, bottom tracking, radiometric distortion correction, slant range correction and geocoding must be conducted to form the SSS geocoded images with clear seabed targets and coordinates [1][2][3]12].After these preprocessing steps, radiation aberration and geometric distortion can be weakened and the image pixel coordinates are obtained.Based on the geocoded images, the image mosaic region can be determined and the gray values are attributed to every pixel in this region based on the known SSS image gray values.
For nonoverlapping areas, only one SSS image gray value exists and can directly serve as the candidate for every pixel in this region.For the overlapped areas, there may exist more than one gray values from different SSS images.Because of the variation of the towfish heights, ocean environment and time varying gain parameters when surveying along different track lines, the gray levels may be different between adjacent strips.To ensure the smooth transition of grayscales in the overlapping areas, a multiresolution fusing method is adopted [24], which is depicted in Figure 1.In this procedure, each image is first decomposed into a set of band-pass filtered component images using a Laplacian pyramid.Each band is then multiplied by a smooth weighting function whose extent is proportional in size to the wave lengths represented in the band.Finally, these band-pass mosaic images are summed to obtain the desired image mosaic.After the coarsely mosaicking process, a large-scale seabed image is formed and the spatial location of the mosaicked image is confirmed.However, this mosaicked image is affected heavily by the inaccurate image coordinates caused by the factors [10][11][12][13] as mentioned and analyzed in Section 1.Thus, the problem of target ghosting and dislocation may exist in the overlapping areas of formed image and the following processing must be conducted.

Determination and Segment of the Overlapping Area
Target ghosting and dislocation resulted from the coarsely mosaicking process only exists in the overlapping areas of SSS images.Thus, the refined adjustment can just be carried out in these overlapping areas, which are determined by the image coordinates.From a practical point of view, the port and starboard sides can be considered independent and the corresponding images are treated as if coming from two one-sided systems [2,3].For each side, echoes are received by the transducer to form the image.In Figure 2, four surveying lines obtained in Zhujiang River, Shenzhen are presented.Each surveying line is presented by the pink line and covers the seabed area along a supposedly straight line.The black arrows denote the tracking directions.The orange dots and blue ones separately denote the two sides of a surveying line, one for port side and the other for starboard side.Because the surveying lines are designed in opposite directions and the overlapping ratio is 50%, the port side of one surveying line may completely overlap the starboard side of the other one.The overlapping areas in Figure 2 are these between adjacent surveying lines.In the measuring process, the proper acquisition of the seabed area located below the sonar requires a redundant coverage of about 40%, which can be used for image registration [19,25].When the ratio is higher After the coarsely mosaicking process, a large-scale seabed image is formed and the spatial location of the mosaicked image is confirmed.However, this mosaicked image is affected heavily by the inaccurate image coordinates caused by the factors [10][11][12][13] as mentioned and analyzed in Section 1.Thus, the problem of target ghosting and dislocation may exist in the overlapping areas of formed image and the following processing must be conducted.

Determination and Segment of the Overlapping Area
Target ghosting and dislocation resulted from the coarsely mosaicking process only exists in the overlapping areas of SSS images.Thus, the refined adjustment can just be carried out in these overlapping areas, which are determined by the image coordinates.From a practical point of view, the port and starboard sides can be considered independent and the corresponding images are treated as if coming from two one-sided systems [2,3].For each side, echoes are received by the transducer to form the image.In Figure 2, four surveying lines obtained in Zhujiang River, Shenzhen are presented.Each surveying line is presented by the pink line and covers the seabed area along a supposedly straight line.The black arrows denote the tracking directions.The orange dots and blue ones separately denote the two sides of a surveying line, one for port side and the other for starboard side.Because the surveying lines are designed in opposite directions and the overlapping ratio is 50%, the port side of one surveying line may completely overlap the starboard side of the other one.The overlapping areas in Figure 2 are these between adjacent surveying lines.In the measuring process, the proper acquisition of the seabed area located below the sonar requires a redundant coverage of about 40%, which can be used for image registration [19,25].When the ratio is higher than 40%, more redundant data exist and will be useful for the image registration.When the ratio is lower, the overlapping areas will become smaller and few common features will be detected.
When the overlapping areas are identified, the segment operation can be conducted along the track line to segment the overlapping areas into several parts according to the target distribution or image textures, which is shown as the black dotted rectangles in Figure 2. Affected by the water current, wave, vessel velocity, cable length variation and complicated topography, the local image distortions may be different in different regions.The segment operation considers these local image distortions and the refined adjustment operation is conducted independently in these segmented partial areas.than 40%, more redundant data exist and will be useful for the image registration.When the ratio is lower, the overlapping areas will become smaller and few common features will be detected.
When the overlapping areas are identified, the segment operation can be conducted along the track line to segment the overlapping areas into several parts according to the target distribution or image textures, which is shown as the black dotted rectangles in Figure 2. Affected by the water current, wave, vessel velocity, cable length variation and complicated topography, the local image distortions may be different in different regions.The segment operation considers these local image distortions and the refined adjustment operation is conducted independently in these segmented partial areas.

Image Registration in the Segmented Overlapping Areas
To refine the image fusion in the overlapping areas, image registration must be conducted to match the common features.For SSS images, the areas with features are paid more attention.The seabed targets, topography characteristics, textures or sediment distributions can be reflected in the SSS images as typical features.
The SURF algorithm, which has invariant ability in rotation, scale, brightness and contrast [26][27][28], can locate key points of high variation and focus on the spatial distribution of gradient information.Thus, the SURF algorithm can be used to find the candidate interest points or these typical features in the SSS images.As the same features presented in the common area of adjacent SSS images may be almost identical, the SURF algorithm can be suitable to detect the CFPs of the adjacent SSS images.Next, the RANSAC algorithm is used to remove mismatches [22,29].Because of using the segment operation, the FP detection and the RANSAC operation are conducted separately in different partial images, which saves much time [20,22].
As SSS image is primarily used for understanding seabed relief and searching seabed targets, the areas with features are often paid more attention, whereas the areas with less features can be ignored.Therefore, for these segmented overlapping areas with no or few CFPs, nothing is done and the previous mosaicked result is preserved, while, for these with abundant CFPs, the refined adjustment must be conducted in the following steps to solve the problem of the target ghosting and dislocation resulted from the coarsely mosaicking process.

The Constraint of Track Line Positions
The detected CFPs can be used to establish a geometric transformation model between adjacent

Image Registration in the Segmented Overlapping Areas
To refine the image fusion in the overlapping areas, image registration must be conducted to match the common features.For SSS images, the areas with features are paid more attention.The seabed targets, topography characteristics, textures or sediment distributions can be reflected in the SSS images as typical features.
The SURF algorithm, which has invariant ability in rotation, scale, brightness and contrast [26][27][28], can locate key points of high variation and focus on the spatial distribution of gradient information.Thus, the SURF algorithm can be used to find the candidate interest points or these typical features in the SSS images.As the same features presented in the common area of adjacent SSS images may be almost identical, the SURF algorithm can be suitable to detect the CFPs of the adjacent SSS images.Next, the RANSAC algorithm is used to remove mismatches [22,29].Because of using the segment operation, the FP detection and the RANSAC operation are conducted separately in different partial images, which saves much time [20,22].
As SSS image is primarily used for understanding seabed relief and searching seabed targets, the areas with features are often paid more attention, whereas the areas with less features can be ignored.Therefore, for these segmented overlapping areas with no or few CFPs, nothing is done and the previous mosaicked result is preserved, while, for these with abundant CFPs, the refined adjustment must be conducted in the following steps to solve the problem of the target ghosting and dislocation resulted from the coarsely mosaicking process.

The Constraint of Track Line Positions
The detected CFPs can be used to establish a geometric transformation model between adjacent SSS images in the overlapping areas.In this process, one SSS image is selected as the reference and the other is consequently as the sensed one.Through the established model, registration from sensed image to reference image can be achieved to ensure the target location uniqueness, which solves the problem of target ghosting or dislocation resulted from the coarse image mosaicking.Because of the local shifts in the SSS image, a rigid transformation model cannot consider the local distortions around every CFP.To solve this problem, an elastic transformation model can be adopted to obtain a commendable relationship between every CFPs.However, only using the CFPs to establish the adjustment model will result in distortions in other remainder image areas.
SSS image coordinates are calculated using the towfish positions, heading data and attitude data in the track line.No external location information exists in the mosaic process.The track line positions can define the measurement area boundary.Thus, in overlapping areas, evenly distributed points in the track line can be adopted as a constraint to ensure the image range invariance on the whole when conducting the adjustment operation.Combining the obtained CFPs and selected track line positions, an elastic model which is used to refine the coarsely mosaicked image can be established to both accommodate the local distortions between CFPs and keep the global stability in the overlapping areas.

Adjustment Model Using CFPs and Track Line Positions
The adjustment model must consider both the local distortions and global stability.The thin plate spline (TPS) function can model complex local distortions without undue effects being propagated throughout the remainder of the image [17,[30][31][32].When n points (x i , y i ) are given in the plane, a TPS function can be defined in Equation (1).
, where a i and b i are the transformation parameters to be solved by this equation.
Combining the detected CPFs (Xf, Yf ) and selected track line positions (Xc, Yc), the adjustment model based on TPS function is established in Equation (2). where Remote Sens. 2018, 10, 953 6 of 22 where (X, y) and (X , y ) separately denote reference and sensed image coordinates, (X , y ) has the similar expression formula as (X, y); (Xf, Yf ) and (Xc, Yc) separately denote the CFP and track line position coordinates; and Ax, Ay, Bx, and By are the transformation parameters to be calculated.
To solve this model, Equation (3) is obtained.
According to the obtained transformation parameters in Equation ( 3), the sensed image in the overlapping area can be transformed into the reference frame and the CFP coordinates in both images would become accordant with each other.In Equation ( 2), the coordinate transformation model between the reference and sensed image is established by taking full advantage of both the CFPs and track line positions.Such an adjustment model has the following three advantages: (1) The local coordinate dislocations in segmented overlapping areas are considered by using the CFPs when establishing the adjustment model.(2) The transformed sensed image will remain stable globally because of using the track line positions as the constraint.
(3) The established model based on TPS function can satisfy the above two requirements.

Gap Filling and Image Fusion
When transforming the sensed image into the geographic frame defined by the reference image, there may exist gaps in the transformed image.To fill these gaps, the adjacent nonblank pixel values can be used and the pixel filling method displayed as Equation ( 4) is adopted to acquire gray values in the gap areas.
where (x, y) is the gap position in the transformed sensed image and (x , y ) is its corresponding pixel position in the raw sensed image.g(x, y) is the interpolated gray value in the gap position (x, y), and f(x , y + 1) are the nonblank gray values of four adjacent point that are closest to the point (x , y ).
In overlapping areas with abundant CFPs, the reference image and the transformed sensed image are fused to obtain the finally mosaicked one using a multiresolution fusing method, as shown in Figure 1.The new obtained image is to replace the previous one which is coarsely mosaicked.In overlapping areas with no or few CFPs and in nonoverlapping areas, the problem of target ghosting and dislocation does not exist.Thus, no operation is conducted and the coarsely mosaicked result is preserved.

The SSS Image Mosaic Process
The above mosaic process can be described in Figure 3.In the proposed method, the first step determines the spatial position of the mosaicked image and guarantees that the following adjustment step will be conducted in the valid measurement area.Next, the refined adjustment operation is only conducted in overlapping areas and other remaining areas keep invariant, which ensures that the original image coordinate information can be utilized to the maximum.Moreover, the combination of the CFPs and the constraint of corresponding track line positions can both accommodate local dislocations between CFPs of both images as well as ensure global stability of the overlapping area in the mosaic process.As a result, the position errors between adjacent images are disintegrated inside each overlapping area and no accumulated errors will occur.This proposed method takes full advantage of both the geocoding method and feature-based method, which can not only avoid the appearance of accumulated errors, but also solve the problem of target ghosting and dislocation.

Experimental Data
To verify the proposed method, an experiment was carried out in Zhu Jiang River with water depth ranging from 10 m to 15 m.In this experiment, EdgeTech 4125 with the operating frequency of 400 KHz and towed by a cable of 5 m behind the vessel was adopted.

Experimental Data
To verify the proposed method, an experiment was carried out in Zhu Jiang River with water depth ranging from 10 m to 15 m.In this experiment, EdgeTech 4125 with the operating frequency of 400 KHz and towed by a cable of 5 m behind the vessel was adopted.Four surveying lines with the lengths of 2015 m, 2390 m, 2474 m and 1936 m were, respectively, set up with the 200 m swath width and 50% overlapping rate.Raw surveying results were recorded in ".xtf" files.Through data decoding, radiation correction, slant-range correction and geocoding, the SSS geocoded images with the pixel resolution of 0.6 m and clear targets were obtained, as shown in Figure 4, which can be used to reflect the seabed surface.(1) Geocoding mosaic.Use the SSS geocoded images to obtain the mosaicked one in the same geographic frame using the method described in Section 2.1.The mosaicked image in Figure 5a shows there exist target ghosting and dislocation which can be clearly seen as L1-L4.It can be found in Figure 5 that common features exist in the overlapping areas of adjacent SSS images, which are named TL1-TL4 in one SSS image and TR1-TR4 in the other one.The mosaicked image is supposed to show the four corresponding features in the same positions.However, affected by the wave and currents, the vessel velocity, the cable length variation and the topography changes as mentioned and analyzed in Section 1, the obtained SSS image feature coordinates may contain errors, which will result in the problem of target dislocation.As a result, the same feature may locate in more than one positions in the mosaicked image and the phenomenon is called target ghosting.In positions L1-L4 of Figure 5a, ghosts accompany the four targets, which means that target ghosting and dislocation exist in these positions.(1) Geocoding mosaic.Use the SSS geocoded images to obtain the mosaicked one in the same geographic frame using the method described in Section 2.1.The mosaicked image in Figure 5a shows there exist target ghosting and dislocation which can be clearly seen as L1-L4.It can be found in Figure 5 that common features exist in the overlapping areas of adjacent SSS images, which are named TL1-TL4 in one SSS image and TR1-TR4 in the other one.The mosaicked image is supposed to show the four corresponding features in the same positions.However, affected by the wave and currents, the vessel velocity, the cable length variation and the topography changes as mentioned and analyzed in Section 1, the obtained SSS image feature coordinates may contain errors, which will result in the problem of target dislocation.As a result, the same feature may locate in more than one positions in the mosaicked image and the phenomenon is called target ghosting.In positions L1-L4 of Figure 5a, ghosts accompany the four targets, which means that target ghosting and dislocation exist in these positions.
(2) Conducting the overlapping area determination and segment.Determine the overlapping areas according to the SSS image coordinates and segment these areas into several equal regions along the track line based on the target distribution and image textures.Considering that enough independent targets should be contained in each segmented part and split lines should lie on the edge of these distinct targets, this experiment separated the overlapping area along the track line into several parts with the length of 420 m, which are obviously shown as the black rectangles (A)-(D) in Figure 5.
the four corresponding features in the same positions.However, affected by the wave and currents, the vessel velocity, the cable length variation and the topography changes as mentioned and analyzed in Section 1, the obtained SSS image feature coordinates may contain errors, which will result in the problem of target dislocation.As a result, the same feature may locate in more than one positions in the mosaicked image and the phenomenon is called target ghosting.In positions L1-L4 of Figure 5a, ghosts accompany the four targets, which means that target ghosting and dislocation exist in these positions.(3) Detecting CFPs.When two surveying lines are designed in the same area, the overlapping areas exist.Even measured along different surveying lines using SSS, the obtained SSS images could reflect the same seabed features after data preprocessing, such as the radiometric correction, slant range correction, geocoding, etc.The features and its neighborhoods in SSS images can be detected and described by the SURF algorithm, which will generate feature vectors, namely descriptors.By matching the descriptors of the features in different SSS images, similar features can be obtained.The matching is based on a distance between the feature vectors, e.g., the Euclidean distance [26].To eliminate the incorrect matched FPs, RANSAC algorithm is then adopted.After this step, the CFPs are finally paired.For the areas with no or few CFPs detected, no further processes are conducted.While for these with abundant CFPs, evenly distributed couple ones are selected for the subsequent image registration and mosaicking.The detected CFPs in the areas (A)-(C) of Figure 5 are shown in Figure 6.When the SURF and RANSAC algorithms are conducted for the adjacent SSS images in area (D), no distinctive CFPs are detected because no significant features exist.Thus, no image registration will be conducted in area (D).(a (4) The constraint of track line positions.For every segmented part, the corresponding track line positions are selected as the constraint to ensure the sensed image not deviate from the actual measurement region after geometric transformation.In this experiment, points with the interval of about 30 m in the track line are selected, which are shown in Figure 7 and denoted by the " " symbols.When establishing the geometric transformation model between the reference image and the sensed one, both the selected CFPs and the constraint of track line positions are used to accommodate the local distortions and global stability.(5) Adjustment model.Some evenly distributed CFPs and track line positions are combined together to construct the adjustment model similar to Equation ( 2), which can accommodate local distortions to match the CFPs as well as be globally stable with the constraint of track line positions.For the separated overlapping areas (A)-(C) shown in Figure 5, the adjustment models are respectively established.After the sensed image has been transformed by the adjustment model, a gap filling method, as displayed as Equation ( 4), should be conducted to obtain a completed image.
(6) Mosaicking images.Fusing the transformed sensed image and the reference one in the overlapping areas (A)-(C), a new mosaicked image can be obtained.Replacing the coarsely mosaicked image with these newly obtained ones, a better mosaicked image can be formed.
Through above processing, the finally mosaicked image is shown in Figure 5b.The CFP coordinates in R1-R4 are more accordant with each other than these in L1-L4 as shown in Figure 5a, which is more helpful for image interpretation and application.Because segment operation is conducted in the overlapping area in Figure 5b, the CFPs can be detected and selected separately in these segmented parts as described in Step (3).In this experiment, topography changes exist in the measurement area and are reflected by the image intensities, which means that the positions with obvious image intensity gradients correspond to these with topography undulations.As SURF algorithm is adopted to generate feature vector to describe the intensity gradients, the CFPs detected in this measurement area are mainly located in the positions with distinctive topography changes.

Evaluation of the Mosaic Method
According to the relationship of selected CFP coordinates, the problem of target ghosting and dislocation can be solved by the adjustment model as displayed in Equation (2).Meanwhile, the refined adjustment operation is only conducted independently in segmented overlapping areas and the track line constraint is used.Thus, the position error can be disintegrated inside these individual segmented areas and will not propagate into other areas.To assess the performance of the image mosaic method, the CFPs and track line points not adopted to establish the adjustment model are chosen to calculate the position deviations before and after the mosaic process, respectively, which is described in detail below.

The Consistency of CFP Coordinates
Choosing the CFPs which are not involved in establishing the adjustment model and computing their coordinate deviations between the reference image and sensed image before and after the mosaicking process.The distributions of their coordinate deviations are shown in Figure 8, where "☆" and "▷" symbols respectively denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○" and"◇" symbols denote these between the reference image and raw sensed image.The corresponding statistical parameters of these coordinate deviations are shown in Table 1. Figure 8 and Table 1 manifest that the former deviations distribute evenly around the mean deviation of 0.0 m, whose standard deviation are, respectively, 1.06 m and 1.76 m in both directions.Meanwhile, the maximum deviations in the east- The relative positions of CFPs of both SSS images in area (a) are shown in Figure 7, which are, respectively, denoted by "+" and " " symbols.It can be seen that the coordinate deviation of CFPs is variational in different positions.In some positions, the FPs in reference image lie to the east of the ones in sensed image, while, in other positions, the situation becomes reversed.In the areas (b,c), a similar situation also exists and implies that the local distortions in the SSS geocoded image are diverse.
(5) Adjustment model.Some evenly distributed CFPs and track line positions are combined together to construct the adjustment model similar to Equation (2), which can accommodate local distortions to match the CFPs as well as be globally stable with the constraint of track line positions.For the separated overlapping areas (A)-(C) shown in Figure 5, the adjustment models are respectively established.After the sensed image has been transformed by the adjustment model, a gap filling method, as displayed as Equation ( 4), should be conducted to obtain a completed image.
(6) Mosaicking images.Fusing the transformed sensed image and the reference one in the overlapping areas (A)-(C), a new mosaicked image can be obtained.Replacing the coarsely mosaicked image with these newly obtained ones, a better mosaicked image can be formed.
Through above processing, the finally mosaicked image is shown in Figure 5b.The CFP coordinates in R1-R4 are more accordant with each other than these in L1-L4 as shown in Figure 5a, which is more helpful for image interpretation and application.Because segment operation is conducted in the overlapping area in Figure 5b, the CFPs can be detected and selected separately in these segmented parts as described in Step (3).In this experiment, topography changes exist in the measurement area and are reflected by the image intensities, which means that the positions with obvious image intensity gradients correspond to these with topography undulations.As SURF algorithm is adopted to generate feature vector to describe the intensity gradients, the CFPs detected in this measurement area are mainly located in the positions with distinctive topography changes.

Evaluation of the Mosaic Method
According to the relationship of selected CFP coordinates, the problem of target ghosting and dislocation can be solved by the adjustment model as displayed in Equation ( 2).Meanwhile, the refined adjustment operation is only conducted independently in segmented overlapping areas and the track line constraint is used.Thus, the position error can be disintegrated inside these individual segmented areas and will not propagate into other areas.To assess the performance of the image mosaic method, the CFPs and track line points not adopted to establish the adjustment model are chosen to calculate the position deviations before and after the mosaic process, respectively, which is described in detail below.

The Consistency of CFP Coordinates
Choosing the CFPs which are not involved in establishing the adjustment model and computing their coordinate deviations between the reference image and sensed image before and after the mosaicking process.The distributions of their coordinate deviations are shown in Figure 8, where " " and " " symbols respectively denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; " " and" " symbols denote these between the reference image and raw sensed image.The corresponding statistical parameters of these coordinate deviations are shown in Table 1. Figure 8 and Table 1 manifest that the former deviations distribute evenly around the mean deviation of 0.0 m, whose standard deviation are, respectively, 1.06 m and 1.76 m in both directions.Meanwhile, the maximum deviations in the east-west and north-south directions are, respectively, less than 3.10 m and 4.31 m, which are smaller than those of the raw coordinate deviations with the maximum deviations of 16.58 m and 16.94 m.As the pixel resolution is set as 0.6 m, the coordinate deviations presented in the refined mosaicked image would be less than 8 pixels.Compared with those between the reference image and raw sensed image, which are more than 30 pixels, the proposed method is testified to have rectified the local distortions around the CFPs.These statistics prove that the coordinates of CFPs in the mosaicked image become consistent.west and north-south directions are, respectively, less than 3.10 m and 4.31 m, which are smaller than those of the raw coordinate deviations with the maximum deviations of 16.58 m and 16.94 m.As the pixel resolution is set as 0.6 m, the coordinate deviations presented in the refined mosaicked image would be less than 8 pixels.Compared with those between the reference image and raw sensed image, which are more than 30 pixels, the proposed method is testified to have rectified the local distortions around the CFPs.These statistics prove that the coordinates of CFPs in the mosaicked image become consistent.As the track line positions are selected as the control points to construct adjustment model, the coordinates of these positions are supposed to be stable or change very little when conducting the refined adjustment operation.Choosing the points in the track line not used to establish the adjustment model and computing their coordinate deviations before and after the adjustment.The distributions of these coordinate deviations are shown in Figure 9, where "----" and "-" symbols, respectively, denote the coordinate deviations in east-west and north-south directions.The corresponding statistical parameters of these coordinate deviations are shown in Table 2. Figure 9 and Table 2 show that the maximum deviations in the east-west and north-south directions are less than 0.2 m.Meanwhile, the deviations distributed evenly around the mean deviation of 0.0 m, whose standard deviations are, respectively, 0.02 m and 0.03 m.Such statistics indicate that the track line positions remain almost unchanged after the adjustment operation and the coordinate variations presented in the mosaicked image are less than one pixel, which proves that the track line positions in the current measurement area are unchanged and the mosaicked image remains in the valid measurement region.These statistics and analyses show the validity of using the constraint of track line positions to ensure global stability in constructing the adjustment model.

The Variation of Track Line Positions
As the track line positions are selected as the control points to construct adjustment model, the coordinates of these positions are supposed to be stable or change very little when conducting the refined adjustment operation.Choosing the points in the track line not used to establish the adjustment model and computing their coordinate deviations before and after the adjustment.The distributions of these coordinate deviations are shown in Figure 9, where "--" and "-" symbols, respectively, denote the coordinate deviations in east-west and north-south directions.The corresponding statistical parameters of these coordinate deviations are shown in Table 2. Figure 9 and Table 2 show that the maximum deviations in the east-west and north-south directions are less than 0.2 m.Meanwhile, the deviations distributed evenly around the mean deviation of 0.0 m, whose standard deviations are, respectively, 0.02 m and 0.03 m.Such statistics indicate that the track line positions remain almost unchanged after the adjustment operation and the coordinate variations presented in the mosaicked image are less than one pixel, which proves that the track line positions in the current measurement area are unchanged and the mosaicked image remains in the valid measurement region.These statistics and analyses show the validity of using the constraint of track line positions to ensure global stability in constructing the adjustment model.It can also be seen that coordinate deviations in Figure 8 and Table 1 are much bigger than those shown in Figure 9 and Table 2 by almost one order of magnitude.The CFPs in different SSS images are entitled different coordinates because of local shifts deduced from the factors [13][14][15][16] as mentioned and analyzed in Section 1.Even after the geometric transformation using the adjustment model, the coordinate deviations of CFPs can only be weakened and not be eliminated completely, whereas the position coordinates in the track line are unique and served as the constraint in constructing the adjustment model.Their role is to ensure the global stability of the overlapping area in the image mosaic process and naturally their positions will remain almost unchanged.Thus, the coordinate deviations of the track line positions will be much smaller than those of the CFPs.

Comparison of the Proposed Method to Others
To further evaluate the reliability of the proposed method, two comparison experiments are conducted below.Experiment 1 is to compare the computational time between the feature-based image matching algorithm used in the proposed method and the area-based matching algorithm.Experiment 2 is to compare the performance of the mosaicking methods with and without constraint of track line positions.

Comparison Experiment 1
Image matching must be conducted before the image mosaic process and common features in overlapping areas are automatically selected.Image matching methods are generally classified into two categories: feature-based and area-based methods.The former method is to use similarities of features extracted from images to achieve CFPs between images and SURF and RANSAC algorithms are used in the proposed method.Area-based method is to use similarity metrics to detect CFPs between images through a template matching strategy.Cyril et al. analyzed different kinds of similarity metrics and concluded that the correlation ratio and mutual information (CR and MI) perform best of all and can determine robust matching blocks between SSS images [19].
To compare the computing time by using the two different image matching methods, the experiment below is performed using the same hardware with the i7, 3.40 GHz Intel Core and 8.00 GB RAM.These two methods are used for SSS image matching operation in areas (A)-(C) and the matching results are shown in Figure 10.The computing time by selecting the same number of CFPs between adjacent images are compared in Table 3. Figure 10 and Table 3 indicate that both methods can achieve abundant CFPs, but using CR and MI to match the feature points will cost more time.The reason may be that the computation of CR and MI of the FPs between different SSS images  It can also be seen that coordinate deviations in Figure 8 and Table 1 are much bigger than those shown in Figure 9 and Table 2 by almost one order of magnitude.The CFPs in different SSS images are entitled different coordinates because of local shifts deduced from the factors [13][14][15][16] as mentioned and analyzed in Section 1.Even after the geometric transformation using the adjustment model, the coordinate deviations of CFPs can only be weakened and not be eliminated completely, whereas the position coordinates in the track line are unique and served as the constraint in constructing the adjustment model.Their role is to ensure the global stability of the overlapping area in the image mosaic process and naturally their positions will remain almost unchanged.Thus, the coordinate deviations of the track line positions will be much smaller than those of the CFPs.

Comparison of the Proposed Method to Others
To further evaluate the reliability of the proposed method, two comparison experiments are conducted below.Experiment 1 is to compare the computational time between the feature-based image matching algorithm used in the proposed method and the area-based matching algorithm.Experiment 2 is to compare the performance of the mosaicking methods with and without constraint of track line positions.

Comparison Experiment 1
Image matching must be conducted before the image mosaic process and common features in overlapping areas are automatically selected.Image matching methods are generally classified into two categories: feature-based and area-based methods.The former method is to use similarities of features extracted from images to achieve CFPs between images and SURF and RANSAC algorithms are used in the proposed method.Area-based method is to use similarity metrics to detect CFPs between images through a template matching strategy.Cyril et al. analyzed different kinds of similarity metrics and concluded that the correlation ratio and mutual information (CR and MI) perform best of all and can determine robust matching blocks between SSS images [19].
To compare the computing time by using the two different image matching methods, the experiment below is performed using the same hardware with the i7, 3.40 GHz Intel Core and 8.00 GB RAM.These two methods are used for SSS image matching operation in areas (A)-(C) and the matching results are shown in Figure 10.The computing time by selecting the same number of CFPs between adjacent images are compared in Table 3. Figure 10 and Table 3 indicate that both methods can achieve abundant CFPs, but using CR and MI to match the feature points will cost more time.The reason may be that the computation of CR and MI of the FPs between different SSS images requires the surrounding image pixels within a small search region.Thus, SURF and RANSAC algorithms are recommend for SSS image matching operation.4.   In Figure 11 and Table 4, the transformed positions with constraint remain closer to the raw ones and change less compared to those without constraint.Besides, the transformed sensed image  4.   In Figure 11 and Table 4, the transformed positions with constraint remain closer to the raw ones and change less compared to those without constraint.Besides, the transformed sensed image without constraint have several other distortions in the edge positions.The reason may be that, when  In Figure 11 and Table 4, the transformed positions with constraint remain closer to the raw ones and change less compared to those without constraint.Besides, the transformed sensed without constraint have several other distortions in the edge positions.The reason may be that, when constructing the adjustment model without track line constraint, the coordinate relationships of the CFPs are only considered and the global stability of the whole area is ignored, whereas the established adjustment model with constraint uses the track line positions to guarantee the global stability as well as accommodate the local distortions between CFPs when conducting geometric transformation.Both Figure 11 and Table 4 show that the transformed sensed image with constraint will not deviate from the current track line and remain in the valid measurement area.

Multi-Strip Image Mosaics and Analysis
Constructing four strip image mosaics together using the proposed method and the result is in Figure 12.As no obviously common targets present in Strip III and IV, only geocoding method is used when mosaicking the two images.While there exist distinct objectives in both Strips I and II, and Strips II and III, the refined adjustment operation must be performed and the complete mosaic process must be conducted.The mosaicked image in Figure 12 shows clear target distribution and ensures the coordinate consistency of CFPs in the overlapping areas.

Multi-Strip Image Mosaics and Analysis
Constructing four strip image mosaics together using the proposed method and the result is shown in Figure 12.As no obviously common targets present in Strip III and IV, only geocoding method is used when mosaicking the two images.While there exist distinct objectives in both Strips I and II, and Strips II and III, the refined adjustment operation must be performed and the complete mosaic process must be conducted.The mosaicked image in Figure 12 shows clear target distribution and ensures the coordinate consistency of CFPs in the overlapping areas.In the mosaic process using the four strip images, only geocoding method is used in the nonoverlapping areas or overlapping areas with few CFPs.While for overlapping areas with abundant CFPs, the refined adjustment operation is subsequently performed.When conducting the adjustment operation, the port and starboard of SSS images are considered independent and the overlapping area is segmented.For every segmented part, the corresponding track line positions are jointed with the CFPs to serve as the constraint to construct an adjustment model to refine the mosaicked result.Using the obtained adjustment model, the position error can be disintegrated inside individual segmented areas and will not propagate over other regions.In such a situation, there only exist position errors resulted from the adjacent image mosaic operation even when multi-strip images are processed and the problem of positional accumulated errors is solved.Thus, the proposed method ensures the global stability of the sensed image even after the geometric transform and raw measurement information can be reserved and utilized to the maximum.
The In the mosaic process using the four strip images, only geocoding method is used in the nonoverlapping areas or overlapping areas with few CFPs.While for overlapping areas with abundant CFPs, the refined adjustment operation is subsequently performed.When conducting the adjustment operation, the port and starboard of SSS images are considered independent and the overlapping area is segmented.For every segmented part, the corresponding track line positions are jointed with the CFPs to serve as the constraint to construct an adjustment model to refine the mosaicked result.Using the obtained adjustment model, the position error can be disintegrated inside individual segmented areas and will not propagate over other regions.In such a situation, there only exist position errors resulted from the adjacent image mosaic operation even when multi-strip images are processed and the problem of positional accumulated errors is solved.Thus, the proposed method ensures the global stability of the sensed image even after the geometric and raw measurement information can be reserved and utilized to the maximum.
The last two SSS images (Strip II and Strip III) are used to evaluate the mosaic method when conducting multi-strip image mosaics.Strip II is used as the reference and consequently Strip III is used as the sensed one.The coordinate deviations of the CFPs between the reference and sensed image before and after the mosaicking process are analyzed.These coordinate deviations are shown in Figure 13, where " " and " " symbols denote coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; " "and " " symbols denote these between the reference image and raw sensed image.The corresponding statistical parameters of these coordinate deviations are shown in Table 5.It can be seen in Figure 13 and Table 5 that the deviations after the refined adjustment become much smaller than the raw deviations and the maximum deviations in east-west and north-south directions are separately 3.90 m and 4.07 m.Meanwhile, compared to the standard deviations of raw coordinates in both directions (3.85 m and 4.67 m), the ones after adjustment become much smaller (1.48 m and 1.59 m).In the experiment, the pixel resolution is set as 0.6 m, the coordinate deviations presented in the refined mosaicked image would be less than 7 pixels after adjustment.These statistics prove that the adjustment model have rectified the local distortions around the CFPs and the coordinates of CFPs in the mosaicked image show good consistency.The last two SSS images (Strip II and Strip III) are used to evaluate the mosaic method when conducting multi-strip image mosaics.Strip II is used as the reference and consequently Strip III is used as the sensed one.The coordinate deviations of the CFPs between the reference and sensed image before and after the mosaicking process are analyzed.These coordinate deviations are shown in Figure 13, where "☆" and "▷" symbols denote coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○"and "◇" symbols denote these between the reference image and raw sensed image.The corresponding statistical parameters of these coordinate deviations are shown in Table 5.It can be seen in Figure 13 and Table 5 that the deviations after the refined adjustment become much smaller than the raw deviations and the maximum deviations in east-west and north-south directions are separately 3.90 m and 4.07 m.Meanwhile, compared to the standard deviations of raw coordinates in both directions (3.85 m and 4.67 m), the ones after adjustment become much smaller (1.48 m and 1.59 m).In the experiment, the pixel resolution is set as 0.6 m, the coordinate deviations presented in the refined mosaicked image would be less than 7 pixels after adjustment.These statistics prove that the adjustment model have rectified the local distortions around the CFPs and the coordinates of CFPs in the mosaicked image show good consistency.
Besides, the track line positions in the mosaicked image are also analyzed.As Strip III and IV are mosaicked only by the geocoding information, the track line positions of Strip IV maintain unchanged, whereas Strips I, II and III are mosaicked using the complete process of the proposed method.Choosing the track line points of Strip III not used in establishing the adjustment model, their coordinate deviations before and after the mosaicking process are computed and shown in Figure 14, where dotted and solid lines which can be seen as "----" and "-" symbols, respectively, denote the coordinate deviations in east-west and north-south directions.The corresponding statistical parameters of these coordinate deviations are shown in Table 5. Figure 14 and Table 5 show that the maximum deviations in both the east-west and north-south directions are less than 0.2 m.Meanwhile, the deviations distributed evenly around the mean deviation of 0.0 m, whose standard deviations are respectively 0.02 m and 0.04 m.These statistics prove that the mosaicked result is not deviated from the original measurement range.The last two SSS images (Strip II and Strip III) are used to evaluate the mosaic method when conducting multi-strip image mosaics.Strip II is used as the reference and consequently Strip III is used as the sensed one.The coordinate deviations of the CFPs between the reference and sensed image before and after the mosaicking process are analyzed.These coordinate deviations are shown in Figure 13, where "☆" and "▷" symbols denote coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○"and "◇" symbols denote these between the reference image and raw sensed image.The corresponding statistical parameters of these coordinate deviations are shown in Table 5.It can be seen in Figure 13 and Table 5 that the deviations after the refined adjustment become much smaller than the raw deviations and the maximum deviations in east-west and north-south directions are separately 3.90 m and 4.07 m.Meanwhile, compared to the standard deviations of raw coordinates in both directions (3.85 m and 4.67 m), the ones after adjustment become much smaller (1.48 m and 1.59 m).In the experiment, the pixel resolution is set as 0.6 m, the coordinate deviations presented in the refined mosaicked image would be less than 7 pixels after adjustment.These statistics prove that the adjustment model have rectified the local distortions around the CFPs and the coordinates of CFPs in the mosaicked image show good consistency.
Besides, the track line positions in the mosaicked image are also analyzed.As Strip III and IV are mosaicked only by the geocoding information, the track line positions of Strip IV maintain unchanged, whereas Strips I, II and III are mosaicked using the complete process of the proposed method.Choosing the track line points of Strip III not used in establishing the adjustment model, their coordinate deviations before and after the mosaicking process are computed and shown in Figure 14, where dotted and solid lines which can be seen as "----" and "-" symbols, respectively, denote the coordinate deviations in east-west and north-south directions.The corresponding statistical parameters of these coordinate deviations are shown in Table 5. Figure 14 and Table 5 show that the maximum deviations in both the east-west and north-south directions are less than 0.2 m.Meanwhile, the deviations distributed evenly around the mean deviation of 0.0 m, whose standard deviations are respectively 0.02 m and 0.04 m.These statistics prove that the mosaicked result is not deviated from the original measurement range.Besides, the track line positions in the mosaicked image are also analyzed.As Strip III and IV are mosaicked only by the geocoding information, the track line positions of Strip IV maintain unchanged, whereas Strips I, II and III are mosaicked using the complete process of the proposed method.Choosing the track line points of Strip III not used in establishing the adjustment model, their coordinate before and after the mosaicking process are computed and shown in Figure 14, where dotted and solid lines which can be seen as "--" and "-" symbols, respectively, denote the coordinate deviations in east-west and north-south directions.The corresponding statistical parameters of these coordinate deviations are shown in Table 5. Figure 14 and Table 5 show that the maximum deviations in both the east-west and north-south directions are less than 0.2 m.Meanwhile, the deviations distributed evenly around the mean deviation of 0.0 m, whose standard deviations are respectively 0.02 m and 0.04 m.These statistics prove that the mosaicked result is not deviated from the original measurement range.
In Figures 8, 9, 13 and 14 and Tables 1, 2 and 5, it can be seen that the coordinate deviations of CFPs and track line positions in the multi-strip image mosaics have the same order of magnitude with that in the adjacent image mosaics.These statistics and the above theoretical analysis prove that the proposed method can not only accommodate the local distortions of CFPs, but also ensure the mosaicked image in the valid measurement area.
Comparing the proposed method to others, the SSS image mosaic results, the computing time and the accuracy statistics are separately shown in Figure 15 and Table 6.The geocoding method performs faster of all, but obvious target ghosting can be seen in positions (A)-(E) in Figure 15a.The other three mosaic results avoid the target ghosting problem even though more computing time is used.Considering the mosaic performance, computing time and the accuracy statistics, the image mosaic method using the CFPs with constraint of track line positions is recommended.During SSS measurement, the SSS transducer emits the acoustic waves into the sea and receives the backscattered ones at fixed time intervals.Such a mechanism could generate one ping for each time and the SSS image can be finally formed by stacking all these ping scanning lines together [1].The distortions may be present in these pings and are different due to varying underwater environment, vessel velocities, towfish depth and heading as well as positioning accuracies.The distortions are different from those in satellite or drone [33].Therefore, the mosaicking process used in satellite or drone images is unsuitable for SSS images.
In this paper, the method in consideration of seabed sediment variation is used to reduce SSS image distortions [34].The image registration operation using SURF and RANSAC algorithm is respectively performed for the raw SSS images and the radiation corrected ones.It can be seen in Figure 16 that few CFPs are detected from the raw SSS images and the detected ones are visibly incorrect and not lie in the distinct feature positions.While for the radiation corrected images, clearly correct CFPs are detected and matched, which will be helpful to establish an adjustment model.The comparison in Figure 16 illustrates the importance of the SSS image quality and suggests the SSS image preprocessing must be conducted before the mosaicking process.

Position Accuracy of the Mosaicked Image
The accuracy of track line positions is important in the proposed method.Generally, the positions of SSS towfish or tracking lines can be determined by two methods.One is the location estimation method and the other is Ultra-short Base Line (USBL) acoustic positioning method.The principle of the first method is shown in Figure 17.The positions of towfish can be determined by ( )

cos( ) sin( )
where Dlayback is the projected horizontal length of the cable, L is the cable length, h is the height of the GPS receiver above the sea surface, fd is the towfish depth, A is the heading data that can be recorded by the compass, and (xtowfish, ytowfish) and (xvessel, yvessel) are separately the coordinates of the towfish and the vessel.

Position Accuracy of the Mosaicked Image
The accuracy of track line positions is important in the proposed method.Generally, the positions of SSS towfish or tracking lines can be determined by two methods.One is the location estimation method and the other is Ultra-short Base Line (USBL) acoustic positioning method.The principle of the first method is shown in Figure 17.The positions of towfish can be determined by (5) where D layback is the projected horizontal length of the cable, L is the cable length, h is the height of the GPS receiver above the sea surface, f d is the towfish depth, A is the heading data that can be recorded by the compass, and (x towfish , y towfish ) and (x vessel , y vessel ) are separately the coordinates of the towfish and the vessel.
principle of the first method is shown in Figure 17.The positions of towfish can be determined by ( ) where Dlayback is the projected horizontal length of the cable, L is the cable length, h is the height of the GPS receiver above the sea surface, fd is the towfish depth, A is the heading data that can be recorded by the compass, and (xtowfish, ytowfish) and (xvessel, yvessel) are separately the coordinates of the towfish and the vessel.The location estimation method is often adopted in SSS measurement, and has various positioning accuracy with vessel velocity and heading.To get accurate positions of towfish, accurate vessel-mounted GPS positioning and stable vessel velocity and heading are needed in SSS measurement.For the second method, the relative coordinates between towfish and GPS can be obtained by towfish-mounted USBL, and the absolute position of towfish can be calculated by combining GPS positioning solution.The accuracies of the two positioning methods are not high and lead to the target dislocation in a geocoded image and the target ghosting appears in the geocoding mosaicked images.To remove the problems, the common features and the constraint of track line are adopted in the proposed mosaic method.The proposed method removes target distortion, target ghosting and accumulated positioning error in the mosaicked image, but can only achieve the position accuracy of the image which is consistent with that of the track lines.To further improve the The location estimation method is often adopted in SSS measurement, and has various positioning accuracy with vessel velocity and heading.To get accurate positions of towfish, accurate vessel-mounted GPS positioning and stable vessel velocity and heading are needed in SSS measurement.For the second method, the relative coordinates between towfish and GPS can be obtained by towfish-mounted USBL, and the absolute position of towfish can be calculated by combining GPS positioning solution.The accuracies of the two positioning methods are not high and lead to the target dislocation in a geocoded image and the target ghosting appears in the geocoding mosaicked images.To remove the problems, the common features and the constraint of track line are adopted in the proposed mosaic method.The proposed method removes target distortion, target ghosting and accumulated positioning error in the mosaicked image, but can only achieve the position accuracy of the image which is consistent with that of the track lines.To further improve the position accuracy of the mosaicked image, the external reference, such as high-accuracy underwater landmarks, can be adopted as a constraint in constructing the adjustment model, if these underwater landmarks exist in the measurement area.

Impact of Overlapping Ratios
The proposed method contains two steps: coarse mosaic and refined adjustments.The latter is only conducted in overlapping areas.In the above experiments, the overlapping ratio between adjacent SSS images is about 50%, which guarantees that enough CFPs in overlapping areas could be detected.When the overlapping ratio becomes smaller, less or even no CFPs could be found, which will affect the subsequent adjustment operation.To analyze the impact of different overlapping ratios, three SSS data sets are used.The first set includes three SSS images (A)-(C) marked in Figure 5 and the overlapping ratios are reduced to 40% and 30%, separately.The other two data sets (D and E) are from the SSS images measured in Bohai Sea, which have the original overlapping ratio 40%.The overlapping ratio is also reduced to 30% in the following experiment.
The image matching results with different overlapping ratios are shown in Figure 18, which indicates that when the overlapping ratio is 40%, enough CFPs can be detected by the SURF and RANSAN algorithms.However, when the ratio is reduced to 30%, almost no correct CFPs could be detected even when common features exist in the overlapping areas.The reason may be that, in addition to the correct CFPs in the overlapping areas, there also exist more wrong ones.The RANSAC algorithm is sensitive to the high outlier rate and unable to eliminate these wrong CFPs when they account for a large proportion [35,36].Thus, 40% is set as a typical overlapping ratio for the SSS measurement task [19,25].E) are from the SSS images measured in Bohai Sea, which have the original overlapping ratio 40%.The overlapping ratio is also reduced to 30% in the following experiment.
The image matching results with different overlapping ratios are shown in Figure 18, which indicates that when the overlapping ratio is 40%, enough CFPs can be detected by the SURF and RANSAN algorithms.However, when the ratio is reduced to 30%, almost no correct CFPs could be detected even when common features exist in the overlapping areas.The reason may be that, in addition to the correct CFPs in the overlapping areas, there also exist more wrong ones.The RANSAC algorithm is sensitive to the high outlier rate and unable to eliminate these wrong CFPs when they account for a large proportion [35,36].Thus, 40% is set as a typical overlapping ratio for the SSS measurement task [19,25].

Interpolation to Fill Gap
To fill gaps derived from the geometric transformation, the mean value of four closest nonblank pixels is used.To analyze the performance of the adopted interpolation method and the impact of gap size on the interpolation, an experiment is conducted.Meanwhile, the nearest interpolation method, which fills the gap with the gray level of the nearest pixel, and the linear interpolation method, which fills the gap by linear interpolation using the gay levels of its surrounding pixels, are introduced for the comparison with the adopted method.In Figure 19, four different gaps with the pixel size of 5 × 5 pixels, 10 × 10 pixels, 15 × 15 pixels, 20 × 20 pixels are shown in the 2nd column and the interpolation results using the three methods are respectively shown in the 3rd, 4th and 5th columns.The red rectangles denote the gap area.Comparing the gap areas in raw image with those in the interpolated ones, the three methods can obtain satisfying visual results when the gap size is small (e.g., 5 × 5 pixels and 10 × 10 pixels).When the gap size (e.g., 15 × 15 pixels and 20 × 20 pixels) becomes larger, vague visual effects appear in these interpolated gap areas.To compute the differences between the interpolated pixel values and the raw ones, the standard deviations are

Interpolation to Fill Gap
To fill gaps derived from the geometric transformation, the mean value of four closest nonblank pixels is used.To analyze the performance of the adopted interpolation method and the impact of gap size on the interpolation, an experiment is conducted.Meanwhile, the nearest interpolation method, which fills the gap with the gray level of the nearest pixel, and the linear interpolation method, which fills the gap by linear interpolation using the gay levels of its surrounding pixels, are introduced for the comparison with the adopted method.In Figure 19, four different gaps with the pixel size of 5 × 5 pixels, 10 × 10 pixels, 15 × 15 pixels, 20 × 20 pixels are shown in the 2nd column and the interpolation results using the three methods are respectively shown in the 3rd, 4th and 5th columns.The red rectangles denote the gap area.Comparing the gap areas in raw image with those in the interpolated ones, the three methods can obtain satisfying visual results when the gap size is small (e.g., 5 × 5 pixels and 10 × 10 pixels).When the gap size (e.g., 15 × 15 pixels and 20 × 20 pixels) becomes larger, vague visual effects appear in these interpolated gap areas.To compute the differences between the interpolated pixel values and the raw ones, the standard deviations are shown in Table 7, which indicates that larger standard deviation appears when the gap size increases.Meanwhile, the linear interpolation method and the adopted method perform better than the nearest interpolation method.In experiments, the gap is generally less than 10 pixels or 6 m, and the adopted interpolation method can get accurate value to reflect the gray level of the gap.
pixels is used.To analyze the performance of the adopted interpolation method and the impact of gap size on the interpolation, an experiment is conducted.Meanwhile, the nearest interpolation method, which fills the gap with the gray level of the nearest pixel, and the linear interpolation method, which fills the gap by linear interpolation using the gay levels of its surrounding pixels, are introduced for the comparison with the adopted method.In Figure 19, four different gaps with the pixel size of 5 × 5 pixels, 10 × 10 pixels, 15 × 15 pixels, 20 × 20 pixels are shown in the 2nd column and the interpolation results using the three methods are respectively shown in the 3rd, 4th and 5th columns.The red rectangles denote the gap area.Comparing the gap areas in raw image with those in the interpolated ones, the three methods can obtain satisfying visual results when the gap size is small (e.g., 5 × 5 pixels and 10 × 10 pixels).When the gap size (e.g., 15 × 15 pixels and 20 × 20 pixels) becomes larger, vague visual effects appear in these interpolated gap areas.To compute the differences between the interpolated pixel values and the raw ones, the standard deviations are shown in Table 7, which indicates that larger standard deviation appears when the gap size increases.Meanwhile, the linear interpolation method and the adopted method perform better than the nearest interpolation method.In experiments, the gap is generally less than 10 pixels or 6 m, and the adopted interpolation method can get accurate value to reflect the gray level of the gap.

Conclusions
To overcome the drawbacks of existing mosaic methods to obtain a large-scale SSS image, this paper proposes a novel method using CFPs with constraint of track line positions.The geocoding mosaic operation is firstly conducted using the SSS geocoded image to form the first-step image.Next, the overlapping areas between different strips are defined and segmented.In the segmented areas, SURF and RANSAC algorithm are used to detect CFPs.When abundant CFPs are detected, an adjustment model is established using these detected ones with the constraint of corresponding track line positions in every segmented part.The adjustment model takes full advantage of CFPs and track line positions to accommodate both local distortion and global stability, which can not only weaken or eliminate local distortions around CFPs, but also can ensure that the transformed sensed image does not deviate from the current track line.
Experimental results proved the effectiveness of the proposed method.The adjacent image mosaic result showed the CFP coordinate deviations were less than 8 pixels and the track line position variations were less than 1 pixel.The statistical parameters of the multi-strip mosaicked image had the same accuracy with that of the adjacent image mosaics in the order of amplitude.Image drifting will not happen and positional accumulated errors is effectively avoided.Such a mosaicked image can

Figure 1 .
Figure 1.The procedure of image fusion using a multiresolution fusing method.

Figure 1 .
Figure 1.The procedure of image fusion using a multiresolution fusing method.

Figure 2 .
Figure 2. Determination and segment of the overlapping area.

Figure 2 .
Figure 2. Determination and segment of the overlapping area.

22 Figure 3 .
Figure 3.The flow diagram of SSS image mosaic process.

Figure 3 .
Figure 3.The flow diagram of SSS image mosaic process.
Remote Sens. 2018, 10, x FOR PEER REVIEW 8 of 22 decoding, radiation correction, slant-range correction and geocoding, the SSS geocoded images with the pixel resolution of 0.6 m and clear targets were obtained, as shown in Figure4, which can be used to reflect the seabed surface.

Figure 4 .
Figure 4. Four SSS geocoded images.(I) (II) (III) and (IV) are four SSS images with 50% overlapping 3.2.Mosaic Process Using Adjacent Strips Strip I and Strip II are used as the examples to illustrate the mosaic method according to the flow diagram shown in Figure 3 in detail.(1)Geocoding mosaic.Use the SSS geocoded images to obtain the mosaicked one in the same geographic frame using the method described in Section 2.1.The mosaicked image in Figure5ashows there exist target ghosting and dislocation which can be clearly seen as L1-L4.It can be found in Figure5that common features exist in the overlapping areas of adjacent SSS images, which are named TL1-TL4 in one SSS image and TR1-TR4 in the other one.The mosaicked image is supposed to show the four corresponding features in the same positions.However, affected by the wave and currents, the vessel velocity, the cable length variation and the topography changes as mentioned and analyzed in Section 1, the obtained SSS image feature coordinates may contain errors, which will result in the problem of target dislocation.As a result, the same feature may locate in more than one positions in the mosaicked image and the phenomenon is called target ghosting.In positions L1-L4 of Figure5a, ghosts accompany the four targets, which means that target ghosting and dislocation exist in these positions.

Figure 6 .
Figure 6.The CFPs detected in three segmented overlapping areas (a,b,c); no CFPs are detected in area (d); the S-L means the starboard side of left SSS image; the S-R means the starboard side of the right SSS image; the blue arrows mean the track line directions.

Figure 6 .
Figure 6.The CFPs detected in three segmented overlapping areas (a-c); no CFPs are detected in area (d); the S-L means the starboard side of left SSS image; the S-R means the starboard side of the right SSS image; the blue arrows mean the track line directions.

Figure 7 .
Figure 7.The distribution of CFPs and corresponding track line positions in area (A).

Figure 7 .
Figure 7.The distribution of CFPs and corresponding track line positions in area (A).

Figure 8 .
Figure 8. Coordinate deviations of CFPs."☆" and "▷" symbols denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○" and "◇" symbols denote these between the reference image and raw sensed image

Figure 8 .
Figure 8. Coordinate deviations of CFPs." " and " " symbols denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; " " and " " symbols denote these between the reference image and raw sensed image.

Figure 9 .
Figure 9. Coordinate deviations of the track line positions after adjacent image mosaics.

Figure 9 .
Figure 9. Coordinate deviations of the track line positions after adjacent image mosaics.

Figure 10 .
Figure 10.The comparison of image matching results using the MI and CR and the proposed method.(I) and (II) are three image matching results using the MI and CR and the proposed method.

Figure 11 .
Figure 11.Transformed sensed image with (a) and without (b) constraint of track line positions.The black lines denote raw positions and the red lines denote transformed positions.

Figure 10 . 22 Figure 10 .
Figure 10.The comparison of image matching results using the MI and CR and the proposed method.(I) and (II) are three image matching results using the MI and CR and the proposed method.

Figure 11 .
Figure 11.Transformed sensed image with (a) and without (b) constraint of track line positions.The black lines denote raw positions and the red lines denote transformed positions.

Figure 11 .
Figure 11.Transformed sensed image with (a) and without (b) constraint of track line positions.The black lines denote raw positions and the red lines denote transformed positions.

Figure 12 .
Figure 12.The mosaicked image using the proposed method.

Figure 12 .
Figure 12.The mosaicked image using the proposed method.

Figure 13 .
Figure13.Coordinate deviations of CFPs."☆" and "▷" symbols denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○" and "◇" symbols denote these between the reference image and raw sensed image.

Figure 14 .
Figure 14.Coordinate deviations of the track line positions after multi-strip image mosaics.

Figure 13 .
Figure13.Coordinate deviations of CFPs." " and " " symbols denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; " " and " " symbols denote these between the reference image and raw sensed image.

Figure 13 .
Figure13.Coordinate deviations of CFPs."☆" and "▷" symbols denote the coordinate deviations in east-west and north-south directions between the reference image and transformed sensed image; "○" and "◇" symbols denote these between the reference image and raw sensed image.

Figure 14 .
Figure 14.Coordinate deviations of the track line positions after multi-strip image mosaics.

Figure 14 .
Figure 14.Coordinate deviations of the track line positions after multi-strip image mosaics.

Figure 16 .
Figure 16.Image registration operation for the raw SSS images (a) and radiation corrected ones (b).

Figure 16 .
Figure 16.Image registration operation for the raw SSS images (a) and radiation corrected ones (b).

Figure 17 .
Figure 17.The schematic diagram to obtain the track line positions of towfish.

Figure 17 .
Figure 17.The schematic diagram to obtain the track line positions of towfish.

Figure 18 .
Figure18.The image matching results with different overlapping ratios: (A-C) are from Figure5; and (D,E) are from the SSS images measured in Bohai Sea.

Figure 18 .
Figure 18.The image matching results with different overlapping ratios: (A-C) are from Figure 5; and (D,E) are from the SSS images measured in Bohai Sea.

Figure 19 .
Figure 19.The performances of three interpolation methods to fill gap with different sizes.Figure 19.The performances of three interpolation methods to fill gap with different sizes.

Figure 19 .Table 7 .
Figure 19.The performances of three interpolation methods to fill gap with different sizes.Figure 19.The performances of three interpolation methods to fill gap with different sizes.

Table 1 .
Statistical parameters of coordinate deviations of 60 CFPs.

Table 1 .
Statistical parameters of coordinate deviations of 60 CFPs.

Table 2 .
Statistical parameters of coordinate deviations of 1200 track line positions after adjacent image mosaics.

Table 2 .
Statistical parameters of coordinate deviations of 1200 track line positions after adjacent image mosaics.

Table 3 .
The computing time of different kinds of image matching algorithms (Unit: second).

Table 4 .
Statistical parameters of coordinate deviations between raw and transformed track line positions with and without constraint.

Table 4 .
Statistical parameters of coordinate deviations between raw and transformed track line positions with and without constraint.

Table 4 .
Statistical parameters of coordinate deviations between raw and transformed track line positions with and without constraint.

Table 5 .
Statistical parameters of coordinate deviations of 45 CFPs and 1200 track line positions after multi-strip image mosaics.

Table 6 .
The computing time and statistical parameters of coordinate deviations of 1200 track line positions after multi-strip image mosaics using different mosaicking methods.