Interpolated Plane Wave Imaging (IPWI): Overcoming Spatial Aliasing

This paper introduces an approach to addressing the challenges posed by sub-Nyquist spatial sampling in ultrasound image generation, particularly in the context of Plane Wave Imaging (PWI). Sub-Nyquist sampling can lead to artifacts that compromise image accuracy. To overcome this, we present an interpolation technique for PWI, which resolves the limitations imposed by the sampling theorem. This technique not only improves image quality but also enhances acquisition speed, offering a valuable solution for real-time inspections.


Introduction
The quest for high-resolution ultrasound imaging at inspection speeds comparable to conventional methods drives the need for innovative solutions. While Full Matrix Capture (FMC) combined with Total Focusing Method (TFM) delivers superior image quality, its realtime performance falls short in demanding applications. Plane Wave Imaging (PWI) offers a promising alternative, enabling faster imaging while maximizing energy transmission. However, spatial sampling constraints and aliasing issues can hinder its performance.
The FMC-TFM technique, integrating data acquisition and processing stages, achieves optimal image quality by utilizing all potential combinations of transmitter-receiver pairs, thereby maximizing information content [1]. When employing N elements, every possible firing combination results in the collection of N^2 A-scans. FMC speeds depend on three critical parameters: the number of firing elements, the pulse repetition frequency (PRF), and its upper limit determined by the total wave travel time through the material. This fundamental physical threshold restricts firing speed. These conditions collectively contribute to FMC's relatively restrained performance in practical scenarios, often yielding frame rates in the order of tens of frames per second (fps), equivalent to scanning speed in mm/s with a 1 mm inspection step [2].
Unlike the Full Matrix Capture (FMC) method, Plane Wave Imaging (PWI) employs a distinct approach where only a single firing cycle is required to simultaneously activate all transducer elements. This novel strategy offers a dual benefit: frame rates are determined solely by wave travel time, and there is a significant increase in energy transmission into the material. In a study by Le Jeune et al. (2016) using CIVA software [3], it was observed that at an imaging depth of approximately 56 mm, the field amplitude of FMC is approximately 24 dB weaker compared to PWI. This difference arises because the field of FMC diminishes proportionally to the reciprocal square root of depth due to cylindrical spreading, whereas the uniform plane wave field of PWI remains consistent [4].
The self-healing attributes of the plane wave, as described by Jon Claerbout [5], underlie this uniformity and constancy. This concept underscores the plane wave's exceptional ability to maintain its coherent structure even when faced with obstacles or perturbations during propagation. Such inherent consistency contributes to the unique advantages of Plane Wave More info about this article: https://www.ndt.net/?id=28691 e-Journal of Nondestructive Testing -ISSN 1435-4934 -www.ndt.net Imaging (PWI), ensuring that the wavefront remains intact and preserves its properties, ultimately enhancing the quality and reliability of the imaging process. The Plane Wave model depicted in Figure 1, simulated using CIVA [3], demonstrates a remarkable phenomenon: the wavefront regenerates itself following its passage through three side-drilled holes. While the combination of PWI with TFM yields improvements in both speed and Signal-to-Noise Ratio (SNR), the resulting image quality does not attain the heights achieved by FMC-TFM. This discrepancy can be attributed to the diminished dataset employed in PWI. Furthermore, the process of generating plane waves within PWI hinges on Huygens' principle, necessitating the presence of closely spaced array elements to effectively create plane wavefronts [6].
Consequently, the imperative for closely positioned array elements, defined by the pitch, serves the purpose of ensuring constructive interference of the plane wave and fulfilling the spatial Nyquist sampling criterion [7]: where ℎ is the maximum pitch, is the slowest velocity, is the maximum frequency component and represents the transducer angle with respect to the normal. For a zero-degree test ( = 90°), the pitch corresponds to half-wavelength.
When evaluating a 64-element, 5 MHz phased array (PA) probe with a 0.6 mm element pitch in steel, wavelength (λ) is approximately 1.2 mm. From Equation 1, and considering the center frequency, we confirm that the maximum pitch of 0.6 mm is satisfied. However, it's crucial to note that the center frequency does not represent the maximum frequency of the transducer.
To illustrate, employing a customary 60% bandwidth yields frequencies around 6.5 MHz, accompanied by a corresponding wavelength of 0.9 mm. Given that the pitch surpasses halfwavelength under these circumstances, the potential for aliasing emerges, leading to the inception of artifacts within the images. This phenomenon primarily stems from the inclusion of higher frequency components within the same phased array (PA) probe.
Contrarily, when examining thicker structures, larger elements become necessary. However, this requirement does not pose a challenge for regular pitch dimensions. The amplitudes of the plane waves remain constant within the near-field depth, which can be estimated approximately as 2 /4 , where D represents the PA aperture and lambda is the wavelength [8]. For the 0.6 mm pitch, the near field extends beyond 300 mm in steel.
To surpass the challenges posed by spatial sampling rules and enable the use of bigger element spacing without causing unwanted distortions in the images, we adopt an interpolation technique in our study. Our strategy is influenced by proven methods widely employed in seismic processing to handle situations where the sampling falls short (sub-Nyquist sampling), frequently leading to the appearance of unwanted patterns like aliasing and grating lobes.
Earlier investigations have addressed similar challenges; for instance, [5] employed Least Square Migration to interpolate PWI images through sub-sample array data. However, our approach follows a distinct trajectory by leveraging the Fourier domain technique. This approach allows for the use of subsampled data or the creation of a smaller virtual element pitch that satisfies the sampling criteria. This paper begins by explaining the methods employed for interpolation and reconstruction. Following that, we delve into the results section featuring various datasets. Finally, we conclude with our findings and outline future research prospects.

Methods: Interpolation and Reconstruction
A plane wave is characterized by its consistent wavefront that extends infinitely in a particular direction. This uniformity enables accurate prediction of its propagation attributes, encompassing velocity, direction, and amplitude. This inherent predictability can be leveraged for interpolation through the application of linear prediction filters in the Fourier domain [9]. To tackle the issue of spatial aliasing, this section introduces Spitz interpolation [10], followed by Fourier domain reconstruction using Stolt f-k migration [11][12].

Interpolated PWI
In the realm of numerical analysis, interpolation is conventionally defined as the process of creating new data points within the span of a discrete collection of established data points. The concept of mathematics can be likened to reading between the lines. In the context of Non-Destructive Testing (NDT), this analogy can be extended to interpreting data by reading between the A-scans. Spitz (1991) introduced a technique aimed at conducting interpolation of seismic traces, designed specifically to manage spatially aliased events [10]. Spitz's approach leverages the insight that linear events within a segment constituted of uniformly spaced traces can be accurately interpolated, regardless of the initial spatial interval. This methodology is particularly effective for 2D (two-dimensional) interpolation, aligning with the context of the PWI data we are addressing.
As derived from the original formulation [10], it is demonstrated that N evenly spaced A-scans containing L linear events can be precisely interpolated when L is smaller than N. If we consider a section in the time-space domain (t-x) composed of N equally spaced A-scans, each containing L linear events, transforming each trace via Fourier transformation yields a section of complex data in the frequency-space domain (f-x). Let's denote the data in the f-x domain as g(f) = (g1(f), g2(f), …, gN(f)) T . It's important to note that in the context of PWI, the quantity of traces N corresponds to the number of firing elements in the transducer. Each A-scan gk may be represented in the frequency domain as: where ( ) is the Fourier transform of the ultrasound pulse corresponding to the jth event. And ( ) = 2 is the event phase shift at the frequency f, corresponding to the time shift between adjacent A-scans.
Our objective is to perform interpolation between adjacent traces by utilizing a prediction filter that operates over the length of the prediction filter between the traces, resulting in a halved pitch-a form of first-order interpolation. In this scenario, the forward-backward one-step prediction filter is denoted as ( ) , and can be estimated in the least squares sense by equations 3 and 4: where * denotes the complex conjugate.
As previously demonstrated, the interpolation process relies on the application of predictionerror filters within the Fourier domain, effectively interpolating plane-wave data. This methodology is realized through the adaptation of a Fourier-based beamforming algorithm. The interpolation procedure is employed on PWI data after the initial phase of beamforming, which results in the creation of a plane-wave image.
Interpolation is essential in two-dimensional (2D) data due to the applicability of the Nyquist principle, which extends to 2D data. In this context, sample rates are specified for each axis, and each axis is associated with a Nyquist frequency that acts as a constraint on the frequency composition of the data [13].
In our current context, the application of Spitz's findings from 1991 reveals that the unit step prediction filter, calculated at a frequency of f/2 using data with a specific pitch, can be utilized at a frequency f to estimate data with a halved pitch. Porsani (1999) has further validated and enhanced this approach, which will be the focus of forthcoming research [14].
Spitz's technique was originally devised for linear events, yet it demonstrates the ability to effectively address events with moderate curvatures that exhibit some lateral amplitude variations. Notably, this adaptability does not involve compromising vertical resolution in favor of lateral resolution [10].

Reconstruction by Stolt f-k Migration
After the interpolation phase, image reconstruction is executed through Fourier domain migration [10] [11]. The frequency-wavenumber (f-k) migration, initially introduced by Stolt in 1978, continues to be the most effective strategy for uncomplicated velocity models [15]. Although the lateral and axial velocity variations of the Earth present a notable restriction in seismic applications, this constraint does not extend to non-destructive evaluation situations. In Plane Wave Imaging (PWI), the wavefront maintains a near-constant velocity, harmonizing with the foundational assumption of the Stolt method.

Results and Discussion
PWI data was collected using a 64-element transducer with a 1 mm pitch and a frequency of 5 MHz (with 60% bandwidth). The ASTM E2491 PA (Type B) assessment block was used, which included 18 Side-Drilled Holes (SDHs) with a 2 mm diameter and an additional 18 SDHs with a 1 mm diameter. It is a 75 mm thick 1018 steel block.  Figure 3 shows the raw PWI data in the left plot, with a 1 mm separation between A-scans matching the transducer pitch. The right plot displays the interpolated PWI data, with the pitch halved to 0.5 mm. Importantly, all reflections showcase a uniform curvature arising from their interaction with the plane wave. This interpolated image demonstrates a notably smoother transition compared to the unprocessed data. The efficacy of the Spitz interpolation method becomes evident when applied to PWI data, where reflectors exhibit quasi-linear behavior or moderate curvatures, and the ultrasound speed maintains a relatively consistent profile. The adoption of this technique facilitates a reduction in the original transducer pitch by half, resulting in improved reconstruction quality. Additionally, the method's adaptability is highlighted in firing configurations involving alternating elements. For instance, selectively activating solely the odd elements within a 64element transducer allows for the generation of a comprehensive array through interpolation, thereby expanding the method's range of capabilities and adaptability.
Furthermore, the application of the Stolt migration technique in the reconstruction process offers notable advantages when employed on a densely populated matrix containing a larger number of A-scans. Executed within the Fourier domain, this migration approach reduces the necessity for extensive zero-padding to accomplish lateral interpolation, thus enhancing both the efficiency and accuracy of the reconstruction procedure. Figure 4 and Figure 5 depict the outcomes of Stolt migration before and after interpolation, respectively.

Conclusions
In summary, our study underscores the significance of interpolation as a rapid and efficient technique that maintains high frame rates in ultrasound imaging. While Plane Wave Imaging (PWI) might not fully achieve the image quality of Full Matrix Capture (FMC), our approach brings it closer to this goal. The interpolation process effectively transforms a 64-element probe into a virtual 128-element probe, enhancing the imaging capabilities.
Maintaining a small element pitch is crucial in ultrasound imaging to achieve constructive interference and meet spatial Nyquist sampling criteria. However, challenges arise when pursuing larger element pitches. To overcome this, an interpolation technique inspired by seismic processing is proposed. Using the Fourier domain technique, a smaller virtual element pitch is created to meet sampling criteria, enabling larger element pitches without introducing aliasing. This not only enhances ultrasound image quality but also provides flexibility in choosing PA probes.