Clutter suppression and classification using twin inverted pulse sonar (TWIPS)

T. G. Leighton, D. C. Finfer, P. R. White, G.-H. Chua, J. K. Dix


This paper describes the detection and classification of targets against clutter by distinguishing between linear and nonlinear scatterers and, further, by distinguishing those nonlinear targets that scatter energy at the even-powered harmonics from those that scatter in the odd-powered harmonics. This is done using twin inverted pulse sonar (TWIPS), which can also, in some manifestations, require no range correction (and therefore does not require the a priori knowledge of the environment needed for most remote detection technologies). The method applies, in principle, to a range of sensor technologies, including the use of radar to distinguish between circuitry, metal and soil; Light Detection and Ranging (LIDAR) to detect combustion products; and Magnetic Resonance Imaging (MRI). A sonar application is demonstrated, detecting objects in bubbly water (including in the wake of a ship of 3953 gross register tonnage). A man-made sonar that can operate in bubbly water is relevant: Cold War sonar is not optimized for the shallow coastal waters that typify many current operations. The US Navy use dolphins in such waters. TWIPS arose as a demonstration that echolocation was possible in bubbly water in response to a video showing dolphins generating bubble nets when hunting: if echolocation were impossible in these nets, then during this hunt, the dolphins would have blinded their sonar.

1. Introduction

The decades of sonar experience built up for deep-water applications during the Cold War are not ideally suited for many of the shallow-water scenarios of interest today (Leighton 2008). Although active sonar has seen enormous development, in particular with respect to high-frequency imaging, these advances are confounded by the clutter produced by oceanic bubble clouds (as generated, for example, by ship wakes, breaking waves or biogenic and geophysical activity). This paper demonstrates a new form of sonar, twin inverted pulse sonar (TWIPS), and shows that it can be superior to the standard sonar in the detection of targets in bubbly water (using laboratory tank tests), and in distinguishing between bubbles and other scatterers (using sea trials). Such enhancements of both detection and classification of targets in bubbly waters are key goals of shallow-water sonar. Furthermore, the underlying physics could be exploited with a range of other radiations (e.g. radar for the detection of improvised explosive devices or covert circuitry).

2. Theory

TWIPS exploits the nonlinearity in the response of a bubble to a sound field that drives it to pulsate, such that the bubble radius R(t)=R0+Rε(t) is characterized by a displacement Rε(t) about an equilibrium position R0. The bubble responds as a lightly damped oscillator, the inertia being invested mainly in the liquid, which must move as the bubble wall pulsates, and the stiffness being imparted by the gas contained within the bubble. One source of the nonlinearity can be seen by assuming that the bubble contains only a volume V of gas at pressure pg, which behaves polytropically, such thatEmbedded Image 2.1 where κ is the so-called polytropic index. This varies between γ (the ratio of the specific heat of the gas at constant pressure to that at constant volume) and unity, depending on whether the gas is behaving adiabatically, isothermally, or in some intermediate manner. It can account for the non-reversible heat transfer (i.e. net losses) across the bubble wall only if its value varies throughout the acoustic cycle, and the constant value used most frequently therefore only describes reversible heat transfer. Neglecting, for a moment, surface tension and vapour effects, in the absence of an insonifying field, the net force on the bubble wall arises from the difference between pg and the constant static pressure in the liquid outside the bubble, p0 (Leighton 1994). If a positive wall displacement occurs (Rε>0), the restoring force Fr is towards the origin (the bubble centre), i.e. Fr=4πR2(pgp0)<0. The stiffness of the gas is therefore given byEmbedded Image 2.2 (using (2.1)). Even when the small-amplitude approximation |8πR(p0pg)|≪|12πκpg(t)R(t)| is applied in equation (2.2), the stiffness of the gas in the bubble is not constant over the oscillation cycle, as it depends on the instantaneous values (pg(t),R(t)) of the internal gas pressure and the bubble radius (and of κ if it varies). The dynamics are therefore inherently nonlinear, and the well-known formulations for the linear model of a bubble only appear if the time-varying quantities (pg(t), R(t)) are replaced by the ‘equilibrium’ values, i.e. the equilibrium internal gas pressure pg,e, which occurs in the bubble when the radius takes the equilibrium value R0 (and the magnitude of the wall speed, and the kinetic energy of the liquid, are both maximal at Embedded Image and Embedded Image respectively), with the common assumption of constant κ. This approximation allows the stiffness to take a constant value in the model, by substitution of equation (2.1) into equation (2.2) as follows:Embedded Image 2.3 If the bubble is assumed to pulsate with small-amplitude simple harmonic wall motion in an infinite body of liquid of equilibrium density ρ0, the maximum kinetic energy Embedded Image associated with that pulsation can be expressed in terms of Embedded Image. This is done by integrating up (from the bubble wall to radius Embedded Image the kinetic energies of concentric spherical shells of liquid, noting that if the liquid is incompressible, then Embedded Image,Embedded Image 2.4 which allows a statement of the effective inertia of the pulsation (mradRF) for such linear oscillations. The undamped natural circular frequency of bubble pulsation is found from the ratio of the stiffness (2.3) to the inertia (2.4) as follows:Embedded Image 2.5 where the so-called ‘Minnaert’ frequency is found by assuming adiabatic conditions (κ=γ) and neglecting the vapour pressure (pv) and surface tension (σ) when replacing pg,e with Embedded Image (the static pressure in the liquid far from the bubble, usually equal to the sum of the atmospheric and hydrostatic pressures).

However, with increasing pulsation amplitude, the small-amplitude assumptions that underpin equations (2.3) to (2.5) become less germane, and TWIPS exploits the resulting nonlinearities. The discrepancies between the time-varying quantities pg(t) and R(t) in equation (2.2) and their constant counterparts in equation (2.3) reveal the amplitude dependence of the stiffness term. Consequently, the nonlinear nature of the bubble pulsation becomes increasingly apparent as the amplitude of pulsation grows, usually as a result of a high-amplitude driving pressure signal p(t) or closeness to resonance. These nonlinearities are commonly expressed in the following Keller–Miksis equation (Keller & Miksis 1980):Embedded Image 2.6 where c0 is the speed of sound in the liquid for waves of infinitesimal amplitude, μ is the shear viscosity of the liquid and pL(t) is the liquid pressure on the external side of the bubble wall, which is related to the bubble gas pressure pg(t) via the vapour pressure pv (Leighton 1994)Embedded Image 2.7 In the family of equations of dynamics for the pulsating bubble, of which equation (2.6) is a member, terms resembling Embedded Image are associated with the inertia (as can be seen by differentiating the term Embedded Image in equation (2.4) with respect to R, noting that Embedded Image. In equation (2.6), the terms on the right-hand side incorporate the stiffness, as in equation (2.3), through the gas pressure (2.7), while the viscous and radiation damping arise, respectively, through the terms involving μ and c0. Thermal damping is not included.

The experiments in this paper are supported by simulations based on this theory. Building on earlier TWIPS simulations (Leighton et al. 2005a,b), equation (2.6) is used to generate the bubble radial response, from which the pressure radiated by an individual bubble is calculated at range r as equal to Embedded Image, where tr is the retarded time t−(r/c0) (Leighton 1994). The simulations are conducted by calculating the pressure radiated to the hydrophone position by representative bubbles, each of which characterizes the bubbles in the size bin encompassing bubbles of similar radii. The three-dimensional volume of liquid is divided into spatial cells into which various bubble populations, and the target, may be placed. The emission from all the bubbles from that volume is calculated by convolving the bin-representative emissions with the bubble-size distribution assumed for that volume. The bubble responses are assumed to be uncoupled (valid for the concentrations modelled here; Leighton et al. 2004), and the bubble cloud is assumed not to evolve between the two pulses in a pair (but will evolve between pulse pairs, for example under buoyancy, turbulence or source motion). The losses during transmission to and from the bubble cloud are imposed by applying the attenuation that would be given by linear bubble pulsations (Commander & Prosperetti 1989), the attenuation of water/seawater (as appropriate; Francois & Garrison 1982a,b) and geometrical losses.

These simulations are used to predict the echo scattered by a bubble cloud when subjected to a driving pressure p(t), which uses two consecutive pulses, the second pulse being a phase-inverted replica of the first. The principle by which this TWIPS operates is shown schematically in figure 1 (Leighton 2004a,b). Consider the following problem scenario: conventional sonar fails to detect a linearly scattering body (the ‘target’, e.g. a fish, mine or seabed) because the returned sonar signal is cluttered by the scatter from bubble clouds. If the insonifying field had sufficient amplitude to generate a nonlinear response from the bubbles, it might be possible to enhance scatter from the target while simultaneously suppressing it from the bubbles, allowing the sonar operator to distinguish between the two (‘classification’). The insonifying TWIPS field p(t) has the following time series (figure 1(i)):Embedded Image 2.8 that is, a pulse containing two components based on a pressure function Γ(t); the second component starting at time Δ after the first and having opposite polarity to it (figure 1(i)). It is assumed that the pressure echo at the receiver can be regarded as the superposition of responses from individual scatterers (valid for the less than 10−4 void fractions modelled here (Kargl 2002; Leighton et al. 2004)). The received signal from an object that scatters linearly (plin(t); figure 1a(ii)) can be related to the emitted pulse p(t) via a convolution integral, i.e.Embedded Image 2.9 Here, the kernel function (or impulse response) hlin(t) incorporates the effects of propagation (assumed to be linear) and scattering from the object. The symbol Hlin[] is used to represent the operator that computes, from its argument, the signal received from the linear scatterer.

Figure 1.

Schematic of the formation of the twin inverted pulse sonar (TWIPS) signal in ideal circumstances (i.e. in the absence of noise, pulse overlap, etc.). The use of a Volterra series expansion allows incorporation of dynamic features (such as bubble motion once the driving pressure has ceased) that simpler expansions (e.g. the power-series expansion used in Leighton et al. 2008a) will not.

The signal pnl(t) received from an object that scatters nonlinearly, such as a bubble, cannot be expressed in a simple universal form. One approximation to a general nonlinear function is provided by the Volterra series (Shetzen 1980). Such a representation can be regarded as a generalization of the familiar Taylor series expansion of a memory-less nonlinear function. A Volterra series representation is not completely general, as there are classes of nonlinear functions that it cannot represent, but previously it has been successfully used to model single bubble scattering (White et al. 1997). In this representation, the signal from the nonlinearly scattering signal (figure 1b(ii)) is given byEmbedded Image 2.10 where τ, τ1, τ2 and τ3 are dummy variables within the integrals. Here, hk(.,.,⋯ ,.) and Hk[] are referred to as the kth-order Volterra kernels and functionals, respectively. The Volterra kernels can also be related to harmonic generation in the nonlinear system. In particular if k is odd, then the kth-order functional gives rise only to odd harmonics of order less than or equal to k. Similarly, when k is even, only even harmonics of order less than or equal to k are generated. The functionals also satisfy symmetry properties that reflect the oddness and evenness of their argument k. Specifically, it is evident that Embedded Image 2.11and Embedded Image 2.12 A nonlinear scatterer insonified by a pulse p(t) of the form shown in (2.8) produces a received signal (figure 1b(ii)) that can be represented as follows:Embedded Image 2.13 The basis of TWIPS processing is to form linear combinations of delayed versions of the received signal. Specifically, one creates p(t) and p+(t) defined asEmbedded Image 2.14 in which pRx(t) is the signal at the receiver. In practice, p(t) and p+(t) are computed by dividing the received signal into two sections, delaying the first by Δ seconds, and then computing the time histories generated through the sums and differences of these sections. For linear scatterers, such processing leads to Embedded Image 2.15and Embedded Image 2.16 whereas the equivalent expressions for the nonlinear scatterers, based on a Volterra representation, are as follows: Embedded Image 2.17and Embedded Image 2.18 Failure to satisfy these ideal predictions (especially equation (2.15)) could be used with the appropriate radiation to diagnose significant disturbance or evolution of the scatterers between the two pulses. Even without nonlinearity, this could be used to detect seismic, sedimentary or biomedical displacements, or be used for fault detection. The quantity p+(t) is already used to increase the effectiveness of biomedical ultrasonic contrast agents (UCAs), microscopic-stabilized gas bubbles that are injected into the blood stream to enhance the scatter from the blood flow as opposed to the surrounding tissue (which scatters the ultrasound with far less nonlinearity; Simpson et al. 1999). It is particularly effective in such an application because the range of bubble sizes is narrowly constrained, so that most of the bubbles can be driven close to resonance to elicit nonlinearities from them.

For oceanic applications, the bubble population is not so constrained: the bubble-size distribution usually covers many decades of radius, and the problem is consequently more difficult. TWIPS makes use of the signal p(t), for which the even-powered components scattered from the bubbles are suppressed (figure 1b(iv)), while the reflections from the solid body (and the odd-powered components in the scatter from the bubbles) are enhanced (figure 1a(iv)).

In deriving equation (2.17), it has been implicitly assumed that the delay Δ is large enough such that the response to Γ(t) has decayed sufficiently so as to be negligible prior to the start of Γ(tΔ). This imposes a lower limit on the selection of Δ because one has to use an inter-pulse delay large enough so that all echoes from the first pulse have decayed before the second pulse is emitted. Conversely, one can regard the inter-pulse delay as defining a maximum range over which a TWIPS-based system can operate, the maximum range being given by c0Δ/2 (which corresponds roughly to 37.5 m for the 50 ms inter-pulse delay most frequently used here). For a given environment, there is also a maximum value of Δ that can be tolerated by a TWIPS system. This derives from the fact that TWIPS processing assumes that the field of scatterers, for example a bubble cloud, does not evolve significantly during the interval between the start time of the two pulses (note that interval is the sum of the inter-pulse delay and the pulse duration). The signals scattered from the nonlinear scatterers should be coherent. The exact value for this maximum time strongly depends upon the prevailing environmental conditions. Critical factors that need to be considered include the motion of the scatterer, e.g. bubbles rising through buoyancy, motion of the sonar platform and motion of the surrounding medium, e.g. through turbulence. In order to understand how robust TWIPS processing is to variations in the environment, a series of simulations was conducted using a bubble cloud in which the individual bubbles were allowed to move; details of these simulations are provided in Finfer (2009). In these simulations, the performance of TWIPS was reduced by 50 per cent when the path-length differences corresponded to 0.37λ, where λ is the acoustic wavelength.

To create representations for display purposes, smoothed envelopes of the basic signals p(t), p+(t) and p(t) are computed. These envelopes are evaluated by band-pass filtering the signals, then computing their envelope (exploiting the Hilbert transform) and finally smoothing the result by averaging over the duration of the outgoing pulse. These smoothed envelopes are denoted here using capital P notation, so that the envelopes of p(t), p+(t) and p(t) are denoted by P, P+ and P−, respectively. The band-pass filters used in the initial stage of this processing can be tuned to enhance particular harmonics. When such processing is performed, then the notation is augmented so that Pm+ denotes the envelope of p+(t) after applying a band-pass filter whose centre frequency is m times the centre frequency of the insonifying pulse. Logically, p(t) should generally be filtered about odd harmonics (m odd) and p+(t) should be filtered about the even harmonics (m even). As an addendum to this notation, the absence of m denotes the use of a band-pass filter at the frequency of the outgoing sonar pulse, i.e. m=1, and on occasion, we shall use notation such as P1+ to emphasize the character of the filter used. It is useful to use, as a minimum requirement, checks such as the fact that the function P1+ should be zero at all times and places when the bubble-free data are processed if the system has sufficient integrity for TWIPS.

The TWIPS performance will be compared with so-called ‘standard sonar’ processing, which is achieved as follows. First, the returned signal is band-pass filtered about the centre frequency of the outgoing pulse. Second, the energy of the return is computed by temporally averaging the envelope using a period that corresponds to the duration of the original pulse. The final step is to average the results from both TWIPS pulses; both pulses are exploited so that the standard technique is not inherently disadvantaged relative to TWIPS processing.

The enhancement of bubbles relative to linear scatterers through P+ is effective, particularly where the bubble-size distribution is so narrow that all bubbles are insonified close to resonance, as occurs with UCAs. However, the calculation of P does not completely remove the influence of nonlinear scatterers (e.g. bubbles). Instead, it results in the enhancement of the linear and odd-powered nonlinearities and the suppression of even-powered nonlinearities. Further enhancement in contrast can be obtained through forming functions based on the ratio P/P+. A particularly important feature of the ratios P/P+ and P+/P is that they do not require range correction through time-varying gain (Leighton et al. 2009a,b). This has been independently noted for use of contrast agents (Tang et al. 2008). The ratio is immune to the requirement of a priori knowledge on the beam and environment, since it self-corrects for range, as both the denominator and numerator fall off with range at the same rate (although at great range, when the components P and P+ are lost in the noise, this self-correcting property is lost).

Conversely, one disadvantage of using ratios is that they can be prone to instability. Small changes in the denominator of the ratio can lead to large fluctuations in the output. The fluctuations in the denominator can be reduced through averaging of values of P+ over sets of adjacent returns (arithmetic, geometric and reciprocal averaging have all been found effective, though for conciseness, only geometric averaging of the denominator is used, as noted in the figure captions). For each image plot shown, a noise-level thresholding is also applied to eliminate spurious variation owing to very low echo returns and noise. The noise-level thresholding is taken to 1 per cent of the maximum value within the image plot.

The ability to calculate two variables, such as P+ and P/P+, from a single echo pair, provides discrimination between linear and nonlinear scatterers (more particularly between those objects that scatter in the even harmonics and those that scatter in the odd ones). Comparison of the outputs of the two functions provides a tool to distinguish clutter from targets if their respective scattering can be forced into different harmonic bands.

In order for TWIPS to be successful, the fidelity of the outgoing transducers and amplifiers is important, since any hardware-induced second harmonic distortion in the outgoing waveform would be interpreted as bubble-mediated by the TWIPS algorithm. These issues are discussed in Leighton et al. (2008b).

3. Method

Both tank experiments and sea trials were conducted to test the efficacy of TWIPS.

(a) Method for tank tests

The A. B. Wood tank at the Institute of Sound and Vibration Research, University of Southampton, UK, contains a body of fresh water measuring 8×8×5 m, with a water temperature of approximately 11°C and a sound speed (in bubble-free conditions) of approximately 1450 m s−1 (figure 2). TWIPS was compared against standard sonar in detecting a target (a steel disc of 230 mm radius and 50 mm thickness) with and without bubbles present. The centre of the disc was on the acoustic axis at a depth of 2.5 m, and placed 2.6 m in front of the source faceplate. This target was used both face-on to the source (when its target strength is −10 dB) and when it was angled (after rotation about a vertical axis) to produce a target strength of −15 dB. The sonar source consisted of four GeoAcoustics T135D transducers in a Maltese-cross configuration. The resulting directivity of this configuration (figure 3a,b), and of that used in the sea trials (figure 3c), has side-lobes, echoes from which cannot be distinguished from on-axis returns using only a single omnidirectional hydrophone. The system emits pairs of pulses, each with a Gaussian envelope and a centre frequency of 6 kHz (figure 4). Only the pulse amplitude and the inter-pulse time, Δ, are within the user’s control. While the direction of an arrival could have been inferred from the use of two or more hydrophones, only a single hydrophone was available (Blacknor Technology D140, serial number 18938 with built-in preamplifier, calibrated by the National Physical Laboratory). It was mounted 0.1 m in front of the source, and 0.3 m above the acoustic axis of the source to reduce the reception of backscatter of the target echoes from the source. Hydrophone data were acquired onto a PC using a four-channel National Instruments sound card acquiring data at 200 kHz on each channel. One channel acquired a trigger signal from a trigger box and the second channel acquired the acoustic data. Acoustic signals were passed through a pair of Krohn-Hite model 3203 filter banks. The high pass was set at 0.2 kHz to eliminate mains contamination, and the low pass was set at 100 kHz to avoid any frequency-folding effects. Data were taken with and without the target in position, with and without bubble clouds in place.

Figure 2.

Schematic of the apparatus for the experimental verification of TWIPS. The shaded plane corresponds to the floor of the laboratory, below which is an underground water tank (8×8×5 m).

Figure 3.

Normalized amplitude far-field directivity patterns of the four transducers at 6 kHz in: (a) the horizontal plane and (b) the vertical plane, where the acoustic axis is horizontal at 0°, as used in the tank (i.e. in the Maltese-cross configuration). (c) One-sided normalized directivity pattern of the sources at 6 kHz in the 2×2 array used in the sea trials, where 0° corresponds to the direction in which the pistons face (i.e. downwards in the sea trial).

Figure 4.

Hydrophone signals as recorded in the A. B. Wood tank in the absence of bubbles and the target. (a) First half of outgoing signal at a distance of 1 m from the source. (b) Second half of outgoing signal at a distance of 1 m from the source. (c) First half of outgoing pulse as measured at target location. (d) Second half of outgoing pulse as measured at target location. Axes labelling and numbering run across all axes ticks aligned with them.

Large bubbles can degrade the TWIPS performance since their slower response times hinder a strong nonlinear response (Leighton et al. 2005c, 2008a; Leighton 2007). In a typical oceanic bubble-size distribution, such bubbles occur in relatively small, but not insignificant, numbers. The objective in the tank was to generate a bubble cloud which resembles that found under oceanic breaking waves (figure 5). Direct injection of gas to the base of a tank does not result in such a bubble-size distribution because all injected large bubbles reach the base of the water column, an occurrence which buoyancy opposes under breaking waves. After trying a range of candidate techniques (e.g. electrolysis, catalytic conversion of hydrogen peroxide; Leighton et al. 2007b), the bubble cloud in the tank was generated by first producing a ‘milk’ of small bubbles in a settling tank using a Venturi system (figure 6), and then pumping this ‘milk’ to the base of the A. B. Wood tank. The bubble-size distribution so produced (figure 5) was confirmed as resembling historical at-sea bubble-size distributions by applying the technique of Leighton et al. (2004) to invert the attenuation of a train of 14 pulses, repeated every second (figure 7).

Figure 5.

Bubble-size distributions (the number of bubbles per cubic metre per micrometre increment in bubble radius range). The line connecting the open-circle data points is the average of six measurements taken approximately 1 second apart in the water tank for the bubble cloud used in this paper. Historic measurements are included for comparison (noting that the environmental conditions (wind speed, measurement depth, etc.) reported in the original sources vary). Open symbols are used to indicate open ocean measurements, specifically of Phelps & Leighton (1998) (+); Breitz & Medwin (1989) (×); Farmer & Vagle (1989) (*); Johnson & Cooke (1979) (•) (noting that the photographic technique of the latter might have undercounted the smaller bubbles; Walsh & Mulhearn 1987; Breitz & Medwin 1989; Farmer et al. 1993; Leighton 1994). Closed symbols are used to indicate the four surf zone datasets of Leighton et al. (2004) (Embedded Image; Deane & Stokes (1999) (Embedded Image); Phelps et al. (1997) (Embedded Image; Meers et al. (2001) (◃).

Figure 6.

The apparatus used to generate the bubble population.

Figure 7.

For the bubble population shown in figure 5, the plot shows the additional attenuation owing to bubbles as measured between two hydrophones spaced 0.62 m apart in the bubble cloud. Each pulse had a duration of 1 ms, an inter-pulse interval of 20 ms, and between them spanned from 3 to 197 kHz (Coles & Leighton 2007). The figure shows six separate readings, spaced approximately 1 second apart, giving an impression of the variability of the attenuation as the bubbles rise.

Although, in principle, TWIPS would work with any pulse type (chirp, pseudo-random sequence, etc.), apparatus limitations meant that only a single pulse type (shown in figure 4) could be tested in the experiments described in this paper. The general limitations (introduced in §2) with respect to pulse duration and the inter-pulse time, Δ, need to be quantified for the specific deployment used. Details of this exercise are described by Finfer (2009), and are summarized as follows. The 0.37λ allowable path-length difference (described in §2) corresponds to about 10 cm at the centre frequency of 6 kHz used in this experiment. The maximum inter-pulse delay for the tank experiments can be estimated by comparing this distance to the motion of the bubbles through buoyancy in the tank (there is little turbulence in this environment and the sensor system is rigidly mounted). Estimates based on bubble rise times in the tank population suggested a maximum usable inter-pulse delay in the tank of 100 ms, which was experimentally confirmed (Finfer 2009). The minimum inter-pulse delay was determined by considering the reverberation time of the tank experiment. Experimentation revealed that, for the A. B. Wood tank, the reverberation time, measured at the location of the Blacknor hydrophone was 240 ms, and correspondingly the reverberation was found to take 65 ms to decay by 15 dB. In practice, a minimum inter-pulse time of 50 ms was found to produce acceptable results.

(b) Method for sea trials

On 27 February 2008, the TWIPS sources were mounted in a downwards-looking 2 × 2 configuration (giving the beam pattern of figure 3c) and towed at 4 knots (2.06 m s−1), 2–4 m behind the stern of the RV Bill Conway at approximately 1.5 m depth. The route was from the National Oceanography Centre (Southampton) (50°53′33′′ N, 1°23′38′′ W) to Calshot Castle (50°49′11.53′′ N, 1°18′23.17′′ W). These busy waters handle 7 per cent of the UK’s entire seaborne trade. Data were taken as the sonar source, which was continually in the wake of the RV Bill Conway, was towed through the centre of the additional wake generated by several vessels of opportunity, and TWIPS was compared with standard sonar for its ability to discern the seabed (which here varies in depth from 10 to 20 m). Assessment of effectiveness is through visibility/invisibility of the seabed and wake in the sonar images, with uniform plotting conditions: every plot provides a colour scale linear representation of the output and is normalized to the maximum return in the field; noise-level thresholding was set at 1 per cent of the maximum value; and where TWIPS ratios are calculated, geometric averaging of the denominator of the ratios was carried out for each group of 10 lines in the plot. The more formal measures of performance (receiver operating characteristic (ROC) curves), which are used to examine the tank test results, are not feasible for the sea-trial data because one cannot conduct pairs of trials with the target (in this case the seabed) being absent and then present. In these data, a time-varying gain (proportional to r2(t), where r(t) is the penetration depth at time t) was applied to all the echoes before processing to allow fair comparison between the conventional sonar and TWIPS results (noting that any such corrections cancel out in the TWIPS functions P+/P and P/P+). As with the tank tests, an inter-pulse delay, Δ, of 50 ms was used in the sea trials (at a towing speed of 4 knots (2.06 m s−1), the source moves approximately 10 cm in 50 ms, in keeping with the 0.37λ criterion of §2 for the 6 kHz used here). In these sea trials, the effect of motion of the towing vessel and the turbulence induced by the recent passage of a large vessel make it desirable to use as small an inter-pulse delay as can be achieved without introducing the ambiguity of making it so short that there are two pulses in the water at once. This choice of inter-pulse delay allowed detection of targets to a range of approximately 37.5 m, which was greater than the water depth in this area. The results for all vessels were similar. For conciseness, results are shown for the wake of the Southampton to East Cowes Raptor Class Red Funnel car ferry MV Red Osprey (3953 gross register tonnage; 93.22 m in length and 17.5 m beam; with a capacity for 895 passengers plus 220 cars).

4. Results

(a) Results for tank tests

Figure 8 shows two sequences of consecutive hydrophone records from the laboratory water-tank tests, arbitrarily chosen, to demonstrate the effect that the presence of bubbles has on the detectability of the target. In 10 consecutive hydrophone time histories, the target, which was clearly visible when no bubbles were injected (figure 8a), is no longer visible (figure 8b) under the bubbly conditions shown in figure 5 (the label T indicates roughly where one would expect to detect energy from the disc target echo). Although it should be recalled that sonar processing would employ correlation and filtering in interpreting such signals, this is a useful preliminary indication not only of the degree to which the target scatter is reduced, but also the increase in clutter caused by scattering from the bubbles. The vertical axes of the graphs are normalized to a common acoustic pressure.

Figure 8.

The effect of bubble scattering on target detectability in the A. B. Wood tank hydrophone time histories (all normalized to the same acoustic pressure). (a) Ten consecutive echolocation traces taken 1 s apart in the absence of bubbles. The target is clearly visible as a high-amplitude region at 3.5–4.5 ms (labelled T). The signal that precedes this consists of reverberation of the outgoing pulse. For the data in (b), the same experimental conditions are maintained, except that bubbles have now been added (the population typified in figure 5), with the bubble cloud filling most of the space between the source and target. Axes labelling and numbering run across all axes ticks aligned with them.

When such returns are processed and stacked in figure 9 (with amplitude represented by colour, as defined in the scale bar), then an image is formed, showing the repeatability of the test as the environment changes with bubble cloud evolution. Note that in these tank experiments, the basic TWIPS quantities P+ and P have not had time-varying gain applied to them owing to the short ranges involved in this set-up.

Figure 9.

Vertical plots of the same 100 time series (which includes the data in figure 8) in order to compare the variable processing operators: (i) standard sonar; (ii) P1−/P2+; (iii) Embedded Image; (iv) P2+; (v) P2+/P1−; (vi) P1− for tank tests with a bubble cloud always located between 1.5 and 3 ms (the time axis is as for figure 8). In (a) the target is absent, and in (b) the target is located between 3 and 4 ms (target strength TS=−15 dB), all other conditions remaining unchanged, except for the natural variability of the bubble population (figure 7). Each panel is produced by stacking the result of processing the same 100 consecutive echo time histories (inter-pulse time of 50 ms). Each colour scale is normalized to a maximum value in each plot, which for (a) is (i) 0.16, (ii) 4.8×104, (iii) 1.0×1010, (iv) 4.4×10−3, (v) 0.044, (vi) 0.65 and for (b) is (i) 0.22, (ii) 3.9×104, (iii) 1.1×1010, (iv) 3.8×10−3, (v) 0.060, (vi) 0.86. The echo from the back wall of the tank occurs at around 6.5 ms and of course can also be treated as a secondary target for TWIPS to enhance. Smaller echoes from other walls also occur (e.g. at around 4.5, 5.5 ms). For all representations, noise-level thresholding was set at 1% of the maximum value and (for (ii), (iii) and (v)) geometric averaging of the denominator was carried out for each group of 10 runs.

Figure 9 uses the filtered versions of the basic functions (i.e. P2+ and P1− in place of P+ and P). With no target present, when processed through the standard sonar method (figure 9a(i)), the echo shows reflections from the bubble cloud (located in front of the source at 1.5–3.0 ms), multi-path echoes at 4.5 and 5.5 ms (avoidance of these was not possible for tests at approx. 6 kHz in this tank) and a strong reflection from the back wall of the tank at approximately 6.8 ms (the latter can serve as a useful secondary target). The main target (the disc of figure 2) is added in figure 9b. Let us suppose that we do not have this a priori information as to the location of the various items in the tank. TWIPS allows easy identification of the items compared with the standard sonar plot (figure 9b(i)) because the feature at 1.5–3.0 ms, which disappears in P1−/P2+ (figure 9(ii)) and Embedded Image (figure 9(iii)) but is strong in P2+ (figure 9(iv)), must be a bubble cloud, since it scatters nonlinearly. Similarly, the two strong returns in P1−/P2+ (figure 9b(ii)) and Embedded Image (figure 9b(iii)) must be linear scatterers (and we know that the one between 3.5 and 4 ms is the target, and the one between 6.5 and 7 ms is the back wall, since only the latter remains in rows (ii) and (iii) of figure 9a). As expected from earlier comments, P2+/P1− (figure 9(v)) is less reliable than P2+ (figure 9(iv)) in identifying bubbles (although P2+/P1− retains the benefit of not requiring time-varying gain). This is because P2+ is sufficiently potent without having to resort to the inherently problematical benefits of dividing by the noisy function P1− (which, as shown in figure 9(vi), is not powerful enough to identify the linear scatterers without use of the ratio in figure 9(ii)). The results of the simulation (which contained only bubbles and target, and so did not feature reverberation from tank walls) show the same trends (figure 10).

Figure 10.

Simulation of conditions for tank tests of figure 9 for the region close to the source. In (a), the target is absent, and in (b) the target (TS= −15 dB) is located 2.6 m from the source. Each colour scale is normalized to a maximum value in each plot, which for (a) is (i) 9.3×1012, (ii) 9.6×101, (iii) 1.1×10−9, (iv) 1.3×1013, (v) 2.6, (vi) 3.7×1013 and for (b) is (i) 1.3×1013, (ii) 2.0×102, (iii) 2.0×10−9, (iv) 1.3×1013, (v) 2.6, (vi) 5.0×1013.

From such data as these, ROC curves are constructed in the established manner for assessing the performance of target detectors (Fawcett 2003). Figure 11 uses the ROC curves of the data of figure 9 to compare the performance of the P1−/P2+ and Embedded Image TWIPS functions against the performance of conventional sonar. The data are analysed to determine the accuracy of a binary-classification decision (target present/absent) based on these data as the threshold for detection is varied (the signal must exceed the threshold to indicate ‘target present’).

Figure 11.

ROC curves for the data of figure 9 (TS=−15 dB). Standard sonar compared with the TWIPS functions: (a) P1−/P2+ and (b) Embedded Image, where the black circles are the ROC curve of standard sonar, the blue triangles are that of TWIPS and the solid line indicates the 50 : 50 line. Axes labelling and numbering run across all axes ticks aligned with them.

In plotting the probability of true positives as a function of the probability of false positives, as derived from experimental data, a viable sonar detection system should produce an ROC curve that lies above the 50:50 line. The shape of the curve is also very important. For example, if a sonar system correctly detects every target but also generates a high rate of false alarms, the sensor may well be useless for practical purposes. An example of this would be if the presence of a false alarm triggered a response in the user that was costly (in terms of finance, time, personnel, public relations, etc.): for a sonar in shallow-water, deployment of a dolphin or diver in response to the false alarm to investigate the contact, or destruction of every alarm (negative as well as positive) by the use of ordnance, represent a costly response. A ‘perfect classifier’ would occupy the point in the upper left of the graph, giving 100 per cent true positives and no false positives.

The performances of TWIPS functions are observed to be better compared with the performance of a standard sonar when the target strength (TS) of the target is low (TS=−15 dB). Moving up the mantissas of figure 11, as the threshold is reduced, the two TWIPS functions achieve approximately 40 per cent true positives before they give a false positive. However, the conventional sonar only achieves less than 10 per cent of true positives before it starts to give false positives.

In tests where the target strength was higher (TS=−10 dB), the performance of the standard sonar improves, as expected: the conventional sonar will be able to detect a very strong target in bubbly water, and the capabilities of TWIPS show little if any advantage with respect to target detection. However, it is very important to note that ROC curves only assess the reliability of a sensor for target detection—they do not reflect the ability of TWIPS to provide discrimination, that is, classifying targets into linear/nonlinear scatterers, and TWIPS retains this advantage over the conventional sonar, even with strong targets. The sea-trial results will focus on this aspect.

(b) Results for sea trials

ROC curves of the type shown in figure 11 could not be constructed for the sea-trial data since it is not possible to remove the ‘target’ (the seabed in this case) while keeping all other conditions equal. Therefore, in the sea trial, the objective was to compare the classification abilities of TWIPS and standard sonar. When the RV Bill Conway entered at 4 knots (2.06 m s−1) the wake of the MV Red Osprey (figure 12), P1−/P2+ offered enhanced detection of the seabed over the conventional sonar and, when coupled with P2+, the ability to discriminate between seabed and wake since P1−/P2+ emphasizes the seabed (figure 12b) and P2+ identifies the bubbles in the wake (figure 12c). Conventional sonar shows returns from both the wake and the seabed (figure 12a) and without the a priori knowledge that we have of the wake at the top of the water column and the seabed at the base, we would not be able to discriminate between the two. Similar observations can be made when the simulation is run for the wake tests (figure 13) for a seabed at 10 m depth and an ideally stable platform, and wake statistics (for which there is little choice in the open literature) based on the bubble fields of Culver & Trujillo (2007).

Figure 12.

Data taken in the wake of the MV Red Osprey (figure 12). Comparison of three processing types for the same set of raw data (taken with an inter-pulse time of 50 ms and presented using linear colour scale having a maximum value shown in {}): (a) standard sonar Embedded Image; (b) Embedded Image; and (c) Embedded Image. For all three representations, noise-level thresholding was set at 1% of the maximum value and (for (b)) geometric averaging of the denominator was carried out for each group of 10 sweeps. Data are not shown for the commercial depth sounder (Wheel house unit Simrad CR50, with transducer Simrad combi C50/200 dual 50 kHz/200 kHz operating at 200 kHz), which was fitted to the RV Bill Conway and which was inoperable in this wake.

Figure 13.

Simulation of conditions for the ship wake test of figure 12. The results of the three processing types are presented using a linear colour scale having a maximum value shown in {}: (a) standard sonar Embedded Image; (b) Embedded Image; and (c) Embedded Image. For all three representations, noise-level thresholding was set to 1% of the maximum value and, in (b), geometric averaging of the denominator was carried out for each group of 10 sweeps.

It is possible to quantify the clutter reductions seen in the data of figure 12. The clutter reduction ratio is calculated as Embedded Image. Here, ETWIPS is the ‘energy’ associated with the TWIPS function in question (P1−/P2+ or Embedded Image, and this is computed by summing the TWIPS values over a region where there is wake. Its value is compared with ESTAND, the ‘energy’ associated with the standard sonar processing over the same location. To prevent simple gains in the system from affecting the measure, ETWIPS and ESTAND are both normalized with the respective average value of the echo returns from the sea floor (which is the target here). The resulting clutter reductions in figure 12 are −22.7 dB for P1−/P2+ and −46.7 dB for Embedded Image.

5. Discussions

The tank tests demonstrated the ability of TWIPS to distinguish between linear and nonlinear scatterers, by suppressing and enhancing each in turn depending on the TWIPS function used. Classification of scatterers as either linear (wall reflections, target disc) or nonlinear (bubbles) was demonstrated in the tank tests by comparing the echoes when processed by the TWIPS functions P/P+, P/P2+ and P+/P. In addition to this classification enhancement, the TWIPS functions (P/P+) and (P/P2+) were shown by the use of ROC curves to be better than the conventional sonar technique at detecting the target in the TS=−15 dB tests of figure 11 (although as the target strength increases, all the sonar methods tested here proved to be effective). The wake tests demonstrated the ability of TWIPS to provide a superior classification compared with a standard sonar.

The techniques described in §2 are not specific to acoustics, but could be applied to any radiation and any pulse type, within application-specific limitations. The tank tests and sea trials both placed constraints on the inter-pulse time, for example, which must not be so long as to allow the scattering environment to change significantly, nor too short (and the degree of amplifier/transmitter distortion must not be so great) that the two pulses, when reflected from a linear scatterer, would not be the inverse of each other.

An additional constraint is that the sound field must have sufficient amplitude to drive the bubbles nonlinearly. High amplitudes at distant bubbles are easier to achieve if there is bubble-free water between the sonar source and the cloud in question (as occurred in the tank test, but not the wake test) because bubbles closer to the source will absorb the driving field. For the monostatic deployments considered here, the bubbles to be suppressed or enhanced by nonlinearity must be close to the sonar source because the field emitted by the source decays geometrically with increasing range. However, if the bubbles are distant, bistatic arrangements can be used that place the source close to the potentially nonlinear scatterers, while the observer remains distant (Leighton et al. 2008c).

Current projects are testing TWIPS for a range of sonar applications. These include harbour protection, the detection of bubbles in marine sediments (Leighton & Robb 2008) and manufacturing (Yim & Leighton 2010). Biomedical applications, in addition to UCA studies, include monitoring biomedical shunts and participating in the current debate on distinguishing nonlinearly scattering bubbles from large, linearly scattering ones in tumour treatment using ultrasound (ter Haar 1995; Kennedy 2005; Rabkin et al. 2006; Coussios et al. 2007; Farny et al. 2008; Leighton et al. 2008c; McLaughlan et al. in press). A twin inverted pulse radar might be used to distinguish between sediment (which scatters linearly), rusty metal (which predominantly scatters odd harmonics) and semiconductors (which scatter all harmonics), with application to the detection of improvised explosive devices and in-wall surveillance equipment. Such schemes could be exploited by any radiation (Light Detection and Ranging (LIDAR), Magnetic Resonance Imaging (MRI), etc.), which can be made to scatter nonlinearly by a relevant target with viable excitation amplitudes (Leighton et al. 2008c).

From the earliest days of this study, the impetus in finding a sonar solution for shallow bubbly water had come from the dilemma that some species of echolocating odontocetes (toothed whales) not only inhabit shallow coastal bubbly waters, but at times, also make bubble nets (Leighton 2004a,b). The previous material has shown that such environments are not conducive to the use of standard sonar techniques, but given that a man-made sonar has been developed for such waters, it is interesting to speculate whether any odontocetes use TWIPS. There is evidence that all but one of the key ingredients of a TWIPS system appear in separate species, but no evidence has been found of all appearing in a single species (Leighton et al. , 2009a,b). Multiple pulses have been reported from various species that inhabit the shallow waters where TWIPS would be advantageous. While, for some, the criterion that the pulses be of equal amplitude has not been established (e.g. for the measurements on Phocoena by Awbrey et al. (1979), Evans et al. (1988) and Au (1993)), multiple pulses with nearly constant amplitude and separation times have been reported from at least two species of the genus Cephalorhynchus (Kamminga & Wiersma 1981; Dawson 1988; Evans et al. 1988; Goodall et al. 1988). Some observers have found phase inversion in the second pulses, and although at times these might be from surface reflections (Li et al. 2005a,b), various authors outline evidence for the deliberate generation of multiple pulses by some species (Dawson & Thorpe 1990; Aubauer et al. 2000), although phase measurements on such pulses are rare (Leighton et al. 2007a, 2008b). However, the centre frequency is too high, and the output levels so far measured too low, to generate a TWIPS effect in a broad distribution of bubble sizes (Leighton et al. 2005c, 2008a). Some species do generate sufficiently high amplitudes (e.g. up to approx. 126 kPa zero-to-peak amplitude from Tursiops gilli, Tursiops truncatus and Pseudorca crassidens (Au 1993)), but the centre frequency is again too high and there is no evidence of the deliberate generation of multiple pulses in these species. The authors are not aware of evidence to show that odontocetes can process phase in a manner sufficient for TWIPS processing. However, while echolocating bats are unlikely to encounter any scatterers as nonlinear as a resonant bubble, there are suggestions that they might have developed phase sensitivity (J. A. Simmons 2009, personal communication). Crucially, there is currently no evidence that any of the candidate odontocetes are sensitive to energy at the second harmonic of the centre frequency of their echolocation emissions, which would be necessary to exploit TWIPS (Au 1993; Kastelein et al. 2002). Overall, whilst man-made TWIPS is effective, the current evidence suggests that no other species is undertaking TWIPS processing, although the data are too scarce to be categorical, particularly for those odontocetes that are restricted to shallow bubbly waters.

6. Conclusions

TWIPS has been shown to enhance the detection of a metal target in bubbly water and the detection of the seabed through a ship wake, and to provide superior classification compared to standard sonar for weak targets in bubble clouds. Its use could be extended to other pulse types and radiations. Currently, the admittedly scarce evidence suggests that odontocetes do not use TWIPS, begging the question of what purpose would be served by the deliberate generation of multiple pulses by Cephalorhynchus.


This work received no external funding (support was provided via T.G.L.’s consultancy activity, and a contribution by the Institute of Sound and Vibration Research (ISVR) to the studentship that supported D.C.F.). Since these funds could not cover the full cost of this project, the authors wish to thank a number of people for their free donation of time and facilities. These include, for volunteer labour associated with the preparation of the tank and equipment, Dave Coles (who also helped aboard the RV Bill Conway), Geun Tae Yim, Gary Robb, Chris Powles and Paul Doust. Anne–Marie Liszczyk undertook her undergraduate project on TWIPS and provided valuable help during the tank and sea trials. Thanks are also due to the Skipper (Bob Wilkie) and crew (John Davies) of the RV Bill Conway. The authors are also grateful to Prof. Hugh Griffiths and David Daniels for useful discussions on radar. Raw and processed data can be found at

  • Received March 18, 2010.
  • Accepted May 4, 2010.


View Abstract