## 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”[1] - , each satisfying the boundary conditions and propagating down the layer with its own speed, wave constant and attenuation. The coherent modal sum is used in matched field signal processing for localization of an underwater source in range and depth. It can be used for localization of the source also with a horizontal array because the angular modal spread is proportional to the direction of arrival relative to the axis of the array (Klemm, 2002). The matched field processing for source position localization is presented in most of the investigations as a parameter estimation problem. It is assumed that the detection has been carried out and the frequency of the source is known. Some authors consider joint detection and localization problem [Booth et al., 2000]. In the beginning MFP of narrowband signals was used which was later modified with the more realistic for practical applications broadband case where MFP results for several frequencies are averaged. The aim of this chapter is to introduce the reader to the principles of matched field signal processing and present narrowband MFP experimental results. An available on the web normal mode propagation code and vertical sonar array signals data from a controlled source localization experiment are used. A comparison is made between MFP results of three known matched field processors applied to vertical array sonar experimental data.

## 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[2] -
(CSDM) is formed through estimation of the cross-spectral functions of the signal at the output of the sonar array. The cross-spectral function may be defined with a Furier transform of the cross-correlation function or directly via finite Furier transform of the input time domain signal (Bendat, 2010). If sample records are available from a stationary random processes in the time domain[3] -

where

The functions

In practise the record length is always finite and the expectation operation

(5) |

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:

(6) |

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, z_{s} the source depth, Ψ_{m} the mode amplitude, _{m} is the eigenvalue. In the solution a large argument asimptotic approximation of the Hankel function is used. When the number of modes is limited to M the final expression for the sound pressure is:

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

where

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” [6] - (Van Trees, 2002). It is also referred in the literature as a “white-noise-constraint” (WNC). Adding noise to diagonal elements expands the noise space and eliminates small eigenvalues in the inverse CSDM. Usually the CSDM is not well conditioned matrix[7] - and the diagonal loading is a way toward more stable matrix inversion procedure. The MVDR weight with diagonal loading of the CSDM is (Van Trees, 2002):

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 *Λ* is diagonal matrix of the eigen values:

where

The optimal weight for MVDR is formed:

where the inverse matrix[8] - is formed with rank reduction and diagonal loading is applied (adding noise to the inverse matrix through diagonal loading), controlled by the parameter ε. A similar form of this processor weight was used in (Kolev & Georgiev, 2007).

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.