SOURCE IDENTIFICATION/INVERSE PROBLEMS

SOURCE IDENTIFICATION/INVERSE PROBLEMS

SESSIONS

Inverse Problems in Sound Radiation Earl G. Williams Physical Acoustics, Code 7137, Naval Research Laboratory, Washington D.C. We present an introduction and overview of inverse problems in sound radiation. The typical inverse problem involves determining the pressure and/or velocity and/or surface intensity (the reconstructed fields) on a surface of interest (such as a vibrating body) from measurement of the pressure on an imaginary surface located a small distance away. This inverse problem is ill-posed, requiring special care in the solution for the reconstructed fields. Methods now exist which are extremely successful in the solution of this ill-posed problem for complicated geometries; for both exterior problems, such as the radiation from machines, and interior problems, such as interior noise in aircraft and automobiles. With these methods under the belt one can turn to post processing of the reconstructed data with the aim of source identification. A new and important source identification tool is wavenumber-filtered intensity which provides an unambiguous identification of the radiating areas of the vibrating source by filtering out small wavelength (poorly radiating) waves.

INVERSE NAH PROBLEM

Solution of the integral equations

Nearfield acoustical holography (NAH) [1] has become a very powerful tool throughout the world for the study of noise sources for both exterior problems and interior problems. NAH is usually based on inversion of one or more of the following integral equations:

Discretization of the surface using boundary element methods (BEM)[4, 5] leads to a linear matrix equation

p(M) =

p(M) =









GN (M|Q)vν (Q)dSQ ,

(1)

∂GD (M|Q)p(Q)dSQ , ∂n

(2)

where vν is the normal surface velocity, p(Q) the surface pressure, GN and GD are the surface Neumann and Dirichlet Green functions known analytically for separable surfaces, and derived numerically for non-separable surfaces, M is the closed measurement surface taken conformal to the integration (reconstruction) surface Q; p(M) =





 ∂GF GF (M|Q) − i (M|Q))µ(Q) dSQ , (3) ∂n

where µ(Q) is the unknown source density[2, 3] and GF = eikR /4πR (free-field Green function), R = M − Q. The inversion of these integral equations determines the unknown in the integrand, i.e. p(Q), vν (Q) or µ(Q), which can then be used to predict the field on any other surface P by replacing M with P in the associated integral equation and recomputing the Green functions. Furthermore, the velocity vector field v at P can be determined by taking the vector gradient of the integral equation. For example for Eq. (1): iωρv = ∇P p(P) =





∇P GN (P|Q)vν (Q)dSQ .

(4)

pδ = Gvδ .

(5)

where vδ represents a vector of length N of one of the unknowns p, vν , µ in Eqs. (1-3), pδ the vector of length N of measured pressure and G the matrix resulting from the discretization. The superscript δ indicates a quantity with noise resulting from the measurement. One possible solution vα,δ to Eq. (5) is given in terms of the regularized inverse Rα using Tikhonov regularization theory[6]; (6) vα,δ = Rα pδ , where Rα is called the regularized inverse of G and   1 1 1 α,1 Rα = VF diag ,···, ,···, UH . λ1 λi λN The filter factor F α,1 (a diagonal matrix) is   F α,1 ≡ diag · · ·

(7)



|2

|λi  2 · · · .  α 2 |λi | + α α+|λ |2

(8)

i

and U and V result from the singular value decomposition (SVD) of G, G = UΣVH , where Σ = diag[λ1 , · · · , λN ]. When α → 0 the £lter factor is unity, and when α → ∞ the factor goes to zero. The regularized inverse Rα is derived from the minimization of the Tikhonov equation: ||Gvδ − pδ ||2 + α||Lvδ ||2

(9)

(||·|| represents the L2 norm) with L = (I −F α )VH where  λ2 λ2N  1 , · · · , F α = diag α + λ21 α + λ2N

SESSIONS

(10)

α,δ

||Gv

δ

− p || = δ,

√ δ = Nσ,

Q≈σ

q ∈ Ω,

20

30

40

50

(11) 70

(12)

where the norm is taken only over the set Ω containing the last Q eigenvectors of UH with the smallest eigenvalues, λq . One can view the first term in the L2 norm in the discrepancy equation, Eq. (11), as the predicted and smoothed value pα,δ of the pressure (now depending on α); pα,δ ≡ Gvα,δ = GRα pδ = UF α,1 UH pδ ,

10

60

where σ is the variance of the noise (assumed Gaussian). The variance is determined using δ ||UH q p ||/

NAH (b) vs Accel (r)

Mobility (dB re m/s/N)

is the £lter factor for the solution to Eq. (9) when L = I. Note that (I − F α ) is a high-pass £lter so that (I − F α )(VH vα,δ ) passes the noisy high-wavenumber, smallwavelength Fourier components VH vα,δ while blocking the relatively noise-free low wavenumber components. Thus the penalty function represented by the last term in Eq. (9) does not include the lowest wavenumbers and is weighted towards the high wavenumber components. At this point the value of α is still undetermined. The solution of the following discrepancy equation due to Morozov provides an estimate of the value of α:

(13)

derived from the reconstructed velocity. F α,1 acts to smooth pδ , reducing the high spatial frequency oscillations. Thus the best value of α corresponds to when the £ltered pressure differs from the measured pressure by just the noise, a result we would expect if the £ltered pressure was identical to the actual pressure p (with no noise). This procedure provides a computationally robust estimate of α.

EXPERIMENT A point-driven rectangular plate was scanned with a microphone to measure a 37x25 point planar hologram. The surface velocity of the plate was reconstructed at 2500 different frequencies using Equations (6-8). The £gure shows the reconstructed velocity (actually mobility) at a point versus a surface mounted accelerometer at the same point. The agreement is excellent. The Morozov discrepancy principle was successful at determining a value of α at every frequency. The signal to noise ratio 20 log10 (||p||/σ) was estimated by Eq. (12) and varied from close to zero dB for the lowest frequencies to a maximum 45 dB. ———————————————————-

500

1000

1500

2000

2500

3000

Frequency (Hz)

FIGURE 1. Measured (from the impedance head) versus reconstructed velocity (from NAH) at the drive point on a pointdriven vibrating rectangular plate.

ACKNOWLEDGMENTS Work supported by the Office of Naval Research.

REFERENCES 1. Earl G. Williams. Fourier Acoustics: Sound Radiation and Near£eld Acoustical Holography. Academic Press, London, UK, 1999. 2. Yu. I. Bobrovnitskii and T. M. Tomilina. General properties and fundamental errors of the method of equivalent sources. Acoustical Physics, 41:649–660, 1995. 3. Gary H. Koopmann, Limin Song, and John Fahnline. A method for computing acoustic £elds based on the principle of wave superposition. J. Acoust. Soc. Am., 86:2433–2438, 1989. 4. Mingsian R. Bai. Application of bem (boundary element method)-based acoustic holography to radiation analysis of sound sources with arbitrarily shaped geometries. J. Acoust. Soc. Am., 92:533–549, 1992. 5. S. T. Raveendra, S. Sureshkumar, and E. G. Williams. Noise source identi£cation in an aircraft using near£eld acoustical holography. In Proceedings of the 6th AIAA, number AIAA2000-2097, Lahaina, Hawaii, 2000. 6. Per Christian Hansen. Rank-De£cient and Discrete Ill-Posed Problems. Siam, Philadelphia, PA, 1998.

SESSIONS

Optimal Conditioning of Inverse Problems in Acoustic Radiation P. A. Nelson and Y. Kim Institute of Sound and Vibration Research, University of Southampton, Southampton SO17 1BJ, UK Solution of the discrete inverse problem in acoustics can yield estimates of acoustic source strength from measurements of acoustic pressure in the radiated field. In this paper the conditions are identified for which the inverse problem is optimally conditioned. This involves sampling both source and field in a manner which results in the discrete pressures and source strengths constituting a discrete Fourier transform pair.

INTRODUCTION It has been noted in a previous publication [1] that, for a particular 2D acoustic radiation problem, the matrix of Green functions relating the sampled far field acoustic pressure to the strengths of a discrete array of elementary acoustic sources has unit condition number at a certain frequency. The associated inverse problem is thus very easy to solve at this frequency and the strengths of the elementary sources can be reliably estimated from the inversion of the Green function matrix [2, 3]. It is the purpose of this paper to explore this relationship further and extend the result of this potentially useful observation to the 3D case.

THE 2D RADIATION PROBLEM

¥

M -1 jkˆmDx p ( kˆ ) = G c ( kR ) å u z ( mDx )e (2) m =0 Now if the variable kˆ (which lies in the range –k to +k) is also divided into M equal increments Dkˆ = 2k / M , then we may write M -1

ˆ p(uDkˆ) = G c (kR) å u z (nDx)e j (uDk )( mDx) (3) m =0

The far field acoustic pressure radiated by 2D velocity distribution u z (x) can be written as ˆ p(kˆ) = G c (kR) ò u z ( x)e jkx dx

distribution. Now note that if the source distribution is assumed to be bounded and sampled spatially at N discrete points such that the continuous variable u z (x) is replaced by the discrete variable uz(mDx), then we may write

(1)

-¥

where uz is the velocity normal to the x-axis, the wavenumber k = w / c 0 , where w is the angular frequency and c 0 is the sound speed, the wavenumber variable kˆ = k sin f and G c (kR) = wr 0 e - jkR + jp / 4 /(2pkR)1 / 2 where r 0 is the density [4]. The far field pressure p (kˆ ) is

sensed at field points having coordinates (R,f), where f denotes the angle to the z-axis. Equation (1) clearly demonstrates the Fourier transform relationship between the far field and the source

where u denotes the index associated with the sampled pressure. Note that this is accomplished by sampling the far field pressure at M equal increments of sin f [1]. Finally it can be seen that when the elementary sources are each separated by one half wavelength (Dx = l / 2) then DkˆDx = (2k / M )(l / 2) = 2p / M since k = 2p / l . Under these circumstances the relationship between the sampled pressure and the discrete source distribution becomes p(uDkˆ) = G c (kR)

M -1

å u z ( m Dx ) e

j

2p um M

(4)

m =0

which demonstrates that they are a discrete Fourier transform pair [5]. In matrix terms we may write

SESSIONS

p = G c (kR ) Wu

(5)

where p is the vector of acoustic pressures having elements p(uDkˆ) , u is the source strength vector having elements u (mDx) and W is the Fourier matrix.The latter has a ratio of maximum to minimum singular value of unity and the errors resulting from its inversion are therefore at a minimum [2, 3].

THE 3D RADIATION PROBLEM

Dk y = kDY / R .

consisting of an M ´ N array whose projection on the source plane is a rectangle of dimensions Lx and L y then DX = Lx / M and DY = L y / N . It then follows that Dk x Dx = 2p / M for a source spacing Dx = ( R / L x ) l and also that

Dk y D y = 2p / N provided that Dy = ( R / L y )l . Under these conditions M -1 N -1

p(uDk x , vDk y ) = G c (kR) å

where

¥ ¥

ò ò u z ( x , y )e

j (k x x +k y y )

dxdy

uz

æ 2pum 2pvn ö + jç ÷ N ø (mDx, nDy )e è M

M -1 N -1

and the far field pressure is exactly a two dimensional discrete Fourier transform of the source distribution. It can be shown in this case that the matrix of Green functions relating the composite vector of far field pressures at the sensors to the composite vector of source strengths has the form of a block Fourier matrix having unit condition number. This matrix has unit condition number when the sources are spaced such that Dx = (R/Lx)l and Dy = (R/Ly)l and its inversion is consequently least sensitive to noise under these conditions.

m =0 n =0

REFERENCES

-¥ -¥

(6)

G s ( kR ) = jwr 0 e

- jkR

/ 2pR

and

k x = k cos f sin q and k y = k sin f sin q using the usual coordinate system. The far field pressure radiated by a discrete M ´ N planar array of sources can be written as p(k x , k y ) = G c (kR)

å

m =0 n =0

The far field acoustic pressure associated with radiation from a vibrating plane boundary can be written as [6] p(k x , k y ) = G s (kR)

With the sensor geometry

å

å

u z (mDx, nDy )e

j ( k x mDx + k y nDy )

(7)

Sampling the pressure field at equal increments of k x and k y then shows that M -1 N -1

p(uDk x , vDk y ) = G c (kR) å u z (mDx, nDy )e

2. 3.

å

m =0 n =0 j (uDk x mDx + vDk y nDy )

1.

4. (8)

where u and v are integers. Now note that the wavenumbers k x and k y can be expressed in terms of the coordinates X = R cos f sin q and Y = R sin f sin q such that k x = kX / R and k y = kY / R where X and Y define the coordinates of a field point in the source plane. Similarly one may write Dk x = kDX / R and

5.

6.

P.A. Nelson, Proceedings of 6th International Congress on Sound and Vibration, 1, 7-32, (1999). P.A. Nelson and S.H. Yoon, Journal of Sound and Vibration, 233, 643-668, (2000). P.A. Nelson and S.H. Yoon, Journal of Sound and Vibration, 233, 669-705, (2000). A.P. Dowling and J.E. Ffowcs-Williams, Sound and Sources of Sound, Ellis Horwood, Chichester, 1983. L.R. Rabiner and B. Gold Theory and Application of Digital Signal Processing, Prentice Hall, Englewood Cliffs, New Jersey, 1975. P.M. Morse and K.U. Ingard, Theoretical Acoustics, Princeton University Press, New Jersey, 1968.

SESSIONS

(9)

Noise Source Identification Using the Inverse Frequency Response Function I.-Y. Jeon a, J.-G. Ih a, Jay Kim b and Y.-J. Kang c aDepartment of Mechanical Engineering, KAIST, Taejon, South Korea bDepartment of Mechanical, Industrial and Nuclear Engineering, Cincinnati University, USA cDepartment of Mechanical Engineering, Seoul National University, Seoul, South Korea

One of the indirect identification techniques for noise sources using the acoustic array is the inverse frequency response function (IFRF) method. The main problem in applying this technique to interior problems is the reverberation due to the reflective boundary surfaces of the enclosure. In this paper, IFRF method is applied to the interior problem with the use of absorptive body. A simple box and a half-scaled automotive cabin are taken as an example. The field reactivity is reduced and checked through the reactivity index level. As a result, the noise sources distributed on the boundary are identified using IFRF method and it gives precise information on the magnitude and phase of noise sources. H(ω ) + , is determined by using over-determined (m > n)

INTRODUCTION Indirect source identification techniques using multiple-microphones such as near-field acoustic holography (NAH) and inverse frequency response function (IFRF) method have been developed and they are useful when the extended multiple sources are dealt with. The main problem in applying these techniques to interior problems is the reverberation within the enclosure because the sound field is very coherent everywhere inside enclosure. As the presence of reflected wave causes standing wave modes and illconditioned transfer matrix, it is difficult to identify where the actual sources are. In order to reduce the reverberation effect, a sound absorptive body was suggested to be placed within the interior space [1]. When the amount of absorption is increased, the field reactivity is generally decreased. It is checked through the reactivity index level. IFRF method can be applicable to irregular shaped enclosure. In this paper, experimental works were performed to investigate the applicability of IFRF method and demonstrate the usefulness of inserting absorbing material within a cavity.

BACKGROUND IFRF method is based on the path information between inputs on a surface and outputs measured on the microphones described in Fig. 1. If a system is linear, the FRF relationship between inputs and outputs is described as P(ω ) = H(ω ) Q(ω ) ,

least-squared solution approach and singular value decomposition, unknown inputs of the system can be identified from the measured pressure by the pseudoinverse of H(ω ) as follows : ˆ (ω ) = H(ω ) + Pˆ(ω ) . Q

(2)

Here, ‘+’ denotes the pseudo-inverse, Qˆ(ω ) is the actual source strength to be identified, and Pˆ(ω ) is the measured field pressure response. Equation 2 means that the actual noise source distribution can be identified by using a linear combination of calibrate sources.

REACTIVITY IN THE FIELD In order to evaluate the field reactivity, the reactivity index level, Lk, is used in this study. Lk is defined as Lp- LI , where Lp is the sound pressure level and LI is the sound intensity level [2]. Lk is affected by many parameters such as the reactivity of the field, the phase match of microphones, and the phase errors in the measurement. For phase-matched measurement system, Lk is about 0.14 dB in the free-field, but sometimes more than 10 dB in reverberant sound fields. Extended structural source

Absorptive body

Pm Q

n

reflected wave incident wave

(1)

where Q(ω ) is the input vector at n points, P(ω ) is the field pressure at m points, and H(ω ) is the FRF matrix of the system. Once the pseudo-inverse of FRF, i.e.,

External Internal acoustic source acoustic source

Enclosing boundary

FIGURE 1. Schematic geometry of IFRF method

SESSIONS

m m

1

y x

3

2

12

00

A z

1.5

with absorption without absorption

15

Source contribution

6

5

Reactivity index level (dB)

600mm

4

20

10 5

with absorption without absorption

1.0

0.5

0 250

700mm

500

1000

2000

0.0

4000

1

2

3

4

5

1/3 octave band center frequency (Hz)

Speaker location

(b)

(c)

(a)

6

FIGURE 2. Example of a rectangular box model. (a) Locations of loudspeakers and microphones, (b) measured reactivity index levels, (c) estimated source contribution at 154 Hz.

met when the measurement is performed in reverberant condition, and that it was tried to recover sources at 154Hz, which is the first resonant frequency of the box cavity. Fig. 2(c) shows that the noise source is identified as a combination of sources 1 and 4 when the actual source is located in A position. Contrastingly, when the field reactivity is high, it seems that all six positions contribute to the final result with the highest contribution from the position 5.

1 (0,0,0)

3 Louds pe a ke r

2

20

with absorption without absorption

15 10 5 0 250

500

1000

2000

4000

1/3 octave band center frequency (Hz)

(a)

(b)

1.5

source contribution

Microphone a rra y

Reactivity index level (dB)

The source distributed on a face of a simple box depicted in Fig. 2(a) was used as a test example. Six artificial sources and 15 microphones were used. A half of box cavity was filled with absorptive body consisted of polyurethane foam. One can find that the reactivity field can be reduced very much by the use of absorption material as can be seen in Fig. 2(b). As mentioned earlier, the most difficult condition can be

Phase

EXPERIMENTAL OBSERVATION

1.0

Q1* Q2* Q3

0.5 0.0 180 0 -180 100

200

300

400

500

Frequency (Hz)

(c)

FIGURE 3. Example of a half-scaled automotive cabin model. (a) Locations of loudspeakers and microphones, (b) measured reactivity index levels, (c) estimated source contribution over 100–500 Hz frequency range.

For a half-scaled automotive cabin depicted in Fig. 3(a), 3 speaker locations and 15 microphones were used to determine FRF matrix and 25 mm polyurethane foam was put to every boundaries of the enclosure. Again as shown in Fig. 3(b), the field reactivity is dramatically reduced by introducing a small amount of absorbing material. When the actual out-of-phase sources at position 1 and 2 are excited in the same amplitude, the result is plotted in Fig. 3(c) and almost same source contribution and 180° phase difference between sources can be observed.

CONCLUSION In this paper, noise source identification by using the IFRF method is presented for the application to the sources in reverberant interior spaces. It was shown

that the introduction of absorptive body in the interior cavity and the monitoring the field reactivity by Lk are the key points in the success of the method.

ACKNOWLEDGMENTS This work was partially supported by BK21 Project and NRL.

REFERENCES 1.

2.

S. Dumbacher and D. Brown, “Practical considerations of the IFRF technique as applied to noise path analysis and acoustic imaging” in Proceeding of International Modal Analysis Conference, Orlando, 1997, pp. 16771685. F. Jacobson, Noise Control Eng. J. 35, 37-46 (1990).

SESSIONS

Jørgen Hald %UHO .M U6RXQGDQG9LEUDWLRQ0HDVXUHPHQW$66NRGVERUJYHM'.1 UXP'HQPDUN The basic operation performed in planar Near-field Acoustical Holography (NAH) is to calculate the sound pressure in a plane parallel with the measurement plane. This calculation is performed via a 2D spatial FFT transform, involving a discrete representation of the sound field data in the spatial frequency domain, [1]. The spatial window associated with the spatial FFT and the discrete spatial frequency domain representation will introduce serious spatial window effects, unless the measurement area is large compared to the source. The requirement for a measurement area larger than the sound source is a problem in particular in connection with Time Domain Holography, [2], where scanning is not possible. A Statistically Optimal Near-field Acoustical Holography (SONAH) method, which overcomes these problems, was introduced in reference [3]. The present paper gives an alternative and probably simpler derivation of the SONAH algorithm.

,1752'8&7,21 A plane to plane propagation of a sound field away from the source can be described mathematically as a 2D spatial convolution with a propagation kernel. A 2D spatial Fourier transform reduces this convolution to a simple multiplication by a transfer function. In Near-field Acoustical Holography (NAH) the Fourier transform is implemented as a spatial FFT of the pressure data measured over a finite area. The use of spatial FFT and multiplication with a transfer function in the spatial frequency domain is computationally very efficient, but it introduces some errors. Its implicit introduction of a spatial window causes spectral leakage in the spatial frequency domain, which will show up as “window effects” in the calculation plane. The discrete representation in the spatial frequency domain introduces periodic replica in the spatial domain, causing “wrap-around errors” in the calculation plane, [1]. As a consequence, the measurement area must cover the entire source to avoid very disturbing window effects. This is a problem in particular in connection with Time Domain Holography, [2], where scanning is not possible. The new Statistically Optimal Near-field Acoustical Holography (SONAH) method performs the plane to plane transformation directly in the spatial domain rather than going via the spatial frequency domain, [3]. This opens up a possibility to perform acoustical holography measurements with an array that is smaller than the source, and with acceptable errors.

1(:621$+'(5,9$7,21 The derivation of the SONAH algorithm given below is based on the same sound field model as the model

\

Measurement plane ]

6RXUFH UHJLRQ G ),*85( Geometry

used in reference [4] for interpolation of sparse measurement grids in NAH. We consider a complex, time-harmonic sound pressure field S(U ) = S( [, \, ] ) with frequency I and wave number N = ω F = 2π λ where ω = 2πI is the angular frequency, F is the propagation speed of sound and λ is the wavelength. For the following description we shall assume that the half space ] ≥ −G is source free and homogeneous, i.e. the sources of the sound field are for ] < − G as shown in Figure 1. The array measurements are performed in the plane ]= 0. We introduce now the 2D spatial Fourier transform pair of the sound pressure in any plane ] = FRQVW : 3 (N [ , N \ , ] ) ≡ S ( [, \ , ] ) =

∞ ∞

∫ ∫ S ( [, \ , ] ) H

−∞ −∞ ∞ ∞

1

( 2π )

2

M(N [ [ + N \ \)

∫ ∫ 3(N [ , N \ , ] ) H

G[G\

− M (N [ [+N \ \)

(1) GN [ GN \

−∞ −∞

(2) It can be shown, [4], that

SESSIONS

3 (N [ , N \ , ] ) = 3 ( N [ , N \ ,− G ) H − MN ] ( ] + G )

(3)

where  N 2 − (N 2 + N 2 ) [ \  N] ≡  − M (N [2 + N \2 ) − N 2 

for

N [2 + N \2 ≤ N 2

for

N [2 + N \2 > N 2

away from the source. Use of formula (3) in (2) leads to ∞ ∞

1

∫ ∫ 3 (N [ , N \ ,− G ) S . (U ) GN [ GN \ (5) ( 2π ) 2 − ∞ − ∞ which expresses the pressure field as an infinite sum of plane and evanescent waves of the form S . ( [, \ , ] ) ≡ H

− M ( N [ [ + N \ \ + N ] ( ] + G ))

[$ $] +

(4)

The so-called evanescent wave components of the Plane Wave Spectrum 3 outside the Radiation Circle, N [2 + N \2 = N 2 , are seen to be exponentially attenuated

S (U ) =

We now let the number 0 of elementary waves used to determine the estimation coefficients increase towards infinity:

, . ≡ ( N [ , N \ ) (6)

Notice that these elementary waves have identical amplitude equal to one on the source plane ] = −G . They are uniquely specified by the 2D spatial frequency vector ., because N] is given by eq. (4). We assume that the complex sound pressure S(UQ ) has been measured at 1 positions UQ ≡ ( [ Q , \ Q ,0) in the measurement plane. We wish to estimate the pressure S(U ) at an arbitrary position U ≡ ( [, \, ] ) , and we wish to estimate S(U ) as a linear combination of

[$ E] +

1

(8)

Q =1

Solution of this set of linear equations in a least squares sense means that we obtain the estimator (7) which is optimal for sound fields with equally probable elementary wave components. The estimator is optimized for Plane Wave Spectra 3, which are “white” in the source plane. To solve (8) in a least squares sense we arrange the quantities in matrices and vectors: $ ≡ S . P (UQ ) E(U ) ≡ S . P (U ) F(U ) ≡ [F Q (U )] (9)

]

This allows (8) to be written as follows E ≈ $F The regularized least squares solution to (10) is F = ( $ + $ + θ 2 ,) −1 $ + E

S . (UQ′ ) G.

* = ∑ S. (UQ ) S . (U ) P P

P

* S. (UQ )

(13)

S . (U ) G.

Here, * represents complex conjugate and the integration is over the 2D plane wave spectrum domain. Notice that the regularization parameter θ needs a re-scaling when we switch to the integral representation. The integrals in equations (12) and (13) can be reduced analytically by conversion of . to polar coordinates: . = (N [ , N \ ) = ( . cos(ψ ), . sin(ψ )) . We introduce the [\-position vector 5 ≡ ( [, \ ) , and let 5Q be the [\-component of UQ. From (12) and (6) we get

[$ $] +

QQ ′

* = ∫∫ S. (UQ ) S. (UQ′ ) G. *

= ∫∫ H M ( N ] − N ] ) G H M. ( 5 Q − 5 Q′ ) G.

[$ $] +

QQ ′

∞

*

= 2π ∫ H M ( N ] − N ] ) G - 0 ( .5 QQ′ ) .G. 0 N

∞

+ 2π ∫ H − 2 N

(15)

. 2 −N 2 G

- 0 ( .5 QQ′ ) .G.

where 5 QQ′ ≡ 5 Q − 5 Q′ . Equation (13) can be treated in a similar way. We can now estimate the pressure using equation (7). The particle velocity can be obtained in the same way as a linear combination of the measured pressure signals. To derive the required estimation coefficients, we use (8), but with the particle velocity of the elementary waves on the left-hand side.

5()(5(1&(6 1.

2. 3.

(10) 4.

(11)

(14)

and further by polar angle integration and use of (4):

(7)

[

(12)

0

In order to determine the estimation coefficients FQ we require formula (7) to provide good estimation for a set of the elementary waves of equation (6):

]

Q

* S. (UQ )

= 2π ∫ - 0 ( .5 QQ′ ) .G.

Q =1

[

P

P

 → ∫∫ P→∞

1

S . P (U ) ≈ ∑ F Q (U ) ⋅ S . P (UQ ), P = 1...0

* = ∑ S. (UQ ) S . (UQ′ ) P

 → ∫∫ P →∞

the measured sound pressure data S(UQ ) :

S(U ) ≈ ∑ F Q (U ) ⋅ S(UQ )

QQ′

J.D. Maynard, E.G. Williams & Y. Lee, Journal of the Acoustical Society of America, 78(4), 1395–1413, October 1985. J. Hald, 7LPH'RPDLQ$FRXVWLFDO+RORJUDSK\, Proceedings of Inter-Noise 1995. R. Steiner & J. Hald, 1HDUILHOG $FRXVWLFDO +RORJUDSK\ ZLWKRXWWKHHUURUVDQGOLPLWDWLRQVFDXVHGE\WKHXVHRIVSDWLDO ')7, Proceedings of ICSV6, 1999. J. Hald, 2SWLPDO LQWHUSRODWLRQ RI EDG RU QRQH[LVWLQJ PHDVXUHPHQW SRLQWV LQ SODQDU DFRXVWLFDO KRORJUDSK\, Proceedings of Inter-Noise 2000.

+

where $ is the conjugate transpose of $, , is a unit diagonal matrix and θ is the regularization parameter.

SESSIONS

Holographic Identification of Mechanisms for Sonar Backscattering Enhancements: Application to Tilted Elastic Disks P. L. Marston and B. T. Hefner Department of Physics, Washington State University, Pullman, WA 99164, USA High frequency acoustic holography is found to be an effective method for the identification of the causes of enhancements in the backscattering of sound by tilted elastic disks in water. Consider a thin circular disk that is tilted so that compressional waves (thin-plate symmetric Lamb waves) are excited by an incident tone burst of sound in water. The time-evolution of holographic images was used to confirm the importance of mode conversion to (and from) in-plane shear waves at the edge of the disk [1]. For a different range of tilts, holography shows that (higher-order) antisymmetric edge waves propagate around the perimeter of the elastic plate. Both of these mechanisms produce enhancements that are especially important in the late-time behavior of scattering signatures. The frequencies used in the holograms were respectively 250 kHz and 300 kHz and data was acquired by scanning a hydrophone in a plane. In a related study, holography was helpful for identifying wavefields having a phase dislocation along an axis as in the case of helicoidal waves [2].

Introduction Although a circular elastic plate has a simple geometry, there is no exact solution for scattering of sound incident at oblique angles relative to the surface normal. In fact, the scattering mechanisms can be quite complicated in that the circular edge can produce both focusing and mode conversion of the excited leaky Lamb waves. The geometry of the plate, however, makes it relatively easy to study using planar acoustic holography [3]. In the present study, the motion of the plate surface was imaged during the scattering process allowing the subsequent motion of the excited Lamb waves to be followed and compared to the backscattered response of the plate. This imaging was performed for incidence angles corresponding to the excitation of the compressional wave (thin-plate symmetric Lamb wave) and the second antisymmetric Lamb wave, the a1 wave. Acoustic holography was also applied to assist in the development of a transducer capable of generating a beam having a phase dislocation along its axis [2]. The performance of the transducer was assessed by both imaging the active element and imaging the field at different distances from the transducer to examine the evolution of the beam.

backscattering direction. The second mechanism involves those portions of the excited wave which strike the plate edge at oblique angles. In this case, a majority of the energy in the compressional wave is converted to a shear (SH) wave which travels across the plate to convert into a compressional wave at the next reflection [1,4]. The compressional wave generated after the second reflection can then reradiate into the backscattering direction if the geometry of the reflections corresponds to that shown in Figure 1b. For this scattering mechanism, the shear wave does not radiate energy into the surrounding medium and can generate scattered returns that arrive significantly later than the enhancements due to the first mechanism. In order to examine these mechanisms, acoustic holography was performed to image the backside of the plate surface relative to the incident field as in Figure 2. A transducer was aligned such that an incident 3-cycle 250kHz tone-burst struck a γ

(a)

n

a

2h

kp

(b)

Imaging the plate surface during scattering For the compressional wave excitation on the plate, there are two scattering mechanisms, the ray diagrams for which are shown in Figure 1. The first mechanism (Figure 1a) involves those portions of the excited compressional wave which travel close to the diameter of the plate. This wave strikes the plate edge at normal incidence and reflects to reradiate into the

γ ks

n

kp kp

ks

FIGURE 1. Coupling geometry for compressional wave excitation on circular plate: (a) along a diameter for the direct reflection and (b) with mode-conversion of the compressional wave to (and from) a shear wave at the edge of the plate.

SESSIONS

borofloat glass plate at γ = 16° which correspond to the angle necessary to excite compressional waves on the plate. The plate had a diameter of 150 mm and a halfthickness of 3.25 mm. A hydrophone was mounted on an x-y positioner which sampled the plane behind the plate at ∆x = ∆y = 2.5 mm increments over a 20×20cm aperture. Because the leaky waves radiate energy from both sides of the plate, by sampling the plane behind the plate, the leaky wave field could be imaged without interfering with the incident acoustic field. At each point in the sample plane, the field radiated from the plate was recorded for subsequent processing. Using algorithms from planar holography [3], the sampled field was backpropagated to the surface of the plate to determine the normal velocity of the plate surface as a function of time. The imaged field could then be presented as an animation to follow the evolution of the imaged wavefield [1]. A sample frame from this animation is shown in Figure 3. This frame shows the reflection of the leaky compressional wave from the plate edge and the focusing of the wavefield due the curvature of the plate edge. This technique was also applied to scattering associated with the excitation of the second antisymmetric wave. In this case, from the imaged field, the scattering mechanism was a wave which traveled around the circumference of the plate and could be either a whispering gallery mode or a type of flexural edge wave.

(a)

transducer

γ =16 o tank

glass plate

positioning system with hydrophone

4 cm Y

X

hydrophone

FIGURE 3. Imaged field on plate surface during the scattering process. The black circle indicates the location of the plate edge. See the animation in [1].

Assessment of Transducer Performance In order to assess the performance of specialized transducers capable of generating a beam with a dislocation along the beam axis, acoustic holography was applied to both image the radiating surface as well as examine the evolution of the beam with distance [2]. By imaging the surface of the transducer it was possible to isolate those portions of the transducer that were radiating with a stronger intensity and to determine the field distribution across the transducer face [4]. In addition to imaging the transducer surface, by back propagating the measured field to intermediate distances and through forward propagation of the field, it was possible to examine the evolution of the beam without moving the measurement plane. This simplified the experiment and made it possible to study the beam at propagation distances that extended beyond the size of the tank.

(b)

ACKNOWLEDGEMENTS image plane

y

sample plane

This work was supported by the Office of Naval Research.

REFERENCES

x 1. 2. 4 cm

FIGURE 2. (a) Top-view of the experimental set-up to image the plate surface during scattering. (b) Source and hologram plane.

3.

4.

Hefner, B. T. and Marston, P. L., Acoustics Research Letters Online 2, 55-60 (2001). Hefner, B. T. and Marston, P. L., J. Acoust. Soc. Am. 106, 3313-6 (1999). Williams, E. G., Fourier Acoustics, Academic Press, New York, 1999. Hefner, B. T., Ph.D. dissertation, Washington State University, Department of Physics 2000, pp. 245-254.

SESSIONS

The Rayleigh Hypothesis And Near-Field Acoustical Holography D.N. Ghosh Roya and E.G. Williamsb a

Sachs and Freeman, Inc., 9315 Largo Drive West, Largo, MD 20774, USA The US Naval Research Laboratory, 4355 Overlook Avenue, Washington, D.C 20375, USA

b

An important consideration in problems of radiation and scattering from material objects is the manner in which the fields are to be mathematically represented in various regions of the space Rn, n = 2, 3. It depends critically on the shape of the source distribution for the radiation and of the object in scattering problems. Particularly important is the relation between the representation in the far-field and that on the boundary. This constitutes the core of the much debated Rayleigh hypothesis the discussion of which is the primary objective of this article. The understanding of the Rayleigh hypothesis is highly important in acoustical holographic work in which the velocity field on a radiating surface is to be reconstructed from the information of the pressure distribution on a measurement surface which is usually different from the surface of the source.

INTRODUCTION Consider a radiating or scattering surface Ω occupying a region in Rn, n = 2,3, and bounded by a surface Γ∈ C2. It is well-known that outside the smallest ball (disc) B circumscribing Ω (with respect to a chosen origin inside it), the radiated or scattered field can be represented entirely in terms of the outgoing wavefunctions. That is, for x ∉ B,

u( x) = ∑ u H ( 1) ( kr)e iθ   88

(1)

in two space dimensions (x = (r,θ)) and ∞

u ( x) = ∑  =0



∑

m =− 

u m h (1) (kr )Ym (θˆ), (2)

in three space dimensions (x = (r, θˆ ), θˆ = (θ, φ)). The rest of the symbols is standard. We would like to know if representations (1) and (2) can also be applied to the interior of B, especially, on the boundary Γ, which is assumed to be arbitrary in shape. An answer in the affirmative leads to what is known as the Rayleigh hypothesis [1-3] since it was assumed to be so by Rayleigh [4] in his work on scattering from diffraction gratings. It should be pointed out that the coefficients 88u and 88um in the above equations are determined essentially by the Helmholtz representations [5] on Γ. Rayleigh's assumption, however, was called into question decades ago [6] and it is common knowledge at present that the assumption does not in general apply to arbitrary surfaces other than cannonical geometries of circles and spheres.

Therefore, caution must be exercised in using Eqs. (1) and (2) on boundaries which are not cannonical in shape. However, it is well-known [7] that the sets {ζneinθ} and {ζmYmn( θˆ )}, where ζ stands for either cylindrical or spherical Bessel and Hankel function, are complete in L2(Γ). Therefore, any wavefield which is square integrable on an arbitrary Γ can be expressed exactly as in Eqs. (1) and (2) with the crucial exception that the coefficients in these expansions will not, in general, be the same as in those equations and hence cannot be obtained from the Helmholtz representations on Γ. The latter applies if and only if the Rayleigh hypothesis holds. In work on acoustical holography, the radiated field is frequently calculated by the equivalent source technique[8]. Fictitious sources are placed on a suitably chosen surface inside the true source configuration and these sources are determined by matching their radiated field to the actual field on the boundary. The most critical consideration is the location of the equivalent sources. It is in here that the Rayleigh hypothesis becomes important. We discuss this point below. FIELD SINGULARITIES As was mentioned earlier, whether Eqs. (1) or (2) can be written on the boundary or perhaps in the interior of Ω depends upon the shape of the boundary surface. It is commonly thought that the inapplicability of Rayleigh's assumption on noncannonical surfaces is owing to the presence of both incoming and outgoing waves. This is, however, not borne out in practice. The real reason lies in the location and shape of the ˜ of the field singularities inside Ω. By convex hull Γ a field singularity we mean a point in space where representations (1) and (2) fail to be uniformly

SESSIONS

convergent and, therefore, analytic. Recall that if a series is uniformly convergent in a region, then it is analytic there. This is equivalent to stating that the field can be analytically continued (from Ωe=Rn\ Ω ) ˜. inside the nonphysical region Ω up to the hull Γ Therefore, Eqs. (1) and (2) can be used in the entire ˜ . Moreover, if the hypothesis is valid at region Rn\ Γ a point, say, y ∈ Rn, then the series represented by Eqs. (1) and (2) are uniformly convergent and hence analytic at any point beyond y. In summary, the equivalent source technique would result in fields given by Eqs. (1) and (2) provided that the source surface encloses all field singularities, i.e., ˜ which depends does not lie inside the convex hull Γ on the shape of Γ. This result also explains the familiar observation that the accuracy of approximation in the equivalent source method increases as the source boundary approaches the actual boundary Γ. Next consider the smallest ball (disc) B of radius R circumscribing Ω. Outside B, that is, for x ∈ Rn \ B , the field is given by Eqs. (1) or (2). However, if |x| ≤ R, then a singularity may be present and Eqs. (1) or (2) may not apply. However, the field can still be expressed by the outgoing waves

U( x) = ∑ U H ( 1) (kr)e iθ   88

(3)

where we have assumed 2D purely for the sake of illustration. The coefficients u8  , however, are to be determined from the consideration that u − U L ( Γ ) 2

be minimum. It is known that the series (3) is uniformly convergent in every closed subset of Ωe. In other words,

L

( 1)

limit u − ∑ U H (kr)e  L→∞ =1 8

iθ

→ 0.

To summarize. For a surface of arbitrary shape without singularities, the Rayleigh hypothesis may or may not be valid depending upon the distribution of the field singularities. In such a case, the field within B is to be calculated in the sense of L2(Γ) and is uniformly convergent in any closed subset of Ωe. Moreover, the Rayleigh hypothesis is valid anywhere outside the convex hull of these singularities. If, however, the surface Γ is singular, then the Rayleigh hypothesis does not apply at all. Finally, the analyticity of Γ and the boundary data are no guarantee of the applicability of Rayleigh's assumption. This work was supported by ONR.

REFERENCES 1. R. Petit and M. Cadilhac, C.R. Acad. Sci., Ser. B, 267, 468-471 (1966). 2. R.F. Millar, Radio Science, 8, 785-796 (1973). 3. P.M. Vanden Berg and J.T. Fokkema, IEEE Trans. Ant Prop., AP-27, 577-583 (1979). 4. Lord Rayleign (J.W. Strutt), Proc. Roy. Soc., Ser. A, 79, 399-416 (1907). 5. D. Colton and R. Kress, Inverse Acoustic and Electromagnetic Scattering Theory, Springer, Berlin, 1992. 6. L.N. Deriugin, Dokl. Akad. Nauk SSSR, 87, 913-916 (1952). 7. I.N. Vekua, Dokl. Akad. Nauk SSSR, 90, 715-718 (1953). 8. M. Ochman, Acustica 81, 511-527 (1995).

SESSIONS

Transducer Source Reconstruction From a Limited Number of Acoustic Field Measurements J. R. Blakey QTH Services, Woking GU21 1NY, England It is well known that the surface displacement/velocity distribution of a plane piston-like transducer can be determined indirectly from measurements of the radiated pressure field distribution over a plane parallel to the transducer surface. The Angular Spectrum of Plane Waves (ASPW) method, for example, typically requires scans of 100×100 points with a corresponding computational effort. This is sufficient incentive to investigate alternatives that require only a relatively few measurements yet can predict both source distributions and radiated fields to an acceptable accuracy. Initial computer simulations by the author [ref.2] demonstrated that the surface velocity distributions of plane, circular, pistonlike transducers could be reconstructed from a limited number of pressure measurements taken along the transducer axis. The present contribution extends this work to measurements made on a transverse line parallel to the transducer aperture plane. The method is not restricted to uniform sampling of the acoustic field. Some representative source reconstructions are included as well as a discussion of the range of applicability of the method.

INTRODUCTION It is well known that the surface velocity distribution of a plane piston-like transducer can be determined indirectly from measurements of the radiated pressure field distribution over a plane parallel to the transducer surface. In the Angular Spectrum of Plane Waves (ASPW) method, measurements are made by uniformly sampling the field over a sufficiently large square area of the plane. The source distribution is then deduced by applying a spatial transform, or projection operator to the measured data. The transform itself is derived from the angular spectrum of plane concept [1,2]. The ASPW method often requires a large number of time-consuming precision measurements and a corresponding computational effort, scans of 100×100 points being typical for transducers operating in the Megahertz region. Thus there is an incentive for the investigation of alternative approaches which can predict both source distributions and radiated fields to an acceptable accuracy but which require only relatively few measurements.

THEORY General The acoustic field produced by a bounded vibrating surface S can be found from the Rayleigh integral: p ( r) =

i ωρ exp( − ikr )U ( r ) dS o 2 π r ∫s

( 1)

Let the surface S be subdivided into N regions, the strength of the vibration over each region is assumed to

be sensibly constant. The acoustic field can be expressed as p( r ) ≅

iωρ N ∑ ∫ exp(−ikr)U (r )dS o 2πr i=1 S

(2)

If measurements are made at M different field points, then (2) yields a set of M simultaneous linear equations which can be represented in matrix form as [ P ] = [ Q ][U ] j ji i

( 3)

The formal solution to the system (3) is given by −1

[U ] = [Q ] [P ] i ji j

(4)

and may be solved provided M≥N.

Axial reconstruction The initial work [3] used computer simulations based on the circular plane piston to establish the feasibility of the principle. The radiating surface was divided into N regions using N-1 concentric circles. Reconstructions were achieved using N values of the (complex) axial pressure field Pj. The N×N coefficient matrix was calculated from (1), specialized to the case of the axial field of a uniform circular plane piston i.e.

p ( z) = ρ cU [ e o o

− ik ( z 2 + a 2 )1 / 2

−e

− ikz

]

(5)

Transverse reconstruction Following the promising results obtained from axial reconstruction, it was natural to extend the

SESSIONS

investigations to reconstructions from field data on a line parallel to the transducer surface plane. The plane circular piston surface, subdivided by N-1 concentric circles, was used again. The off-axis field expression required for the co-efficient matrix is derived from (1). The particular form employed is due to Fung, Cobbold and Bascom [4].

RESULTS

CONCLUSIONS The current investigations have shown that an acceptable reconstruction of a transducer ‘s surface vibration can be achieved from a limited number of measurements. The simulation experiments have shown that acceptable reconstructions can be obtained over a range of distances 0.05 a2/λ< z < a2/λ from the source. 200

Axial Reconstructions The original simulations included the uniform piston, a linear taper and a Bessel function distribution Jo(x); successful reconstructions were obtained for each source distribution and a range parameter values. A typical result is given in Figure 1.

100

0

-100

1

-200

0.5

0 0

0

Relative radius

1.0

FIGURE 1. A simulated 5 point reconstruction of a Bessel source distribution for a radius a=75mm at 0.25MHz.

Transverse (Radial) Reconstructions The current investigation includes both computer simulations and experimentally derived data. The same source functions used for the axial reconstructions were employed in the simulations. These were carried out for a range of source-measurement plane distances in order to test the range of applicability of the reconstruction algorithm. It was found that the reconstruction quality is relatively insensitive to the measurement distance. A series of axial and transverse field measurements were made using a PC-controlled scanning tank. The source was a 500kHz, 30mm dia. plane transducer (Dupont) and the field probe was a Precision Acoustics 0.5mm dia. Needle hydrophone. A number of frequencies in the range 200-500kHz were used in order to produce different source distribution shapes. Figure 3 compares the results of a 25×25 point ASPW source reconstruction with an 8 point transverse reconstruction for a 200kHz excitation. The transducer exhibits a radial mode excitation, rather than a pistonlike vibration. The 8 point reconstruction compares reasonably well with the ASPW back-projection result. These results on a nominal 500kHz transducer represents an extreme test of the method, since a transducer under normal operation will have a simpler source distribution, therefore requiring fewer points for its adequate representation.

.25 .

.50 Relative radius

.75

1.0

FIGURE 2. Comparison of source reconstructions for a 30mm dia. transducer, driven at 200kHz. The curve is a radial section of the back projection. The line graph is the 8 point radial reconstruction.

Other experiments have shown (not illustrated here) that the predicted radiation patterns for the few-point source reconstructions compare closely to the true source radiation patterns. Although it was anticipated that ill-conditioning of the system of equations (4) might limit the usefulness of the method, in practice this was not a problem.

ACKNOWLEDGMENTS The Author wishes to thank the Directors of Precision Acoustics Ltd and especially Andrew Hurrell for advice and technical support.

REFERENCES 1.

P.R. Stepanishen and K.C. Benjamin, 71, 803-810 (1982).

2.

R. Reibold, Acustica 63, 60-64 (1987).

3.

J.R. Blakey, Acoustic Sensing and Probing (Ed. A.Alippi), World Scientific Publishing Co.Pte.Ltd. Singapore, 1992, pp. 267-274.

4.

C.C.W. Fung, R.S.C. Cobbold and P.A.J. Bascom, J. Acoust. Soc. Am. 92, 2239-2247 (1992).

SESSIONS

The Acoustic Inverse Problem : On the Isoparametric Element Formulation C. Langrenne and A. Garcia Laboratory of Acoustic, Conservatoire National des Arts et Métiers, 292 rue Saint Martin, 75003 Paris, France

The goal of this paper is to solve the inverse problem, i.e., find the velocity distribution on the structure of a sound source by means of hologram measurements around it. For nonsmooth bodies, isoparametric formulation can be used to solve the acoustic direct problem. Applying this method in the present case of an inverse problem, this attempt has failed when the surface geometry and the acoustic variables on the surface of the body are represented by second-order triangular shape functions. The method proposed in this paper provides a computational method for implementing the Helmholtz integral formula for the acoustic inverse problem. The surface of the geometry and the acoustic variables on the surface of the body are represented, respectively, by second-order and first-order triangular shape functions.

INTRODUCTION For industrial cases, it is very important to know the velocity distribution on the structure of a sound source and localise high amplitude of displacements. Near-field Acoustic Holography (NAH) is used when the shaped of the source is simple – planar, cylindrical – because the treatment in the wave number domain is very fast [1]. For arbitrary shaped three dimensional bodies, the Helmholtz formulation needs to be solve to calculate the propagation operator between measurement pressure and structure velocity. An isoparameter element formulation is usually used for implementing the acoustic direct problem [2]. The structure geometry is represented by shape functions of first or secondorder. In the present case, only triangular elements are discuss. An example on a cello illustrates errors on the solution when second-order triangular shape functions are used to represente acoustic variables and a method is given to avoid the problem.

THE HELMHOLTZ INTEGRAL FORMULATION The starting point for the solution method is the classical Helmholtz integral which may be written in the form :

¶G W (r - s ) - ir0ckv(s )G(r - s (1) p(r ) = òò p(s ) G 4p ¶n

where p(r) is the pressure satisfying the Helmholtz equation D + k 2 p (r ) = 0 and G(r-s) is the free space Green function. The coefficient W has the value 0 for r in the interior domain, 1 for r in the exterior domain, and W has a different value which is a function of the local geometry of G at s [2].

(

)

DISCRETIZATION By the use of an isoparametric element formulation, the problem is discretised. The position, pressure and velocity variables are interpolate with the 6

same N m triangular shape given functions [3] : 6

s (x ) = å N m6 (x ).s m m =1 6

p (x ) = å N m6 (x ). p m

(2)

m =1 6

v(x ) = å N m6 (x ).v m m =1

where sm, pm and vm are the values at node m on the element. With two equations of (1), on the exterior domain and on the surface of the body, the pressure on the surface can be eliminated and the operator which relies the velocity on the structure to the radiated pressure is calculated :

p=F.vs

(3)

SESSIONS

We propose an other interpolation of the pressure and the velocity variables with first-order triangular shape functions : 3

3

m =1

m =1

p (x ) = å N m3 (x ). p m and v (x ) = å N m3 (x ).v m

(2)

where the (1,2,3) nodes correspond to the (1,3,5) nodes of the previous 6 node triangular element. The implementation of this new procedure is very simple and can be deduce from the other one, noting that :

()

()

()

()

()

()

()

()

()

()

()

()

N 13 x = N 16 x + 1 N 26 x + 1 N 66 x 2 2 N 23 x = N 36 x + 1 N 26 x + 1 N 46 x 2 2 N 33 x = N 56 x + 1 N 46 x + 1 N 66 x 2 2

FIGURE 1. formulation.

Backpropagation using isoparametric

The procedure of inversion is not describe here. This is an ill-posed problem which is very sensible at noisy data. A regularization method must be used to filtering the evanescent contribution [4].

ILLUSTRATION AND CONCLUSION Measurements are made around a cello on a regular cylindric mesh of 299 points. The structure of the cello is discretized with 128 elements and 258 nodes. Figure 1 shows the result of the backpropagation at 224 Hz with the isoparametric method. Comparing to the figure 2, which is the result of the proposed method, errors appear on nodes (1,3,5) of each triangular element. This is a numerical artefact. The integration on each element is compute with a Gauss method. The result of this integration is, of course, not to throw back in question. But with 6 node triangular element, it depend strongly of the values of the (2,4,6) nodes. The contributions of the (1,3,5) nodes are so poor that the inverse process failed on these points.

ACKNOWLEDGMENTS We thank A. Chaigne, ENST (Paris), for disposal of the cello, M. Bonnet, Ecole Polytechnique (France) for helpful advice.

FIGURE 2. Backpropagation isoparametric formulation.

using

non-

REFERENCES 1. J. D. Maynard.E. G. Williams & Y. Lee, J. Acoust. Soc. Am., 78(4), 1396-1413(1985). 2. A.F Seybert, B. Soenarko, F.J. Rizzo and D.J. Shippy J. Acoust. Soc. Am., 77(2), 362-368 (1985). 3.

M. Bonnet, Boundary equation methods in solids and fluids, John Wiley & Sons, 1999.

4.

C. Langrenne and A. Garcia, Proceedings of InterNoise 97, August 25-27, Hungary, 1997.

SESSIONS

Time Domain Acoustic Holography: Influence of Filtering on Localisation of Non Stationary Source O. de La Rochefoucaulda, M. Melonb and A. Garciac Laboratoire d’Acoustique, Conservatoire National des Arts et Métiers, 5 rue du Vertbois, 75003 Paris, France a rochefou @ cnam.fr, bmelon @ cnam.fr, cgarcia @ cnam.fr An algorithm based on Fourier transforms both in time/frequency and space/wave vector domains is used to forward and backward project the time dependent pressure field radiated by a sound source. Numerical results concerning an ideal baffled planar piston excited by a gaussian pulsed velocity are presented. An adapted filter should be chosen to minimize amplified noise effects without removing data of interest. Several filters have been tested (Veronesi, Tikhonov and Wiener) so that the algorithm can be optimized. Then the time-dependent pressure field radiated by two baffled loudspeakers driven by a sinus signal is measured by a two-dimensional microphone array. Loudspeaker signals can be shifted in time and noise can be added. The forward and backward projection of the measured pressure field from a first plane to a second one are compared to the measured pressure field on the second plane. The agreement between the measured and the projected quantities remains good, for 20 dB of signal to noise ratio, for propagation distance of about ten centimetres and even for backward projection distance of a few centimetres.

INTRODUCTION Nearfield Acoustical Holography associated with methods of regularisation is a particularly powerful technique for the identification and localisation of vibratory sources. However, when dealing with nonstationary sources like engines throttling up, Time Domain Holography must be used. This method requires the simultaneous measurement of time pressure over a microphone array covering the sound source. In this paper, the time domain method to forward and backward project time dependent pressure fields is first described. Then, numerical results concerning the pressure field resulting from a gaussian impulse into a baffled planar piston are backward propagated and then compared to the direct calculation. Those simulations are very useful for testing different filtering processes. The last section contains experimental set-up and results obtained with two baffled loudspeakers.

THEORY Forbes et al. [1] have proposed the method of forward projecting time dependent pressure field used in this paper. One begins by dealing with the acoustic pressure p(x, y, z0,t) , which must satisfy the wave equation. The different transformations, needed to obtain the pressure field in the plane z = z1 from plane

z = z 0 , are summarised in equation (1) where spatial and time Fourier transforms are used. The time dependent pressure can be either forward projected for z1> z0 or backward projected for z0 >z1>0 .

p ( x , y, z 0 , t ) ® P ( x , y, z 0 , w ) ® P ( k x , k y , z 0 , w ) FFTt

FFTx , y

Pr opagation ¯ (1) p( x, y, z 1 , t ) ¬-1 P( x, y, z 1 , w) ¬-1 P(k x , k y , z 1 , w ) FFTt

FFTx , y

The pressure P(k x,k y, z1,w ) can be calculated from the pressure P(k x,k y, z0,w ) as follows: P(k x,k y, z1,w )= P(k x,k y, z0,w )exp(- jk z (z1 - z0) ) (2)

( )

2

where k z 2 = w -k x2 -k y 2 . c0

NUMERICAL RESULTS A projection algorithm based on a discrete approximation of equation (1) has been developed in the Matlab® environment. The sound source used in the simulation is an ideal baffled circular piston, excited by a gaussian pulse. The analytic expression of the time pressure radiated by this piston is widely available in literature [2]. We are mainly interested in the backward propagation of the pressure field since forward propagation has not posed any major problems. For the backward propagation, the presence of evanescent waves requires the use of a k-space filter, W(k x,k y,w ) , which is applied to the field P(k x,k y, z0,w ) . It makes it possible to reject the components of waves numbers masked by the background noise. Several filters already used in stationary acoustical holography [3] were implemented and compared: filters of Veronesi, Wiener and Tikhonov. The Veronesi filter [4] is a lowpass filter which permits to reject the wave vectors, the

SESSIONS

values of which are superior to the value of the cut-off wave vector kc. Two parameters (kc and the slope α) have to be adjusted. The Wiener filter [5] is applied to the measured signal. It considers the noise included in the signal. The noise level is supposed to be constant in the space domain and is estimated for each frequency with an energetic criteria. The Tikhonov filter [6] is calculated from the expression of the propagator and a constraint function. This filter is non supervised. Only one parameter has to be optimised. The comparison between the theoretical and the backward propagated gives good results, especially with the Tikhonov and Veronesi filters. Due to problems of noise estimation for transient field, the Wiener filtering process does not give valuable results.

EXPERIMENTAL RESULTS Measurements were performed by using a squared planar microphone array composed of sixteen Sennheiser electret microphones equally spaced every 3 cm. The microphones are calibrated in amplitude and phase. The source is made up of three baffled loudspeakers. Two of them are alimented by a sinusoidal signal, which can be shifted in time. Additional white noise can be radiated by the third loudspeaker (in top of the figure). A robot allows the displacement of the array, parallel to the baffle, covering the sound source. The electrical signals from the sixteen microphones are amplified and digitally stored on a HPE1432 acquisition card and then passed to a computer running HP VEE. The signal processing of the measured time waves is then done in the Matlab® 6 environment. For the results presented here, each of the two loudspeakers is alimented by two sinusoidal arches (frequency of 1000 Hz). The second channel is 2.5 ms shifted from the first one. Fig. 1 shows the RMS pressure, measured on plane z = 4.4 cm, as function of radial coordinates x and y. The average pressure over time on the source plane does not allow the two sources to be distinguished.

additional information included in the evanescent waves. However, as for the numerical simulations, the sensitivity of the results to the filtering parameters is significant.

FIGURE 2. Pressure as function of radial coordinates x and y at the time t = 3.46 ms.

Backward propagation to the plane source allows a better localisation of the two loudspeakers, according to time, than the direct measure in plane z = 4.4 cm.

CONCLUSION The pressure field radiated by a planar piston was simulated and then backward propagated. Measurements on tests sources have also been done, so that algorithms can be optimised, and different methods of filtering have been tested. The agreement between direct measurement or simulation and the result after backward propagation is fairly good. Tikhonov and Veronesi filters seem to be more adapted for transient field. Other techniques such as regularisation methods, already successfully tested for stationary acoustic sources [7], have to be applied to transient field. Other time frequency methods are also being studied in our laboratory in order to characterise sources requiring a long time of acquisition. The encouraging results obtained by both simulations and measurements, make it possible to consider various industrial applications for our algorithms such as the localisation of the noise generated by different parts of a car engine.

REFERENCES

FIGURE 1. RMS Pressure as function of radial coordinates, x and y.

Fig. 2 shows, as function of x and y, at the time t = 3.46 ms, the pressure backward propagated from plane z = 4.4 cm to the source plane z = 0 cm. At this time, only the right loudspeaker is in use. We can see that, at a given moment, the localisation is fairly good after backward propagation. This is also due to the

1. Forbes, M., Letcher S. and Stepanishen P., J. Acoust. Soc. Am. 90, 2782-92 (1991). 2. Stepanishen P., J. Acoust. Soc. Am. 49, 1629-38 (1971). 3. Der Mathéossian, J. Y., "Problème inverse régularisé pour l'holographie acoustique en champ proche applications en milieu perturbé", Ph. D. Thesis, Univ. de Paris VII, 1994. 4. Veronesi, W. A., Maynard, J. D., J. Acoust. Soc. Am. 81, 1307-22 (1987). 5. Teukolsky, S. A et al., Numerical recipes in Pascal, Cambridge University Press, 1989. 6. Tikhonov, A. and Arsenine, V., Méthodes de résolution de problèmes mal posés, Editions Mir, Moscou, 1974. 7. Der Mathéossian, et al., "Application to regularisation methods to inverse problem: Nearfield acoustical Holography in a reflective medium", in proceedings of the Journées Imagerie Acoustique, Lyon, France, March 1-2, 1994, pp. 4450.

SESSIONS

Estimation of the Shock Wave Parameters of Projectiles via Wigner-type Time-Frequency Signal Analysis. G. C. Gaunaurd and T. Pham U.S. Army Research Laboratory, 2800 Powder Mill Road, Adelphi, MD 20783-1197 U.S.A. A fired supersonic bullet generates a mach cone, with vertex angle arcsin(1/M), where M=v/c is the Mach number. The Nshaped shock wave propagates perpendicular to the surface of the mach cone. To detect such shock waves and to distinguish them from the associated muzzle blasts is of importance when one tries to predict the classification and trajectory of the bullet. This has applications to determine the location of prospective snipers. This latter goal is achievable by correctly estimating the parameters of the wave transient as the shock goes over a nearby array of observing acoustic sensors. Several techniques have been used earlier to estimate such parameters in the presence of contaminating noise levels. We introduce here a novel timefrequency (t-f) approach, which uses some of the distributions belonging to the “Cohen class”. It reveals unique features of the signature, which become more conveniently displayed and evident in this 2-D domain. This facilitates the classification of the transient and its differentiation from the muzzle blast and it may lead to the identification of the bullet’s trajectory.

INTRODUCTION & THEORY

For projectiles of length L < 1 cm, the wavelength λ of the sound fields that they may radiate when they are in motion will be such that λ >> L for frequencies up to 34 kHz, and then the projectiles will behave as point-sources. A number of works [1-3] have studied the pressure fields emitted by point-sources in uniform rectilinear motion in free-space. If the source is stationary and at a height “h” above a flat impedance boundary, several other works have analyzed the emitted and reflected sound fields they generate [4,5]. These papers have included descriptions of the early controversial debate between H. Weyl and A. Sommerfeld, which is summarized in [5]. For a point source in uniform rectilinear but subsonic motion at a height “h” above an impedance plane, other works [6,7] have also given it attention. The present paper generalizes that situation to the case of supersonic motion of such bullets as projected by a distant gun, and as observed by an array of ground sensors, not always directly under the path of the projectile, but rather located at some lateral miss-distance to one side of the path. For a point source of strength q (t ) moving supersonically (i.e., M = v / c > 1 ) above a plane surface of specific acoustic impedance ξ = Z / ρc with uniform velocity v, the pressure field it generates can be shown to be the time derivative of a velocity potential ψ of the form:  e −ikb 2 R1 e +ikb 2 R1    + R1  qo −ikb 2 ( ct − Mx )  R1 ψ= e   (1) −ikb 2 R2 +ikb 2 R2 4π   e e + ] + C R [ R2 R2  

where b 2 = 1 (1 − M 2 ) and the two distances are: R1 = {(vt − x) 2 − ( M 2 − 1)[( y − h) 2 + z 2 ]}

1

2

(2) 1 R2 = {(vt − x) 2 − ( M 2 − 1)[( y + h) 2 + z 2 ]} 2 . The source moves in the x-y plane a distance “h” above the x-axis, and in the x-direction. The method of images is used to obtain Eq. (1), and C R is a reflection coefficient from the ground that can be found from the boundary conditions to be: {ξ ( y + h) − b 2 [ R2 − M (vt − x)]} . (3) CR = {ξ ( y − h) − b 2 [ R2 − M (vt − x)]} In this supersonic case there are two retarded (or emission) times contributing to the field, which are: t1, 2 = [ Mx − ct ± R1, 2 ] c( M 2 − 1) , (4) which are the two allowable roots of the equation of motion, viz., R 2 = c 2 (t − t1, 2 ) 2 (5) = [ x − v(t1, 2 − R1, 2 c)] 2 + ( y m h) 2 + z 2 . In all the above, only the radiation terms, of order 1 R1 (or 1 R2 ) have been kept. Particular cases of interest can be recovered from the above formulation. For example, if there is no ground, or alternatively, if there is no time for a reflected wave from the ground to develop, then C R = 0 , and then only the first two terms in Eq. (1) remain. Then, ψ ( R, t1,2 ) = [q ' (t1 ) + q ' (t 2 )] 4πR1 , (6) with properly from Eq. (4). sensors then differentiation 1 R1 .

simplified retarded times t1 and t 2 The radiated pressure observed by the follows from Eq. (6) by a time where we only retain terms of order

SESSIONS

SIGNAL PROCESSING THE SIGNATURES. The pressure signals of the projectiles are then recorded by acoustic ground sensors. The muzzle blast signatures, also recorded, exhibit the behavior of a simple radiating monopole at the gun-barrel location. The shock wave signatures look like typical N-waves but distorted by the noise and the fact that the bullet (purposely) not always passes exactly on top of the sensors. Here, the gun-to-sensor distance is 200 yards, which is too short for ground reflected waves to have time to develop. The time-record for the muzzle blast and for the shock wave, for the three considered miss-distances of 20, 60 and 110 feet are shown in Figs. 1 and 2. These signatures are also analyzed in the joint t-f domain in Figs. 3 and 4 for the same three miss-distances. In these plots we have used a pseudo-Wigner [8] distribution (PWD), which is a Wigner-Ville distribution (WD) with windowing. PWD is defined as ~ ' ' W f (x, t ) = 2 ∫−+∞ ∞ f (τ + τ ) f (τ − τ ) , (7) ' w f (τ ' ) w f * (−τ ' )e −2ixτ dτ '

Figure 2. Same as Fig. 1 for the shock waves.

where the asterisk denotes complex conjugation, f (τ ) is the input signature, and w f (τ ) is the (Hamming) window function. The speed of the projectiles is M = 2.7 and their length is L = 7.62 mm. Other sizes and speeds were also considered. Ideal signatures (i.e., the perfect N-waves that would result from the theory shown above) are not displayed, but their comparison with Figs. 3 and 4, illustrate the noise effect present in these latter figures, which makes the t-f plots of the PWD to spread both along time and frequency bands. The present approach clearly distinguishes between the muzzle blast and the associated shock wave, and further details, which do not fit here, will be given elsewhere.

Figure 3. WD of muzzle blasts shown in figure 1.

Figure 4. WD of shock waves shown in figure 2.

REFERENCES

Figure 1. Muzzle blasts of 7.56 mm from a range of 200 yards at miss distances of 20, 60, and 110 feet.

1) P. M. Morse and K. U. Ingard, Theoretical Acoustics, McGrawHill, NY (1968). 2) H. Honl, “On the sound field of a point-source in uniform translatory motion”, Ann. der Physik, 43 (5), 437-464, (1943). 3) A. Dowling and J. E. Ffowcs Williams, Sound and Sources of Sound, John Wiley, NY (1983). 4) A. Wenzel, J. Acoust. Soc. Amer., 55, 956-963, (1974). 5) S. I. Thomasson, J. Acoust. Soc. Amer., 59, 780-785, (1976) 6) T. Norum and C. Liu, J. Acoust. Soc. Amer., 63, 1069-1073, (1978). 7) S. Oie and R. Takeuchi, Acustica 48, 123-129, (1981). 8) L. Cohen, Time-Frequency Analysis, Prentice-Hall, NJ (1995).

SESSIONS

Wavelet based Shockwave Muzzle Blast Discriminator B. Mays and T. Pham U.S. Army Research Laboratory, 2800 Powder Mill Road, Adelphi, MD 20783-1197 U.S.A. A projectile’s trajectory can be estimated by measuring the arrival times of the acoustic energy at several locations in space. In the case of a supersonic projectile fired from a gun, both the acoustic shock wave and muzzle blast may be observed. For acoustic sensor systems attempting to determine a projectile's trajectory, the challenge is to first, correctly classify the transient signal as either a shock wave or a muzzle blast and then, calculate the direction-ofarrival via appropriate arrival times across a sensor array. An incorrect classification will result in large estimation errors of the projectile's trajectory. This paper proposes a wavelet-based approach to the shock wave/muzzle blast classification problem. A normalized product of several Discrete Wavelet Transform (DWT) scales will form the basis of the classifier [1]. Experimental results are presented for proper classification over various round types and miss distances from the sensors.

INTRODUCTION Acoustic systems, which seek to exploit energy produced by supersonic projectiles, must correctly estimate the origin of the energy prior to processing. The supersonic projectile produces acoustic energy via a shock wave created by projectile motion and muzzle blast created at projectile launch. Proper discrimination between these two arrival energies at an acoustic system must be achieved if an estimation of the projectiles trajectory is desired. This stems from the propagation pattern of the two different energies. The muzzle blast energy will appear as a far field plane wave originating from the gun while the shock wave will propagate in the form of an acoustic cone trailing the projectile with angle θ = arcsin(1 / M ) , where M is the Mach number [2]. For an acoustic system that incorrectly classifies an acoustic arrival of a shock wave as a muzzle blast, the estimate of shot origin will be perpendicular to the shock cone incorrectly estimating the trajectory of the projectile. The discrimination process is generally trivial in the case of ideal shock wave and muzzle blast arrival energies as the shock wave duration time is microsecond range while the muzzle blast duration is on the order of milliseconds. Under general conditions the shock wave spectral characteristics are greatly affected by propagation and can be severely attenuated over snow covered terrain. Also, many practical systems may have insufficient bandwidth to preserve rise time characteristics of the shock wave. Muzzle blast signatures can also be affected by multipath situations producing apparently faster rise-times than anticipated.

This paper presents a computational efficient method for discrimination that shows promise when presented with non-ideal data.

THEORY & APPLICATION It is well known that wavelets can be used for detecting and characterizing singularities. The approach here is to extend the method developed by Sadler and Pham [1] to discriminate between shock waves and muzzle blasts using product functions of several scales of the Mallot and Zong’s DWT (MZDWT) of the original signal [3]. Each scale of the MZ-DWT is defined as WT j , j = 1, 2, ..., J , where J is the number of dyadic scales. The general discriminating function takes on the form J

p

D( p, J ) = ∏ (WT j ) , j =1

(1)

where p is order of the wavelet scale. The scales, for example, can be cubed first to emphasize the region containing the singularity while preserving the sign of the results. For this particular application, J = 3 and p = 3 provides the best discriminant and it is referred to as cubic discriminant. Initial results have shown that when the signal rise time is preserved, the shockwave produces multiple peaks in D ( p, J ) while the muzzle blast only produces one peak. For cases with limited bandwidth each class of signal only produces one peak but with opposite signs (see figures 3 and 4). In this paper, this will be the discriminating factor. The maximum value is then extracted from the results and its sign noted. A positive value indicates shockwave while a negative peak indicates a muzzle blast.

SESSIONS

Experimental Data

two peaks in the cross product with the dominant one having the opposite polarity than expected.

FIGURE 1: Muzzle blasts for 5.56 and 7.72 mm caliber rounds at various miss distances in feet.

FIGURE 3: Cubic discriminants for muzzle blasts.

FIGURE 2: Shock waves for 5.56 and 7.72 mm caliber rounds at various miss distances in feet.

FIGURE 4: Cubic discriminants for shock waves.

The experimental data was collected for two caliber rounds, 5.56mm and 7.62mm. The rounds were fired along various firing lanes perpendicular to the recording microphone and at a range of 200 yards. Figures 1 and 2 show the muzzle blast and shock wave signatures for the various miss distances. Again, this data was then processed to form the normalized product of the first 3 MZ-DWT scales cubed, as defined in Eq. (1) for J = 3 and p = 3 to form the cubic discriminant. Figures 3 and 4 shows the cubic discriminates for the muzzle blast and shock wave data respectively. Note the strong peak in the cubic discriminant. This peak is unique and has positive polarity for the shock wave and negative polarity for the muzzle blast in all of the data sets except one. The 80 ft miss distance 7.62mm data, which contradicts the pattern, appears to have a strong echo that severely alters characteristics of the muzzle blast. This echo is presumed to have caused

Conclusions The sign of the normalized product of the first 3 DWT scales cubed (cubic discriminant) is a robust discriminator for shock waves and muzzle blasts. This low complexity calculation is applicable to timecritical processing systems. The work shown here is extremely preliminary. Further investigation needs to be done to determine which functions of the DWT scales are the best to use as sample rate and expected pulse widths varied. Performance in the presence of multi-path and transients other than shock wave and muzzle blast needs to be quantified.

REFERENCES [1] B. Sadler, T. Pham, L. Sadler, “Optimal and wavelet-based shockwave detection and estimation,” Journal Acoustical Society of America, Vol. 104 (2), Pt. 1, August 1998. [2] G.B Whitman, “Flow patterns of a supersonic projectile,” Common Pure Appl. Math. 5, (1952).

[3] S. Mallat, S.Zhong, ”Characterization of signals from multiscale edges,” IEEE Trans. Pattern Anal. March. Intell. 14, 1992.

SESSIONS

Direct Evaluation of a Vessel Passive Detection K. Ugrinovica and I. Mateljanb a

Faculty of Natural Sciences, Mathematics and Education, N. Tesle 12/III, HR-21000 Split, Croatia Faculty of Electrical Engineering, Mechanical Engineering and Naval Architecture, R. Boskovica bb, HR-21000 Split, Croatia

b

The paper treats the problem of the statistical direct evaluation of the vessel passive detection structure. The vessel passive detection structure is based on the statistical optimum passive detection of the vessel underwater acoustic noise as a signal in the presence of deep-sea ambient acoustic noise as interference. The vessel signal is supposed to be a zero mean Gaussian stochastic process composed as a sum of some statistical independent parts. These parts are due to cavitation phenomenon of the vessel propeller(s) and to vibration phenomena of the vessel operating machinery, mechanisms and propeller(s). The interference is supposed to be a zero mean Gaussian stochastic process as well. The statistical optimum detection of the signal is considered in the frequency domain and is realized by means of the likelihood ratio statistical test and the Neyman-Pearson criterion. To evaluate the vessel passive detection structure performance we have to consider the probability density functions of the signal, of the interference and of their sum. The direct evaluation measure is now the value of detection probability versus the value of the false alarm probability with signal-to-noise ratio as a parameter.

INTRODUCTION Each vessel has its own characteristic underwater acoustic field, which is a result of the sound emission into the sea by a lot of sound sources on the vessel. These sound sources are mutually statistically independent. This field is the source of the vessel underwater acoustic signal (vessel noise) which can be detected by passive sonar. The vessel underwater acoustic signal propagates trough the sea as an acoustic medium. The sea has its own underwater acoustic noise (sea noise), that is of stochastic nature. So, the passive vessel detection is statistical detection. The sea noise, as a stochastic process, in a deep sea has Gaussian distribution. The sea-water is dispersive acoustic medium and has stochastic parameters. Therefore the original vessel underwater acoustic signal, on the locations very distant of the vessel, degenerates in stochastic manner. The vessel underwater acoustic signal on the receiving location is assumed to be a sum of the mutually independent Gaussian stochastic processes. These stochastic processes are due to the cavitation phenomenon of the vessel propeller(s) and due to the vibration phenomena of the vessel operating machinery, mechanisms and propeller(s).

Detection Algorithm For the vessel underwater acoustic signal, on the basis of likelihood ratio statistical test with the Neyman-Pearson criterion, the optimum detection algorithm for the alternative hypothesis H1 in the frequency domain becomes [1]:

rk

K

H1: å k =1

(

M

)(

N 2pf k 1 + Mr k

for rk =

)

å ~yki

2

i =1

³ lp

( ) N ( 2pf ) S 2pf k

(1)

(2)

k

where S(2pfk), is the spectrum of the signal s(t) on the frequency fk ; N(2pfk), the spectrum of the sea noise n(t) on the frequency fk; K, a finite number of discrete samples in the frequency domain; M, a number of yki , the Fourier complex compounderwater sensors; ~ nent as an element of the receiving signal matrix for the frequency fk and ith sensor. For the specific case, when the vessel underwater acoustic signal is extremely weak, the ratio (2) becomes very small. So we can approximate in (1) that 1+Mrk»1, and (1) becomes now in a following form: K

H1 :

å k =1

(

rk

N 2pf k

M

)

2

å ~yki ³ l p .

(3)

i =1

Evaluation Direct evaluation of the performance of the optimum statistical detection receiver is described by the function of detection probability versus false-alarm probability. The random variable is the reduced likelihood ratio. The value of the detection probability for the real situations has to be at least 0.5. The relations of the detection probability and falsealarm probability we can solve as the probability

SESSIONS

density functions of the reduced likelihood ratio for the alternative and null hypotheses respectively. The reduced likelihood ratio is of the following form

rk

K

( )

l ~ y =å k =1

(

)

å ~y i =1

ki

(4)

and we see that is composed as the sum of 2K squared members with a variances equal to each set of two elements of the same frequency. It means that the statistical distribution of the reduced likelihood ratio (4) is approximately the gamma distribution with degree of freedom equal 2K. The probability density function of the (4) for the null hypothesis is

é ù l l K -1 expê ú 2 ë 2s 0 ( K ) û p0 l( ~ y) = 2 K s 02 K ( K ) G ( K )

[ ]

(5)

where

s 02 ( K ) =

M 2K

K

år

k

(6)

k =1

é ù l l exp êú 2 ( ) ë 2s 1 K û ~ p1 l y = 2 K s 12 K ( K )G ( K ) K -1

[ ( )]

(7)

where M 2K

K

å ( rk + Mrk2 )

[ ]

D = ò p1 l( ~ y ) dl( ~ y) lp

D = 1 - F ( lD ) for l D =

s1

where F(.) is Laplace integral in the following form F (t ) =

1 2p

t

òe

1 - z2 2

dz .

(13)

-¥

Finally, we can write complete relation for detection probability versus false-alarm probability with rk as a parameter, or

ù é ú ê ú ê ú ê ú ê -1 F (1 - a ) + 2 K ê D = 1-F -2 Kú. K ú ê 2 ú ê rk å ú ê 1 + M k =1 K ú ê rk ú ê å úû êë k =1

(14)

(9)

-2 K ,

The last relation (14) is the complete form of the detection probability function and we can see that the detection probability is a function of a chosen falsealarm probability a, discrete samples K and number of sensors M. At the same time there is parameter rk that tells us about vessel underwater acoustic signal and about underwater ambient sea-noise.

ACKNOWLEDGEMENT Financial support by the Ministry of Science and Technology of the Republic of Croatia (Project No. 177080: Research in Passive Detection Theory of Extremely Weak Underwater Signal) is gratefully acknowledged.

REFERENCES

or for a great value of K approximately 2l p

(12)

(8)

k =1

and where G(.) is complete gamma function. Now we can define the detection probability as ¥

-2 K ,

Conclusion

and for the alternative hypothesis is

s 12 ( K ) =

s0

2

M

N 2pf k

2l p

a = 1 - F (la ) for la =

(10)

1. K. Ugrinovic, Optimum Passive Sonar Structures. Doctoral Thesis in Croatian, University of Zagreb, Zagreb, 1990.

and we can define as well the false-alarm probability as ¥

[ ]

a = ò p0 l ( ~ y ) dl ( ~ y) lp

(11)

or for a great value of K approximately

SESSIONS

Acoustic Imaging Using Synthetic Aperture and High Resolution Arrays Li Jun

a ,b , c

c

b

, Hu Kexin , Peng Jianming , Gong Xianyi

a,b

a

State Key Laboratory of Ocean Acoustics, Hangzhou, China Hangzhou Applied Acoustics Research Institute,Hangzhou, China c Institute of Acoustic, Chinese Academy of Sciences, Beijing, China b

Abstract: Acoustical Imaging is an important technology for detection, classification and localization of small objects which are buried in or proud of the bottom. Among the sensing techniques used for acoustic imaging, the low frequency synthetic aperture and the high frequency high-resolution array processing are promising methods that are given more attention by many researchers. The paper introduces a HF high resolution array processing and a LF synthetic aperture sonar(SAS), and we propose a method that is based on combined time frequency and space-time(CTFSD) processing. We show some results of both HF and LF imaging.

INTRODUCTION Acoustic image is achieved by high-resolution sonar. HF high-resolution sonar is researched and developed for many years, and is widely used to marine exploration and mine detection [1]. There is another method that uses the synthetic aperture technique to achieve the high-resolution through the movement of a small physical aperture. Comparing to HF high resolution sonar, the range of SAS can be much longer, and the field of view can be much more wide [2].

HF ARRAY PROCESSING The study of detection and imaging of short distance object has put emphases on the integrated near field-dynamic focusing and far field-beamforming. We developed a dynamic focusing method to fulfill near field-dynamic focusing. According to the tolerance for loss of unfocusing, we divided the slant range into several focusing-points. The phase term at every focusing-point is stored in memory, and then is read out for compensation. The dynamic focusing system is developed using hardware circuit. The linear array beamforming in far field can be accomplished by spatial FFT transform. We used the operational amplifiers to fulfill the FFT transform beamforming. In order to minimize the sidelobe of the array,

weighting method was adopted.

Results of Lake-test Experiment The lake-test of the HF system was operated in 1995 at Mu-Gan lake with the water depth 15m. The wet-end unit and the targets are deployed below the surface 1.5m and 5m, respectively. The slant distance between the targets and the array is about 5m. The result for imaging the “日”-form target without any post image processing is shown in Fig.1(a). The result for imaging “O”-form target with post image and multi-looking processing is shown Fig.1(b). The results of the lake-test demonstrated that the HF high resolution array processing system can detect and image the targets moored in the water.

FIGURE 1. Lake-test results of imaging the targets, (a) “日” form target, (b) the circle form target

SESSIONS

LF SAS PROCESSING HAARI has studied SAS coorperated with IOA since 1997. The system used in the SAS lake test consists of three major elements: 1)Wet-end unit (array, TCM-20), 2)Towed body and winch, 3)Shipboard electronics (a real time processing equipment, PC). The equipment has the faculty for real-time processing and display. The result shows that two gas tanks proud of the bottom can be detected and imaged clearly by conventional SAS [2]. But it was difficult to detect and image the partially buried targets. Although motion compensation has been adopted to process the data of lake-test, the partially buried objects can not still be detected and imaged. So a new method is presented, called CTFSD algorithm. The CTFSD algorithm [2] is a synthetic aperture array processing with the high order time-frequency analysis and extraction & waveform reconstruction as preprocessing.

The Results of CTFSD Algorithm

The lake-test of the LF SAS system was operated in recently at Mu-Gan lake for many times. The water depth is 17m, and the sound speed profile is of weekly negative gradient. We used 1/4 sonar array with the physics size 0.2m and the synthetic aperture is 12m, provided the platform speed 1.5m/s and the pulse period 0.1s. The cross-form target made up of four gas tanks was partially buried in the bottom. Using the same synthetic aperture, the conventional SAS and CTFSD algorithm were adopted to process the lake test data. The results show that CTFSD algorithm is able to detect the target partially buried in the bottom more clearly than conventional SAS. The range of image is in wide-area (30m 40m). The target is about 24m at slant range.

CONCLUSIONS The HF high resolution array processing has been characterized by the integrated near field dynamic focusing and far field beam forming. The LF synthetic aperture processing has been investigated based on the emphases of the joint time-frequency analysis/waveform reconstruction and space-time (synthetic aperture) processing. The several of lake tests demonstrated that the analysis and extraction of the combined temporal, frequency and spatial features, associated with HF high resolution and LF penetrability, are most important in order to detect and image the objects proud of and buried in the bottom against reverberation and other interferences.

REFERENCES 1. Robert C. Spindel,

High Frequency Acoustics at the

Applied Physics Laboratory”, Proceeding of the China/US workshop on High Frequency Underwater Acoustics, 1998. 2. Li Jun,etc, A Combined Time-Frequency-Space Detecting Method ,Chinese Journal of Acoustics, to be published. FIGURE 2. The lake-test results of Synthetic aperture processing, (a) conventional SAS processing, (b) CTFSD processing.

SESSIONS

Near-field Source Localization using Bottom-mounted Linear Sensor Array in Multipath Environment S. H. Leea , C. S. Ryub and K. K. Leeb a Department of Computer Engineering, Uiduk University, Gyeongju, Korea b Department of Electronic Engineering, Kyungpook Nat’l University, Daegu, Korea

In this paper, we propose a 3-D near-field source localization algorithm using a bottom-mounted linear sensor array in multipath environment. This algorithm utilizes the signals from differnt paths and uses a simple linear sensor array to estimate 3-D location of the source, while conventional 3-D localization algorithms require 2-D array of sensors. The conic angles of the source signal through each path are different, and the position of the source can be estimated using these conic angles and the time difference of these signals. We derive equations of the time difference and the conic angles, and estimate the position of source by simultaneously solving these equations.

INTRODUCTION

x

Surface

Image Target

The passive localization of an underwater acoustic source has received considerable attention in literature. Source localization algorithms extract information of source position from time difference of source signal received at spatially distributed sensors[1-3]. Due to geometry, most of 3-D source localization algorithms utilize 2-D planar sensor array. In this paper, we propose 3-D near-field source localization algorithm using a linear sensor array in multipath environment. It estimates conic angles and time difference of each path signals. Source location can be estimated from the geometric relation of conic angles and time difference.

LOCALIZATION ALGORITHM Consider the source and sensor array in multipath environment as illustrated in Fig. 1. A uniform linear array consisting of M sensors is receiving signals from a source at bearing angle β, range R, and depth z through direct path and surface reflected path. Time difference between signals through these paths is given by (1) from geometric relation in Fig. 1. 1 (RS − RD ) (1) c where, c is sound propagation speed. Path signals are impinging on the sensor array with conic angles θD and θS , given by (2) and (3). τ=

cos θD

= cos αD · cos β

(2)

cos θS

= cos αS · cos β

(3)

The proposed algorithm consists of 3 steps as follows:

P+z

R

Target p = (β, R, z ) RS

0

RD P-z

surface reflected path

y

direct path

αD M ...

...

θS

P : depth

θD

β

Sensor Array

1

αS

Bottom

z

FIGURE 1. Sensor and source geometry.

Step 1: Estimate conic angles (θˆ D , θˆS ) and time difference (ˆτ). The proposed algorithm estimates the conic angles and extracts signals from their directions using source localization algorithm proposed by Q. Zhang, et al.[4] Time difference τ between the signals can be estimated using cross-correlation of the extracted signals. Step 2: Estimate source location. From the geometry, τ can be rewritten as 1 τ= c

q

!

q R2 + (P + z)2 −

R2 + (P − z)2

.

(4)

From (2)-(3), θD and θS can be represented by the functions of source position (β, R, z). cos θD

=

cos θS

=

R · cos β p 2 R + (P − z)2 R · cos β p R2 + (P + z)2

SESSIONS

(5) (6)

150

z [m]

Then the source position can be estimated without search by simultaneously solving (4)-(6), and given by   c2 τ2 cos θˆ D + cos θˆS   z = (7) 4P cos θˆ D − cos θˆS r 1 16P2 z2 R = c2 τˆ 2 + 2 2 − 4(P2 + z2 ) (8) 2 c τˆ p R2 + (P − z)2 cos β = cos θˆ D . (9) R

200

48 49 50

Step 3: Local search. The estimation accuracy is improved using the maximum likelihood(ML) algorithm[5].

250 800

51

850

900

950

1000

1050

1100

1150

β [deg]

52

1200

R [m]

FIGURE 2. Localization result.

SIMULATION

CONCLUSION In this paper, we proposed 3-D localization algorithm using a linear array of sensors. The algorithm utilizes the facts that the conic angles of the signals through different paths are different. We derived the relation of the source position and the conic angles and time difference between the signals from different paths and proposed the equation that can estimate source position. Computer simulation shows the proposed algorithm gives a reliable estimate of the source position.

RMSE of R estimate [m]

The performance of the algorithm is assessed by the following computer simulation. We assume the sensor array of M = 32 half-wavelength spaced sensors located at the bottom P = 400 m. The carrier frequency of the source is 100 Hz, and sampling frequency is assumed 256 Hz. Thus the length of sensor array is 232.5 m. The algorithm uses 150 samples of the signal and sensor noise is assumed to have independent white Gaussian distribution. The source is located at β = 50◦ , R = 1000 m, z = 200 m and signal to noise ratio is taken to be 0dB. From 100 Monte-Carlo runs, the localization result is shown in Fig. 2. The perpendicular line to each plane indicates the true position of source and ellipses on each plane represent 99% confidence region to two corresponding parameters. Bearing estimate errors are very small and range and depth estimates are less then 40m and 10m, respectively. Fig. 3. shows root mean square error(RMSE) in range estimate with respect to β. As shown in the figure, RMSEs are small except the region near β = 90◦ at which the conic angles of two path signals are coincident. RMSEs in other location coordinates are also small and the estimation performance of the proposed algorithm is reliable.

250 200 150 100 50 0 20

30

40

50

60

β [deg]

70

80

90

100

110

FIGURE 3. RMSE in range with respect to bearing angle.

ACKNOWLEDGMENTS This work was supported in part by the Underwater Acoustic Research Center, SNU, Korea and in part by the BK21 project.

REFERENCES 1. M. Hamilton and P. M. Schultheiss, “Passive Ranging in Multipath Dominamt Environments, Part I: Known Multipath Parameters,” IEEE Trans. Sig. Processing, 40, 1, 1-12, Jan. 1992. 2. Y. X. Yuan, G. C. Carter and J. E. Salt, “Near-Optimal Range and Depth Estimation Using a Vertical Array in a Correlated Multipath Environment,” IEEE Trans. Sig. Processing, 48, 2, 317-330, Feb. 2000. 3. C. M. Lee, K. S. Yoon, J. H. Lee, and K. K. Lee, “Efficient algorithm for localising 3-D narrowband multiple sources,” IEE proc. Radar, Sonar, and Nav., 148, 1, 23-26, Feb. 2001. 4. Q. Zhang and J. Huang, “Joint Estimation of DOA and Time-Delay in Underwater Localization,” ICASSP 1999, 2817-2820, 1999. 5. Y. Bresler and A. Macovski, “Exact Maximum Likelihood Parameter estimation of Superimposed Exponential Signals in Noise,” IEEE Trans. Acoust. Speech Sig. Processing, ASSP-34, 1081-1089, Oct. 1986.

SESSIONS

Matched Field Processing for Ocean Acoustic Inversions A. Tolstoy ATolstoy Sciences,8610 Battailles Ct., Annandale VA 22003, USA, [email protected] Matched Field Processing (MFP) has become an important tool in the search for underwater acoustic inverse solutions. MFP typically combines acoustic array data with predictions for those data via an ocean propagation model. It is now critical that those models be validated in phase and amplitude at a range of frequencies in the complicated shallow water regions now of interest. This presentation will address the nature of MFP, applications to deep and shallow water inversions (with particular emphasis on shallow water geoacoustic tomography), and recent efforts to benchmark propagation models used in MFP inverse methods.

INTRODUCTION

HOW GOOD IS MFP?

The speciality of ocean acoustics has long been interested in inverse problems. These problems range from efforts in the localization of a source of sound (such as from a marine mammal or submarine) to the estimation of environmental properties (such as ocean sound−speeds over large volumes also known as ‘‘ocean tomography’’). Most recently efforts have been directed toward the calculation of geoacoustic parameters such as bottom depth, sediment thicknesses, sound−speeds, densities, and attenuations. These environmental parameters are difficult, if not impossible, to measure over a region of any size or variability. Thus, methods which can non−invasively ‘‘measure’’ such values over substantial areas are highly desirable. Most approaches under investigation today employ some form of ‘‘Matched Field Processing (MFP)’’.

MFP has been studied extensively over the past two decades or so. A significant portion of the effort has been to understand the advantages and disadvantages of a variety of MFP techniques and to understand common pitfalls. Some key details include: 1) Sampling. Is the search space fine enough to find the proper peak? 2) Mismatch. Is the propagation model appropriate? Are the ‘‘true’’ input parameters allowed for? 3) Noise. Is the signal too weak to be extracted? These issues were considered in a worskshop aimed at evaluating the wide variety of MFP processors avaliable back in 1993 [3] and continue to be of interest.

WHAT IS MFP? In its simplest form MFP is a signal processing technique based upon the comparison (or correlation) of ‘‘data’’ versus ‘‘model predictions’’. Of course, the ‘‘data’’ are usually measured acoustic field parameters (such as Transmission Loss), but they can also be simulated data as well. Additionally, the ‘‘model predictions’’ can also be measured data, although they usually are not. MFP is essentially an extension of the usual plane wave beamforming (examining energy assuming a variety of incident plane waves coming from different directions or ‘‘beams’’) but which allows for greater sophistication in the modeling of the acoustic field. The predicted fields are usually complex, i.e., consist of amplitude and phase, are computed by highly realistic propagation models [1], and are usually made at more than a single receiver. In particular, predictions are made on an array of receivers. Thus, MFP typically involves the cross−correlation between measured and predicted complex acoustic fields over an array [2].

APPLICATIONS As mentioned in the abstract MFP has now been applied to a number of inverse problems. Each type of problem has issues with regard to: the non−uniqueness of solutions (also referred to as ‘‘sidelobes’’), the robustness of a solution to errors in the modeling or in search intervals, and interference from other sources (or noise). Additionally, the size of a search space can vary dramatically for the different applications leading to many different types of optimizations. New processors are constantly being developed to focus on particular aspects of an application.

Source Localization It is generally accepted that MFP can successfuly find a loud source in deep water using low frequencies and long arrays (usually vertical but also towed horizontal arrays). In general the search space is relatively well defined with only the range and depth of the source to be determined, and these are typically ‘‘findable’’ even in the presence of inaccurate environmental information [4].

SESSIONS

However, sidelobes can still be troublesome. It seems that the low frequencies offer robustness with regard to errors but not much resolution. On the other hand, the higher frequencies offer resolution but at the price of more sensitivity to uncertainties in environmental parameters and to propagation modeling accuracy. Efforts to deal with sidelobes have led to improved array designs, to searches for optimal frequencies, to environmental studies to find more accurate information, and to tracking for a moving source.

Deep Water Tomography Global warming has motivated much recent research into the study of large ocean volumes of water and their temperature character, particularly as a function of time. Since water temperature and acoustic sound−speed are directly related, it is natural to study the behavior of ocean sound−speed profiles via acoustic tomography. It has been suggested that by using MFP plus multiple broadband low frequency sources (such as explosives) plus a number of long vertical arrays carefully distributed throughout an ocean region it may be possible to estimate temperatures to well within 0.5 degrees (C) very rapidly [5]. This still needs to be demonstrated, and major difficulties include the expense of developing tand deploying long vertical arrays. There is additional concern for undetermined environmental impacts from explosive sources dropped into the ocean.

Shallow Water Geoacoustics When early efforts were directed at source localization, one of the major difficulties seemed to be the estimation of environmental parameters which are input to propagation models. This seemed to be a major drawback for shallow water efforts where bottom properties were nearly always impossible to determine. However, that difficulty is now regarded as an advantage if one uses MFP to estimate those parameters outside the context of source localization. That is, if a processor is sensitive to a geoacoustic parameter, such as sediment sound−speed, then that processor can be used to estimate that parameter. Recent work has been aimed at developing methods using such optimization methods as simulated annealing, genetic algorithms, simplex methods, gradient methods in hybrid schemes, neural networks and more. These mathematical techniques have come to be necessary as the search spaces which now include numerous geoacoustic parameters become enormous and ill−behaved, i.e., with many local minima [6].

Propagation Model Benchmarking One of the critical elements of MFP is the propagation model used to generate model predictions for comparison with the measured data. These models need to be as accurate as possible Additionally, since the model must typically be used to evaluate many thousands of candidate parameters, these models need to be as efficient as possbile. Recently, it has also become obvious that model accuracy can dramatically affect the quality of an inversion. Efforts to benchmark propagation modeling include a Special Session of the Acoustical Society of America in 1989 [7] and more recently a session in 2001 [8] plus a very recent Inversion Workshop organized by S. Chin−Bing.

Array Element Localization One final area of recent interest is the localization of individual elements of an array. Often arrays are not purely vertical or horizontal and there is some uncertainly with regard to the exact location of any particular receiver. Errors in receiver locations can easily translate into inversion errors. On the other hand, if the processor is sensitive to that uncertainty, then it can address it as an inverse problem [9].

REFERENCES 1. F.B. Jensen, W.A. Kuperman, M.B. Porter, and H. Schmidt, Computational Ocean Acoustics, Am. Inst. Physics, New York, 1994. 2. A. Tolstoy, Matched Field Processing for Underwater Acoustics, World Scientific Pub, Singapore, 1993. 3. M.B. Porter and A. Tolstoy, J. Computat. Acoust. 2, 161−185 (1994). 4. M.D. Collins and W.A. Kuperman, J. Acoust. Soc. Am. 90, 1410−1422 (1991). 5. A. Tolstoy, J. Computat. Acoust. 2 , 1−10 (1994). 6. A. Tolstoy, N.R. Chapman, and G. Brooke, J. Computat. Acoust. 6 , 1−28 (1998). 7. F.B. Jensen and C.M. Ferla, J. Acoust. Soc. Am. 87, 1499−1510 (1990). 8. K.B. Smith and A. Tolstoy, J. Acoust. Soc. Am. 109, 2332−2335 (2001). 9. S.E. Dosso, G.H. Brooke, S.J. Kilistoff, B.J. Sotirin, V.K. McDonald, M.R. Fallat, and N.E. Collison, IEEE J. Oceanic Eng. 23, 365−379 (1998).

SESSIONS

Detection and Classification of Targets Buried in Sandy Seafloor Using a 3d Sediment Sonar D. Brecht, H. Peine Forschungsanstalt der Bundeswehr für Wasserschall und Geophysik Federal Armed Forces Underwater Acoustics and Marine Geophysics Research Institute Low frequency sediment sonars allow the detection of objects which are buried in the seafloor and are invisible for conventional sidescan systems. Due to the considerably lower resolution of these systems, classification becomes a critical key issue. The resolution can be improved significantly by synthetic aperture processing, provided the motion of the ship can be compensated. Multi aspect analysis can further improve the classification by extracting additional information from the sonar data. Data from a new two dimensional multi aspect sonar system are used to investigate and demonstrate the potential of these concepts. The ship mounted System provides a high degree of flexibility with respect to target selection and aspect angle diversity. We present some experimental results from the first sea trials in 2000.

INTRODUCTION Mines buried in the seafloor are a dangerous threat since they can not be detected with conventional minehunting sidescan sonars. Sediment sonars use lower frequencies which can penetrate into the sediment. Resolution is essential in order to discriminate the attenuated signals from the strong surface and volume reverberation. In order to maintain a reasonable resolution, non-parametric sediment sonars need large apertures and, since their size is restricted, new signal processing concepts. The demonstrator system EXSESO (Experimental Sediment Sonar) allows to investigate detection, classification and signal processing concepts under real conditions. It was developed at STN Atlas and the first sea experiments, mainly for test purpose, were completed in September 2000. As a difference to many other systems, EXSESO is a ship mounted sonar which is not designed for experiments under laboratory conditions but has to struggle with many real world problems. A mayor advantage lies in its flexibility to set up experiments in a lot of different areas, environments and conditions. A main research topic is, consequently, detection and classification of buried objects under a variety of different conditions (sediment types, environments, grazing angles,...). Further topics are validation of simulation models, development of display concepts for 3-/4-dimensional seafloor data and, finally, investigation of applicable signal processing algorithms which may improve resolution, image quality and classification under real

conditions.

CONCEPT The antenna consists of a 20kHz projector and a 2dimensional (square) phased array receiver. The receiver has an aperture of about 1.3m length and is divided into 18 elements in across track (vertical) and 6 elements along track (horizontal) direction. At 20 kHz the beam has a width of 3° and can be steered electronically in both directions. The 6 along track elements may also be used to support autofocus or motion compensation algorithms. The whole system is mounted below the research vessel 'Planet'. It can be mechanically tilted up to 90° in starboard direction. Additionally, the receive beam can be steered electronically within a sector of 30° in across track and 18° in along track direction, which is insonified by the transmit unit. While the ship passes a target, the object is seen multiple times from different aspect angles. The resolution of 3° should be sufficient for detection, but may even be further improved by synthetic aperture processing. The system is designed as a real world demonstrator system which is operated from a moving platform and has to cope with the problems of the real world like ship motion, bad weather and environments, false targets,... On the other hand, it provides much experimental flexibility with regard to different areas, environments, sediment types, object types (shape, size, material), aspect/grazing angles.

SESSIONS

FIRST RESULTS For the first experiments we used a concrete filled metal cylinder which was buried into a sandy sediment. At first, the cylinder was flush buried by divers. The cylinder was placed within the sediment, but it was not completely covered by sand because a small crater had been created by the burying process. Nevertheless, the object could not be detected by a sidescan system which was used for comparison. It did not produce a significant highlight in the image, and the shadow which is essential for classification of mines was missing completely. In the image which was produced by the sediment sonar (Fig. 2), the object was clearly visible. Later on, the cylinder was buried deeper into the sediment. This time it was completely covered by a layer of sand and not visible by a TV system any more. As expected, the detection proved to be a lot more difficult than in the first case. Nevertheless, Fig. 3 shows an example where the target was still clearly detectable. In the sediment sonar image, the object produced a strong echo which can, at least in this example, clearly discriminated from the surrounding reverberation.

FIGURE 2. Flush buried concrete cylinder.

FIGURE 1. Experimental setup of the Exseso Demonstrator System. FIGURE 3. Deep buried concrete cylinder

.

SESSIONS