1. Introduction
Passive underwater sonar array systems were subject to development in the past decades for accomplishing tasks in naval surveillance and geo-acoustical investigation of the ocean bottom. Modern sonar passive array systems are using adaptive matched field signal processing techniques. Passive sonar array signal processing based on ocean physical propagation model is known as “Matched Field Processing” (MFP). It is a specific application of the more general case of space-time adaptive signal processing in the area of sonar array signals processing (Klemm, 2002). It is a cross-correlation technique developed for matching the values of the computed with a propagation model sound pressure with the measured values at the output of a passive sonar vertical array. The maximum in the cross-correlation or ambiguity surface gives an estimation of the position of the underwater source. This technique is first proposed by Clay in 1966 and later on theoretically and experimentally developed by Hinich and Bucker in 1972 and 1975 (Sullivan & Middleton, 1993). In 1980 and 1981 Klemm introduced in the matched field signal processing algorithms the Capon’s minimum variance and Burg’s maximum entropy estimators for better resolution in shallow water environments (Sullivan & Middleton, 1993). The acoustic field in an ocean horizontally stratified waveguide may be presented as a sum of “normal waves” or “normal modes” standing waves
2. Model based matched field processing of vertical sonar array signals
In conventional antenna array systems signal processing the usual methods of beamforming include phase control of the signals in receive or transmit mode in order to form and steer antenna beam pattern in given direction (see Van Trees, 2002). In MFP the steering vector is determined by the predicted with a propagation model sound pressure for a range of source coordinates in a waveguide with known depth sound speed profile (SSP) and bottom acoustical parameters. A flow diagram of MFP is given on fig. 1.
From the time domain array sampled signals the cross-spectral density matrix (CSDM) is estimated. As a steering vector for the matched filed processor the calculated sound pressure level with normal mode propagation model is used. An ambiguity surface is formed through matching or cross-correlating of the steering vector with the received sound pressure field. The maximum of the ambiguity surface gives the most likely position of the sound source in range and depth. There are several general signal processing algorithms applicable for MFP.
2.1. Cross-spectral density matrix estimation
Cross spectral density matrix The name covariance matrix is also used For example signals at the output of two hydrophones from the array
where
The functions
In practise the record length is always finite and the expectation operation
The row in (5) represents the frequency spectrum of the signal from a given hydrophone (channel).
The array cross-spectral matrix for the n-th frequency bin is:
Where (*) is the complex conjugate operator. The cross-spectral matrix for an M element array is square MxM matrix for a given frequency
An estimate of the CSDM
2.2. Steering vector calculation
The frequency range and the parameters of the ocean waveguide determine the sound propagation model for the matched field processing steering vector. In the application a low frequency source positioned in a shallow water ocean waveguide is considered which specifies the normal mode propagation model. For sound propagation numerical calculations in a range-independant environment the normal mode representation of the acoustical pressure actuated by a harmonic source in a horizontally stratified medium is given by the following infinite radical sum (Ferla et al., 1993):
Where r is the range, z the depth, ρ the density, zs the source depth, Ψm the mode amplitude,
where dependence of the pressure – p from distance – r and depth z is obvious.
Each of the modes has unique wave number, attenuation, phase speed and group speed. The sum of the modes for different source positions and distances from the source forms unique vertical distribution of the pressure which is subject to matching with the experimentally measured pressure (see fig.7). The calculated sound pressure is used as a steering vector in the matched field processor.
The steering vector is the sound pressure level for a given possible position of the known source in depth. This forms a sound pressure level data cube presented on fig. 4. For every possible position of the source in depth
2.3. Matched field signal processors
Two matched field processors are subject to investigation in the application. The first one is the classical or Bartlett matched field processor and the second one is the so called minimum variance distortion less processor (MVDR).
2.3.1. Bartlett (conventional) matched field signal processor
In conventional beamforming the sonar array pattern is steered through phase control of the signals at the output of the antenna array (see fig. 1) through steering vector. In matched field processing the steering vector represents the computed pressure in the points of the array elements. The processor output or as in our case a matched field filter output is:
Where
where
This processor is a benchmark for the other more sophisticated processors because it has much power and low probability of no detection.
2.3.2. MVDR matched field signal processor
The MVDR processor is derived with quadratic constraints (Van Trees, 2002). An optimal filter is given with the block diagram on fig. 6.
The signal vector at the input can be written as:
where Or no gain
where In the majority of the literature it is also referred as MVDR (Van Trees, 2002)
Later in the text we are using the name MVDR for weights (14), (15) and specify if the CSDM of the signal or noise is used. This estimator (beamformer) in the literature is known also as Minimum Variance Estimator (MVE) and Maximum Likehood Estimator (Klemm, 2002). The power of MFP-MVDR is given with (Van Trees, 2002):
where
One way to build robust to some degree of mismatch weight vector and improve the numerical stability in matrix inversion is through diagonal loading of the CSDM. This is the so called method of “diagonal loading” Van Trees shows that imposition of quadratic constrain leads to diagonal loading The analysis of the experimental averaged CSDM showed that the ratio between largest singular value to the smallest is more than 1 000 000
Another way to overcome problems with CSDM inversion is to apply singular value decomposition (Van Trees, 2002). With singular value decomposition the square CSDM is presented in the following form:
where
where
The optimal weight for MVDR is formed:
where the inverse matrix Also known in the literature as pseudo-inverse
The rank in the application is determined if the corresponding singular value is under 1/100 from the maximum value. On fig. 8 the singular values of a CSDM are plotted and the rank is determined to be r=5.
3. Vertical sonar array MFP experimental results
Passive array sonar data from SACLANTC 1993 North Elba experiment available on Internet (Saclant sonar data, 1996) is used for processing. Additional information about the experiments and signal processing examples is given in (Gingras & Gerstoft,1995, Gerstoft & Gingrass, 1996). The array data was collected in shallow-water off the Italian west coast by the NATO SACLANT Center in La Spezia, Italy. SACLANT Center has made this data available to the public for the purposes of fostering signal processing research. The original SACLANT time series has been converted to a series of MATLAB files each of which contains a matrix data file from the 48 sensors. Each file represents about 1 minute of data The vertical array consists of 48 hydrophones with spacing 2 m between elements at total aperture length 94 m (18.7 m to 112.7 m in depth). On 27-th of October a source was towed in the surveillance area at a speed approximately 3.5 kn. The source emitted PRN signal with a center frequency of 170 Hz. For this source 10 data files are available. The data for the moving source was subject to MFP. The available on the site sound speed profiles were used for replica computations with KRAKEN normal mode program (Porter, 2007). Technical manual of the software was used from (Porter, 1992). On fig. 9 a datagram of sampled time domain array signals and their spectrum is given.
The CSDM was formed from averaged FFT data. The FFT length was 4096 resulting in 0.2441 Hz FFT binwidth for 1000 Hz sampling data rate. A Hamming window with 50 % overlapping was used for weighting. It was observed from fig. 11 the relatively low level of the power of CSDM in the frequency window of interest with center frequency 170 Hz compared to the low frequency noise (50 – 100 Hz). For estimation of signal CSDM the frequency bin for 169.9 Hz is used.
On fig. 12 the ambiguity surfaces of three matched field processors: Bartlet, MVDR with diagonal loading and MVDR with reduced rank SVD and pseudoiverse matrix are plotted. The position of the source is determined with the maximum in the ambiguity surface. All data was equally treated with rank of the reduced matrix r = 5, and diagonal loading constant ε2 = σ2 = 10. Comparison of the three processors side lobes for an equal depth cut of the ambiguity surfaces is given on fig. 13. The level of the side lobes is lower when introducing spatial constraints with rank reduction of the inverse matrix. In the source localization results on fig. 12 the Bartlett processor is used as a benchmark. Even it has more side lobes in the ambiguity surface there is global maximum corresponding to the source position. During the processing of the available array data it was established that the accuracy of passive localization with the narrow band MFP is very sensitive to environment parameters change. Little inaccuracy of SSP or bottom geo-acoustical parameters leads to substantial errors in localization. The images on fig. 12 are for the moving source data recorded on the 27-th of October. The available environment file on the web site is used in the computations (Saclant Sonar Data, 1996). If in the environment file the SSP is changed with the data from another file - “sspprofiles.txt” for the 27-th of October more than one and halve thousand meters offset in the source range occurs.
4. Conclusion
The basics of matched field signal processing were presented. Three matched field processors were applied on real experimental vertical array sonar data for source localization purposes. The MFP-MVDR processor has a lower side lobe levels especially applied with rank reduction and pseudo-inverse matrix diagonal loading (see fig. 13). For the investigated ranges and shallow water waveguides however the sound structure in range is reverberant and repetitive resulting in multiple ambiguity surface local maximums or false alarms. The current state of the technology leads to development of underwater sensor networks with application of bistatic and multistatic sonar concepts. Results from measurements with a bistatic sonar experimental setup are presented in (Kolev et al. 2009), (Kolev et al. 2010). The interest for developing adaptive beamforming methods is as with uniform but also with non-uniform and random sparse (ad-hoc) arrays with distance between nodes bigger than a half-wavelength (Hodgkiss, 1981, Gerstoft, P & Hodgkiss, 2011). MFP methods are subject to development for application in these underwater networks.
Acknowledgments
The authors are grateful to SACLANT centre scientists who have made experimental sonar array data and normal mode propagation codes available on the web.
References
- 1.
Bendat J. Persol A. 2010 , John Wiley & Sons, 2010.978-0-47024-877-5 - 2.
Booth O. Abawi A. Schey P. Hodgkiss W. 2000 Detectability of Low-Level Broad-Band Signals Using Adaptive Matched-Field Processing with Vertical Array Aperture Arrays. IEEE Journal of Ocean Engineering,25 3 JULY 2000.0364-9059 - 3.
Ferla M. Porter M. Jensen F. 1993 C-SNAP: Coupled SACLANTCEN normal mode propagation loss model. SACLANTCEN Memorandum SM-274, SACLANT Undersea Research Centre, La Spezia, Italy, Dec. 1993. - 4.
Gerstoft P. Hodgkiss W. 2011 Improving beampatterns of 2D random arrays using convex optimization, JASA Express Letters, March 2011. - 5.
Gerstoft P. Hodgkiss W. W. Kuperman W. Song H. 2003 Adaptive Beamforming of a Towed Array During a Turn. IEEE Journal of Ocean Engineering,28 1 January 2003,44 54 0364-9059 - 6.
Gerstoft P. Gingras D. 1996 Parameter estimation using multifrequency range-dependent acoustic data in shallow water. J. Accoust. Soc. Am. 99 (5) May 1996.3589 3596 0001-4966 - 7.
Gingras D. Gerstoft P. 1995 Inversion for geometric and geoacoustic parameters in shallow water: Experimental results. J. Accoust. Soc. Am. 97 (6) June 1995.3589 3596 0001-4966 - 8.
Random Acoustic Arrays.Hodgkiss W. 1981 In: Underwater Acoustics and Signal Processing, Proceedings of the NATO Advanced Study Institute, L. Bjørnø, (Ed.), 327-347, Reidel Publishing Company, Dordrecht, Holland, 1981, ISBN 9789027712554. - 9.
Klemm R. 2002 The Institution of Electrical Engineers, 2002.0-85296-172-3 - 10.
Klemm R. 1993 Interrelations Between Matched-Field Processing and Airborn MTI Radar. IEEE Journal of Ocean Engeneering,18 3 July 1993.,168 180 0364-9059 - 11.
Kolev N. Georgiev G. 2007 Reduced Rank Shallow Water Matched Field Processing for Vertical Sonar Array Source Localization. 15-th International Conference on Digital Signal Processing, Cardiff, UK,87 90 Digital Object Identifier: 10.1109/ICDSP. 2007. 4288525. - 12.
Kolev N. Iliev I. Ivanov P. 2009 Experimental investigation of an underwater sensor network with radio control on the basis of 802.11. Proceedings of the National scientific conference “ACOUSTICS 2009”,81 87 1312-4897 - 13.
Kolev N. Kaloyanchev P. Iliev I. Ivanov P. 2010 A Multistatic system for underwater surveillance, Proceedings 5-th international scientific conference HEMUS 2010, Plovdiv, Bulgaria, 2010,245 254 1312-2916 - 14.
Porter M. 2007 Kraken Normal Mode Program, retrieved from Ocean Acoustic Library http://oalib.hlsresearch.com/Modes/index.html - 15.
Porter M. 1992 The Normal Mode Program, NRLIMRI5120-92-6920, Retrieved from http://stinet.dtic.mil. - 16.
Saclant Sonar Data 1996 retrieved from Signal Processing Information Base (SPIB) http://spib.rice.edu/spib/saclant.html. - 17.
Sullivan J. Middleton D. 1993 Estimation and Detection Issues in Matched Field Processing, IEEE Journal of Ocean Engineering,18 July 1993,156 167 0364-9059 - 18.
Van Trees H. 2002 John Wiley & Sons, 2002.0-47109-390-4
Notes
- standing waves
- The name covariance matrix is also used
- For example signals at the output of two hydrophones from the array
- Or no gain
- In the majority of the literature it is also referred as MVDR (Van Trees, 2002)
- Van Trees shows that imposition of quadratic constrain leads to diagonal loading
- The analysis of the experimental averaged CSDM showed that the ratio between largest singular value to the smallest is more than 1 000 000
- Also known in the literature as pseudo-inverse