A fractional Fourier transform analysis of the scattering of ultrasonic waves

Many safety critical structures, such as those found in nuclear plants, oil pipelines and in the aerospace industry, rely on key components that are constructed from heterogeneous materials. Ultrasonic non-destructive testing (NDT) uses high-frequency mechanical waves to inspect these parts, ensuring they operate reliably without compromising their integrity. It is possible to employ mathematical models to develop a deeper understanding of the acquired ultrasonic data and enhance defect imaging algorithms. In this paper, a model for the scattering of ultrasonic waves by a crack is derived in the time–frequency domain. The fractional Fourier transform (FrFT) is applied to an inhomogeneous wave equation where the forcing function is prescribed as a linear chirp, modulated by a Gaussian envelope. The homogeneous solution is found via the Born approximation which encapsulates information regarding the flaw geometry. The inhomogeneous solution is obtained via the inverse Fourier transform of a Gaussian-windowed linear chirp excitation. It is observed that, although the scattering profile of the flaw does not change, it is amplified. Thus, the theory demonstrates the enhanced signal-to-noise ratio permitted by the use of coded excitation, as well as establishing a time–frequency domain framework to assist in flaw identification and classification.

Many safety critical structures, such as those found in nuclear plants, oil pipelines and in the aerospace industry, rely on key components that are constructed from heterogeneous materials. Ultrasonic nondestructive testing (NDT) uses high-frequency mechanical waves to inspect these parts, ensuring they operate reliably without compromising their integrity. It is possible to employ mathematical models to develop a deeper understanding of the acquired ultrasonic data and enhance defect imaging algorithms. In this paper, a model for the scattering of ultrasonic waves by a crack is derived in the timefrequency domain. The fractional Fourier transform (FrFT) is applied to an inhomogeneous wave equation where the forcing function is prescribed as a linear chirp, modulated by a Gaussian envelope. The homogeneous solution is found via the Born approximation which encapsulates information regarding the flaw geometry. The inhomogeneous solution is obtained via the inverse Fourier transform of a Gaussian-windowed linear chirp excitation. It is observed that, although the scattering profile of the flaw does not change, it is amplified. Thus, the theory demonstrates the enhanced signal-to-noise ratio permitted by the use of coded excitation, as well as establishing a time-frequency domain framework to assist in flaw identification and classification.

(a) Chirp excitations
Coded excitations are an effective way of delivering large amounts of energy using relatively low acoustic pressure amplitudes. Inspiration for coded signal design can be drawn from bioacoustics; bats and dolphins use frequency-modulated sweeps to navigate and hunt [7,8]. The use of coded excitations in signal processing has been shown to improve SNR and lessen trade-offs between sample penetration and image resolution [9,10]. Chirps contain multiple frequencies which vary in time, typically in a linear or exponential manner. Their broad frequency content increases the likelihood of reaching the resonant frequency of a defect, in turn causing stronger vibrations and consequently improving the probability of detection. For the purposes of this work, a Gaussianmodulated linear chirp in time (t) of the form has been used, where m is the gradient of the chirp (the rate at which it sweeps through a prescribed range of frequencies), f 1 is the initial frequency, t 1 is the centre of the Gaussian envelope and σ is its standard deviation. By varying f 1 , m and σ , the bandwidth of the chirp can be altered. Setting m = 0 Hz, the chirp reverts back to a time harmonic signal with frequency f 1 . For a fair comparison of the chirp with a continuous gated waveform (the typical signal emitted by the transducer), it is imperative both are optimized for the same transducer and hence use its full bandwidth. Figure 1 shows one such matched pair. In the time domain, the gated continuous sine wave (in black) spans a far shorter time interval than that of the Gaussian-modulated linear chirp (both signals have the same peak amplitude). Studying the plots of their respective Fourier transforms in figure 1b, it is apparent that they have similar −6 dB bandwidths (approx. 50%) but observe that the chirp contains a far greater amount of energy.

(b) The fractional Fourier transform
An alternative approach to differentiating between flaw scattering and noise can be taken via the analysis of the collected data in the time-frequency domain. In an ideal case, a received ultrasound signal would exhibit signs of scattering at the time interval pertaining to the location of the flaw, facilitating detection and subsequent characterization. However, in practice, scattering by the microstructure of the host media can dominate the signal. To improve identification of defects it is suggested that the scattered signals are analysed for their frequency content; the frequency spectrum of the wave scattered by a flaw should be different from that scattered by a heterogeneity. However, time information must be retained in order to locate the flaw.  achieved via (i) the time windowed Fourier transform [11], where the discrete Fourier transform is applied to short time intervals allowing the frequency content at that specific time to be analysed independently of the rest of the signal, or (ii) the fractional Fourier transform (FrFT) [12][13][14], which enables continuous movement between the time and frequency domains, allowing the simultaneous retention of both frequency and time domain information. In this work, the FrFT is employed to analyse the benefits of chirp excitation over gated continuous wave excitation. As a generalization of the ordinary Fourier transform, the FrFT is more flexible in its applications and hence of potential interest to any area in which the Fourier transform is frequently implemented. Its main advantage is that it allows continuous movement between the time and frequency domains, retaining information from each, thus presenting an alternative to using the time windowed Fourier transform. There exist several conventions for defining the FrFT, each of which gives a slightly different physical interpretation. For the purposes of this work, the FrFT of order a is given as the linear integral transform [12]  and α = aπ/2. When a is an integer, it denotes the number of repeated applications of the ordinary Fourier transform. Hence, setting a = 1 (and, consequently, α = π/2), equation (1.2) simplifies to the ordinary Fourier transform.

Solving the inhomogeneous wave equation in time-frequency space
To build a mathematical framework to allow for analysis of the scattering of an ultrasonic chirp by a flaw, the wave equation with a time-dependent forcing function, q(t), will be solved in the time-frequency space. Note that q(t) is spatially independent and hence acts like a body force, affecting the whole flaw domain simultaneously. To justify this, it is assumed that the length of the flaw is smaller than the wavelength. This is in keeping with the low-frequency assumption made in the Born approximation [15] which is used later. A second assumption inherent to the Born approximation is that the transmission and reception of waves takes place at a distance from the flaw which is much larger than the wavelength: the far-field assumption. To begin, consider the non-homogeneous wave equation and c is the wave speed. We demand that the solution is bounded in time and space and satisfies the Sommerfeld radiation condition guaranteeing that the waves are outgoing and decay sufficiently fast so there exist no sources at infinity. It is also assumed that an initial pressure amplitude of h 0 is present at the ultrasonic array. This simplified case will of course cover the detection of an object in a fluid host but is also relevant to certain restricted classes within elastodynamics such as the propagation of horizontal shear waves in an isotropic solid. As mentioned above, the applications of interest will involve heterogeneous materials (with spatially dependent wave speeds) but this is a reasonable model with which to start this investigation and the methodology will extend naturally to these more general cases. By [12], the FrFT, taken with respect to time, of a derivative of a function is given by Hence, taking the FrFT of every term in equation (2.1) gives the non-homogeneous wave equation in time-frequency space.

(a) The homogeneous solution
To solve equation (2.4), we start by finding the solution to the homogeneous differential equation via separation of variables. The solution is written in product form .
for some b ∈ R. This can then be separated into two equations and from which the temporal and spatial components of the homogeneous solution can be derived.
To solve equation (2.7), a chirp-like ansatz of the form is chosen. Substituting this into equation (2.7) gives Equating coefficients of powers of u, one can calculate γ and β as γ = −π i tan α and β = ± −b 2 sec 2 α + 2π i tan α, (2.12) where the square root is taken so as the real part is positive. Hence However, to ensure that the solution is bounded in u, d 2 is set equal to zero. Turning our attention now to the spatially dependent component of the homogeneous solution, equation (2.8) can be recognized as the homogeneous Helmholtz equation wherek = b/c (c is the plane wave speed and b is analogous to the circular frequency). An explicit approximation of the scattered wave at a specified b, h(x, b), is derived via the Born approximation [15], giving given by 2.17) and the unit vectors u 1 , u 2 and u 3 lie along the axes of the flaw. Thus, a solution to the homogeneous equation can be written as However, to ensure that this integral is finite, d 1 must be chosen as some function of b with which to bound the solution. For this work, the function has been chosen, where A, b 1 and σ are chosen to mimic the Gaussian envelope observed in the plot of the power spectrum of the forcing function ( figure 1b). Thus, the solution is effectively summed over the range of frequencies received by the transducer array.

(b) The inhomogeneous solution
To find a time-dependent inhomogeneous solution f p a (u) to equation (2.1), an ansatz of the form This form is chosen since, using the identities tan α = − cot(α + π/2) and sec α = csc(α + π/2), one can rewrite s a (u, b) as and Substituting these into equation (2.4) shows that the left-hand side can be written as In the work below, q 0 (u) will be written as q(t) and q 1 (u) as q(b), to allow for easier physical interpretation. By letting D(b) = (−4b 2 π 2 − iπ sin 2α)v a (b) (and remembering s a (u, b) = K a+1 (u, b)) equation (2.27) can be written as which holds for all orders a. Hence, by the index additivity rules it is shown that By taking the FrFT of order a = −1 (analogous with the inverse Fourier transform) of both sides, it follows that  where E 0 = exp(−t 2 1 /σ 2 ), p = 1/σ 2 , r = 2π f 1 m and w 0 = 2π f 1 + 2t 1 i/σ 2 . Applying the inverse Fourier transform gives (2.32) The standard integral for a general quadratic exponent is calculated using assuming that the real part of a is positive. It follows that Note that the spatial component of the solution, which gives rise to the scattering profile of the defect and encapsulates the geometry of the flaw, is amplified by the chirp excitation and this suggests that an increased SNR will result.  be assumed that the integral is zero outside this interval (the received signal will not contain any frequencies outside of the bandwidth). This approach also resolves any issues around the singularity in the integrand when b = 0 and α = 0. The inhomogeneous solution is plotted as a time-order plot in figure 2 for the cases where the forcing function is set as (i) a gated continuous wave and (ii) a Gaussian-modulated linear chirp. Each point within the plot is the additional amplitude at that point in time-frequency space obtained by using the particular excitation. It can be seen that, since the linear chirp contains more energy, it provides a marked increase in the scattering amplitude over the entire time-frequency space. In the case of the gated continuous wave in figure 2a,  (a) Choosing the optimal order a It is shown in [16] that the FrFT of a Gaussian function has the form of a Gaussian for all orders a. The standard deviation of these Gaussian functions, σ a , varies in time-order space, with the narrowing of the function being synonymous with an increase in maximum amplitude. Hence, it can be concluded that the optimal value of a at which to employ the FrFT is the value at which the minimum σ a occurs. As the chirp rate m increases, this maximum peak (which occurs at min σ a ) moves further away from the frequency domain. This is explained schematically in figure 4. The long-dashed horizontal line in figure 4a represents a continuous wave (where m = 0) and hence results in a single value on the frequency axis. As the gradient increases (see the dotted and solid lines), the breadth of the frequency spectrum increases. The curves in figure 4b demonstrate how the width of the Gaussian changes in the fractional Fourier domain (of course, the Gaussian is infinite but here the width is approximated by 6σ a as 99.73% of the signal lies within this interval). The fractional order which exhibits the widest Gaussian is orthogonal to the order with the narrowest distribution (i.e. there is a difference of 1 between the orders at which these extremes occur). It then follows that, as the bandwidth of the chirp decreases, the order at which the narrowest Gaussian form arises approaches the frequency domain. Now, plotting the rate of frequency change (with respect to time) of the linear chirp results in the plot as seen in figure 5. The angle made with the frequency axis can thus be calculated as α = tan −1 ( 1 2 f 1 m). This provides the optimum angle at which to take the FrFT [17] and translates to order a = 2 tan −1 ( 1 2 f 1 m)/π . To assess the formula's success in predicting this order in regards to the inhomogeneous solution derived above, f p a (u) (as defined in equation (2.36)) is plotted over orders −2 ≤ a ≤ 2 in figure 6. As the Gaussian function is not centred at zero, the plot appears skewed; however, the narrowing phenomena can still be observed. Owing to the low gradient of the chirp excitation, the location of the narrowest distribution approaches the frequency domain as predicted. The simple algebraic formula derived for the optimal order of the FrFT of a Gaussian-windowed linear chirp exhibits an error in application to the inhomogeneous solution and does not incorporate the maximum amplitude (which is circled in black). However, the error is small (within 0.1 of the order at which the maximum does occur) and the formula could potentially guide the implementation of the discrete FrFT [18], effectively reducing the neighbourhood (and, subsequently, the computational expense) over which the FrFT is taken. It is hoped that the general solution as derived above (equation (2.37)) will act as a basis for further work on improving the extraction of the optimal order a at which to implement the FrFT for signals which have encountered a defect and been subsequently scattered, thus eventually reducing the time-order space to one dimension for numerical implementations.

Conclusion
A wealth of information regarding the detection, imaging and sizing of flaws is contained within the scattering matrices as constructed by data arising from ultrasonic phased array inspections.  by choosing an ansatz that, once substituted into the inhomogeneous wave equation, resulted in a linear integral equation which could be solved by formulating a Fourier transform pair for a Gaussian-modulated linear chirp. Since the excitation parameters are exclusively contained in the additive term provided by the inhomogeneous solution, the comparison between gated continuous wave excitation and chirp excitation can be drawn by focusing solely on that term. It was plotted in a time-order plot for the cases where the forcing function was set as (i) a gated continuous wave and (ii) a Gaussian-modulated linear chirp and it was shown that, since the linear chirp contained more energy, there was a marked increase in the scattering amplitude. This was reinforced by plotting and comparing the corresponding scattering matrices for a chosen peak in the time-order plot, which further demonstrated the increased amplification provided by the chirp. Thus it is anticipated that an improved SNR will result when applied to experimental data. Recent studies on the use of the FrFT with regard to chirp signals [16,17] suggest a potential methodology for studying the wave scattering problem considered in this paper. Indeed, the findings from that work were assessed for their applicability in this paper and, whilst the discrepancies were clear, they did provide a reasonable estimate of the optimal order at which to implement the FrFT. It is envisaged that the analytical formulation of the general solution in §2d will allow for improved extraction of the optimal a for the more complicated case of a chirp which has been scattered by a defect. The benefit in doing so would be the effective reduction of the time-frequency space to one dimension, thus enabling reduced computational cost in NDT applications.