WO2013055236A1

WO2013055236A1 - System and method for modelling multibeam backscatter from a seafloor

Info

Publication number: WO2013055236A1
Application number: PCT/NZ2012/000185
Authority: WO
Inventors: Jiashun YU; Stuart Anthony HENRYS; Colin Robert BROWN; Ivor Francis MARSH
Original assignee: Institute Of Geological & Nuclear Sciences Limited; National University Of Ireland Galway
Priority date: 2011-10-13
Filing date: 2012-10-12
Publication date: 2013-04-18

Abstract

The invention provides a method of modelling multibeam backscatter from a seafloor. The method comprises determining, with a computer processor, an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least partly by a local beam incidence angle and a beam amplitude value; calculating, with a computer processor, an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculating, with a computer processor, one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor. The invention further provides a system configured to model multibeam backscatter from a seafloor.

Description

SYSTEM AND METHOD FOR MODELLING MULTIBEAM BACKSCATTER

FROM A SEAFLOOR

FIELD OF INVENTION

The invention relates to a method and system of modelling multi-beam backscatter from the seafloor, using a combination of a boundary integral and Lambert's law.

BACKGROUND TO INVENTION

In ship-mounted multi-beam surveys the array or swath of beams, orthogonal to the vessel's direction, are generated by a single ping. The data measured includes co-registered bathymetric and acoustic backscatter data at many different locations along a strip of the seafloor ensonified by the ping.

With appropriate calibration, multi-beam echo sounders provide accurate bathymetry and a measure of the seafloor acoustic backscatter strength as a function of grazing angle. The data have been used inter alia to:

^■ interpret seafloor topography

^■ understand surficial sediment processes

^■ classify benthic habitats

^■ identify geotechnical hazards and sub-seafloor gas seeps

^■ assess environmental impacts.

Much effort has been made to transform swath acoustic multi-beam data into maps that provide quantitative geological and biological information. However, compared to existing sophisticated algorithms used for multi-beam data processing, there is a paucity of practical tools for the inference of seafloor physical properties.

The beams comprising a multi-beam ping propagate at different angles through the water column until they intersect the seafloor where the propagating wavefront is subject to scattering. The - - scattering redistributes the incident acoustic energy in multiple directions. The energy returned to the transducer from each beam carries important information about the seafloor morphology and physical properties.

The backscatter angular response provides the optimum data set for the inference of the seafloor acoustic properties.

Backscatter is caused by the complex interaction of a wavefront with roughness and contrasts in velocity and density at the water-seafloor interface. Backscatter is also caused by inhomogeneities due to velocity and density contrasts within the sediment volume. The seafloor might be modelled as fluid, elastic or poroelastic media and might include the effects of attenuation, stratification, anisotropy, density and velocity gradients, interface waves, shear waves,

compressional and shear wave conversions, Biot slow waves, gas bubbles and sub-surface volume heterogeneities.

The scattering mechanisms are coupled. The contribution of each mechanism to the observed backscatter is frequency-dependent and multiple scattering may take place. The combined effects of all scattering mechanisms can be accurately simulated using finite-difference methods (FDM). However, such models are computationally intensive for modelling the multi-beam response in the typical frequency range of swath systems. This typical frequency range is 12-450 kHz.

Detailed modelling taking into account all of these factors is also impractical because multi-beam surveys are not normally accompanied by geophysical and geological surveys to determine poro- and visco- elastic rock properties necessary to model adequately the data or interpret the results.

It is also unclear whether scattering observed over a wide range of grazing angles and frequencies is dominated by the roughness of the sediment-water interface alone or whether volume scattering is also important. The contributions of seafloor roughness and volume scattering to the backscatter response may vary from one location to another. - -

Furthermore, seafloor roughness that is on a length scale below that of a sonar wavelength is at the limit of resolution but may still cause significant wave scattering.

It is an object of the invention to provide a modelling technique that enables a practical inversion of remotely acquired backscatter data for the inference of seafloor acoustic properties, or that will at least provide the public with a useful choice.

SUMMARY OF THE INVENTION

In broad terms in one form the invention comprises a method of modelling multibeam backscatter from a seafloor, the method comprising determining, with a computer processor, an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least partly by a local beam incidence angle and a beam amplitude value; calculating, with a computer processor, an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculating, with a computer processor, one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor.

Preferably the method further comprises applying a weighting function, with a computer processor, to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam-form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

Preferably the method further comprises calculating the approximate beam amplitude at the source location for the two or more beams within the incident beam wavefield.

Preferably the incident beam wavefield comprises two or more beams directed at different azimuths for a single ping.

Preferably the method further comprises calculating a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor. - -

In broad terms in another form the invention comprises a system configured to model multibeam backscatter from a seafloor, the system comprising a memory; and a processor programmed to determine an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least pardy by a local beam incidence angle and a beam amplitude value; calculate an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculate one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor.

Preferably the processor is further programmed to apply a weighting function to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam- form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

Preferably the processor is further programmed to calculate the approximate beam amplitude at the source location for the two or more beams within the incident beam wavefield.

Preferably the processor is further programmed to calculate a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor.

In broad terms in another form the invention comprises a computer-readable medium having computer-executable instructions that, when executed by a processor, cause the processor to perform a method of modelling multibeam backscatter from a seafloor, the method comprising determining an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least pardy by a local beam incidence angle and a beam amplitude value;

calculating an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculating one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor. Preferably the method further comprises applying a weighting function to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam-form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

Preferably the method further comprises calculating a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor.

BRIEF DESCRIPTION OF THE DRAWINGS

Preferred forms of a system and method for modelling multi-beam backscatter from a seafloor will now be described with reference to the accompanying figures in which:

Figure 1 shows the process of emitting beams toward a seafloor.

Figure 2 shows a preferred form environment in which the invention is intended to operate.

Figure 3 shows an overview of a preferred form method of the invention.

Figure 4 shows a more detailed view of the beams of Figure 1.

Figure 5 shows a model of a water layer and a seafloor.

Figure 6 shows the power spectrum corresponding to seafloor models.

Figure 7 shows examples of synthetic seismograms. - -

Figure 8 shows examples of directionally beam forming FDM time domain backscatter seismograms.

Figure 9 shows a preferred processor in accordance with the invention. DETAILED DESCRIPTION

In multi-beam surveys, as shown in Figure 1, a vessel 100 floating on the surface 110 of a body of water 120 directs an array 130 of beams toward the seafloor 140. This array 130 of beams generated by a single ping ensonifies a strip of seafloor 140. This allows the simultaneous measurement of high-resolution depth data and co-registered acoustic backscatter at many different locations of the seafloor.

The swath is wide in the across-track direction and is approximately three times the water depth. The swath is narrow in the along-track direction. This technique is therefore suitable to a two- dimensional seafloor model where the bathymetry is invariant in a direction perpendicular to the swath.

The swath of ensonified seafloor 140 can be tens of metres wide or thousands of metres wide. The width of the swath depends upon the multi-beam sonar system and the height of the sonar system above the seafloor surface 140. The seafloor surface 140 has topographical variations at multiple scales. These topographic variations range from sub-wavelength (mm scale) to a large number of multiple wavelengths (km scale).

For multi-beam surveys using sonars operating in the frequency range 12-450 kHz, this results in an extremely large parameter model defining the seafloor surface 140 beneath the swath. It is not computationally practical to solve the wave equation system for such a model. The concept of a footprint combining micro-scale roughness with macro-scale topographic changes in the seafloor model has the potential to provide a solution. The proposed solution is not limited to the typical frequency range 12-450 kHz. Figure 2 shows a preferred form environment 200 in which the invention is intended to operate. The vessel 100 includes a multibeam echo sounder unit 205. The echo sounder unit 205 is typically mounted on the hull of the vessel or over the side, or is towed some distance behind the vessel. Echo sounder 205 typically includes a transducer that is long in the along ship direction and short in the athwart ship direction. This enables the unit 205 to generate a swath of sound that is wide in athwart ship and narrow in the along ship direction. Unit 205 also includes a receiver transducer 215. The receiver transducer 2 5 is wide in the athwart ship direction and narrow in the along ship. This enables the unit 205 to receive the wide swath produced by the transmitter , but still have narrow beams in the along track direction.

The sound generated by the transmit array 210 is reflected by the seafloor at different angles and is received by the receiver transducer 215 at slightly different times. The signals are processed by a computer processor 220 programmed and operating as a transceiver unit. A data processing module 225 controlling operation of the processor 220 converts the depth values and generates a visualisation in the form of a bathymetric map on a display 230.

Computer processor 220 is further configured to export processed data 235. These processed data 235 could be transmitted over a network in real time, could be transmitted over a network in a batch format, or could be stored on a computer readable medium which is then physically transmitted. In one preferred form the exported data 235 represents a series of acoustic beam values, each beam value having an associated amplitude and angle.

Computer processor 240 is configured to operate under control of data modelling module 245. The data modelling module 245 causes the computer processor 240 to perform the data modelling techniques described below. The modelled data 250 is exported using one or more of the methods described above in relation to processed data 235.

The modelled data is further processed by computer processor 255. The computer processor 255 operates under the control of visualisation module 260. The computer processor 255 and visualisation module 260 together transform the modelled data 250 to generate a modelled map of the seafloor suitable for display on computer display 265. - -

The invention provides a preferred form alternative to the prior art, which involves the transfer of such data 235 directly from computer processor 220 to computer processor 255. An owner of the preferred form method 300 is illustrated in Figure 3.

Defining a footprint model 305

The preferred form technique includes calculating a set of coordinates defining a footprint comprising the segment of the seafloor 140 ensonified by an acoustic beam directed as a function of azimuth from a source point 100.

Referring to Figure 4, the echo sounder 205 from vessel 100 directs a beam array 130 toward the seafloor 140. A footprint is defined as a segment 400 of the seafloor 140 with effective coverage of a beam 130. The midpoint 410 of the footprint 200 is the intersection of the seafloor 140 with the central axis of the beam 130.

The continuous seafloor is thus separated into a set of footprints. These footprints can be overlapping or contiguous. The techniques below define the multi-beam amplitude at a given beam grazing angle as the backscatter from the associated footprint. The modelling of a multi- beam ping is therefore reduced to the computation of the backscattered wave field from the footprint associated with each beam within a ping. As the size of a footprint is normally approximately 100^th of the size of the seafloor ensonified by a ping, the computational requirements have the potential to be significantly reduced.

The model is characterised by survey parameters, geometrical parameters and physical parameters.

The survey parameters comprise the frequency f of the signal generated by the ping, the beam local incidence angle Θ between the central axis of the beam cone and the nadir direction, and the angular width φ of the beam. This angular width is twice the angle between the effective edge and the central axis of the beam cone. - -

The geometrical parameters comprise the two-dimensional spatial coordinates (x_c , z_c ), which represent the horizontal and vertical distances of the midpoint of the footprint from the source, the dip angle μ of the footprint, and the roughness parameter η characterising the small scale random variation of the footprint. With regard to the dip angle, μ>0 is defined for a seafloor dip in the direction of x. For example, μ<0 for the geometry in Figure 1.

The physical parameters comprise the density and acoustic velocity of the water ( _w , v_w ) above the footprint, and the underlying material forming the seafloor (p_s , v_s).

Modelling topographic variations within the footprint 310

A preferred form embodiment includes calculating a stochastic model of the seafloor characterising fine-scale topographic variations (referred to as roughness) within the footprint.

The footprint 400 is the superposition of a sloping linear segment of the seafloor 140 and fine scale topographic variations referred to as roughness. The set of linear segments approximates the topographic variations of the seafloor at the scale commensurate with the horizontal resolution of the swath bathymetry. The model parameters may vary from footprint to footprint. The footprint morphology Ω is therefore defined as

z = z_c + x s μ + z (x) cos μ; x_x < x≤ x_N

where (3c, z) are the coordinates relative to the localised coordinate system. 3c is the distance from the midpoint (x_c , z_c ) along the linear segment of the point of interest within the footprint.

Figure 4 shows 3c, and x_N as the intersection points of the linear segment with the effective edges of the beam, given as

_ _ z_c t (0 - <p/ 2) - x_c

——————— _——_____ ^ (2) cos μ - sin // tan(# - ^ / 2) - -

z_c tan(0 + (p/2)-x_c

cos μ - sin μ tan(# + ω/2)

The coordinate z ( ) , normal to the sloping linear segment, is given by the random roughness amplitude of the footprint ζ =η . (4)

The variable η is defined as a random number from a random series Y≡ (η_ι, η₂ , ·^■■ , η_Ν )^Γ given by the expression:

Y=LW, (5) where W≡ (ς_ι , g₂ ,^■ ·^■ , g_N )^T is a standard white noise series, and L is a N x N matrix such that

LL^T =R, (6) where the Nx N matrix, R is determined by the correlation coefficient function, r x_i— .), i,j = \,···,Ν of the amplitude of the random footprint.

It is preferable to use two parameters (fl,/) to characterise the footprint, where a is the standard deviation of the roughness amplitude within a footprint, and / is a distance between two random amplitudes where their correlation reduces to half. The parameters a and / are presumed to be indicative of the average height and proportional to the width of scatterers, respectively, on the seafloor surface, such that

rfljc. -Xj \=l)=a^z/2. (8) - -

The correlation function can take many forms which satisfy eqs. (7) and (8).

One preferred form definition is:

r(x_j - x_J )≡2 ¹ a² (9)

A further description, relative to the localised coordinate system (x , z) , of the same footprint is: Q≡{x,z \ z = z (x ), x_x≤ x≤ x_N }. (10)

Then Ω is a projection of Ω on the universal coordinate domain.

Calculating scattered wavefield on the footprint 315

The preferred form technique includes calculating the incident beam wavefield on the footprint for a given frequency and calculating the scattered wavefield on the footprint using a combination of a boundary integral and Lambert's law approximation.

The amplitude p , of a wave at location r > 0 in an infinite, homogenous, isotropic medium that is generated by a point source at r = 0 , is given by the solution of the wave equation for a point source. This is the Green's function for the wave equation ikr

p(r) = S— , (11)

r where the wave amplitude S , at a location one unit distance away ( r = 1 ) from the source, represents the strength of the source. - -

Huygens' Principle states that each point of an advancing wavefront is the centre of a fresh disturbance and the source of a new train of waves, and that the advancing wave as a whole may be regarded as the sum of all the secondary waves arising from points in the medium already traversed.

One preferred form involves approximating the amplitude at r = 0 of a wave scattered back from the footprint Ω as coming from many point sources located on Ω . The point sources are then summed together using a surface boundary integral

where p^sc ( ) , x e Ω, is the amplitude of a backscattered wave on the surface at . The path distance is r (x, z) e Ω . (13)

The wave number in the water is

a>— Ί , and ds is an area element, which is a length element in the two-dimensional model

Applying a weightingfunction 320

The preferred form technique involves optionally applying a weighting function to the scattered wavefield to beam-form the backscatter wavefield to the source point. - -

The recorded backscatter is sampled by a beam over the footprint. To simulate the energy distribution of the beam, the integral need to be weighted so that

where the incident angle 9 is given as

,9 =arctan(x/z) (x,z)e Q . (17)

Within the footprint, i.e. (x, z) e Ω , the incident angle ranges within the effective flare angle of the beam, i.e. θ - φΙ 2≤3≤θ + φ Ι 2 . The weighting function w( ff) is determined by the configuration of the recording instrument.

In general, it would be weighted strongly at the centre, and weaker towards the edge of the beam. It is assumed to follow a Gaussian distribution because it has a simple analytic expression for the beam divergence and a predictable minimum width for a given propagation distance:

The factor (φ I n)¹ is the variance of the Gaussian distribution. The larger n is, the stronger the central part of the beam is weighted. Hence n>1 can be used to shape the weighting function for the beam.

On the other hand, the beam width φ and n should be constrained so that the truncated tail beyond the two ends of the effective edge of the beam is negligible. For example, n—2 leaves the truncated tails accounting for ~ 4% of the distribution. The truncated tails are negligible. -

Calculating backscatter wavefield 325

The preferred form technique further includes calculating the backscatter wavefield over a range of acoustic beam azimuths and frequencies to form backscatter response functions.

It is assumed that the roughness of the footprint is such that the scattering follows Lambert's Law. Therefore, for the backscattering case, of which the incident and scattering are at the same angle, the following applies: p^sc (x) = Α οο8² (3 + μ) ρ^ίη (χ) G Q , (19) where the incident wave on the footprint is ik_w r

pⁱⁿ (x) = S-— ϊ ε Ω , (20)

r

Inserting (19) and (20) into (16), yields

2ik_w r

p^sc (0) = AS f ^— - cos² (5 + i)w(5) ifc . (21)

r

The term 1 _s = A cos² (i + μ) is the backscattering intensity of the rough seafloor. For normal incidence, i.e., i + = 0 , I _s = A . Hence, Λ can be interpreted as the strength of backscattering at normal incidence. To use eq. (21) to calculate the backscatter, the value for Λ needs to be determined.

Observations show that Λ varies over a very wide range of seafloor types (-20.2 to -5.6 dB) and survey signal frequencies used. However, the fact that Λ generally increases with seafloor impedance encourages use of a proportional relationship to approximate the change of Λ with respect to the normal incidence reflection coefficient ξ between the sea water and the sub seafloor medium, i.e.,

Λ= οξ (22) where

and where cis a constant independent of medium properties.

The finite difference benchmark modelling discussed below has found that using c—2.83 (or 101og₁₀ c = 4.5 dB) is a very accurate approximation.

Numerical solution of the integral

The variable [x_{ ,x_N ] is digitised with an interval of δ

7 = 1,···, N, (24) where the number of samples N is defined as

N≡ x_N -χ^/δ + Ι . (25)

The following equations apply

x ≡ x_c+ Xj cos μ - Zj sin μ, j = h-;N, (27) z_j≡z_c+x_js μ + z_Jcosμ, j = !,^■■■, N , (28)

Ζ ≡ ζ(χ_})≡η j = l,— ,N . (29)

Equation (21) can therefore be changed into a discrete form

where

7 = 1

dS_j =^(dx)² +(dz_j)² — S² ₊(^^)² j = 2,.; ,NN--ll , (31)

j = N

with

9 _j = arctan(x / z_j ) j = !,^■··, N. (33)

Computational calculations

The following algorithm describes the procedures to model the backscatter of a single multibeam ping.

1. For an arbitrary footprint on a (two-dimensional) seafloor ensonified by a ping in a multibeam survey, the end coordinates X, and X_N of the sloping footprint segment are - - calculated according to eqs. (2) and (3). Lambert's coefficient Λ is determined using eqs. (22) and (23).

2. For a given interval of δ and number of samples N determined by eq. (25), the discrete coordinate values, X - , j = 1, · · · , N ate calculated according to eq. (24).

3. A random series, η j = l,- - -,N characterising the roughness of the footprint, is assigned to the Zj , j = 1, · - ·, Ν coordinate according to eq. (29).

4. Coordinates Xj and Zj , j = \, - - -, N are calculated using eqs. (27) and (28). The path distance r, , , j— 1, · · · , N , incident angle = !,·^■·, N and flare angle difference d j , j = 1, · · · , N are calculated using eqs. (26), (33), and (31).

5. The weighting function Wj , j = 1,· · · ,Ν is calculated according to eq. (32) and the backscatter at the receiver point p_sc (0) is calculated according to eq. (30) .

6. Steps 1 to 5 are repeated for all the footprints over the angular range of the two-dimensional model to complete the synthesis of the backscatter record for a single ping.

RESULTS

Example calculations using the techniques described above for the problem of a rough seafloor are given in this section. The results of these examples are validated by comparison with calculations made using a fourth-order time-domain acoustic finite difference method (FDM).

The FDM method can be used to yield a full acoustic wave solution to the surface scattering problem over a broad frequency spectrum and without assumptions of the underlying physical geometry. This allows the direct analysis of the rough seafloor scattering physics and permits a benchmark of the techniques described above. The standard model in Figure 5 consists of a water layer 500 that is 30 m deep overlying a rough seafloor. The acoustic wave velocities are 1500 m/s, for the seawater and 1600 m/s for the seafloor sediment (reflection coefficient, ξ- 0.032). The seafloor is horizontal, but four variations of random roughness are introduced characterised by root mean square (RMS) amplitudes (a) ranging from near smooth 0.5 cm for Model A, 1.0 cm and 2.0 cm for Models B and C, respectively, to rough 4.0 cm for Model D. Roughness correlation length ( ) is 40 cm for all models.

The power spectra corresponding to the seafloor models are shown in Figure 6. The slope of each curve is the same (-3.6) corresponding to the constant correlation length, but the intercept value (at spatial frequency equal to 1) increases with RMS amplitude. Although power spectra of seafloor microtopography together with grain size distribution, physical properties data, and photographs are few, there are no reliable correlations of seafloor statistics and geology. In general, however, fine grained sediments have lower values of RMS roughness and models described here with reference to the above table are loosely representative of a shallow seafloor comprising siltstones and sands to course grained sediments.

The techniques described above represent the monostatic response as a function of grazing angle. These solutions were calculated for beams of grazing angles ranging from -56° to 56°.

For the acoustic FDM, the model space of interest is limited within the triangular area between the two outer rays (grazing angle -56 ° to 56 °) and the seafloor as shown in Figure 5. Since the backscatter amplitude from the rough seafloor is extremely weak, compared with the direct arrivals or reflections, the modelling needs to be able to produce very clean waveforms without contamination from boundary reflections.

The models used for the FDM were 76 m by 76 m (30401 by 30401 grid points), with the point source located in the centre and 30 m above the seafloor and 8 m thick layer of sediments as shown in Figure 5. To adequately represent the rough seafloor boundary in this model, with minimal grid dispersion, it requires a grid interval of 0.25 cm, corresponding to 50 grid points per wavelength for a 12 kHz sonar signal represented by a Ricker wavelet. This spacing interval, which presents the roughness amplitude as fine as 0.125 cm, is 1/4 to 1/32 of the RMS amplitudes of Models A to D.

A Ricker wavelet was used as source signal, with the amplitude equal to 1 at 12 kHz. The synthetic seismograms shown in Figure 7 are recorded at a level of 38 m and traces are selected across 3.58 m (1432 seismograms) at the source point to simulate the multibeam receiver array.

Sixty thousand (60000) time steps (interval of 9.0x10^"7 s) were calculated, representing a time history of 0.054 s for the wavefield.

To directionally beam-form FDM time domain backscatter seismograms, first a Fast Fourier Transform (FFT) was applied to each of the 1432 traces within the time window of 0.0390≤/≤0.0483 s, as shown in Figure 8. Then the 12 kHz frequency components were extracted and the multibeam backscatter response from grazing angles between -56° to 56° was calculated using the FFT in the space domain. This allows a direct comparison of the response to the FDM results for the same model as shown in Figure 8. In practice, the port and starboard backscatter responses of the FDM models were calculated separately and then averaged in Figure 8. Similarly, for the responses for two random models, with the same seafloor statistics, were computed and averaged. The average stacked curve for 100 random responses is also shown for each model in Figure 8.

Ping values for all four models were calculated using a constant c (eq. 23) of 2.83 (4.5 dB). The RMS values of the difference between the two methods are given in Table 1.

TABLE 1. Comparisons between BILL and finite difference calculated ping values. 2 and 100 represent a stack of 2 model calculations and a stack of 100. The average difference and root mean square (RMS) difference were calculated by first taking the difference of the backscatter values of each beam, and then average or RMS of the differences of all beams. - -

Average difference was calculated by first measuring the difference between the values of the two methods, for each beam, and then determining the average of the differences over all beams in the ping. In summary, modelling matches very well the FDM benchmark models of different roughness by no more than an RMS error ±2.1 dB.

The techniques described above provide a combined method using boundary integral and the Lambert's Law (BILL) to forward model multibeam acoustic swath backscatter from a two- dimensional model of the seafloor with roughness and physical property contrasts across the ocean/ sea-floor interface.

The examples discussed above demonstrate the potential to calculate highly accurate solutions to seafloor wave scattering problems that give excellent agreement with finite-difference calculations made over a range of seafloor types at 12 Hz. The techniques have been developed to simulate monostatic multibeam seafloor backscatter and are able to calculate full angular response curves across a range of frequencies. That is, the BILL solution is not limited to the typical frequency range of 12-450 kHz or to the narrow grazing angles (between -56 ° to 56 °) of the model described in Figure 5. Examples presented here include the same two-dimensional seafloor statistics for each beam within a single ping. However, different seafloor parameters can be accommodated in each beam footprint.

Fast and accurate solutions to the forward problem of seafloor scattering have the potential to enable ping calibrated backscattered multibeam data acquired across a range of frequencies and angular responses to be directiy inverted for seafloor properties. Figure 9 shows a simplified block diagram of a device forming at least part of computer system 220, 240 and/ or 255 in the example form of a computing device 900. Sets of computer executable instructions are executed within device 900 that cause the device 900 to perform the methods described above. Preferably the computing device 900 is connected to other devices. Where the device is networked to other devices, the device is configured to operate in the capacity of a server or a client machine in a server-client network environment. Alternatively the device can operate as a peer machine in a peer-to-peer or distributed network environment. The device may also include any other machine capable of executing a set of instructions that specify actions to be taken by that machine. These instructions can be sequential or otherwise.

A single device 900 is shown in Figure 9. The term "computing device" also includes any collection of machines that individually or joindy execute a set or multiple sets of instructions to perform any one or more of the methods described above.

The example computing device 900 includes a processor 905. One example of a processor is a central processing unit or CPU. The device further includes main system memory 910 and static memory 915. The processor 905, main memory 910 and static memory 9 5 communicate with each other via data bus 920.

Computing device 900 further includes a data input device 925. In one embodiment a display device 930 is touch sensitive. The display with the touch sensitive display device 930 and any icons and indicia displayed on the device 930 together form the data input device 925. The user is able to select various icons and indicia on the display device 930 by depressing or touching certain areas of the device 930.

It will be appreciated that the data input device 925 could take other forms for example a set of physical keys that either form part of the computing device or are in communication with the computing device. Alternatively the data input device 925 includes voice recognition to receive and act on spoken commands and other sounds.

Computing device 900 optionally also includes reader unit 935, network interface device 940, optical media drive 945, cursor control device 950, and signal generation device 955. - -

Reader unit 935 is able to receive a computer readable medium 960 on which is stored one or more sets of instructions and data structures, for example computer software 965. The software 965 uses one or more of the methods or functions described above. Reader unit 935 includes a disc drive and/ or a USB port. In these cases the computer readable medium includes a floppy disc and a static storage device such as a thumb drive. Where the optical media drive 945 is used, the computer readable medium includes a CD Rom.

Software 965 may also reside completely or at least partially within main system memory 910 and/ or within processor 905 during execution by the computing device 900. In this case main memory 910 and processor 905 constitute computer-readable tangible storage media. Software 965 may further be transmitted or received over network 970 via network interface device 940. The data transfer uses any one of a number of well known transfer protocols. One example is hypertext transfer protocol (http).

Computer-readable medium 960 is shown in an example embodiment to be a single medium. This term should however be taken to include a single medium or multiple media. Examples of multiple media include a centralised or distributed database and/ or associated caches. These multiple media store the one or more sets of computer executable instructions. The term

"computer readable medium" should also be taken to include any medium that is capable of storing, encoding or carrying a set of instructions for execution by a processor and that cause the processor to perform any one or more of the methods described above. The computer-readable medium is also capable of storing, encoding or carrying data structures used by or associated with these sets of instructions. The term "computer-readable medium" includes solid-state memories, optical media and magnetic media.

The foregoing describes the invention including preferred forms thereof. Modifications and improvements as would be obvious to those skilled in the art are intended to be incorporated in the scope hereof, as defined by the accompanying claims.

Claims

1. A method of modelling multibeam backscatter from a seafloor, the method comprising: determining, with a computer processor, an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least partly by a local beam incidence angle and a beam amplitude value; calculating, with a computer processor, an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculating, with a computer processor, one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor.

2. The method of claim 1 further comprising: applying a weighting function, with a computer processor, to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam-form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

3. The method of claim 2 further comprising calculating the approximate beam amplitude at the source location for the two or more beams within the incident beam wavefield.

4. The method of claim 1 wherein the incident beam wavefield comprises two or more beams directed at different azimuths for a single ping.

5. The method of claim 1 further comprising calculating a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor.

6. A system configured to model multibeam backscatter from a seafloor, the system comprising:

a memory; and

a processor programmed to:

determine an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least partly by a local beam incidence angle and a beam amplitude value;

calculate an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law

approximation; and

calculate one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor.

7. The system of claim 6 wherein the processor is further programmed to apply a weighting function to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam-form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

8. The system of claim 7 wherein the processor is further programmed to calculate the approximate beam amplitude at the source location for the two or more beams within the incident beam wavefield.

9. The system of claim 6 wherein the incident beam wavefield comprises two or more beams directed at different azimuths for a single ping.

10. The system of claim 6 wherein the processor is further programmed to calculate a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor.

11. A computer-readable medium having computer-executable instructions that, when executed by a processor, cause the processor to perform a method of modelling multibeam backscatter from a seafloor, the method comprising: determining an area of the seafloor ensonified by an incident beam wavefield, the beam wavefield represented at least partly by a local beam incidence angle and a beam amplitude value; calculating an approximate incident beam wavefield, scattered back from the ensonified area of the seafloor, using a boundary integral and Lambert's Law approximation; and calculating one or more beam incidence angles in the beam wavefield to model the effect of surface roughness of the seafloor.

12. The computer-readable medium of claim 11 wherein the method further comprises applying a weighting function to a backscatter wavefield scattered back from the ensonified area of the seafloor to beam-form the backscatter wavefield and to calculate the approximate beam amplitude at a source location.

13. The computer-readable medium of claim 12 wherein the method further comprises calculating the approximate beam amplitude at the source location for the two or more beams within the incident beam wavefield.

14. The computer-readable medium of claim 11 wherein the incident beam wavefield comprises two or more beams directed at different azimuths for a single ping.

15. The computer-readable medium of claim 11 wherein the method further comprises calculating a stochastic model of the seafloor to represent the approximate surface roughness of the seafloor.