SURFACE TEXTURE DETERMINATION METHOD AND APPARATUS
The present invention relates to the analysis of signals that are scattered from a. surface in order to determine the type of surface. This may be used for example, for geological mapping of seabeds, biological monitoring and environmental control, underwater archeology, pipeline surveys and underwater cabling and mineral resource mapping.
In particular, the invention is applicable to sonar apparatus, such as sidescan (including phased and and bathymetric) or multibeam sonar apparatus, and may used for determining the texture of seabeds which may, for example, be composed of mud, clay, sand, gravel, rock, stones or other materials.
Sidescan sonar apparatus is typically mounted on a ship or on a "fish" (a float that may be controlled to maintain a substantially constant depth) towed behind the ship. A pulse (or ping) is transmitted towards the surface and a receiving device records the signal that is returned. As shown in Figure 1 of the accompanying drawings, the signal is scattered from the surface and therefore the sidescan sonar device SS monitors the backscattered signal from a range of distances across the surface in a direction perpendicular to the direction B in which the ship is travelling. The range of distances d correspond to a range of angles of elevation θ (known as the grazing angle) and a range of azimuth angles α. Since signals from more distant sections of the surface take longer to return, increasing time within a scan sample (the section of the received signal that corresponds to a single sonar pulse) is, for slightly rough surfaces, approximately equivalent to increasing distance. Therefore in the following description reference to "time" may be considered to be analogous to the azimuth angle and therefore is. related to the grazing angle θ.
After a given time, in which the ship will have progressed a certain distance, another ping is transmitted towards the seabed, allowing a further scan of the seabed to be performed. By combining a plurality of successive scans of the seabed, a strip of the seabed is surveyed as the ship passes.
The amplitude of the received signal, referred to as a backscatter profile, is a function of the attenuation caused by the depth of water, the angle of reflection, interference caused by the reflection from a plurality of independent scatterers which is dependent on their spatial distribution and the directivity of the transducers used. Consequently it is not possible to classify the texture of the surface by simple visual analysis of the received signal.
Previously, methods of analysing the received signal have been known in which a two dimensional image of the strip of the surface which has been scanned is produced and subjected to standard image processing techniques. This typically requires a large amount of computation which may prevent it from being performed at the same time as the signals are received (i.e., in real-time). Furthermore, in the processing prior to the image processing a lot of potentially useful information is lost. Other methods have relied on extracting spectral parameters from an entire scan sample. However these methods are not suitable if the strip of surface being scanned contains more than one texture since the scan sample may contain a combination of two or more different textures. Furthermore, these methods typically require preprocessing of the signal.
It is therefore an object of the present invention to provide a method and apparatus for determining the distribution of the type of surface from a signal scattered from the surface which requires less computation than previous methods and apparatus, which can determine the presence of more than one type of surface in a signal and which does not necessarily require processing of the backscatter profile to correct for the depth of the water, the angle of incidence, the directivity of the transducers or the alteration of the signal in the water, for example. The present invention achieves this and other objectives by providing a method of determining the distribution of a type of surface, comprising: conducting wavelet analysis of a plurality of signals scattered from a surface; and
using a plurality of energy coefficients from one or more scale bands of the wavelet analysis to determine the type of surface at a plurality of discrete locations on said surface; and using the determined types of surface at said plurality of discrete locations to determine the spatial distribution of one or more said types of surface on said surface.
This method is advantageous because it can be used to map areas with different types of surface using much less computation than previous methods and can determine the different types of surface in different sections of a single signal. A further advantageous feature of this method is that it may be performed either at the same time as the signals are being received or alternatively the signals may be stored on a storing medium and the method may then subsequently be performed. The method may therefore also be used to analyse old signals that pre-date the present invention. A further advantage of this method is that it does not require preprocessing of the signal to correct for factors such as alteration of the signal, depth of the water, the angle of incidence or the directivity of the transducer.
In a preferred embodiment, the invention further comprises a step of calculating a moment of the energy coefficients from each of said one or more scale bands. The plurality of coefficients correspond to a section of the backscattered signal and therefore a section of the surface. One or more of the moment values may be used to characterise the seabed. This is advantageous because it reduces the effect of any anomalous sections of the signal. The moments used may, for example, be the mean or the variance of the energy coefficients.
In a further preferred embodiment, coefficients from two scale bands are used to characterise the surface. This is advantageous because it means that the surface may be characterised by a minimal amount of information, which may in consequence reduce the time take to analyse the information.
In a yet further preferred embodiment, the coefficients used to characterise the surface are taken from scale bands corresponding to high frequencies. The use of
these scale bands is advantageous since they more clearly distinguished different types of surfaces.
In a yet further preferred embodiment, the invention comprises calculating for each of one or more scale bands one or more moments of a plurality of the energy coefficients. Each of the moments of the energy coefficients within a scale band corresponds to a different part of the received signal and therefore to a different part of the surface. This is advantageous as by this means it is possible to determine the type of surface at a plurality of locations within a single scan of the surface.
In a yet further preferred embodiment, a plurality of successive received signals are processed and the results combined to produce a representation of the arrangement of the different types of surfaces on a portion of the surface. Advantageously this may then be presented on a display to enable a user to view the arrangement of the different types of surfaces.
According to another aspect of the present invention, there is provided a computer program comprising computer code means for performing the method of determining the distribution of a type of surface.
According to a yet further aspect of the present invention, there is provided an apparatus for determining the distribution of a type of surface.
The invention iu nυw ue described by way of non-limitative examples with reference to the accompanying drawings, in which:
Figure 1 shows the geometry of a sidescan sonar apparatus;
Figure 2 shows a backscatter profile and its continuous wavelet transform; Figure 3 shows a dyadic wavelet transform of a similar backscatter profile to the one shown in Figure 2;
Figures 4a, 4b and 4c show pairs of time-averaged squared coefficients, calculated from Pace's data, plotted against, each other;
Figure 4a shows the time-averaged wavelet squared coefficients, calculated with the Dyadic Wavelet Transform, from the first and third high frequency bands with scales 21 and 23.
Figure 4b shows the time-averaged wavelet squared coefficients calculated with the Discrete Wavelet Transform using the Daubichies wavelet of order 4, from the first and third high frequency scale bands;
Figure 4c shows the time-averaged wavelet squared coefficients calculated with the wavelet packet transform, from the high frequency scale bands.
Figure 5 shows a backscatter profile that was artificially created to simulate a surface that is half stone and half sand;
Figure 6a shows time-averaged squared coefficients for a single scale band calculated from the profile of Figure 5 using the Dyadic Wavelet Transform;
Figure 6b shows the time-averaged dyadic wavelet coefficients from a single scale band for the backscatter profile of Figure 5 using a moving window with a size of 256 samples;
Figure 7 shows the dyadic wavelet coefficients from a single scale band for a sequence of backscatter profiles;
Figure 8 shows pairs of time-averaged wavelet squared coefficicina, calculated from airborne sonar data plotted against each other, calculated with the Dyadic Wavelet Transform from the first and third high frequency bands with scales 21 and 23; and
Figure 9 shows the region of the plot of Figure 8 that is close to the origin.
In a first example, a sonar arrangement such as that shown in Figure 1 is used to provide backscatter profiles of a seabed. Each profile corresponds to the backscatter from a sonar pulse transmitted towards the surface, in this case the seabed. The backscatter. profile is the amplitude of the received signal at a range of angles (corresponding to the distance from the ship or fish and received over a sample of time). Each profile therefore provides information about the seabed on a line perpendicular to the direction of motion of the boat. By performing successive scans, information about a strip of the seabed can be determined. The width of the
strip (i.e. the length of each scan) is dependent on the depth of the water and the design of the sonar but typically may, for example, be between a few tens of metres and several hundred metres.
The backscatter profile may be recorded for future analysis or may be analysed in realtime depending on the intended use of the resulting information. The resulting information may include the type of surface within each sample, the location of different types of surface within each scan sample and the boundaries between the different types of surface. The information for a plurality of scans can then be combined to produce a representation of the arrangement of the types of surface on the seabed. The location of the sections of the seabed corresponding to each backscatter profile or part of each profile may be determined using, for example, a global positioning system and recorded with the profile so that the representation of the arrangement of the types of surface on the seabed can be overlaid on charts of the seabed or combined with additional data such as, for example, bathymetric data.
In the first example (using Pace's data set, described in N. G. Pace and H Gao "Swathe Seabed Classification" IEEE Journal of Oceanic Eng. 13(2): 83-90, April 1988) the sidescan sonar operates at a centre frequency of 48 kHz, transmitting a pulse of 1 ms duration and a nominal bandwidth of 2 kHz. The backscatter profile was recorded with a sampling period of 125 μs and each profile consists of 1024 samples. The data were pre-processed to account for attenuation and the angle of reflection. However the present invention may also be used with data that has had little or no pre-processing.
In order to analyse the data, a wavelet transform is taken of each backscatter profile (corresponding to a single scan). Many wavelet transforms may be used, such as a dyadic wavelet transform, a discrete wavelet transform or wavelet packets, using wavelet bases such as the Daubechies wavelet, the Symmlet wavelet or the Coiflet wavelet. Each of these is computationally efficient.
Figure 2 shows a backscatter profile for a single surface type~and its continuous wavelet transform. The dyadic wavelet transform (which is discrete in
scale) of a similar profile is shown in Figure 3. On the vertical axis are shown the coefficients for each scale (analogous to frequency), decreasing towards the top such that those scale bands representing the lowest frequencies are at the bottom. The horizontal axis represents distance across the scan (which is also time received). In order to determine the type of surface, a section of the coefficients
(referred to as a window) are squared and averaged. For example, the coefficients relating to the first quarter of the distance across the scan may be squared and then the average of the squares calculated. In the present case a window size of 512 samples was used. Since this relates to a section of the duration of the received signal, these are referred to as time-averaged squared coefficients. Two or more , scale bands are selected and the corresponding time-averaged squared coefficients are studied for each sample. Plotting the time-averaged squared coefficients for two or three scale bands enables the differences between distinct types of surface to be observed, as further described below. Figures 4a and 4b show such plots for Pace's data using pairs of scale bands, in this case the first and third scale bands. The Figures correspond to squared coefficients calculated using Dyadic and Discrete Transforms, respectively. Figure 4c shows an equivalent plot for which the coefficients were calculated using a Wavelet Packet Transform. The different shapes of the points relate to different types of surface, in particular to clay, mud, sand, stone, gravel and rock (cross, circle, diamond, asterisk; star, dot, respectively). As shown in these Figures, the different types of surface are associated with different regions on the plot. As seen from the Figures, in this case the results from using the Dyadic Wavelet Transform are especially effective in segregating different types of surface into distinct regions of the plot. For an area ofunknown surface type, the relevant time-averaged squared coefficients can be calculated from the wavelet transform of the backscatter profile. Then, by comparing them to a selected pair of squared coefficients for known surface types, it is possible either approximately by eye or by means of an automatic classifier to determine the type of surface that is producing the backscatter profile. An example of a suitable automatic classifier is a multi-layer perceptron, although
other automatic neural networks, such as radial basis functions and support vector machines may also be used. Squared coefficients from scale bands corresponding to frequencies that are close to the Nyquist frequency (half the sampling frequency) are especially good for distmguishing different types of surfaces. Although in the example shown squared coefficients from two scale bands are used, in alternative embodiments of the present invention squared coefficients from a single or three or more scale bands may be used. For example use of squared coefficients from several scale bands may be used to not only determine the type of surface, such as gravel, for example, but may also permit a determination of the particle size, for example the grade of the gravel or sand.
The extent of the sample that is averaged to produce the time-averaged squared coefficients depends on the situation. Where a surface is known to be relatively homogeneous a relatively large time average may be used. However, in cases where there is likely to be a variety of surface types within a single scan a smaller time average is advantageous. In a development of the first example a moving window is used to produce averages of the squared coefficients from a single scale band. This enables the determination of the location of different types of surface within a single scan.
Figure 5 shows a combined backscatter profile that was created to simulate a surface with two types of surface. Half of the signal is from a stone seabed and the other half is from a sand seabed. Figure 6a shows the squared coefficients for a single scale band calculated from this backscatter profile using the Dyadic Wavelet Transform and Figure 6b shows a moving average of the selected squared coefficients. At the halfway point there is a marked drop, indicating the transition from one type of surface to another. Once the moving averages of the squared coefficients have been calculated, the determination of the boundary may be approximately judged by eye from a display, may be determined using, an edge detection technique or the values may be analysed using an automatic classifier, as . described above. Additionally, moving averages of squared coefficients from two or more scale bands may be combined to determine the location of the boundary.
The information from a. plurality of scans may be combined to show the boundary between two, or more, regions of differing surface types within the complete strip of seabed scanned. Such a plot is shown in Figure 7, which shows the non-averaged squared coefficients from a high frequency scale band of the dyadic wavelet transform on the vertical axis and the distance across the scan and the distance along the strip of the surface on the horizontal axes. The boundary between the two types of surface is clear. Combining the information from a plurality of scans and/or a plurality of strips of scans, the distribution of the types of surface over any given area can be determined. In a second example ofthe present invention, airborne sonar is used. Data to produce Figures 8 and 9 was collected using a Continuous Tone Frequency Modulation (or CTFM) system. This is provided with a single receiver and a transmitter fixed at a given angle of elevation. The transmission frequency has a sawtooth pattern, the frequency varying between 45 kHz and 90 kHz (the wavelength is approximately 3 to 7 mm). The sampling frequency is 25 kHz, each signal being comprised of 4014 samples representing amplitude versus range.
As in the first example, dyadic wavelet transforms were used and coefficients from two selected scale bands are squared, time-averaged and plotted. For the results shown in Figures 8 and 9, the first and third scale bands were selected and the size of the windows was 512 samples. Figures 8 and 9 show the pairs of coefficients for carpet (diamond), plastic (dot), grass (asterisk), gravel (cross) and asphalt (circle) surfaces. Figure 9 shows the plot zoomed in towards the origin. Again the different types of surface are easily distinguished by a single pair of time-averaged squared coefficients. The pairs of time-averaged squared coefficients could therefore be displayed graphically to allow approximate judgement ofthe type of surface by eye or an automatic classifier could be used, as described above.
Although in the above description and examples, the signals used have been sonar, the present invention is not limited to such signals. The invention can be performed using any form of signal creating a backscatter profile. .