Strategies for data acquisition using ultrasonic phased arrays.

Ultrasonic phased arrays have produced major benefits in a range of fields, from medical imaging to non-destructive evaluation. The maximum information, which can be measured by an array, corresponds to the Full Matrix Capture (FMC) data acquisition technique and contains all possible combinations of transmitter-receiver signals. However, this method is not fast enough for some applications and can result in a very large volume of data. In this paper, the problem of optimal array data acquisition strategy is considered, that is, how to make the minimum number of array measurements without loss of information. The main result is that under the single scattering assumption the FMC dataset has a specific sparse structure, and this property can be used to design an optimal data acquisition method. An analytical relationship between the minimum number of array firings, maximum steering angle and signal-to-noise ratio is derived, and validated experimentally. An important conclusion is that the optimal number of emissions decreases when the angular aperture of the array increases. It is also shown that plane wave imaging data are equivalent to the FMC dataset, but requires up to an order of magnitude fewer array firings.

Introduction

Ultrasonic arrays are now routinely used in industry for non-destructive testing and assessment of the structural integrity of materials. Modern computational power makes it possible to collect raw array data first and then perform image reconstruction in post-processing. One common data collection method, which is now available in many commercial ultrasonic array systems is Full Matrix Capture (FMC) [1-4]. In this case, the ultrasonic wave is consecutively transmitted by each array element and the reflected signal is measured by all array elements simultaneously. The advantage of this method of array data acquisition is that an FMC dataset contains the maximum possible information that can be measured by an array from a given position. It means that data corresponding to any array acquisition method can be synthesized from the FMC dataset and various imaging techniques can be applied [1,5-10]. Moreover, advanced signal processing becomes possible. For example, scattering matrices of defects can be extracted and used for their characterization [11-15].

However, the FMC data acquisition method also has some disadvantages. The first important issue is related to the volume of measured array data. If an ultrasonic array has N elements, then the corresponding FMC dataset consists of N2 time signals. In many industrial cases, an inspection is carried out by physically scanning a phased array to cover a large area. Although the data volume for a single FMC set is relatively small, the total amount of inspection data in this case can be extremely large, containing terabytes of data.

Secondly, in order to measure the FMC data N transmissions are required, which can take too long for some applications. Note, that the technological advances in the last few years have made it possible to perform ultrasonic imaging using the FMC data acquisition method in near ‘real time’, where this means frame rate of tens frames per second [16]. This is close to the fundamental physical limit, which is related to the finite speed of ultrasonic wave propagation and their decay in materials. Although this imaging rate is fast enough in some industrial cases, there are multiple applications which require a much faster data acquisition rate of several hundred (or, in some cases, even thousands) frames per second. Regarding industrial non-destructive testing applications, well-known examples are pipeline inspection using pigging [17] and railway track inspection [18]. In these cases, the inspection vehicle moves at tens of metres per second so an ultrasonic frame rate of 10 frames per second is wholly insufficient. In addition, a very high frame rate opens the possibility for many other applications, such as real-time monitoring of fatigue crack growth.

One possible solution to reduce the amount of measured data is based on compressive sensing theory [19], which provides a method to reconstruct data from far fewer samples than required by Nyquist sampling theorem. The main condition under which this reconstruction is possible is that the data are sparse. In the application of this theory to ultrasonic imaging the data to be reconstructed represents an image, and sparsity is usually interpreted as a small number of scatterers present in the material [20,21]. This approach is not generalizable and lacks a formal definition to assess its performance in different inspections.

The problem of reduction of the number of ultrasonic array firings has been intensively studied by many researchers. The common source method uses all elements to transmit simultaneously, so that a plane wave propagates into the structure and an image is created from the scattered waves [22-25]. Another possibility is to use just a few array elements as transmitters in order to produce the best image for a given number of firings [26,27]. However, in both cases the reduction in the data acquisition time comes at the expense of significantly reduced image quality. Alternatively, iterative imaging methods have been suggested in order to suppress the aliasing noise [28,29]. These techniques help to reduce imaging artefacts, but the signal processing time increases in proportion to the number of iterations. The main progress in this area has been achieved in the medical imaging field using plane wave excitations [30,31]. It has been shown that it is possible to achieve good quality imaging with a relatively small number of transmitted plane waves. The image quality in plane wave imaging is affected by a number of parameters, specifically the number of transmitted plane waves and the maximum steering angle. A method for optimizing these parameters based on Pareto optimality has been developed in [32]. It has also been demonstrated that plane wave imaging techniques give good results in non-destructive testing applications with much smaller number of emissions compared to the FMC method [33,34]. However, the physical reason for this has not been explained. Also, it remains unclear if there is a fundamental relationship between maximum steering angle and the number of plane wave transmissions.

The key question investigated in this paper is how to optimize the array data acquisition process so it measures the same amount of information as the FMC method, but using fewer measurements. The crucial observation is that under the single scattering assumption an FMC dataset has a specific sparse structure, which is independent of the number of scatterers. This property can be used to construct the optimal array data measurement procedure. It is shown that plane wave imaging represents one possible implementation of the optimal data acquisition approach. Also, the angular sampling criterion, which relates the number of plane wave to the maximum steering angle, is derived. The proposed plane wave sampling theory is experimentally validated on an aluminium sample with volumetric defects. Finally, the key applications of this approach to a range of fields are discussed.

Structure of ultrasonic array data

Formulation of the problem

In this paper, for simplicity, two-dimensional imaging using a one-dimensional linear array is considered. However, all results can be extended to three-dimensional imaging using a two-dimensional array. The system is illustrated schematically in figure 1. A two-dimensional elastic isotropic half-space is considered with Cartesian coordinate axes (x, z) defined with the z-axis normal to the surface of the half-space. In the simulations shown below the material is aluminium (Young's modulus: 69 GPa, Poisson's ratio: 0.334, density: 2700 kg m−3). A 5 MHz, 64 element array is used as an example throughout this paper, and the element pitch of the array is 0.6 mm (i.e. 0.5λ at the centre frequency).

Figure 1.

Array measurement geometry.

In general, the ultrasonic array data, g, can be measured in many different ways, for example, by applying various time-delays to the array elements or using different subsets of array elements in transmission and reception. However, the maximum possible amount of information that can be collected by the array in a particular position is contained in the full matrix of transmitter–receiver array data, measured using the FMC technique. The FMC data depend on three variables: position of transmitter array element, xT, position of receiver array element, xR and time, t, and can be visualized as a three-dimensional cube of data, gfmc(t, xT, xR) (figure 2a). One example of an alternative array data capture technique is the plane wave imaging method [30,33]. In this case, the measured array data are given by a dataset gpw(t, γT, xR), which depends on the angle of transmitted plane wave, γT, position of receiver array element, xR and time, t.

Figure 2.

Modelling results for five point scatterers: (a) FMC data, −40 dB isosurface image; (b) back-propagation two-dimensional image; (c) generalized back-propagation image showing −40 dB isosurfaces; (d) slice of the generalized image at xR = 0. Dashed line indicates the direction of side-lobes (angle θ) given by (2.7), angle φ0 is the array half angular aperture.

The main question, investigated in this paper, is if there is an optimal array data measurement technique, which allows the same amount of information to be collected as FMC, but using a smaller number of measurements. Mathematically, it means, that the data g has a smaller number of components, than the FMC data, gfmc, and there exists a one-to-one relationship between them: where T is the mapping operator from the new measurement space g into the FMC measurement space gfmc.

Expression (2.1) means that measurements gfmc and g are equivalent to each other. One important practical implication is that in this case all data processing methods available for the FMC data, gfmc, can also be applied to the data g. For example, scattering matrices for different defects can be extracted from the raw measured data and then used for defect characterization [15]. Another important observation from relationship (2.1) is, that if the number of measurements in g is smaller than the dimension of the FMC data gfmc, then different measurements in the FMC dataset are not independent. Therefore, the main goal of this paper is to find possible dependencies in the FMC data and then use them to design the optimal measurement procedure. It is important to note that these sparsity relationships will be based on fundamental properties of the FMC data and array and not rely on the assumption of a small number of scatterers or other qualitative comparisons.

Concept of imaging as a data compression method

Conventionally, the amount of information in the array datasets gfmc and g, is assessed by comparing the ultrasonic images, Ifmc and I, obtained using data gfmc and g, respectively. If Lfmc and L are imaging operators for the gfmc and g datasets, then

Traditionally, it was assumed, that if images have similar quality (for example, based on the signal-to-noise ratio or contrast metrics), then there is no loss of information. However, in general, the equality Ifmc = I does not necessarily imply that the data gfmc and g are equivalent in the sense of relationship (2.1), key if more advanced processing is to be carried out. A sufficient condition for the existence of a mapping operator T between two datasets is that the corresponding imaging operators are reversible. In this case, The construction of an inverse imaging operator requires the range of the imaging operator to be known, which in turn requires a specific model of the scattering mechanism. For example, if scatterers are considered as a collection of frequency-independent point scatterers, then an image represents a reconstruction of the scatterer distribution function. In this case, the image-to-data operator is represented by the forward model of wave scattering from a known distribution of point scatterers. Alternatively, the imaging operator Lfmc can be considered globally [12], as an operator in the vector space gfmc. The advantage of this approach is that the reversibility of the imaging operator becomes model-independent. For example, it was shown that the back-propagation imaging method is reversible in such a global sense [12].

The analysis in this paper is based on the assumption that multiple scattering between different scatterers is small compared with the single scattering contribution, so a single scattering model is valid. It is important to stress that no other assumptions about the physical nature of scattering are made.

The modelling example of FMC data for five-point scatterers, located at a depth of 20 mm from the array is shown in figure 2a. It can be seen that the scattering data are spread over all combinations of transmit–receive signals, and it is difficult to see any general dependencies between different components of the FMC dataset. However, the main characteristic of the single scattering data are that for every scatterer all scattered signals come from a single location, corresponding to the position of the scatterer. The majority of the array imaging methods are based on this property and effectively focus the scattered signals back at the scatterer location. If the imaging method is reversible, then it is possible to convert the image back to the array data. Therefore, in this case the image represents an alternative form of the array data. Moreover, because scattered signals from different scatterers are localized in the imaging space, the imaging algorithm can be considered as an array data compression operation.

Reversible back-propagation imaging method

In this section, the analysis of the reversible back-propagation imaging operator is performed which helps to reveal data redundancies in the FMC dataset.

In general, the transmitter and receiver elements are designed to be sensitive to the longitudinal wave mode only, so shear waves and mode conversion effects are not considered. The back-propagation imaging method can be represented by a linear operator, B, which converts FMC data, gfmc(t, xT, xR), into the generalized image, b(z, xT, xR) [12], Here F is a two-dimensional Fourier transform with respect to the spatial coordinates xT, xR, F−1 is the inverse Fourier transform and H is the back-propagation of the angular spectrum operator. Note that the physical meaning of the generalized image, b(z, xT, xR), is transmitter–receiver array data, measured at time t = 0 by an array located at depth z [12]. Alternatively, it can be considered as beamforming with different transmit, (xT, z), and receive, (xR, z), focusing.

The operator H converts the time data into a function of propagation distance and can be written in the form of a one-dimensional Fourier transform (see appendix Aa). Therefore, each operator in expression (2.4) for the back-propagation operator B is reversible and the generalized image b(z, xT, xR) can be converted back into the array data gfmc(t, xT, xR): where B−1 is the inverse imaging operator.

The conventional two-dimensional image of scatterer position I(x, z) is given by the pulse-echo data xT = xR of the generalized image Therefore, the generalized image b(z, xT, xR) contains more information than is necessary for localization of the scatterers. However, the extra information corresponding to the non-diagonal data xT≠xR is crucial for the inverse imaging [12].

Point spread function in the generalized image domain

Figure 2b shows the conventional two-dimensional back-propagation image, i.e. in the pulse echo xT = xR plane, for the modelling example of five-point scatterers. The corresponding generalized image is presented in figure 2c. It is seen that the back-propagation operator focuses the energy from each scatterer into the vicinity of its location, and, therefore, the data in the generalized image domain are localized around this pulse-echo plane xT = xR. Note that the generalized image can be considered as a different representation of the original FMC dataset shown in figure 2a, as they can be transformed between the two representations using the back-propagation and inverse back-propagation operators. Comparison of figure 2a,c reveals the fundamental sparse nature of the FMC dataset. It is clear that the degree of sparseness is defined by the point spread function (PSF) of the array in the generalized image domain. The structure of the PSF is investigated in this section.

For simplicity of the initial analysis, it is assumed that the array elements are omnidirectional and a point scatterer is located at depth z0 directly below the geometrical centre of the array (x0 = 0). An example of the PSF in the generalized image domain is shown in figure 2d. The amplitude of the side lobes can be estimated analytically using asymptotic analysis of the back-propagation imaging operator. The detailed derivation is given in appendix Ab, and only final results are summarized below.

The main side lobes are located in the (xT, z), (xR, z) planes in the directions θ = θ (figure 2d): where the angle θ is an elevation angle in the (xT − x0, xR − x0, z − z0) space, and θ = 0 corresponds to the direction of the positive z-axis. The angle φ0 corresponds to the maximum scattered angle measured by the array (half angular aperture of the array), where L is the array aperture.

The quantity of interest is the relative side lobe level as a function of distance r and array angular aperture φ0: where P(r, θ) is the point spread function and r is the radial distance in the (xT − x0, xR − x0, z − z0) space with the centre at the location of the scatterer.

It is convenient to introduce dimensionless distance , where λ is the ultrasonic wavelength at the centre frequency f0. Then the asymptotic behaviour of the PSF in the far field is given by

if

if

Qualitatively, expression (2.10) is applicable when the half array angular aperture φ0 is relatively small. This corresponds to the situations when the array size is small or the scatterer is located far from the array. On the other hand, when the array angular aperture φ0 increases, the array focusing capability becomes better. Consequently, the decay rate of the PSF also increases according to expression (2.11).

Figure 2d shows the xT − z plane (xR = 0) of the generalized image of the PSF for the scatterer located in the centre of the array. The dashed line corresponds to the theoretically predicted side-lobe direction. Note that the half angular aperture of the full array is 40°. However, large angles contribute mainly to the noise and, therefore, an additional angular filter of φ0 = 35° was applied in the wavenumber domain. In this case, according to expressions (2.7), (2.11) the side-lobes direction is θ = 72.5° and the side-lobes length at −40 dB in xT direction is 9 mm. It can be seen that these theoretical estimations are in a good agreement with figure 2d.

Optimal data acquisition strategies

General approach

In this section, the sparse nature of the FMC dataset in the generalized image domain is used to design the optimal data acquisition procedure. For the analysis, it is convenient to rewrite the back-propagation imaging operator (2.4) as an operator acting on the angular spectrum G(t, k, k) = F[gfmc(t, xT, xR)]: Note, that the Fourier transform operator F acts on the spatial coordinates xT, xR, and the back-propagation of angular spectrum operator H acts on the variable z only.

If the array data are collected using the FMC technique, then the sampling in the wavenumber domain (k, k) is defined by the array aperture L as In the previous section, it was shown that the generalized image b(z, xT, xR) has a specific sparsity, when all data are localized around the pulse-echo plane xT = xR. This means that the sampling rule (3.2) is too conservative and results in oversampling. The concept of the optimal data acquisition method is to design the optimal sampling in the wavenumber domain, so that it reflects the sparse property of the generalized image.

From the practical point of view, an important question is how to implement the sampling scheme in the wavenumber domain using array measurements. There are multiple possible sampling schemes which are consistent with the sparsity of the generalized image, however, not all of them allow practical implementation using conventional ultrasonic array instrumentation. Therefore, in this paper an additional optimization criterion of minimal number of array firings is used. In this case, it is convenient to represent the Fourier transform operator as F = FTFR, where FT and FR are Fourier transform operators with respect to transmitter and receiver coordinates. Then the back-propagation imaging operator can be written in terms of the transmit angular spectrum GT(t, k, xR) = FT[gfmc(t, xT, xR)] as From this expression and using conventional sampling arguments, it follows that the optimal sampling in the wavenumber domain k is where W is the maximum linear size of the PSF in the generalized image domain in the xT, xR directions. It should be noted that, strictly speaking, the size of the PSF in the generalized image domain (points where the PSF is not equal to zero) is unbounded. If the parameter W is finite, then this results in additional noise in the image due to aliasing. The level of this noise can be controlled by using some threshold, and then the corresponding value of W can be calculated from expressions (2.10) or (2.11), depending on the array size and scatterer location.

Angular spectrum imaging

One possible practical implementation of the optimal data acquisition approach described above is to transmit angular spectrum components, corresponding to the sampling rule (3.4), and then use all array elements on reception simultaneously similar to the FMC technique. The maximum measured wavenumber is given by . Then the number of array firings is given by , and using (3.4) can be written as Then the image is produced by applying the back-propagation operator in form (3.3). In this case, all information is preserved and the FMC data can be reconstructed from the generalized image using inverse imaging (2.5).

However, there are some difficulties associated with the angular spectrum excitation. Each angular spectrum component represents a one-dimensional wave propagating in the z-direction. For high wavenumbers k, this wave is highly dispersive and has very small amplitude compared with the plane wave excitation at k = 0, which makes it very sensitive to experimental noise. Therefore, a practical implementation of the angular spectrum approach is very challenging.

Plane wave imaging

The problems associated with the angular spectrum method can be overcome by using an alternative approach—plane wave imaging, which is described in this section. Figure 3a schematically shows the sampling in the wavenumber-frequency domain, k − ω, corresponding to the angular spectrum excitation. Each angular spectrum component is associated with the fixed wavenumber k, which is independent of frequency. It means that the phase velocity of the angular spectrum wave is frequency-dependent. In order to illuminate the dispersion effect the wavenumber k in the sampling scheme can be taken to be proportional to frequency This corresponds to non-uniform sampling and is schematically shown in figure 3b. Each component (3.6) represents the plane wave propagating in the direction φ relative to the array. Note, that practically each plane wave can be generated by applying appropriate time-delays to the array elements. The reception is performed simultaneously on all array elements as before.

Figure 3.

Sampling schemes in wavenumber-frequency, k − ω, domain, corresponding to (a) angular spectrum excitation; (b) plane wave excitation.

Because the sampling in the wavenumber domain is non-uniform, the back-propagation image cannot be calculated using the Fast Fourier Transform algorithm. Therefore, a numerical evaluation of Fourier integrals on a non-uniform grid must be performed instead. However, a computationally more efficient approach is to use the asymptotic form of the back-propagation method [35]. In this case, the imaging algorithm has a form similar to the total focusing method, and, hence, very efficient numerical implementation is possible. The detailed derivation is given in appendix Ac.

The angular interval, Δφ, of the plane wave excitation can be approximately calculated from expressions (3.4) and (3.6) as

For a scatterer located in the centre of the array (x0 = 0), the PSF is symmetrical in the xT, xR plane. In this case, W can be taken as the size of the PSF in the xT direction, and can be written as W = 2rsinθ = 2rcos0.5φ0, where r is the side-lobe size in the generalized image domain (figure 2d). Here relationship (2.7) between size-lobe direction θ and the angular array aperture φ0 was used. Then the angular sampling interval can be expressed as

The number of plane waves is N = 2φ0/Δφ. Finally, using expressions (2.10), (2.11) for asymptotic behaviour of the PSF, it is possible to obtain the following expression for the optimal number of firings: and Here the variable δ describes the relative noise level due to aliasing.

Expression (3.9) allows us to draw some important practical conclusions. Note, that if the location of a scatterer is fixed, then the angle φ0 monotonically increases when the array aperture increases. Then from (3.9) it follows that for a relatively small angular range φ0 the required number of firing is constant and independent of φ0, and, consequently, independent of the array aperture. Moreover, when the number of array elements increases (hence, angular range φ0 increases), the optimal number of firing decreases as (φ0tan0.5φ0)−1/3. This behaviour can be explained by the fact that increased angular range results in better focusing, and, hence, smaller size of the PSF.

It should be stressed that the number of transmitted plane waves given by expression (3.9) preserves the PSF in the generalized image domain (to the accuracy given by the relative noise level δ). It means that in this case the plane wave array data is equivalent to the FMC data and can be converted to the FMC data using inverse imaging. On the other hand, the conventional two-dimensional image corresponds to the pulse-echo xT = xR plane in the generalized image domain. Therefore, if it is required to reconstruct the two-dimensional image only, then the number of plane wave excitations, N, is half that given by expression (3.9): In this case, the two-dimensional image is preserved (up to the relative noise level δ), but the part of the generalized image outside of the pulse-echo plane is corrupted by aliasing noise.

It is noted that the results (3.9), (3.11) for the optimal number of plane wave transmissions do not explicitly depend on the position of the scatterer relative to the array and are entirely defined by the angular aperture only. However, it is assumed that the angular aperture is symmetrical, [ − φ0, φ0]. In practice, it corresponds to the imaging area directly below the array centre. If a point scatterer is located at the position (x0, z0), x0≠0, that is offset from the array centre, then the angular range, [φ1, φ2], is asymmetrical: , . However, the important result obtained above is that the maximum angular step, Δφ, of the plane wave transmissions is inversely proportional to the size of the point spread function, W, in the generalized image (equation (3.7)). This condition is true for any location of a scatterer, and, therefore, can be used to estimate the minimum number of plane waves. Qualitatively, if a scatterer is located with an offset relative to the array centre, then the focusing performance of the array worsens and the relative side lobe amplitude of the PSF increases. Therefore, for a given relative side lobe level, δ, the size of the PSF increases. Consequently, the plane wave angular step becomes smaller.

The structure of the PSF in the generalized image domain for an offset scatterer can be investigated similar to the centred scatterer (see appendix Ab). However, the analysis in this case becomes more complex, and, in general, the result cannot be represented in a compact form as in (3.9). Nevertheless, in a particular case of a scatterer located close to the edge of the array, x0 ≈ L/2 (L is the array aperture), it is possible to estimate the minimum number of plane waves as (see (3.11)) and where N1, N2 are given by (3.10). Note that in this expression it is assumed that the angular range is [0, φ0], and corresponds to the imaging area around the scatterer location at x0 = L/2. If it is required to perform imaging in the entire area below the array, −0.5L≤x0≤0.5L, then the angular range is [ − φ0, φ0], and the minimum number of plane waves is given by Therefore, in this case the optimal number of array firings is between 22/3 ≈ 1.6 and 2 times larger than for a centred scatterer.

Experiments

Scatterer in the centre of the array

The theory developed in the previous sections was experimentally validated on an aluminium specimen with a 1 mm diameter side drilled hole located at 20 mm depth. A linear array with the same parameters (5 MHz, 64 elements, 0.6 mm element pitch) as in the modelling example was used. Firstly, the measurements were performed by positioning the array exactly above the defect. The case of a scatterer located with the offset relative to the array centre is considered later. Figure 4 shows the back-propagation imaging results obtained using the FMC dataset. Note that similar to §2d the half angle filter of 35° was applied in the wavenumber domain. Figure 4a,b shows the two-dimensional back-propagation image and the generalized image, respectively. Figure 2d shows the xT − z plane of the generalized image at xR = 0. The dashed line corresponds to the theoretically predicted side-lobes direction. It can be seen that in general the experimental PSF structure agrees well with the theoretical estimations (see also figure 2d).

Figure 4.

Experimental images for 1 mm side drilled hole obtained from the FMC data: (a) back-propagation two-dimensional image; (b) generalized back-propagation image, −40 dB isosurfaces; (c) slice of the generalized image at xR = 0. Dashed line indicates the direction of side-lobes (angle θ) given by (2.7), angle φ0 is the array half angular aperture.

The plane wave imaging results are presented in figure 5 for 9, 17 and 25 transmitted plane waves. The maximum steering angle is 35° and the angular filter of 35° is also applied on reception in the wavenumber domain. Note that according to expression (3.11) only nine plane waves are required for the maximum noise level of −40 dB in the two-dimensional image. It can be seen that in all cases the two-dimensional images are practically identical and the noise level is below −40 dB (see also figure 4a for the two-dimensional image obtained from the FMC data).

Figure 5.

Experimental images for 1 mm side drilled hole obtained with 9 (a,b), 17 (c,d) and 25 (e,f) emitted plane waves and aperture half angle limit of 35°: (a,c,e) back-propagation two-dimensional image; (b,d,f) slice of the generalized image at xR = 0.

Figure 5b,d,f illustrates how the aliasing noise in the generalized image domain changes with the number of emitted plane waves. For example, if only nine plane waves are used then the generalized image is completely corrupted with noise in the region outside the pulse-echo plane |xT − xR|≥2 mm. However, if the number of emissions is increased to 17 plane waves, then the aliasing noise in the generalized image domain does not affect the data as predicted by expression (3.9).

Once the generalized image is constructed, it can be converted to the FMC data using the inverse imaging process. Preliminarily, the part of the generalized image free from the aliasing noise was filtered out by applying the spatial filter |x|≤Δx, z1≤z≤z2. Based on figure 5, it can be determined that z1 = 18 mm, z2 = 23 mm and Δx = 2, 8 and 10 mm for 9, 17 and 25 plane wave emissions, respectively. An example of the reconstructed time-trace corresponding to transmitter element 40 and receiver element 20 is shown in figure 6, where for comparison the same time-trace from the experimentally captured FMC dataset is also shown. It can be seen that in all three cases the reconstruction is very good, even in the case of only nine plane wave transmissions. This surprising result is explained by noting that the majority of the PSF data in the generalized image domain, which is needed for the inverse imaging, is concentrated in a relatively small area around the defect location. In order to quantitatively illustrate this fact the relative side lobe amplitude of the PSF as a function of array half angular aperture, φ0, and the distance in xT,R direction is shown in figure 7. The calculations were performed using expression Here δ, δ are given by formulae (2.10) and (2.11) and the distance in the side lobe direction was replaced by , where is the distance in xT,R direction. For example, it is seen that for φ0 = 35° the PSF data with a relative amplitude above −20 dB is contained in the region |x|≤1.4λ, where the ultrasonic wavelength λ = 1.2 mm at 5 MHz.

Figure 6.

Time trace corresponding to transmitter element 40 and receiver element 20 (red lines) and the same time trace reconstructed using the inverse imaging method from the plane wave dataset with (a) 9, (b) 17 and (c) 25 emissions (blue lines).

Figure 7.

Relative side lobe amplitude of the PSF in the generalized image domain as a function of distance r in the xT,R direction and angular aperture. (Online version in colour.)

The possibility to reconstruct the FMC dataset means that all signal processing operations available for the FMC data can also be applied to the plane wave dataset. For example, the scattering matrices extracted from the FMC data and the plane wave data with 9, 17 and 25 emissions are shown in figure 8. For the plane wave data, the FMC data were reconstructed first using the inverse imaging as described above and then the scattering matrix was extracted from the FMC data. It is seen that the results are very similar confirming that the number of plane waves given by expression (3.11) is enough to preserve the majority of FMC information. Note that the scattering matrices extracted from the FMC data and nine plane wave data (figure 8a,b) differ mainly at the maximum angles ±φ0 = ±35°. This is consistent with the asymptotic analysis performed in appendix Ab, where it has been shown that the main side lobe contribution is given by the maximum steering angles and the stationary point (0°). However, this information is partly removed from the plane wave data during the preliminary filtering in the generalized image domain before the application of the inverse imaging operator.

Figure 8.

Experimental scattering matrices for 1 mm diameter hole extracted from the FMC data (a) and the plane wave data with (b) 9, (c) 17 and (d) 25 emissions and aperture half angle limit of 35°. Scattering matrices are normalized to its maximum values. (Online version in colour.)

Finally, the signal-to-noise ratio (SNR) in the two-dimensional image as a function of number of plane waves and angular aperture was calculated. The noise amplitude was estimated as the maximum image amplitude in the area outside the rectangular box containing the PSF. The half x-size of the box was taken equal to , and the size in the z-direction was 6 mm. The result is presented in figure 9a. Also figure 9b shows the SNR based on the theoretical formulae (3.9) and (3.11). It can be seen that there is a reasonable agreement, however, for angular apertures smaller than 20° the experimental images show slightly higher SNR compared with the theoretical predictions. This can be explained by the fact that asymptotic expression (3.9) is based on the assumption of regular sampling in the wavenumber domain. In this case, the PSF pattern is periodically repeated in the xT direction. Plane wave imaging however corresponds to irregular sampling, which is also frequency-dependent (figure 3). Consequently, the grating lobe structure in the generalized image domain is different (figure 5). It is important to note that this modelling approach therefore represents a good conservative approach.

Figure 9.

Signal-to-noise ratio in the two-dimensional image as a function of the number of plane waves and angular aperture. (a) Experimental results for 1 mm diameter hole; (b) Theoretical predictions (expressions (3.9) and (3.11)).

Scatterer offset from the array centre

In this section, a scatterer located with an offset to the array centre is considered. In order to illustrate the plane wave imaging performance in this case a modelling example is presented first. The linear array aperture for the array parameters given in §2a is L = 37.8 mm, the point scatterer is located at a depth of z0 = 20 mm and its lateral position x0 varies from 0 (array centre) to x0 = L/2 (array edge). The maximum plane wave steering angle is φ0 = 35° as in the case of the centred scatterer. Therefore, the angular aperture of the scatterer is [ − φ, φ0], 0≤φ≤φ0, and φ = 0° corresponds to the scatterer located at the array edge and φ = φ0 corresponds to the scatterer at the array centre.

Initially, the FMC array data were modelled using the far-field ray tracing model [11,12]. Then the FMC data were converted into the plane wave array data by applying appropriate time-delays to the array elements on transmission. Figure 10a shows the signal to aliasing noise ratio in the conventional two-dimensional image as a function of the number of emitted plane waves and the scatterer location. Similar to the previous section, the noise amplitude was estimated as the maximum image amplitude in the area outside the rectangular box containing the peak of the PSF. The size of the box was chosen to satisfy the condition that the amplitude of the PSF in the back-propagation image (FMC dataset) is less than −40 dB in the area away from the selected region.

Figure 10.

(a) Modelling results of the signal-to-noise ratio in the two-dimensional image as a function of the number of plane waves and lateral position of the scatterer. The scatterer is located at the depth of z0 = 20 below the array, the maximum plane wave steering angle is φ0 = 35°. Angular aperture of the scatterer is [ − φ, φ0], 0≤φ≤φ0. The number of plane waves corresponds to the full angular aperture [ − φ0, φ0]. (b) The case φ = φ0 for the scatterer at the array centre; (c) the case φ = 0° for the scatterer located at the array edge. Dashed lines in (b) and (c) correspond to theoretical prediction (3.14).

It is important to note, that in figure 10a the number of plane waves for a given angle φ is related to the full angular aperture [ − 35°, 35°], and, consequently, corresponds to the case when the scatterers are located in the whole area −x0≤x≤x0. It can be seen that for a given number of plane waves the relative level of aliasing noise gradually increases as the position of the scatterer approaches the edge of the array, φ → 0. Additionally, figure 10b,c separately shows the SNR for the centred scatterer and the scatterer located at the array edge together with the theoretical predictions (3.11) and (3.14). It can be seen that the agreement is very good. In particular, the number of plane waves required to image the entire area below the array, −L/2≤x0≤L/2, with the SNR greater than 40 dB is 2−1/3N2 ≈ 15 as predicted by the lower bound of the theoretical condition (3.14).

So far only a single scatterer case has been considered. However, because imaging operations are linear, all results remain valid even when more than one scatterer is present, providing that the single scattering contribution is dominant. This is supported by the experimental measurements shown in figure 11. In this case, an aluminium specimen was used with a row of 1 mm diameter side drilled holes with a pitch of 10 mm at a depth of 20 mm. Figure 11a shows the back-propagation imaging results obtained using the FMC dataset. The half angle filter of 35° was applied in the wavenumber domain. The plane wave images for 9 and 15 transmitted plane waves are presented in figure 11b,c. The steering angle is [ − 35°, 35°] and the angular filter of 35° was also applied on reception. The aliasing noise is visible in the case of nine plane waves predominantly between and above the scatterers. However, the noise amplitude is at approximately −36 dB, which is half that predicted by the theoretical estimation of −30 dB (figure 10c). This is explained by the fact that the highest aliasing noise is related to the side scatterers, but their absolute amplitude is smaller compared with the scatterer at the array centre due to the reduced angular range. In the case of 15 plane waves, it can be seen that the image is practically identical to the back-propagation image and the noise level is below −40 dB as predicted theoretically.

Figure 11.

Experimental images for 1 mm side drilled holes with a pitch of 10 mm at a depth of 20 mm: (a) back-propagation image obtained from the FMC array data; (b) plane wave image with nine transmitted plane waves; (c) plane wave image with 15 transmitted plane waves. In all cases, the angular aperture on transmission and reception is [ − 35°, 35°].

Other factors affecting plane wave imaging

One of the main results of this paper is the fact that the optimum number of plane waves is determined by the array PSF in the generalized image domain. Therefore, the structure of the PSF can be used to assess the influence of different parameters on the plane wave imaging performance. As an example, two factors which are important from the practical point of view are briefly considered in this section: array element pitch and multiple scattering. Note, that only a qualitative illustration is given below, and the detailed quantitative analysis of these cases will be the subject of future work.

The developed theory implicitly assumes that the array element pitch satisfies to the Nyquist sampling criterion (λ/2). However in many fields, including medical ultrasonic imaging it is common to use arrays with λ-pitch. Figure 12 shows the modelling example of a conventional two-dimensional image and the generalized image of the PSF for the same 5 MHz 64 element array as in §2, but with λ element pitch. Additionally, an angular filter with the half angular aperture of 20° was applied on transmission and reception [36]. It can be seen that grating lobes are successfully suppressed in the two-dimensional image. However, there are significant grating lobes in the generalized image domain outside the pulse-echo plane. These grating lobes do not appear on the two-dimensional image, but effectively result in an increased size of the PSF, and, therefore, affect the minimum number of array firings required to preserve all information. This provides the quantitative framework to explain why images can look fine with a reduced number of firings and specifies the number required to achieve this, something previously not reported.

Figure 12.

Point spread function for the 5 MHz, 64 element array with 1λ element pitch and angular filter of φ0 = 20°: (a) back-propagation two-dimensional image; (b) generalized back-propagation image, −40 dB isosurfaces; (c) slice of the generalized image at xR = 0. Dashed line indicates the direction of side-lobes (angle θ) given by (2.7).

Another important aspect is related to the single scattering assumption. Strictly speaking, a multiple scattering contribution always exists, although in many cases it is small compared with single scattering [37]. In the generalized image, domain multiple scattering results in additional noise distributed outside the conventional two-dimensional image plane xT = xR. This ‘multiple scattering’ part of the generalized image will result in aliasing noise in the two-dimensional image during the plane wave imaging with reduced number of firings. Essentially, this will degrade the performance of the two-dimensional image and the magnitude and nature of this effect must be taken into account when designing the reduced plane wave firings.

Conclusion

The problem of optimization of ultrasonic array data acquisition has been investigated. The conventional FMC technique involves collecting all possible transmitter–receiver signals and, therefore, contains the maximum information which can be measured by an array. However, the volume of measured data is very large and the data acquisition is not fast enough for some applications.

It has been shown that under the single scattering assumption the FMC dataset has a specific sparse structure. This property can be used to design an optimal data acquisition method. Two possible strategies, corresponding to different sampling schemes in the wavenumber-frequency domain, have been considered: the angular spectrum and plane wave imaging techniques. The analytical relationship between the number of emissions, angular array aperture and signal-to-noise ratio has been derived. An important counterintuitive conclusion from this relationship, is that the optimal number of emissions decreases when the angular aperture of the array increases.

The analysis shows that plane wave imaging is more suitable for practical implementation than the angular spectrum method. The theory has been validated experimentally using 5 MHz, 64 element array measurements on an aluminium sample with a 1 mm diameter holes. In particular, it has been demonstrated that the theoretical prediction of the optimal number of plane waves is in a good agreement with the experiment. Also, it has been shown that the FMC dataset can be reconstructed from the plane wave data with a reduced number of emissions. This result means that all signal processing operations available for the FMC data can also be applied to the plane wave dataset. For example, the scattering matrices of the defects can be extracted and used for defect characterization.

The possible applications of the developed theory are twofold. Firstly, the precise knowledge of the PSF structure gives an opportunity to construct an efficient spatial filter in the generalized image domain for the accurate extraction of defect scattering matrices. Secondly, the optimal number of array firings for various practical configurations can be predicted and used for fast imaging without loss of scattering information, or with detailed knowledge of what information has been discarded.