Multi-perspective bistatic ultrasound imaging and elastography of the ex vivo abdominal aorta

Knowledge of aneurysm geometry and local mechanical wall parameters using ultrasound (US) can contribute to a better prediction of rupture risk in abdominal aortic aneurysms (AAAs). However, aortic strain imaging using conventional US is limited by the lateral lumen–wall contrast and resolution. In this study, ultrafast multiperspective bistatic (MP BS) imaging is used to improve aortic US, in which two curved array transducers receive simultaneously on each transmit event. The advantage of such bistatic US imaging on both image quality and strain estimations was investigated by comparing it to single-perspective monostatic (SP MS) and MP monostatic (MP MS) imaging, i.e., alternately transmitting and receiving with either transducer. Experimental strain imaging was performed in US simulations and in an experimental study on porcine aortas. Different compounding strategies were tested to retrieve the most useful information from each received US signal. Finally, apart from the conventional sector grid in curved array US imaging, a polar grid with respect to the vessel’s local coordinate system is introduced. This new reconstruction method demonstrated improved displacement estimations in aortic US. The US simulations showed increased strain estimation accuracy using MP BS imaging bistatic imaging compared to MP MS imaging, with a decrease in the average relative error between 41% and 84% in vessel wall regions between transducers. In the experimental results, the mean image contrast-to-noise ratio was improved by up to 8 dB in the vessel wall regions between transducers. This resulted in an increased mean elastographic signalto-noise ratio by about 15 dB in radial strain and 6 dB in circumferential strain.


Multiperspective Bistatic Ultrasound Imaging and
Elastography of the Ex Vivo Abdominal Aorta

I. INTRODUCTION
A BDOMINAL aortic aneurysms (AAAs) are local dilations of the abdominal aortic wall with an aortic diameter of 3.0 cm or more [1], [2]. The prevalence of AAAs is reported to be between 1.9% and 18.5% in men older than 60 years and 0% and 4.2% in women older than 60 years, and there is a high risk of rupture as the diameter of the aorta increases [2]. The overall mortality rate of AAA rupture is reported up to 80% [3]. Surgical intervention is considered if the risk of rupture exceeds the risk associated with intervention. Current thresholds for surgical intervention are an aortic diameter of 5.0 cm for women and 5.5 cm for men [1]. However, the results of autopsy studies highlight the limitation of simply using the aortic diameter for risk assessment of AAAs. Some large diameter aneurysms (7-10 cm) remain stable, while other aneurysms can rupture below 5 cm [4]. This suggests that additional variables contribute to aneurysm rupture other than aneurysm diameter alone.
From a biomechanical perspective, aneurysm rupture is mechanical failure of the vascular wall. An aneurysm ruptures where the wall stress becomes higher than the tensile strength of the vascular tissue. Knowledge of local mechanical aortic wall parameters is, therefore, necessary to correctly analyze rupture risk and could contribute to a better prediction of rupture risk than the AAA diameter [5], [6]. Previous studies have demonstrated that the peak wall stress is an accurate indicator for rupture risk of AAAs, which can be estimated from imaging-derived geometries using a finite element analysis [7], [8]. Ultrasound (US) speckle tracking and strain imaging can be used for model personalization, yielding mechanical properties such as the wall stiffness [6], [9]. Moreover, strain can also be an important parameter in itself, as it may indicate local weakening of the AAA wall [10].
US strain imaging or elastography is a technique that was first described by Ophir et al. [11]. Since then, it has shown to be an accurate method for mechanical characterization and strain estimation in the vessel wall [12]. However, noninvasive aortic US imaging is challenging considering the typical imaging depth [5] and poor lateral resolution [13]. Furthermore, tissue interfaces perpendicular to the propagation direction of the US wave reflect well, whereas the interfaces parallel to the wave may not reflect the wave at all [14]. As a result, the images lack anatomical information for the reconstruction of the geometry of AAAs. This also has implications for This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see ht. tps://creativecommons.org/licenses/by/4.0/ the strain estimation in the vessel wall. To improve US image quality and strain estimations, the image area can be insonified from multiple beam-steering angles after which images or displacement fields are averaged by means of spatial compounding [15], [16]. Despite the improvements in both image quality, and motion and strain estimation, this technique cannot be applied to the abdominal aorta or deeperlying structures in general, considering the limited overlapping region in compounding at large depths [14]. To overcome this insufficient overlap, multiple probes can be used to achieve accurate motion and strain estimation over larger circumferential segments of the aortic wall. By minimizing the transmit time between transducers, image decorrelation can be limited and allows for image and displacement compounding of multiperspective (MP) acquisitions. This can be achieved using ultrafast imaging, which can reconstruct the entire US image using a single transmit and receive event with all elements in the US transducer [17].
We previously introduced ultrafast MP strain imaging, in which the displacements from the two perspectives were compounded after displacement estimation on single probe, ultrafast, compounded US data [18]. Ultrafast MP strain imaging showed improved motion tracking and circumferential strain estimates in the aortic wall in ex vivo experiments with porcine aortas [18]. However, by alternately transmitting and receiving with all elements of one probe in a MP monostatic (MP MS) or interleaved scanning sequence, accurate motion and strain estimation can only be performed in about 50% of the wall circumference. The full potential of high-framerate MP imaging can only be leveraged by using all available data in a dual-receive acquisition scheme. It allows for the use of MP bistatic (MP BS) imaging, where both probes receive simultaneously on each transmit event. Four signals can be reconstructed using an MP BS scanning sequence, instead of only two signals using MP MS imaging, without increasing the number of acquisitions. In Fig. 1, the single-probe (T 1 R 1 and T 2 R 2 ) and trans-probe signals (T 1 R 2 and T 2 R 1 ) in MP BS imaging are illustrated, in which T indicates the transmitting transducer, that is, either 1 or 2 and R indicates the receiving transducer, that is, either 1 or 2.
The addition of the trans-probe data allows for the reconstruction of an additional reflection on the aortic wall that occurs in between the transducers. Hence, accurate motion and strain estimation can be performed across almost the entire circumference of the aorta, which is also visualized in Fig. 1. However, the use of trans-probe data implies the need to develop completely new methods to derive MP compounded images, motion estimations, and strain in the aortic wall, which are presented in this article. The performance of single-perspective monostatic (SP MS), MP MS, and MP BS ultrafast image acquisition configurations using two curved array transducers is compared. Experimental strain imaging was performed in US simulations, as well as ex vivo in a study involving four porcine aortas. The US simulations were used to investigate the benefits of aortic MP BS strain imaging and to quantify the strain accuracy between image acquisition configurations. In the ex vivo study, different image and displacement compounding strategies were tested to retrieve Fig. 1. Illustration of the dual-receive acquisition sequence and all four signals in MP BS imaging (T 1 R 1 , T 1 R 2 , T 2 R 1 , and T 2 R 2 ). The 2-D imaging plane of the aorta is shown, together with the curved array US transducers. The blue and red circumferential segments indicate regions in which accurate estimation of vessel wall geometry and motion can be performed using each signal. the most useful information of each received US signal, and these strategies were also individually compared.
Finally, apart from the conventional sector image grid in curved array US imaging, a polar grid with respect to the vessel's local coordinate system is introduced on which the US RF signals were reconstructed. Instead of combining displacement estimations from the two coordinate systems of each receiving transducer, the polar grid allows for a more direct displacement estimation in the local radial direction of the aorta using coherently compounded MP image data. This allows for improved motion estimation of the aortic wall since motion is predominantly radial.

A. Simulation Design
US simulations were performed using the open-source MATLAB toolbox k-Wave [19], [20]. A medium was simulated to mimic the mechanical properties of the abdomen using the parameters in Table I. The grid dimensions were 11.6 cm in the x-direction and 8.7 cm in the z-direction, with a grid spacing of 0.0755 mm in both directions. The simulated aorta was implemented at a depth of 4.5 cm below the transducer surfaces and had an inner diameter of 7.5 mm with a uniform wall thickness of 1.7 mm [ Fig. 2(a)]. Speckle was simulated in the soft tissue and the vessel wall by adding Gaussian white noise to both the speed-of-sound and mass density maps, with a standard deviation of 2.5 m/s and 3.5 kg/m 3 , respectively.
To accurately simulate the same curved array transducers as used in the experiments, each of the 128 elements was described individually, including the exact element positions, an element width of 1.12 mm, and the tilt angle. Moreover, the simulated curved array transducers had the same center frequency of 3.7 MHz, with a transmit and receive bandwidth of 60%. Since a binary definition of such tilted source and sensor pixels would lead to errors in the element width and to staircasing on the discretized grid, an analytical description of the element surfaces was band-limited and then discretized [25]. This resulted in weights covering the entire image grid since a band-limited element of a certain dimension corresponds to a convolution of a rect function with a sinc function. However, since the weights far from the element location were very small, only the weights inside a binary source mask could be used to spread the input pressure signal of the transducer element over the corresponding pixels. The weight masks were applied within a radius of 10 pixels around each element location. This limited the computation time and the memory required to perform the simulation but did not have a large influence on the image quality, compared to not using a binary source mask. In case of pixels belonging to multiple elements, the assigned pressure input signals for those pixels were superposed.
The relative transducer positions and interprobe angle in the simulations mimicked the geometry used in the experiments. However, since relative probe positions were slightly different in each experimental dataset, they were not exactly the same. The exact transmit sequence from the experiments was directly used in the simulations as well, with the same frame rate of 180 Hz. The sampling frequency was 30 MHz.
Simulating deformation of the aortic wall required a displacement map to be applied to the input medium parameters in the mass density map, attenuation map, and the speed-ofsound map. An equal expansion in all directions (isotropic expansion) was simulated that could easily be compared to a ground truth (GT) to quantify accuracy of the strain estimates. The simulated radial displacements at the vessel lumen-wall border mimicked the systolic phase of aortic inflation that was modeled in 30 frames. The radial displacements were derived from the ex vivo study, by multiplying an input circumferential strain ε circ , with the radius of the vessel at the lumen-wall border r b . In this way, expansion of the aorta ex vivo could be mimicked since radial displacement dominated the deformation pattern, whereas circumferential displacement was negligible. Due to the simulated axisymmetric shape of the aorta, the radial displacement u r decreased toward the outer border with 1/r when the radius r of the position in the grid was larger than r b [26]. When the radius of the position in the grid was smaller than r b , the displacements decreased proportional to r . The total radial displacement field [ Fig. 2(b)] was described by (1), where C is a constant equal to r 2 b ε circ .
It should be noted that in the definition of the displacement fields, the assumption was made that the vessel wall and the surroundings were of the same material. In reality, the differences in stiffness between the vessel wall and the surroundings as well as the presence of the spine and the axiasymmetric shape of the aorta cause the deformation to be slightly different.
The local displacements were achieved by band-limited interpolation of the k-Wave grids, which could be efficiently implemented using zero padding, after which linear interpolation was used to look up the correct field values at the displaced positions. After the simulation was performed, the pressure field was multiplied with the weights allocated to each pixel within every element mask. The pixel signals within an element mask were summed, to obtain the resulting RF signal for each element. No additional electrical noise was added to the RF data. The channel RF data were reconstructed according to the method described in Section II-C1.

B. Experimental Setup
In this study, four porcine thoracoabdominal aortas were used, originating from healthy young pigs from the local slaughterhouse. These aortas were 'residual waste' that could be used for research purposes. They were individually stored in a phosphate-buffered saline solution (PBS) in a −20 • freezer until usage. The excessive tissue was removed, after which side branches were closed firmly using sutures. After preparing the aorta, it was transferred to an in vitro mock-circulation setup to mimic physiological conditions of aortic inflation, which is described in more detail in [18], [27], and [28].
To mimic the abdomen, the container in which the aorta was embeddded was filled with an 8 wt% gelatin solution (300 bloom, FormX, Amsterdam, The Netherlands) that also contained 2 wt% silicon carbide (∼400 mesh particle size, Sigma-Aldrich Chemie GmbH, Steinheim, Germany) as an acoustic scatterer. A 3-D printed solid spine (Polyamide 12, Shapeways Inc., Eindhoven, The Netherlands) was positioned close to the posterior side of the aorta. The presence of the spine influenced the deformation of the aorta, as it does in an in vivo situation. The aortas were imaged using a Vantage 256 system (Verasonics, Kirkland, WA, USA) and two C5-2v curved array Verasonics probes, with a center frequency of 3.7 MHz. The distance between the transducers' surfaces at its closest point was about 1.5 cm, with a relative angle between the beam axes of each transducer of approximately 85 • . Care was taken with the positioning of the two US probes, to make sure that they were properly aligned in the same imaging plane. The position of the aortic center was located approximately 4 cm below the transducer surface. An interleaved acquisition scheme was designed ( Fig. 1), in which each transducer consecutively performed 11 transmits over different steering angles between −10 • and 10 • with respect to the element normals. The high frame rates required for ultrafast acquisition were realized by transmitting diverging waves using all elements in the transducer. The time between the transmit events was set to 0.25 ms to prevent interference between the acquisitions. As a result, the frame rate for each transducer was equal to 180 Hz. The backscattered RF signals were then received by all the elements in both probes, therefore creating a dual-receive setup. The raw channel data were extracted and stored, from which the individual four signals T 1 R 1 , T 2 R 1 , T 1 R 2 , and T 2 R 2 were reconstructed offline. SP MS, MP MS, and MP BS image acquisition configurations could thus be obtained from the same acquisition scheme.

C. Multiperspective Bistatic Strain Imaging
MP BS strain imaging was performed on both the simulated and experimentally acquired US data. An overview of the workflow is presented in Fig. 3.
1) Reconstruction: The backscattered RF signals were individually reconstructed by a commonly used delay-and-sum (DAS) beamformer. However, to allow for the reconstruction of the trans-probe data, the beamforming process needed to be extended to the current dual-transducer setup, which is explained by the methods in [29]. The transmit events were coherently combined by summation over the transmit angles to generate each final beamformed signal that was stored as in-phase quadrature (IQ) data for further processing in MATLAB 2019a (The MathWorks, Natick, MA, USA).
As a final step, a coordinate system had to be chosen to reconstruct the beamformed IQ data on. For this study, two different coordinate systems were used for the reconstruction: the conventional sector grid in curved array US imaging and a polar grid with respect to the vessel's local coordinate system (Fig. 4). The conventional sector grid allowed for the estimation of displacements along the US beam direction (axial direction) of each transducer, after which the axial displacements from the two different coordinate systems could be combined by means of a radial projection, which will be explained later. The sector grid consisted of 256 lines in the lateral direction (two lines per pitch) and 1961 samples in the axial direction with a pixel spacing of 0.051 mm [(1/8)λ] and comprised the entire recorded field of view. The US signals were reconstructed on the sector grid corresponding to its receiving transducer. Apodization by means of a Hanning window was applied to the beamformed US signals that were reconstructed on the conventional sector grid to reduce the effect of artifacts on the displacement estimation.
The use of the polar reconstruction grid came from the idea to determine aortic wall motion in the direction of the presumed displacement, which is radial. The radial displacements could be estimated in a single step using coherently compounded image data, without the need for a radial projection, interpolation, and displacement compounding (Fig. 3). The polar grid was defined by a center point that was placed in the middle of the aorta [ Fig. 4(b)]. In the ex vivo study, this center point was estimated for every frame individually by applying a Star-Kalman algorithm for automatic detection of the lumen-wall border [30]. For each degree in the circumferential direction, 15 mm of RF data in the radial direction was reconstructed, with a radial pixel spacing equal to the axial pixel spacing in the conventional sector grid. Updating the midpoint for every frame was necessary for the local radial direction to stay aligned with the expansion of the aorta. After all, the aorta did not expand equally in all directions because of the 3-D printed spine inserted below the aorta. The simulated US data were also reconstructed on the polar grid, however, updating the midpoint was not required because of the isotropic expansion of the simulated vessel.
For accurate calculation of the transmit and receive delays in the experimental trans-probe data (T 1 R 2 and T 2 R 1 ), it was important to know the exact location of the two probes. Since the surface of the other probe was visible in each reconstructed image, it was possible to define a rigid transformation between the two probes using the T 1 R 1 and T 2 R 2 signal. This transformation was found by using the circular Hough transform for automatic detection of the transducer surfaces, a method adapted from [18]. For this purpose, the T 1 R 1 and T 2 R 2 signals were reconstructed on a sector grid containing one line per pitch, after which a Canny edge detection was applied to the envelope detected signals. Connected structures containing less than 200 true pixels in the edge detected images were removed, making sure that the final image contained two main edges corresponding to the pixel intensity changes between the visible transducer surface and the surroundings. After applying the circular Hough transform to these images, this resulted in two point locations of the highest intensity. The algorithm then automatically selected the midpoint of the curvature corresponding to the inner side of the edge-detected reflection of the transducer surface since the reflection originated from the acoustic lens. Hence, the element positions were located closely behind it. The resulting rigid transformation, defining the relative transducer element positions, was then used to reconstruct the T 1 R 2 and T 2 R 1 signal. The mean registration error over all four datasets, equal to the difference in Euclidean distance between the detected origins of the transducers in the T 1 R 1 and T 2 R 2 signal, was equal to 0.17 mm. Although this registration error does provide an indication of the mismatch between transducer positions, the actual mismatch between image features in the image is much smaller since a slight offset in transducer element location will still lead to almost the same path lengths to a certain pixel in the image.
2) Image Fusion: To combine the image data from the two coordinate systems using the conventional sector grid, interpolation on a Cartesian grid was required, this was performed using a scattered interpolant function. The image data received by transducer 2 were transformed before interpolation using the rotation and translation found by the circular Hough transform since the Cartesian grid was aligned with transducer 1. The isotropic pixel spacing in the Cartesian grid was equal to the pixel spacing in the axial direction of (1/8)λ. The actual fusion was performed by incoherent compounding, which does not take the phase of the signal into account to avoid errors introduced after interpolation. To retrieve the most useful information from each US signal, masks were designed, which were aligned with the US beam direction of each receiving transducer. The masks M were created using a cosine function [18]: Here, θ i, j represents the angle between each pixel (i, j ) in the interpolated images of each transducer and the line between the receiving transducer origin and the midpoint of the aorta [ Fig. 5(a)]. The masks were normalized by dividing each mask with the sum of all used masks before being used for image compounding. Note that this could slightly shift the mask shown in Fig. 5(a), depending on the degree of overlap, since normalization will make sure that the sum of all masks is equal to 1. Two methods were created to perform weighted incoherent image compounding, which were defined as two-part weighted compounding (2PWC) and fourpart weighted compounding (4PWC). In 2PWC, the singleand trans-probe data of each receiving transducer were first coherently combined on the grid of the receiving transducers, before being weighted and summed incoherently using two masks. In 4PWC, all four signals in MP BS imaging were weighted individually. For this purpose, a third mask was created that was aligned with the location of the trans-probe reflections on the aortic wall that occurred in the wall segments between the transducers. The way to obtain the fused image I , using 2PWC or 4PWC, is shown by the following equations: whereM P1 andM P2 are, respectively, the normalized masks aligned with probes 1 and 2. Moreover, the signs c + and i + indicate a coherent or incoherent summation, respectively.
In (4) for 4PWC,M TP corresponds to the mask used to weigh the trans-probe signals. The factor (1/2) was used to make sure that the brightness of the reflections across the circumference of the aorta remained constant since the reflections of the T 1 R 2 and T 2 R 1 signal on the aortic wall will always appear at the same location for arbitrary probe positions, because of reciprocity. For reference, an MP MS image was also obtained, for which 2PWC was applied, but only the single-probe signals (T 1 R 1 and T 2 R 2 ) were taken. The experimental and simulated US data reconstructed on the polar grid could be combined by means of coherent compounding since all four signals in MP BS imaging were reconstructed in the same coordinate system. Hence, no interpolation was required. Coherent 2PWC and 4PWC could then be reformulated by: In the simulations, the relative transducer positions were exactly known, and therefore, the coherently compounded image was created by summation of the individual signals to fully benefit from the increased aperture.
After compounding, the resulting B-mode images were obtained after log compression to a decibel (dB) scale and normalization by subtraction of the maximum pixel intensity in the dataset, similar to conventional US imaging.
3) 2-D Displacement Estimation: A vessel segmentation was made on the fused image by manually selecting points on the lumen-wall border in a clockwise manner. The corresponding points on the outer border were found by extrapolating the selected points with a uniform wall thickness of 1.7 mm. The wall thickness was determined from measurements in the animal samples and estimated to be 1.7 mm ± 0.2 mm across the four vessels. Linear interpolation between the inner and the outer walls resulted in a mesh of 5 × 101 points (radial × circumferential direction).
The Cartesian coordinates of the wall segmentation obtained from the incoherent compounded image were transformed back to the coordinate systems of each individual transducer since displacement estimation works better with the availability of RF image data in each signal, before fusion. In the polar grid, displacement estimation could directly be performed on the fused image since it was obtained by coherent compounding, which preserved phase information.
Displacement estimation was performed by speckle tracking, using a 2-D coarse-to-fine displacement estimation algorithm, developed by Lopata et al. [31]. The coarse displacement estimation was performed on the envelope signal, which is the absolute value of the IQ data, after which the fine displacement estimation was performed on the RF signal, which is the real part of the IQ data for the axial reconstruction spacing used [(1/8)λ]. Displacement estimation in the conventional sector grid was performed using a kernel size of 2.6 mm × 4.5-5.9 mm (axial × lateral) for the coarse displacements and a kernel size of 0.8 mm × 2.1-2.7 mm for the fine displacements. The search area was defined to be 3.0 mm × 4.9-6.5 mm for the coarse displacements and 1.0 mm × 2.5-3.2 mm for the fine displacements. The lateral kernel sizes were larger than the axial kernel sizes because of the limited resolution in this direction. In the axial direction, image resolution was much better and the fine kernel size could be set to approximately half the wall thickness to be able to measure radial strain. The change in wall thickness over one cardiac cycle only covered approximately four RF data points, with a spacing of (1/8)λ, in the axial direction. However, displacement estimation is typically performed using overlapping windows. As a result, sufficient data were available for radial displacement estimation. Subresolution displacement estimations were improved by interpolation of the cross correlation function, as described in [31] and [32]. As a last step, the axial and lateral displacements were filtered with a median filter of 0.6 mm × 4.5-5.9 mm (11 × 11 pixels).
After displacement estimation, the axial displacements u ax were radially projected according to u ax / cos(θ ) to obtain the radial displacement field u rad [33]. Here, the same definition of θ is used, as shown in Fig. 5(a). The lateral displacements were not considered because of the poor lateral resolution and the lack of phase information in this direction.
The radial displacement fields from the two different coordinate systems were then interpolated on the Cartesian grid, after which displacement compounding was performed similar to the way the images were combined, by using angular weighted masks. The angular weighted normalized masksM * were also defined using a cosine function similar to (2), and however, they contained cut-out regions between the 70 • and 110 • angle with respect to each receiving transducer in order to prevent infinite radial displacements [18]. The compounded radial displacement field was obtained by applying either 2PWC or 4PWC based on the methodology chosen to do the incoherent image compounding. Equation (7) explains 2PWC for the radial displacements, in which the radial displacement fields for each probe u rad,P1 and u rad,P2 were combined using two masks. Similarly, 4PWC for the radial displacements is summarized in (8). Because of the cut-out regions, different masks were now obtained to weigh the displacements from the T 1 R 2 and T 2 R 1 signals (M * T 1 R 2 andM * T 2 R 1 , respectively).
In the polar grid, the total radial displacement field could be obtained using a single displacement estimation step on the coherently compounded image data (Fig. 3). The kernel sizes used for estimating the displacements were 2.1 mm × 41 • (radial × circumferential) for estimating the coarse displacements and 1.3 mm × 25 • for estimating the fine displacements. The corresponding search area in the radial direction was equal to 2.7 mm for the coarse displacements and 1.6 mm for the fine displacements. In the circumferential direction, the maximum allowed displacement was set to 0, thereby not estimating the displacements in this direction, Fig. 6. Illustration of the aorta at end-systole. In the polar grid with respect to the vessel's local coordinate system, a and b represent the vectors between a point in the vessel wall segmentation, and the reference midpoint (frame 1) or current midpoint, respectively.
which mainly caused a drift in the displacement estimation results, and motion is predominantly radial. Note that setting circumferential displacement to 0 did not mean that resulting radial and circumferential strain maps were the same since the displacement estimation did not enforce incompressibility. After obtaining the total displacement field, the previously made segmentation of the aortic wall was used to perform motion tracking of the aorta.

4) Strain Calculation:
The radial and circumferential strains were estimated from the coordinates of the vessel wall segmentation after motion tracking. For each frame, the displacements of the tracked coordinates were calculated with respect to the first frame, after which a 2-D least-squares strain estimator described in [34] was applied. A strain kernel with a fixed size of five radial × nine circumferential mesh points (1.7 mm × 4.6 mm) was used in the SP ultrafast acquisition. In the MP acquisitions, a smaller strain kernel of five radial × five circumferential mesh points (1.7 mm × 2.6 mm) was used to increase the strain resolution, at the possible cost of accuracy.
In the experimental results using the polar coordinate system, a correction was performed on the tracked x-and zcoordinates of the vessel wall segmentation (x c and z c ) to take the displacement of the midpoint into account, that occurred between each reconstructed frame: This was necessary to obtain the correct relative displacements with respect to the first frame for the calculation of the strains. In (9) and (10), (| a| − | b|) represents the difference in vector length between each tracked point in the vessel wall segmentation and the midpoint in the current frame, compared to the reference midpoint of the first frame. Moreover, θ 0 represents the angle between each point in the vessel wall segmentation and the location of the reference midpoint. The vectors a and b as well as the angle θ 0 are shown in Fig. 6.

5) Analysis of Results:
The analysis of the results was performed in a regional manner as a function of the vessel circumference. The aorta was divided into eight regions (R1-R8), where each region had a span of 45 • [ Fig. 5(b)].
The accuracy of the simulations was quantified in the last simulated frame, corresponding to the end-systolic phase of aortic inflation, by comparing the strain estimation after speckle tracking with the GT strains that were calculated analytically. The GT circumferential strain ε circ,GT was equal to u r /r [35], with u r defined in (1), and because of the isotropic expansion, the GT radial strain ε rad,GT was equal to -ε circ,GT . Using these GT strains, the average relative error (ARE) was calculated over all n strain estimates in the inner layer of the vessel wall segmentation in R1-R8: It expresses the ratio between the error of the strain estimation and the total GT strain.
The results of the ex vivo study were analyzed by means of image quality, tracking performance, and strain estimation precision. As shown in Fig. 1, the added benefit of the transprobe data in image quality and strain estimation is visible in the wall regions between transducers that correspond to R4 and R8 in Fig. 5(b). Hence, each measure was locally quantified in these regions. The global results were also reported, in which the entire circumference of the aorta was considered. In this analysis, the results at sides of the aorta with a span of 40 • were excluded in the MP image acquisition configurations. These regions benefit the least from the current dual-transducer setup, due to the interprobe angle and luminal artifacts at these positions.
To assess the change in image quality as a result of weighted image compounding, the contrast-to-noise ratio (CNR) between the vessel wall and the lumen was computed. The CNR is defined by where μ w and μ l denote the mean envelope pixel intensities in the aortic wall and the lumen, respectively, and σ 2 w and σ 2 l their variance. For the analysis of the tracking performance, the mean drift error (ME) was used. This error was estimated between the n estimated positions of the middle wall layer in the vessel segmentation before and after inflation. With a higher precision of the estimated displacements, the ME decreases since the aortic wall comes back to the same position after inflation. The ME is defined by where x i and y i denote the positions in the aortic wall before systole (bs) and end-diastole (ed).
To be able to assess the accuracy of the tracking, the Dice similarity coefficient (DSC) [36] was calculated at end-systole. For this purpose, the aorta was manually segmented at the endsystolic phase of aortic inflation with a smaller wall thickness of 1.5 mm, and the segmentation was compared against the tracking result. The wall segmentation and the tracking result were converted into binary masks to calculate the DSC.
The strain estimation precision was quantified by the elastographic signal-to-noise ratio (SNRe) using the mean (μ ε ) and  the standard deviation (σ ε ) of the radial and circumferential strain estimates as measured in the middle layer of the segmented wall: SNR e = 20 log 10 μ ε σ ε . (14) This value is expected to be larger in the regional analysis of the results, compared to the global results since a certain variation is also present in a completely correct estimation due to the spine's influence, which will lower the SNRe.

A. Simulations
In Fig. 7, the strain results obtained from SP MS, MP MS, and MP BS simulations are shown, as well as the GT strain. It is shown that the strains in the MP BS simulation were more homogeneous compared to the SP MS and MP MS simulation and resembled the GT strains most out of the three image configurations. Overestimations of radial strain along with underestimations of circumferential strain remained visible at the 3 and 9 o'clock positions since almost no backscatter could be received from these wall locations.
In Fig. 8, the regional accuracy of the strain estimates from the simulations are compared between the image acquisition configurations. It gives a good overview of how each signal helps to make the strain estimation more robust, and it can directly be related to the illustration in Fig. 1. Comparing SP imaging to MP imaging, the addition of the second transducer resulted in a better accuracy in both radial and circumferential strain in all regions except R6. Using both transducers in dualreceive mode by MP BS imaging, the ARE in R4 was further decreased by 41% and 78% in the radial and circumferential directions, respectively, compared to MP MS imaging. In R8, the radial ARE decreased by 84% and the circumferential ARE by 78%.
B. Experiments 1) Image Quality: In Fig. 9, the appearance of the aorta between the image acquisition configurations in the experimental study is shown. It can be seen that by the addition of the trans-probe data in MP BS imaging, an additional reflection on the aortic wall appeared, which occurred in between the transducers. This improved the visibility of the wall-lumen interface and the wall-gelatin interface at these positions.
The mean and standard deviation of the CNR of every aorta (N = 4) are listed in Table II. Here, the differences between 2PWC and 4PWC strategies are also shown. In both grids, it can be seen that by using MP BS imaging, the mean CNR in vessel wall regions between transducers increased compared to MP MS imaging, by about 8 dB in R4 and 5 dB in R8. In the global results, an increase of about 2 dB was observed. The small differences that exist between the two grids can be explained by the apodization, which was applied to the beamformed US data reconstructed on the sector grid, but not to the beamformed US data reconstructed on the polar grid.
2) Tracking Performance: The tracking performance was quantified by the ME and DSC in Table II. Comparing the SP MS results from the sector grid with those from the polar grid, it is shown that displacement estimation on the polar grid resulted in less lateral drift across the circumference of the aorta, with a reduction in global ME from 1.1 to 0.2 mm on average. This improvement in 2-D displacement estimation accuracy in SP images, in general, also leads to an improvement of the tracking performance in MP images on the polar grid, compared to the results on the conventional sector grid. An exception was observed in MP MS images in R8, where the mean DSC decreased from 0.86 to 0.78 and the ME increased from 0.1 to 0.2 mm.
In the comparison between MP MS and MP BS imaging, different outcomes were observed within the sector grid and the polar grid. Within the sector grid, the tracking performance of MP BS imaging, according to the DSC and ME, decreased in R8 compared to MP MS imaging. Here, the ME increased by 0.1 mm on average, and the mean DSC decreased from 0.86 to 0.54 or 0.68, using 2PWC and 4PWC, respectively. This was not reflected by the results from the polar grid, where the ME in R8 decreased by 0.1 mm on average and the mean DSC increased from 0.78 to 0.83 using both MP BS strategies.
An increase in tracking performance was also observed in R4 on the polar grid, with a decrease in ME by 0.1 mm using MP BS imaging compared to MP MS imaging.
3) Strain Imaging: In Fig. 10, the strain results between image acquisition configurations on the sector and polar grid  are shown. It can be seen that these experimental results in Fig. 10 are in agreement with the simulation results in Fig. 7.
In the MP MS image acquisition configuration, there was an overestimation in radial strain in the vessel wall regions between transducers. Using MP BS imaging, this overestimation did not occur anymore. Moreover, the circumferential strain became more homogeneous in these regions. Between sector and polar grids, results on the polar grid were also observed to be more homogeneous in the posterior wall. Strain estimation precision was quantified by the SNRe circ and SNRe rad for every aorta and the means and standard deviations are summarized in Table II. Between the sector and polar grid, global SNRe circ using MP BS imaging was observed to be between 2 and 4 dB larger on the polar grid. Moreover, a large increase of mean SNRe circ was observed in R8 by 10.6 dB using 2PWC and 5.2 dB using 4PWC. MP results in SNRe rad on the polar grid were also generally higher than the results on the sector grid, although no large differences were observed in the global SNRe rad . An exception occurred in R8 using 2PWC where the mean SNRe rad was 1.7 dB smaller on the polar grid, compared to the sector grid.
In the comparison between MP MS and MP BS imaging within the sector grid, it was observed that the mean SNRe circ decreased using MP BS imaging in R8 by 12.3 dB using 2PWC and by 6.7 dB using 4PWC. This was not observed in R4, where the mean SNRe circ increased using MP BS imaging by 7.5 dB using 2PWC and 7.7 dB using 4PWC. There was no increase observed in the average global SNRe circ using MP BS imaging compared to MP MS imaging.
Looking at the comparison between MP image acquisition configurations within the polar grid, an increase in mean SNRe circ was observed in both R4 and R8. In R4, this increase with respect to MP MS imaging was equal to 6.2 and 5.5 dB using 2PWC and 4PWC, respectively. In R8, an increase of 3.8 dB using 2PWC and 4.0 dB using 4PWC was observed. A minor decrease in mean global SNRe circ of 0.9 dB was observed using 2PWC, but using 4PWC, the mean global SNRe circ was increased by 0.3 dB.
Comparing the SNRe rad between MP MS and MP BS imaging within each grid, it was observed that the mean global SNRe rad of MP BS imaging was decreased, compared to MP MS imaging in both grids. Moreover, there was a large variance observed between the four datasets, which was related to the varying effect of luminal artifacts on the estimated strain. However, looking at the vessel wall regions between the transducers specifically, an increase in the mean SNRe rad was observed. In the sector grid, the mean SNRe rad increased by 18.8 and 17.2 dB using 2PWC and 4PWC, respectively, in R4. In R8, this increase was equal to 29.7 dB using 2PWC and 18.5 dB using 4PWC. In the polar grid, it was observed that the increase in mean SNRe rad was slightly lower at 10.3 dB using 2PWC and 8.9 dB using 4PWC in R4. In R8, the mean SNRe rad increased by 13.9 and 13.4 dB using 2PWC and 4PWC, respectively.

IV. DISCUSSION
In this study, ultrafast MP BS imaging was applied to aortic strain imaging in both US simulations and ex vivo experiments involving porcine aortas. It was shown that using two curved array transducers in an MP BS image acquisition configuration, better images can be made of the abdominal aorta, visualizing a larger part of the vessel circumference, compared to SP MS or MP MS imaging. This large increase in the wall-lumen contrast enabled the estimation of more accurate wall motion and strain in almost the entire aortic wall. Moreover, the introduced polar reconstruction allowed for a more tailored and direct displacement estimation of the aortic wall which, to the best of our knowledge, has not been done before. It resulted in a high strain estimation precision despite limitations such as the low center frequency used and the small wall thickness.
In the US simulations, it was found that the accuracy of both radial and circumferential strains was substantially increased using MP BS imaging, compared to SP MS and MP MS imaging in vessel wall regions between transducers. Comparing the accuracy in radial and circumferential strain, it was found that the ARE of the radial strain estimates was much larger. This can be explained by the relatively small wall thickness since small distortions in the displacement field can lead to larger differences in the radial strain. In the US simulations, these erroneous displacements were mainly caused by the reverberation artifacts in the upper part of the vessel lumen. The reverberations occur specifically in this phantom since they were caused by the strong reflections on the inner and outer aortic walls. Because these artifacts occur in the axial direction of each US probe, it explains why a lower error in R6 was observed using SP MS imaging, compared to MP MS and MP BS imaging. This was an interesting finding since it suggests that the strain results at the sides of the aorta in MP BS imaging could be improved by only using the signal information from the SP MS images not containing this artifact at these locations. In future research, the simulations could be made more realistic by means of finite element modeling that can incorporate other organs, e.g., a spine, to mimic the more inhomogeneous deformation that occurs in an aorta in vivo. Moreover, electrical noise could be added to the channel RF data. While this could have contributed to make the US simulation more realistic, we chose to omit noise to fully isolate the effect of multiple views from the impact of noise.
In the experimental results, differences were observed between the sector and the polar grid. It was shown by comparing Fig. 10(c) and (f) that strain results on the polar grid were more homogeneous in the posterior wall. Inaccuracies in the strain result from the sector grid could be caused by badly estimated displacements in individual single-probe or transprobe US data on positions that lack good image information. Even though an angular-based weighted compounding method was implemented, inaccurately estimated displacements could still end up in the final compounded radial displacement field. In the polar grid, radial displacements could directly be estimated from the MP BS US data, which makes the displacement estimation less prone to inaccuracies.
The experimental results were limited by the lack of a GT, a common limitation, and thus, the differences in the strain results could only be analyzed in terms of strain estimation precision. It is important to realize that an improvement in strain estimation precision does not necessarily imply an improvement in the accuracy as well. This is especially visible in R8 on the conventional sector grid, where the mean SNRe circ decreased using MP BS imaging, compared to MP MS imaging. At the same time, the increase in mean SNRe rad in MP BS imaging was quite high in that region at 29.7 dB for 2PWC and 18.5 dB for 4PWC. This anomaly was caused by the beamwidth artifacts next to the reflections on the aortic wall, introduced by the large pitch of the transducer used. These artifacts affected the block matching algorithm, causing both the inner and outer walls of the vessel segmentation to expand more than the actual trans-probe reflection from the aortic wall in R8. The beamwidth artifacts were more reduced using 4PWC compared to 2PWC since 4PWC was a sum of four images instead of two and was, therefore, more robust to individual distortions. In the polar grid, SNRe circ was improved using MP BS imaging in R8 compared to MP MS imaging, which was also expected. Moreover, no differences between 2PWC and 4PWC were observed. Because of the improved image quality of the coherent compounded MP BS US image compared to that of the individual signals, displacement estimation in the polar grid also seems to be less affected by image artifacts. In future research, beamwidth artifacts could be further reduced by using more transmit angles, at the cost of frame rate, which always remains a tradeoff. In vivo, their effect is expected to be less since a stronger speckle background is expected together with more gradual changes in acoustic impedance.
Comparing MP MS and MP BS imaging in the polar grid, a local increase in both SNRe circ and SNRe rad in R4 and R8 was found, which corresponded to the expectation based on the increase in wall visibility. However, despite these local improvements, global SNRe rad was not improved. This could be explained by the artifacts in the lumen, which were also present in the single-probe data, but more predominant and partially overlapping the vascular wall in the trans-probe data, which affected the radial strain results in R3 and R5. Again, this effect was larger using 2PWC compared to 4PWC. These artifacts were also linked to the large pitch of the transducer used and are explained in [37]. In future research, adaptive weighting based on the image content may be investigated to improve both the image quality and the resulting wall strains even more since one can freely choose which region to use where to obtain the final MP BS result.
The accuracy of the radial and circumferential strain in the polar grid relies on a good midpoint estimation, which could be determined automatically for every frame using the Star-Kalman algorithm. The location of the reference midpoint is especially important to obtain the right strains. A small deviation from the right position will lead to a slight overestimation of circumferential strain in vessel wall locations that are closer to the reference midpoint than in reality and vice versa. A deviation in any other midpoint from the right location could influence the vessel wall tracking since the direction of expansion is then not aligned with the grid. Polar beamforming may also be interesting for use in conventional single-probe techniques. In the future application of strain estimation in AAAs, the polar grid may be adapted according to the vessel's local coordinate system since AAAs come in all forms and shapes and are not necessarily circular. Fitting the ellipsoid shape of an AAA could be done using a prolate spheroidal coordinate system.
This study was performed in simulations and an ex vivo setup. In vivo studies need to be performed to investigate the feasibility of dual-transducer acquisition in the clinic. The use of multiple probes will potentially increase the operational difficulty for the person performing the scan, particularly because of the importance of probe alignment in the same imaging plane. However, this may be solved by the use of probe holders that will allow the operator to hold multiple probes in one hand or 3-D ultrasonography. The performance of the method in vivo will be dependent on the interprobe angle that can be reached, as well as the depth of the aorta, that may be deeper than the depth that could currently be reached in our experimental setup (approximately 4 cm). These factors influence how much of the wall circumference will be well visible by specular reflections. Nevertheless, image and displacement compounding of MP data can still improve image quality and strain estimation compared to conventional single-probe techniques, even with relatively small interprobe angles.
Coherent compounding of multitransducer US data may also become more challenging in vivo in the presence of aberrations, caused by anisotropic speed-of-sound. This problem would particularly present itself in transabdominal US imaging of obese patients. However, in a study by Peralta et al. [38], it was shown that coherent bistatic imaging offers opportunities to correct for phase aberration errors in case of aberration clutter by the cross correlation of backscattered echoes from individual receive elements. These, and other methods that can correct for aberrations, will be investigated in future research.
The current registration algorithm was based on the detection of the probe surfaces and enabled an accurate reconstruction of trans-probe data, despite a small mismatch between detected probe positions. Since the visibility of the probe surface of the other transducer may be limited in vivo, the current registration algorithm may be adapted in future research to include other well visible in vivo structures, such as the spine as shown in [14].
Finally, risk management of AAAs could strongly benefit from a translation of the dual-receive acquisition into 3-D US imaging to retrieve a more complete assessment of the AAA geometry and its mechanical behavior. The 3-D US imaging captures out-of-plane motion and removes the requirement of imaging the exact same plane with both transducers. Therefore, the fixed transducer setup may likely be changed in the future to allow for a more free-hand operation. Previous research has shown that the 3-D speckle tracking of AAAs is feasible [5], as well as 3-D ultrafast imaging in vivo [39]. Future research could improve these acquisition methods even more by implementing 3-D MP BS imaging that will contribute to accurate 3-D biomechanical characterization of AAAs.

V. CONCLUSION
This study contributes toward ongoing efforts in improving the characterization and risk management of AAAs. Ultrafast aortic strain imaging was performed using US simulations and ex vivo experiments involving porcine aortas and two curved array transducers. The use of MP BS imaging demonstrated an improvement in image quality with respect to SP MS and MP MS imaging by providing more information of the vessel wall circumference in aortic wall regions between transducers. This improved the accuracy of radial and circumferential strain estimates in US simulations and resulted in more homogeneous strain results in ex vivo experiments. Moreover, this study introduced a polar grid with respect to the vessel's local coordinate system that has shown to return accurate displacement estimations in aortas, despite limitations such as the low center frequency used and the relatively small wall thickness. Future research will investigate the applicability of these methods in vivo and using 3-D ultrafast US imaging.