High-Resolution Ultrasound Imaging Enabled by Random Interference and Joint Image Reconstruction

In ultrasound, wave interference is an undesirable effect that degrades the resolution of the images. We have recently shown that a wavefront of random interference can be used to reconstruct high-resolution ultrasound images. In this study, we further improve the resolution of interference-based ultrasound imaging by proposing a joint image reconstruction scheme. The proposed reconstruction scheme utilizes radio frequency (RF) signals from all elements of the sensor array in a joint optimization problem to directly reconstruct the final high-resolution image. By jointly processing array signals, we significantly improved the resolution of interference-based imaging. We compare the proposed joint reconstruction method with popular beamforming techniques and the previously proposed interference-based compound method. The simulation study suggests that, among the different reconstruction methods, the joint reconstruction method has the lowest mean-squared error (MSE), the best peak signal-to-noise ratio (PSNR), and the best signal-to-noise ratio (SNR). Similarly, the joint reconstruction method has an exceptional structural similarity index (SSIM) of 0.998. Experimental studies showed that the quality of images significantly improved when compared to other image reconstruction methods. Furthermore, we share our simulation codes as an open-source repository in support of reproducible research.


Introduction
Medical imaging refers to different tools and techniques used in hospitals for diagnostic purposes. Noninvasive imaging modalities, particularly ultrasound imaging, play an important role in the early identification of abnormalities. In this study, we are interested in improving the spatial and contrast resolution of ultrasound imaging modalities. Conventional ultrasound imaging methods use beamforming to spatially divide the region of interest (ROI) into multiple scanlines [1,2]. Each scanline is then imaged separately using a focused ultrasound pulse. Conventional ultrasound techniques that utilize beamforming are capable of achieving an image resolution that is approximately equivalent to the width of the transmitted ultrasound beam [3]. The smallest beam width is achieved at a focal point. For a central frequency of 3 MHz, the best spatial resolution is approximately equal to 1 mm [3,4]. However, owing to wave diffraction, it is impossible to achieve an infinitely narrow ultrasound beam that travels precisely along a single scanline [4,5]. As a result of the diffraction effect, the focused ultrasound pulse is reflected not only from the tissue within the desired scanline but also from the adjunct areas. The echo signals reflected from adjunct areas degrade the image resolution by corrupting it with speckle noise [6][7][8].
In the past, different techniques have been used to enhance the spatial and contrast resolution of ultrasound systems [9,10]. For instance, in [11,12], speckle reduction techniques using statistical models have been proposed. Speckle reduction techniques use mathematically modeled filters to suppress speckle noise. To achieve significant improvements in resolution, we believe that it is insufficient to modify only the signal processing (software) part without modifications to the imaging protocol (hardware). An excellent example of symbiosis between software and hardware is the introduction of contrast agents, where two parts complement each other to produce significant improvements in resolution [13][14][15]. Contrast microbubbles are used to create a spatially localized sound source that significantly improves the ultrasound resolution.
In ultrasound, coded excitation signals are a well-known technique used to improve the quality of images [16,17]. However, to the best of our knowledge, using excitation signals to generate an ultrasound field with a spatially varying property such as one presented in this paper is a relatively new research direction. In the past few years, there is a growing interest in the development of imaging systems that use unfocused ultrasound fields with spatially variant properties. In [18], an object much smaller than the wavelength of the signal is reconstructed using a priori measurements of an incident ultrasound field with the spatially variant property. In [19], a 3D ultrasound imaging setup is shown that uses only one sensor. The sensor is equipped with an acoustic lens that yields individual objects in the ROI to be uniquely identifiable by the echo signals. In [20,21], specially tailored optoacoustic lenses are used to generate arbitrary acoustic fields and holograms. It was demonstrated that imaging with random ultrasound fields has the potential to improve the spatial resolution of ultrasound scanners beyond present limits.
In [22], we proposed a novel interference-based ultrasound imaging method, where instead of a traditional focused ultrasound beam, we use the transmission of an unfocused ultrasound wavefront of random interference. Here, we further improve the spatial and contrast resolutions of the interference-based imaging method by introducing a joint image reconstruction scheme. The new joint reconstruction scheme allows us to directly reconstruct high-resolution ultrasound images by utilizing radio frequency (RF) signals from all elements of the sensor array in a joint optimization problem. The joint reconstruction scheme provides superior contrast and image resolutions. According to the simulation study, the proposed joint reconstruction method leads to significant resolution improvements compared to different beamforming-based methods. The proposed method achieves a mean-squared error (MSE) of 0.00469, a peak signal-to-noise ratio (PSNR) of 33.3 dB, a signal-to-noise ratio (SNR) of 27.6 dB, and structural similarity (SSIM) of 0.998. This constitutes to a two-time resolution improvement over the previously proposed interference-based compound method. We share our simulation codes as an open-source repository. The link to the code repository can be found in Supplementary Materials. The remainder of this paper is organized as follows. In Section 2, we describe the basic concept of ultrasound imaging using random interference and compressive sensing theory. In Section 3, we explain the proposed joint image reconstruction scheme. In Section 4, we provide the simulation and experimental results obtained. In Section 5, we discuss areas of future study. Finally, in Section 6 we summarize the contributions of this research.

Background
Traditionally, ultrasound interference has been treated as an undesired effect that degrades the quality of ultrasound images. We have recently shown that an ultrasound wavefront of random interference can be used as a means to differentiate between individual point scatterers in the ROI [23]. In this section, we review the effects of random interference and compressive sensing theory.

Effect of Random Interference
In conventional ultrasound imaging, beamforming is used to spatially divide the ROI into several scanlines. An ultrasound image is then acquired one scanline at a time by using beamforming to focus the transmit and receive ultrasound pulses. The image is reconstructed on the basis of the assumption that the received ultrasound signals consist of echoes reflected only from inhomogeneity within the given scanline. However, in practice, owing to the diffraction effect, the incident wavefront reflects from the medium within the scanline and its adjunct areas [23]. Owing to the imperfect coherence of echo signals, the beamformed signal is corrupted by interference patterns known as speckle noise [6][7][8]. Thus, in conventional focused methods, ultrasound interference is a highly undesired effect that degrades the image resolution.
Unlike conventional beamforming-based methods, we intentionally create an unfocused ultrasound wavefront of random interference, which yields spatial impulse responses of individual point scatterers in the ROI to be mutually incoherent [22]. Thus spatial resolution can be achieved by identifying individual point scatterers on the basis of their spatial impulse responses.
Here, we demonstrate the effect of random interference, where a small change in the ROI incurs a significant difference in the received RF signals. In Figure 1, we show simulated RF signals obtained when a wavefront of random interference reflects from a group of scatterers [24][25][26]. In Figure 1, we show the echo signals of the incident wavefront of random interference reflected from the scatterers as grayscale images. In Figure 1a, the ROI includes a group of scatterers consisting of four points located at an axial distance of 45 mm and lateral distances of −5, −1, 1, and 5 mm, respectively. In Figure 1b, a group of point scatterers consists of four points located at lateral distances of −4.75, −1, 1, and 5 mm. In Figure 1b, the position of the first scatterer was changed from −5 mm to −4.75 mm. Figure 1c depicts the absolute difference between the images shown in Figure 1a,b. From Figure 1c, we can see that the effect of random interference yields a significant difference in the received echo signals even when we make a small change to the scatterer map. Using random interference and its effect on the scatterers in the ROI, we can design an imaging scheme where an ultrasound image is reconstructed using a priori information about wave propagation [18,22]. We represent the entire ROI as a set of individual spatial points. Then, the RF signals received at the sensor array can be decomposed into a set of echo signals reflected from individual point scatterers in the ROI. In a simulation study, a priori information can be acquired by generating RF signals for every point scatterer in the ROI. In the experimental study, a priori information can be generated from measurements of a plastic plane submerged into a water tank [22]. We summarize the steps used to acquire a priori information in Algorithm A1.

Compressive Sensing
With the proposed imaging method, we chose to use the compressive sensing (CS) theory because it provides a robust solution to underdetermined systems of linear equations. Originally, compressive sensing theory was introduced to reconstruct a signal using fewer measurements than that proposed by the Shannon-Nyquist sampling theorem [27][28][29]. However, later studies have shown that CS theory can be used to significantly improve the resolution of imaging modalities. During the past few decades, CS theory has been successfully used to improve the resolution of spectrometers, microscopes, and many other imaging modalities [30][31][32][33][34][35][36].
The concept of classical CS theory can be explained using a system of linear equations y = Ax = ΦΨx, where x ∈ R N is a signal to be reconstructed, A ∈ R M×N is a specially constructed matrix consisting of a sensing matrix Φ and a sparsifying basis Ψ, and y ∈ R M is a vector representing the signal acquired using a measuring instrument. The goal of CS is to reconstruct x from y and G when M << N. In such a case, y = Ax is an underdetermined system because we have more unknowns than equations. The matrix A must characterize the transmission properties of the imaging system, therefore, it is called the transmission matrix. In the case of interference-based imaging, the transmission matrix characterizes the propagation properties of the wavefront of random interference. The CS theory guarantees the successful reconstruction of x under the following two conditions. The first condition requires the signal x to be sparse in some domain. The signal is called sparse if K non-zero values of the signal satisfy K < M << N. The second condition requires a sensing matrix to be incoherent. The matrix is called incoherent when the cross-correlation of its columns is small. Many good algorithms were recently introduced to solve the l 1 -norm minimization problem and reliably reconstruct the signal of interest. In particular, we use the YALL1 MATLAB package solver [37]. decomposed into a set of echo signals reflected from individual point scatterers in the ROI. In a simulation study, a priori information can be acquired by generating RF signals for every point scatterer in the ROI. In the experimental study, a priori information can be generated from measurements of a plastic plane submerged into a water tank [22]. We summarize the steps used to acquire a priori information in Algorithm A1.

Method
In [22], we proposed a new ultrasound imaging method based on random interference. In this paper, we further improve our results by proposing a joint image reconstruction scheme that combines RF signals from all elements of the sensor array in a single l 1 -norm minimization problem. As a result, the new proposed imaging scheme is capable of reconstructing ultrasound images with higher resolution.

Observation Model
Let us consider a linear transducer array with 128 piezoelectric elements. We let the vector r i for each i ∈ {1, 2, . . . , N Rx } indicate the position of the array elements. In the proposed method, an acoustic wave is generated by simultaneously exciting all array elements using unique random signals. We use random excitation signals to drive the elements of the sensor array, which yields an ultrasound wavefront of random interference. We summarize the steps used to generate random excitation signals in Algorithm A1(Step 1). Let us consider a square image from 35 mm to 55 mm in the axial direction and from −10 mm to 10 mm in the lateral direction. Let us represent a medium within the ROI as a collection of point scatterers that are evenly distributed at equal distances d = 0.25 mm apart. Thus, the image of interest consists of a total of N Sc = 6561 point scatterers. We let the vector r k for each k ∈ {1, 2, . . . , N Sc } indicate the position of the kth scatterer in the ROI. We denote an ultrasound image as a vector f i := [f 1 f 2 . . . f N sc ] for i ∈ {1, 2, . . . , N Rx } and k ∈ {1, 2, . . . , N Sc }, where elements f k represent the inhomogeneity of the medium, which gives rise to the echo signal reflected off the scatterer at the kth spatial location. In Figure 2, we show a diagram of the proposed imaging method utilizing random interference. In Figure 2a, the random excitation signals are used to drive the transducer elements that yield an ultrasound wavefront of random interference. When the wavefront of random interference Sensors 2020, 20, 6434 5 of 15 travels through the medium, its energy is reflected from the scatterers. Here, we consider a simple 3 × 3 image that can be represented as a vector f i : Figure 2b, we show that an echo signal p i received at the ith receiving element can be represented as the sum of individual impulse responses of the scatterer.
. In Figure 2b, we show that an echo signal i p received at the ith receiving element can be represented as the sum of individual impulse responses of the scatterer. We can express an ultrasound echo signal received at the ith transducer element as , is a vectorized image, as observed by the ith array element. According to the principle of superposition [9], the received ultrasound echo signal i p is the sum of individual signals reflected off the scatterers i f . Our goal is to find an estimate of the image f given the RF signals i p and transmission matrices , this can be accomplished using CS theory [27][28][29].

Joint Image Reconstruction
In [22], we proposed an interference-based imaging method that used the compounding of several independently reconstructed images ˆi f given a set of echo signals i p and transmission matrices i G . In this study, we propose a joint image reconstruction scheme that further improves the resolution of interference-based ultrasound imaging by directly reconstructing a high-resolution image.
When an ultrasound wavefront of random interference propagates through the medium, its energy is partially reflected and received at the transducer array. The RF signals carry information about the same object image Object f observed at slightly different angles. In the interference-based joint image reconstruction method, we propose directly reconstructing a high-resolution image by utilizing all signals from the array in a single optimization problem.
First, to estimate the image f , we need to obtain transmission matrices that carry information about the propagation of the proposed ultrasound wavefront of random interference. We use the spatial impulse responses of individual point scatterers as a priori information. We can then reconstruct the ultrasound image using the following optimization problem: : arg min subject to .
In Equation (1), the image i f is reconstructed from an echo signal i p acquired at a single receiving channel i. In a single pulse-echo transmission, we use all elements of the array to receive the reflected We can express an ultrasound echo signal received at the ith transducer element as where p i ∈ R M is the received RF signal, G i ∈ R M×N is the transmission matrix, and f i ∈ R N is a vectorized image, as observed by the ith array element. According to the principle of superposition [9], the received ultrasound echo signal p i is the sum of individual signals reflected off the scatterers f i . Our goal is to find an estimate of the image f given the RF signals p i and transmission matrices G i for i ∈ {1, 2, . . . , N Rx }, this can be accomplished using CS theory [27][28][29].

Joint Image Reconstruction
In [22], we proposed an interference-based imaging method that used the compounding of several independently reconstructed imagesf i given a set of echo signals p i and transmission matrices G i . In this study, we propose a joint image reconstruction scheme that further improves the resolution of interference-based ultrasound imaging by directly reconstructing a high-resolution image.
When an ultrasound wavefront of random interference propagates through the medium, its energy is partially reflected and received at the transducer array. The RF signals carry information about the same object image f Object observed at slightly different angles. In the interference-based joint image reconstruction method, we propose directly reconstructing a high-resolution image by utilizing all signals from the array in a single optimization problem.
First, to estimate the image f, we need to obtain transmission matrices that carry information about the propagation of the proposed ultrasound wavefront of random interference. We use the spatial impulse responses of individual point scatterers as a priori information. We can then reconstruct the ultrasound image using the following optimization problem: Sensors 2020, 20, 6434 6 of 15 In Equation (1), the image f i is reconstructed from an echo signal p i acquired at a single receiving channel i. In a single pulse-echo transmission, we use all elements of the array to receive the reflected echo signals. Therefore, we have 128 different versions of Equation (1), one for each receiving element in the array, as follows: Owing to the effect of random interference, the received ultrasound signals p i and the transmission matrices G i for each i ∈ {1, 2, . . . , N Rx } carry unique information about the image of interest f i . In [22], an ultrasound image is reconstructed by applying the optimization problem in Equation (2) to the individual equations in Equation (3). The obtained imagesf i for i ∈ {1, 2, . . . , N Rx } are combined to form a high-resolution image as follows: However, a better and more accurate approach is to utilize the fact that the image of interest f Object does not change depending on the receiving element number i. Therefore, the RF-signals p i for i ∈ {1, 2, . . . , N Rx } include information about the same object image f Object . In such a case, we can use a joint image reconstruction scheme to further enhance the reconstruction accuracy and image resolution. As a simple solution to the joint reconstruction problem, we will next discuss a matrix inversion approach. First, we rearrange the measurement signals p i and transmission matrices G i in Equation (3) as follows: For simplicity, let us use a subscript T to denote the tall matrix and tall vector in Equation (5) as follows: Then, a matrix inversion solution to Equation (6) can be written as where (·) denotes a transpose. Then, the object image is equal to In such a case, we can directly reconstruct a high-resolution image by utilizing information across the elements of the array. Similarly, we can use the optimization problem in Equation (2) to recover the imagef Object as follows: where · 1 denotes the l 1 -norm and ε ≥ 0. To successfully reconstruct an object imagef Object , the transmission matrix in Equation (9) needs to be incoherent. We achieved incoherency by transmitting an ultrasound wavefront of random interference. We measure the incoherence of the transmission matrix G using the Gram matrix D = G * G. We assume that each column of the matrix G is normalized to one. The absolute value of the off-diagonal elements d i,j of matrix D represents a cross-correlation value between the corresponding pair of columns i and j of the matrix G.

Simulation Study
In our simulation study, we use the Field II ultrasound software [24][25][26]. First, we simulate a linear transducer with 128 identical elements. Second, each array element was excited with a random sequence, as described in Section 2.1 and [22]. The central frequency was set to 3 MHz. The sampling frequency was set to 40 MHz.
To demonstrate the resolution capabilities of the newly proposed joint reconstruction scheme, we simulate a Shepp-Logan phantom. In Figure 3, we provide a side-by-side comparison of different reconstruction methods. The scattering map of the phantom is shown in Figure 3a. In Figure 3b, we show an image reconstructed using a conventional focused B-mode method. The image was reconstructed using 128 scanlines. The focal point was set to 45 mm depth. The image of the phantom appears with large sidelobes and the details of the phantom cannot be properly reconstructed.
In Figure 3c, we show an image reconstructed using synthetic aperture imaging provided in the "Beamformation toolbox" [38]. The image was simulated using 65 emissions and beamformed with 180 lines. In Figure 3d, we show an image reconstructed using the interference-based compound method proposed in [23]. The reconstructed image features a more accurate location and intensity of the scatterers that precisely matches the original scatterer map of the phantom image. However, in the case of the compound method, the image in Figure 3d is corrupted by noise from the compounding operation. The presence of noise affects the visibility of small elements of the phantom at a lateral distance of 0 mm. In Figure 3e, we show an image reconstructed using the proposed joint image reconstruction scheme. The proposed joint reconstruction scheme allows us to successfully reconstruct a high-resolution image of the phantom. Figure 3e shows an image obtained using the proposed joint image reconstruction scheme. The image in Figure 3e has accurate spatial and contrast resolutions. The noise specific to the image in Figure 3d was eliminated, and the small details of the phantom at a lateral distance of 0 mm are accurately reconstructed. The image reconstructed using the newly proposed method exhibits less speckle noise and a more accurate representation of small components of the Shepp-Logan phantom. We compare the results from the Shepp-Logan phantom study in Table 1. The proposed interference-based joint reconstruction method, among all methods, has a lower MSE, higher PSNR and SNR, and an extremely high value for SSIM. Table 1. Mean-squared error (MSE), peak signal-to-noise ratio (PSNR), signal-to-noise ratio (SNR), and structural similarity index (SSIM) values for conventional and proposed imaging methods.

Experimental Study
In addition to the simulation results, we conducted an experimental study using a custom research ultrasound scanner supplied by Alpinion medical systems, Seoul, South Korea and a tissue-mimicking phantom Model 040GSE manufactured by CIRS, Virginia, USA. The ultrasound research equipment was equipped with an arbitrary wave generator (AWG) and dedicated memory to store the array of excitation signals. The AWG can store an array of size 128 × 2048 of eight-bit data, where each row

Experimental Study
In addition to the simulation results, we conducted an experimental study using a custom research ultrasound scanner supplied by Alpinion medical systems, Seoul, South Korea and a tissuemimicking phantom Model 040GSE manufactured by CIRS, Virginia, USA. The ultrasound research equipment was equipped with an arbitrary wave generator (AWG) and dedicated memory to store the array of excitation signals. The AWG can store an array of size 128 × 2048 of eight-bit data, where each row is dedicated to the corresponding element of the sensor array. We used a linear transducer with 128 identical piezoelectric sensors (L3-12 manufactured by Alpinion, Seoul, South Korea). The central frequency was set to 3 MHz. Likewise, the sampling frequency was set to 40 MHz.
Throughout this paper, we use the following parameters for different imaging techniques mentioned in the experimental study. In the case of the conventional focused B-mode, the center Throughout this paper, we use the following parameters for different imaging techniques mentioned in the experimental study. In the case of the conventional focused B-mode, the center frequency was set to 6 MHz, the number of scanlines was set to 128, and the focal point was set at 5 cm depth. In the case of the plane-wave imaging, the center frequency was set to 6 MHz and the number of plane wave emissions was set to 7. Likewise, in the case of synthetic aperture imaging, we used 128 transmissions, one for each element in the sensor array, and the image consisted of 250 lines.

Subarray Imaging
In Equation (5), we suggest arranging echo signals p i and matrices G i into a tall vector p T and a tall matrix G T , respectively. Although this is the best approach in practice, it is feasible only for small size problems. For example, in our experimental study, the dimension of RF signals is 128 × 3300, and the number of spatial points in ROI is 38,801. The dimension of the transmission matrix for a single ith receiving element is 3300 × 38,801. To use Equations (5) and (9), we would have to work with an extremely large tall matrix G T . Here, we propose simplifying the joint image reconstruction method given in Section 3.2 and Equation (5) by reducing it to a set of smaller problems that divide the ROI into smaller imaging zones while still utilizing the joint reconstruction method.
First, we perform 10 pulse-echo cycles each time with a different set of random excitation signals. In Figure 4, we show how RF signals from a single pulse-echo transmission can be rearranged into 18 different sets of subarrays signals. Then, we use subarray signals to reconstruct a subarray image using Equations (5) and (9). and a tall matrix T G , respectively. Although this is the best approach in practice, it is feasible only for small size problems. For example, in our experimental study, the dimension of RF signals is 128 × 3300, and the number of spatial points in ROI is 38,801. The dimension of the transmission matrix for a single ith receiving element is 3300 × 38,801. To use Equations (5) and (9), we would have to work with an extremely large tall matrix T G . Here, we propose simplifying the joint image reconstruction method given in Section 3.2 and Equation (5) by reducing it to a set of smaller problems that divide the ROI into smaller imaging zones while still utilizing the joint reconstruction method. First, we perform 10 pulse-echo cycles each time with a different set of random excitation signals. In Figure 4, we show how RF signals from a single pulse-echo transmission can be rearranged into 18 different sets of subarrays signals. Then, we use subarray signals to reconstruct a subarray image using Equations (5) and (9). In Figure 5a, we show a diagram of a tissue-mimicking phantom that represents the ROI. In Figure 5b, we show three rectangles that divide the ROI into three imaging zones. By reconstructing images separately for different imaging zones we can reduce the computational complexity. Likewise, in Figure 5c, we show the acceptance angle of a single piezo-electric sensor. Due to the sensor's acceptance angle, not all points in the ROI contribute to the received RF signal. Therefore, echo signals reflected from spatial points that lay outside of the acceptance angle do not contribute to the RF signal acquired by the sensor. By removing columns that represent spatial responses of such points we can further reduce the computational complexity and increase image reconstruction accuracy. In Figure 5d-f, we show the acceptance angle for elements in the subarray. In Figure 5g-i, we show low-resolution images reconstructed for the selected subarrays using data acquired during a single pulse-echo transmission. We use data from 10 pulse-echo transmissions, each divided into 18 subarray configurations, and 3 imaging zones. In total, we have 540 subarray images that combined to form a final high-resolution image. In Figure 6, we compare images of a dental floss submerged into a water tank and reconstructed using different methods. In Figure 6a, we show an image reconstructed using a conventional focused (d-f) acceptance angle for subarrays; (g-i) low-resolution images obtained using joint reconstruction and subarrays.
In Figure 5a, we show a diagram of a tissue-mimicking phantom that represents the ROI. In Figure 5b, we show three rectangles that divide the ROI into three imaging zones. By reconstructing images separately for different imaging zones we can reduce the computational complexity. Likewise, in Figure 5c, we show the acceptance angle of a single piezo-electric sensor. Due to the sensor's acceptance angle, not all points in the ROI contribute to the received RF signal. Therefore, echo signals reflected from spatial points that lay outside of the acceptance angle do not contribute to the RF signal acquired by the sensor. By removing columns that represent spatial responses of such points we can further reduce the computational complexity and increase image reconstruction accuracy.
In Figure 5d-f, we show the acceptance angle for elements in the subarray. In Figure 5g-i, we show low-resolution images reconstructed for the selected subarrays using data acquired during a single pulse-echo transmission. We use data from 10 pulse-echo transmissions, each divided into 18 subarray configurations, and 3 imaging zones. In total, we have 540 subarray images that combined to form a final high-resolution image.

Single Wire Study
In Figure 6, we compare images of a dental floss submerged into a water tank and reconstructed using different methods. In Figure 6a, we show an image reconstructed using a conventional focused B-mode method that can be found in all modern ultrasound systems. In Figure 6b, we show an image reconstructed using plane-wave imaging with 30 plane wave transmissions. In Figure 6c, we show an image reconstructed using synthetic aperture imaging. In Figure 6d, we show an image reconstructed using the proposed interference-based joint reconstruction method. The images in Figure 6a-c appear with side lobes caused by strong reflections of ultrasound waves that corrupt the image scanline during the beamforming operation. Meanwhile, the proposed method is capable of reconstructing a more accurate image of a wire without side lobes. The image in Figure 6d was reconstructed using data from a single pulse-echo transmission of a random interference wave.
Sensors 2020, 20, x FOR PEER REVIEW 11 of 16 Figure 6. Experimental results using a nylon wire submerge into a water tank: (a) an image reconstructed using the conventional focused B-mode; (b) an image reconstructed using plane-wave imaging; (c) an image reconstructed using synthetic aperture imaging; (d) an image reconstructed using the interference-based joint reconstruction method using data from a single pulse-echo transmission.

Tissue-Mimicking Phantom Study
In Figure 7, we compare images of a tissue-mimicking phantom reconstructed using different methods. In Figure 7a, we show images of a tissue-mimicking phantom reconstructed using a conventional focused B-mode method. In Figure 7b, we show an image reconstructed using planewave imaging. In Figure 7c, we show an image reconstructed using synthetic aperture imaging. In Figure 7d, we show an image reconstructed using the interference-based compound method. Likewise, in Figure 7e we show an image reconstructed using the proposed interference-based joint reconstruction method. From Figure 7, it can be observed that the proposed method reconstructs cyst and nylon wires with greater accuracy. The diameter of the nylon wires appears much closer to that of the original. Similarly, the speckle noise has been removed. The proposed joint image reconstruction method further improves the resolution of the interference-based ultrasound. The Figure 6. Experimental results using a nylon wire submerge into a water tank: (a) an image reconstructed using the conventional focused B-mode; (b) an image reconstructed using plane-wave imaging; (c) an image reconstructed using synthetic aperture imaging; (d) an image reconstructed using the interference-based joint reconstruction method using data from a single pulse-echo transmission.

Tissue-Mimicking Phantom Study
In Figure 7, we compare images of a tissue-mimicking phantom reconstructed using different methods. In Figure 7a, we show images of a tissue-mimicking phantom reconstructed using a conventional focused B-mode method. In Figure 7b, we show an image reconstructed using plane-wave imaging. In Figure 7c, we show an image reconstructed using synthetic aperture imaging. In Figure 7d, we show an image reconstructed using the interference-based compound method. Likewise, in Figure 7e we show an image reconstructed using the proposed interference-based joint reconstruction method. From Figure 7, it can be observed that the proposed method reconstructs cyst and nylon wires with greater accuracy. The diameter of the nylon wires appears much closer to that of the original. Similarly, the speckle noise has been removed. The proposed joint image reconstruction method further improves the resolution of the interference-based ultrasound. The experimental results shown in Figure 7 are consistent with our simulation study shown in Figure 3. The image in Figure 7e was reconstructed using data from a ten pulse-echo transmission of random interference wave.

Discussion
Here, we share some thoughts on possible future ways to improve the resolution of interferencebased imaging. The contribution of the proposed method is related to the area of image acquisition and reconstruction. Therefore, we did not use any signal/image post-processing techniques such as filtering or smoothing. However, in the future, it would be of great interest to see how modern speckle reduction techniques such as in [39,40] can be adapted to enhance the image quality of the proposed method. Likewise, in the proposed method, it is very important to obtain accurate a priori information. Therefore, it is of great interest to see if advanced beamforming algorithms such as those in [41][42][43] can be used to better estimate the spatial impulse responses.
The theoretical temporal resolution of the proposed method is equal to the time required to perform 10 pulse-echo transmissions. However, a significant amount of time is required to reconstruct high-resolution images. We use a PC equipped with Intel Xeon Gold 6240R, 256 GB RAM, and installed MATLAB environment. The maximum size of the array that can be efficiently processed using available PC is 40 × 3300 × 38,801. A single subarray image is reconstructed in about 110 s. The final high-resolution image is reconstructed in about 16 h. While it is a significant deviation from a common real-time ultrasound imaging, the proposed method achieves 250 μm resolution and has the potential to compete in some applications with 3 T magnetic resonance imaging MRI [44].
Using the proposed method, we were not able to obtain high-resolution in vivo images because

Discussion
Here, we share some thoughts on possible future ways to improve the resolution of interference-based imaging. The contribution of the proposed method is related to the area of image acquisition and reconstruction. Therefore, we did not use any signal/image post-processing techniques such as filtering or smoothing. However, in the future, it would be of great interest to see how modern speckle reduction techniques such as in [39,40] can be adapted to enhance the image quality of the proposed method. Likewise, in the proposed method, it is very important to obtain accurate a priori information. Therefore, it is of great interest to see if advanced beamforming algorithms such as those in [41][42][43] can be used to better estimate the spatial impulse responses.
The theoretical temporal resolution of the proposed method is equal to the time required to perform 10 pulse-echo transmissions. However, a significant amount of time is required to reconstruct high-resolution images. We use a PC equipped with Intel Xeon Gold 6240R, 256 GB RAM, and installed MATLAB environment. The maximum size of the array that can be efficiently processed using available PC is 40 × 3300 × 38,801. A single subarray image is reconstructed in about 110 s. The final high-resolution image is reconstructed in about 16 h. While it is a significant deviation from a common real-time ultrasound imaging, the proposed method achieves 250 µm resolution and has the potential to compete in some applications with 3 T magnetic resonance imaging MRI [44].
Using the proposed method, we were not able to obtain high-resolution in vivo images because biological tissue has a complex multi-layer structure that disperses the incident wavefront in various directions. Therefore, it is hard to generate transmission matrices that fit all possible cases of biological tissue. This is a topic of our ongoing research efforts. One can solve this issue with dynamically adjustable transmission matrices, where the ROI is divided into many thin layers. Then, the image reconstructed for the top layer can be used to estimate the transmission matrices for the next layers. Alternatively, a conventional B-mode image can be used as a reference to better estimate the transmission matrices.

Conclusions
In this work, we claim that by combining a novel interference-based image acquisition method with a joint image reconstruction algorithm we can achieve very accurate ultrasound images at 250 µm resolution with strong SNR and precise SSIM. The proposed joint reconstruction method uses echo signals across all array elements to directly estimate a high-resolution image. The proposed method does not require a focused ultrasound beam, thereby removing the speckle noise caused by diffraction. The interference-based method is different from conventional approaches that use the strength of the echo signals to visualize details of imaged objects. Instead, an intentionally created wavefront of random interference and a l 1 -norm minimization algorithm are used to visualize the object by identifying individual spatial impulse responses of the scatterers. Using simulation and experimental results, we have shown that the joint image reconstruction method improves the resolution of ultrasound images.