Reverberation presents a persistent challenge for active sonar operations in shallow water, primarily due to the generation of target-like clutter within the reverberation signal, which can trigger false alarms and significantly impair sonar system performance (Gruden et al., 2021). The origins of such clutter are varied, including discrete and buried objects on the seabed (Prior, 2005; Ratilal et al., 2005; Holland et al., 2007), non-discrete seabed structures (Holland and Ellis, 2012), fish schools (Weber, 2008), and oceanographic phenomena like nonlinear internal waves (NIWs) (Henyey and Tang, 2013).

Zhou (Zhou et al., 1991) examined the impact of soliton internal waves on shallow-water sound propagation, identifying “acoustic mode coupling” as the cause of abnormal frequency responses. John and Stanley (John et al., 1994) observed broadband fluctuations in a 1000-kilometer acoustic pulse experiment in the Pacific Ocean, suggesting internal waves as a plausible explanation. John (John and Michael, 1998) later proposed a numerical method for simulating sound speed perturbations caused by random internal waves, though this model was limited in its applicability to certain regions. (Dirk et al., 1997) developed a model for sound wave propagation in random shallow-water waveguides, based on experiments off the New Jersey coast, treating sound velocity profile variations as stochastic processes. (Lynch et al., 2010) found that internal wave curvature significantly affects modal amplitudes and arrival angles, but the underlying scattering mechanisms remained elusive.

Traditional methods for measuring the velocity of internal waves include the use of acoustic Doppler current profilers (ADCP), fixed buoys, and satellite remote sensing techniques. Shipborne ADCP (Sun and Shen, 2010) detects internal wave information based on the strength of the sea water echo signal and acoustic Doppler information, requiring a combination of electronic compass and high-precision Global Positioning System (GPS) for correction to calculate the velocity of ocean currents. Marine monitoring buoys are one of the most reliable and efficient platform monitoring methods in routine marine environmental monitoring. Marine monitoring buoys (Chen et al., 2020) can, to a certain extent, be unaffected by adverse weather conditions at sea, continuously acquiring marine environmental data over a long period, and have the advantages of long-term, continuous, and all-weather automatic observation. They are an important, reliable, and stable means in marine observation. Buoys measurements (Chun et al., 2021) rely on the influence of internal waves on marine environmental parameters, coupled with real-time monitoring by sensors such as ADCP, Conductivity-Temperature-Depth profiler (CTD), Temperature-Depth profiler (TD), and current meters. Data collected are then subjected to inversion analysis to deduce the velocity of internal waves. Despite the extensive coverage required for accurate measurements, deploying a network of buoys incurs significant costs. However, they have the disadvantages of high cost and limited observation range. Measuring the velocity of internal waves using synthetic aperture radar (SAR) is also one of the commonly used methods (Chong and Zhou, 2013). Since the 1970s, various bands and polarizations of airborne SAR and spaceborne SAR have obtained a large number of internal wave images, providing extensive 2D information, which has formed a strong supplement to on-site measurements and optical observation methods, providing a rich source of data for internal wave detection and becoming an important remote sensing method for marine internal wave observation. SAR images and other hydrographic data are used to extract hydrodynamics parameters such as the depth, velocity, wavelength, and amplitude of internal waves. Furthermore, the Long Baseline Forward (LBF) method (Zhang et al., 2024) can improve the imaging performance and efficiency issues of traditional imaging algorithms; Zeng (Zeng et al., 2024) by analyzing the internal wave fields simulated using the MIT General Circulation Model (MITgcm), satellite Synthetic Aperture Radar (SAR) observations, and moored temperature-salinity-depth (TSD) chain observations, the source areas and initial propagation times of internal waves can be identified. Discussing the reasons for the generation of different types of internal waves in the region can help improve the understanding of internal wave phenomena in the South China Sea. However, the use of SAR may be affected by weather and has the disadvantage of lower data update frequency, resulting in poor real-time performance.

Henyey and Tang (2013) numerically demonstrated that NIWs can produce significant target-like clutter signals, exceeding the mean reverberation level (RL) by more than 10 dB. They noted that these NIWs deflect acoustic rays to higher grazing angles, leading to strong, slowly varying clutter signals. They also observed a general increase in RL following the arrival of clutter, although the underlying cause was not explored. The dynamic nature of shallow-water waveguides makes the impact of moving NIWs on distant reverberation a compelling area of study.

Conventional monitoring of ocean internal waves has relied on fixed-point measurements from anchored temperature chains. This study experimentally observes the abnormal oscillation of reverberation intensity caused by soliton waves (a single packet form of NIWs) in shallow seas. Theoretical analysis in section 2 supports the notion that soliton waves can induce target-like clutter in reverberation, validated by experimental data. Section 3 delves deeper into the theoretical implications, revealing that moving soliton internal waves trigger a quasi-periodic oscillation in reverberation intensity post-clutter, with the dominant frequency linked to the internal wave’s velocity. This provides a theoretical basis for determining internal wave velocity through active reverberation telemetry, enhancing underwater-acoustic detection of internal waves. However, further research is required to refine the accuracy of internal wave velocity measurements through post-clutter reverberation intensity oscillations. Sections 4 and 5 present the experimental setup and discussion on shallow-water reverberation clutter.

Researchers have advanced multiple reverberation models, such as those by (Ellis, 1995; Grigor’ev et al., 2004; Gao et al., 2010; Tang and Jackson, 2012) with a prevalent approach being the normal mode reverberation model that incorporates mode coupling in forward propagation as detailed by Yang et al (Yang, 2014). However, these models often rely on the simplistic Lambert scattering model for seafloor effects, neglecting the more complex coupled scattering phenomena. For a nuanced examination of seafloor topography’s influence on scattering, innovative methodologies like the coupled reverberation mode, as introduced by Gao et al., are essential. Additionally, Figure 1A is the block diagram of this model.

This study considers bottom reverberation in a two-layer medium with a rough bottom profile characterized by H = H 0 + ζ ( r → ) , where H 0 is horizontal sea bottom depth and rough seabed fluctuation ζ is much smaller than H 0 . The density and sound speed of the upper and lower layers are represented by ρ 0 , c 0 and ρ 1 , c 1 , respectively, with b being the density ratio. k r n is the eigenvalue, k 0 z n , and k 1 z n are the vertical wave numbers in the two media, defined as k 0 z n = k 0 2 − k r n 2 ,   k 1 z n = k 1 2 − k r n 2 . The pressure function φ ( r → ) , originating from a harmonic source at depth z s , is defined as description in Yang (2010)’s Equation 4 and Gao et al (2010)’s Equation 1.

where the local mode (vertical) function   Z n is expressed as

and

Additionally, the element G m n of the coupled matrix G , which is applicable to a rough bottom, is defined as follows:

Given that bottom reverberation arises from first-order perturbations at the seabed interface, the horizontal function ψ n ( r → ) is formulated as

In the aforementioned equations, r   → represents the horizontal vector from the acoustic source to the receiver. The vector r ′ → signifies the horizontal vector from the acoustic source to the seabed scattering point. The vector R → corresponds to the horizontal vector from the seabed scattering point to the receiver. Equation 7 characterizes the horizontal factor of the acoustic field’s potential function. Among them, − 2 k z n 2 B n H 0 ( 1 ) ( k r n r ′ → ) represents forward propagation, ∑ m k r m Γ n G m n H 1 ( 1 ) ( k r m r ′ → ) ( x ′ r ′ → ∂ ζ ∂ x ′ + y ′ r ′ → ∂ ζ ∂ y ′ ) indicates modal coupling due to seabed fluctuations and H 0 ( 1 ) ( k r n R → ) signifies backscattering. In addressing the influence of water column inhomogeneities, such as the presence of NIWs, adjustments must be applied to the vertical mode function Z n . This necessitates the employment of mode coupling, as detailed by Yang (2014), to incorporate the water column’s coupling matrix. For the initial modeling of an NIW packet, we adopt the soliton solution of the Korteweg–de Vries (KdV) equation (Zheng et al., 2001), which is given by

In this context, the hyperbolic secant function is represented by sech, with η ^ signifying the amplitude of the NIW packet. The variable r corresponds to the horizontal coordinate, r x pinpoints the center position of the soliton, and Δ characterizes the width of the soliton.

In the modeling of pressure perturbations when the water column disturbs, the incident mode vector is characterized by Yang (2014)‘s Equation 2:

where T ( r 1 , 0 ) = d i a g { e i k r ( r 1 ) } and T ( R , r p + 1 ) = d i a g { e i k r ( R − r p + 1 ) } represent the mode propagation matrices. Here, A ^ ( R ) denotes the incident mode vector at a specific range R, while A ^ ( 0 ) is the source mode vector, which is associated with a standard eigenvector given by [ φ 1 ( z s ) , φ 2 ( z s ) , … , φ n ( z s ) ] H and is formulated as follows:

Within the mode coupling range from r k to r k + 1 , the entry of the coupling matrix Λ is expressed as:

In the context of mode coupling, the element Λ n m within the range r k to r k + 1 is defined using the Kronecker delta δ n m , with δ c representing the speed perturbation induced by the NIW and c 0 being the ambient sound speed in the absence of NIWs.

The Sound Speed Profile (SSP) alteration due to NIWs is numerically computable, as demonstrated by Henyey and Tang (2013). The simulation employs an environmental model incorporating an NIW, depicted in Figure 1 , to model an NIW packet as coherent structures traversing the acoustic path without deformation. Although the KdV equation-based mathematical model may not precisely describe the leading wave, the primary wave-packet is predominant in causing acoustic clutter, as highlighted by Henyey and Tang (2013), suggesting that the physical mechanisms and observed phenomena are analogous.

Subsequently, the mode coupling within the water column is integrated into the coupled mode reverberation theory. Traditionally, the coupled matrices in seabed reverberation theory are influenced solely by bottom roughness. The horizontal Hankel factor H 0 ( 1 ) and the vertical eigenfunctions for the source s i n ( k 0 z n z s ) and receiver s i n ( k 0 z n z ) are combined and derived into the vertical function Z n ( z , k r n ) and the horizontal function. When accounting for the inhomogeneous water column, these three terms are redefined within the incident mode vector A ^ ( R ) , encapsulating the source’s eigenfunction and the Hankel function’s phase term, the receiver’s eigenfunction, and the spreading factor 1 k r n R , which approximates the Hankel function’s amplitude. This framework facilitates a clear description of the coupling process within the incident mode vector A ^ ( R ) .

Incorporating the NIW-induced incident mode coupling into the propagation of the coupled mode theory for seabed reverberation, the vertical function in Equation 2 is subsequently adjusted to reflect these considerations:

where A ^ n ( r , z s ) denotes the n-th component of the incident mode vector A ^ ( R ) as defined by Equation 9. The horizontal function ψ n ( r → ) retains the form given in Equation 7. Meanwhile, the element G ^ m n of the coupled matrix G ^ , as presented in Equation 8, is transformed to become

Consequently, the pressure function φ ^ ( r → ) , as indicated in Equation 1, can be reformulated as follows:

During the numerical simulation, Equation 7 is utilized to calculate the Reverberation Level (RL). To simplify the scenario, a mono-static condition is assumed, and the partial derivative of the isotropic roughness ζ is taken as zero and transformed the area integration into line integration. This leads to the following expression:

Furthermore, Equation 15 is expanded to its full form as follows:

Thus, an amended expression for the RL that accounts for the presence of NIWs is derived as Equation 18.

A shallow-water reverberation experiment was conducted in the North Yellow Sea from July 8 to 10, 2014. The water depth averaged 45 m, with a low-frequency transducer deployed at a depth of 22 m serving as the sound source. A vertical line array (VLA), comprising 16 element pressure hydrophones, was used to receive the signals. The hydrophones were spaced 3 meters apart, and both the sound source and the VLA were situated on the same vessel to establish a monostatic reverberation testing configuration. Environmental measurements indicated bottom parameters of ρ 2 = 1850 k g / m 3 and c 2 = 1770 m/s. A variety of signals, including CW and Chirp signals, were transmitted starting at noon local time on July 8th.

The purpose of the experiment was to explore the impact of internal waves on shallow-water reverberation. To effectively monitor internal waves, two temperature chains equipped with 16 temperature and depth (TD) sensors each were deployed near the research vessel, with sensors spaced 0.5 m apart. These TD sensors, measuring 150 mm in length and 40 mm in diameter, were capable of measuring temperatures ranging from -2°C to 40°C. The chains were anchored perpendicular to the isobaths, considering that the propagation of NIWs is generally aligned with the isobaths’ vertical direction and influenced by tidal forces. This arrangement ensured that soliton waves, akin to plane waves, would sequentially pass through both temperature chains.

The first temperature chain (Tc-A) was positioned 723 m from the vessel, North Survey I, which has a displacement of 1200 tons, while the second chain (Tc-B) was 1037 m away. The vessel was anchored roughly midway between the two chains. Figure 6A and Figure 6B illustrate a series of NIWs detected by Tc-A and Tc-B. Despite some packets’ shapes changing as they propagated, the initiation point of each packet allowed us to determine that the propagation time across both temperature chains was approximately 1.5 hours (5400 s), with a distance of 1760 m between Tc-A and Tc-B. Consequently, the speed of the soliton wave was calculated to be v r = 0.325 m/s, which aligns with historical data of internal wave speeds in the northern Yellow Sea, ranging from 0.3 m/s to 0.4 m/s (Hu et al., 2020).

A pronounced soliton wave, identified from the sequence of the NIW packets in Figure 7A , was selected for analyzing the clutter effects. This wave was recorded by Temperature Chain-A (Tc-A) at 14.82 hours after noon, as depicted in Figure 7A . Concurrently, at 14.115 hours after noon, a series of CW signals, each 0.2 s in duration, were emitted. The reverberation intensity captured by the 8th, 11th, and 12th hydrophones is displayed in Figure 7B , revealing the emergence of target-like clutter at 2.047s.

The internal wave illustrated in Figure 7A is capable of inducing clutter in the signals received by the test ship. With the soliton wave’s speed established at 0.325 m/s, the distance from the wave to the sound source at 14.115 hours was approximately 1549 m (723 + 826 m). According to the theoretical calculations based on Equation 18, the clutter was expected to appear at a time of 1549 m × 2 1500 m / s = 2.065 s .

As indicated in Figure 7C , the theoretical clutter appearance time of 2.065 seconds closely matches the measured time of 2.047 s. A minor discrepancy between the theoretical and measured clutter intensities is observed, attributable to the lack of a second temperature chain to verify the soliton wave’s precise structure at the 1549 m mark from the ship during the signal transmission at 14.115 hours. The theoretical model’s soliton wave structure is inferred from the Tc-A data recorded at 14.82 hours. Despite the half-hour time difference, the soliton wave’s structure is expected to have undergone slight changes, which is the primary source of the mismatch. Notably, by the time this strong soliton wave reached Tc-B after 1.5 hours, its amplitude had altered significantly, as evident from Figure 6B . Initially, the authors were uncertain about the reliability of estimating the soliton wave structure at 14.115 hours, but the comparison in Figure 7C confirms the estimation’s accuracy. Furthermore, by analyzing data from various time points, it can be deduced that internal waves can induce the generation of clutter, as illustrated in Figure 7D .

It is important to clarify that the clutter is not a result of direct reflection of underwater acoustic wave by the soliton wave. Instead, the clutter effect leads to enhanced mode coupling within the forward-propagating sound field. This implies that energy from lower modes (with lower incident grazing angles) is transferred to higher modes (with higher grazing angles). Notably, another set of solitary internal waves was detected by Tc-A at 16.5 hours, as shown in Figure 8A . At 16.27 hours, the test vessel emitted 30 sets of single-frequency signals, each consisting of ten pulses with widths of 0.2 s, 0.5 s, and 1.0 s, at a frequency of 580 Hz.

At the time of observation, the solitary internal wave was positioned 210.9 m from Tc-A. Based on the established speed of internal waves in the region, the anticipated arrival time of the clutter was determined to be 1.245 s. The testing vessel transmitted a 580 Hz signal, and the data from the 9th hydrophone, located at a depth of 19.5 m, were analyzed, as illustrated in Figure 8B . The clutter was detected at 1.527 s; after accounting for the direct sound travel time of 0.2 s, the clutter’s actual occurrence time was 1.327 s. This aligns well with the predicted clutter time of 1.245 s, resulting from the internal wave’s influence, with a minor discrepancy of 0.082 s. Figure 8B shows that, at this depth, clutter associated with different pulse widths (0.2 s, 0.5 s, 1.0 s) of the CW signal is observable at the expected time intervals, verifying that the clutter was indeed induced by the moving solitary internal wave.

For the analysis of reverberation data with clutter occurrences at different pulse widths, data segments corresponding to 0.2 s, 0.5 s, and 1.0 s pulse durations were examined. As depicted in Figure 9 , for the (A) 0.2 s pulse width, the theoretical prediction for clutter induction by the internal wave was at 1.2074 s, with the actual measurement at 1.2045 s, yielding a negligible error of 0.0029 s. For the (B) 0.5 s pulse width, the calculated clutter time was 1.1289 s, and the measured time was 1.1254 s, with a slight error of 0.0035 s. Lastly, for the (C) 1.0 s pulse width, the expected clutter time was 0.9136 s, against an observed time of 0.9733 s, resulting in a larger discrepancy of 0.0597 s.

When plotting these three curves together, it becomes evident that the clutter position shifts as the internal wave moves, a movement referred to as the “drift” phenomenon, as shown in Figure 9D .

In conjunction with the theoretical framework from Section 2, the reverberation signals from 30 sets of 580 Hz single-frequency transmissions, influenced by the moving internal waves at the specified depth, were averaged at a fixed time. The difference from reverberation signals not affected by the internal waves was calculated to ascertain the change in reverberation intensity at the time of observation. Time-domain interpolation was applied to the intensity change data based on the transmission timing, creating a time-domain reverberation intensity change curve spanning 846 s. A Fourier transform was subsequently applied to convert this into the frequency domain, revealing how the reverberation intensity changes over various time spans.

The relationship between the frequency f m n and the horizontal wavenumber k m n is given by f m n = v r k m n / 2 π , which can be rearranged to k m n = 2 π f m n v r , where v r denotes the velocity of the internal wave. The empirical value of k m n was found to be 19.3 × f m n .

The time-domain and frequency-domain representations of the reverberation signal, processed with a 0.2 s window, were generated using data from 6.0 s to 6.2 s and from 6.3 s to 6.5 s.

The reverberation signals depicted in Figures 10A and 10B were captured between 6.0 and 6.2 s. These signals are associated with the horizontal wave numbers as simulated by the Kraken model, detailed in Table 1 .

Figure 10C and Figure 10D display the reverberation signals recorded in the time interval from 6.3 to 6.5 s. These signals correspond to the horizontal wave numbers as simulated by the Kraken model, with the specific values and comparisons provided in Table 2 .

The reverberation signals presented in Figures 10E and 10F were obtained using a 0.5-second time window, with the data extracted from the time span of 6.0 to 6.5 s. These signals are aligned with the horizontal wave numbers from the Kraken simulation, as detailed in Table 3 . Additionally, the reverberation signals from 6.5 to 7.0 s, processed with the same time window, are also considered as part of the analysis.

Figures 10G and 10H illustrate the reverberation signals acquired between 6.5 and 7.0 s. These signals are matched with the horizontal wave numbers from the Kraken simulation results, which are detailed in Table 4 .

The reverberation signals depicted in Figures 10I and 10J were processed using a 1.0 s window and were extracted from the time interval ranging from 6.0 to 7.0 s. These signals correspond to the horizontal wave numbers obtained from the Kraken simulation, as presented in Table 5 .

The analysis, in conjunction with the simulation experiment, reveals that the principal frequency component of the reverberation intensity oscillation at a constant observation time is correlated with the interference interval between a pair of normal waves.

Thus, it can be demonstrated from the correspondence between data processing and simulation results that internal wave imprint target information on the reverberation curve, and that internal wave velocity can be extracted from the reverberation oscillations in the frequency domain. Compared with traditional methods of measuring internal wave speeds, using reverberation to extract internal wave velocity eliminates the need for high-cost ADCP transect testing using high-frequency signals, as well as traditional long-term fixed-point observations by buoys and extensive SAR data inversion. This approach offers more proactive, intuitive, and convenient features, with lower implementation costs, requiring only a transmit-receive sonar. By utilizing reverberation, which is traditionally regarded as interference in active sonar, the detection of internal waves is achieved.

This research conducted a shallow-water reverberation and clutter experiment in the northern Yellow China Sea in July 2014, focusing on the analysis of clutter caused by NIWs through reverberation data. The findings confirm that the soliton wave’s movement is indicative of target-like clutter, as evidenced by the wave’s temporal progression. Theoretical insights were gained by applying the coupled-mode approach (Yang, 2014) to calculate and interpret the reverberation and clutter dynamics with a moving NIW. The simulations indicated that NIWs, particularly those with substantial amplitude, generate significant and target-like reverberation clutter, corroborating the observations by Henyey and Tang (2013). The study also highlighted the quasi-periodic oscillation in reverberation intensity following the clutter, which is influenced by the soliton’s motion. The Fourier spectrum analysis of these oscillations revealed primary peaks that aligned with the mode-coupling theory’s predictions, suggesting a dependency on the NIW’s position. The main frequency of reverberation oscillation during fixed observation periods was found to be linked to the velocity of the internal wave movement. While this study offers plausible explanations, further validation is necessary using more refined experimental data.