Large elements and advanced beamformers for increased field of view in 2-D ultrasound matrix arrays

Lar ge elements and adv anced beamformers for increased field of vie w in 2-D ultrasound matrix arrays Mick Gardner , Michael L. Oelze Electrical and Computer Engineering, University of Illinois at Urbana-Champaign, 306 N Wright St, Urbana, IL, 61801, USA Abstract Three-dimensional (3D) ultrasound promises various medical applications for abdominal, obstetrics, and cardiov as- cular imaging. Howe v er , ultrasound matrix arrays have extremely high element counts limiting their field of view (FO V). This w ork seeks to demonstrate an increased field-of-view using a reduced element count array design. The approach is to increase the element size and use adv anced beamformers to maintain image quality . The delay and sum (D AS), Null Subtraction Imaging (NSI), directional coherence factor (DCF), and Minimum V ariance (MV) beam- formers were compared. K-wa ve simulations of the 3D point-spread functions (PSF) of NSI, DCF , and MV display reduced side lobes and narrowed main lobes compared to D AS. Experiments were conducted using a multiplexed 1024-element matrix array on a V erasonics 256 system. Elements were electronically coupled to imitate a larger pitch and element size. Then, a virtual large aperture was created by using a positioning system to collect data in sections with the matrix array . High-quality images were obtained using a coupling factor of two, doubling the FO V while maintaining the same element count in the virtual large aperture as the original matrix array . The NSI beamformer demonstrated the best resolution performance in simulations and on the large aperture, maintaining the same reso- lution as uncoupled D AS for coupling factors up to 4. Our results demonstrate how larger matrix arrays could be constructed with larger elements, with resolution maintained by adv anced beamformers. K eywor ds: Beamforming, ultrasound matrix arrays, 3D ultrasound, large elements, null subtraction imaging 1. Introduction Three- dimensional (3D) ultrasound is an emerging and attractive technique in the medical imaging com- munity because it can provide full-volume views of a region of interest while being portable, safe, and capable of real-time imaging. V olume imaging with ultrasound has various medical applications in obstetrics, muscu- loskeletal, and cardio vascular imaging to name a few examples [1]. Howe ver , 3D ultrasound imaging usually requires a fully populated 2D matrix array . Such arrays hav e a prohibitively high number of array elements, e.g. 1024 elements in a 32x32 system. W ith an array pitch typically around one w av elength or less to pre vent grat- ing lobes during beam steering, the field of view (FO V) of such probes can be sev erely limited. Also, most ul- trasound scanners are not equipped with enough input channels to recei ve data from such a high element count, so multiplex ers are often used to switch between sec- tions of the array on di ff erent transmit / receive e vents [2]. Howe ver , in this case, frame rates are limited by the need for multiple transmit / receive ev ents to acquire all aperture data. W ith a 4x1 multiplex er, an ywhere be- tween 4 - 16 transmit ev ents would be needed to acquire, for example, a single plane-wa ve angle [3]. T o handle higher element counts, so called “mi- crobeamformers” have been inte grated into the handles of the transducers [4, 5]. W ith microbeamforming, the delay-and-sum (D AS) operation is partially performed on ASICs inside the transducer handle on patches of elements. Then, only one wire per patch needs to go back to the scanner for the remaining time delay compensation. W ith this approach, probes ha ve been designed with over 9,000 acti ve elements [6]. How- ev er , the micro-beamformers introduce deviations from ideal time delays due to static focusing and quanti- zation [7]. Micro-beamformed patches often hav e a fixed focal depth to simplify the hardware, which is not always aligned with the dynamic receive focus when beamforming across patches [7]. Delay quantization is also introduced by the micro-beamformer architecture, which might use a sample-and-hold method that can only delay signals by discrete time interv als determined by an input clock [8]. These quantizations can reduce the resolution and contrast of high-frame rate images, especially for larger micro-beamformed patches [9]. Other approaches allo w the FO V to extend beyond the footprint of the array . These approaches include con vex and phased arrays [10], as well as panoramas (also called extended field-of-view or EFO V) [11]. W ith con vex and phased arrays, the field of vie w can ex- tend beyond the probe footprint due to array curvature or beam steering, respectively . Howe ver , with this ap- proach, maintaining resolution over the FO V is a chal- lenge because scan lines become more spread out with depth. W ith panoramas, the field of vie w is e xtended by translating the probe to multiple locations, acquiring an image at each location, then applying some image regis- tration method to align the di ff erent views [12, 13]. The issue with panoramas is the y can only acquire static im- ages, nev er real-time video data, because of the need to translate the probe to multiple locations. A few alternative approaches to array design allo w for larger footprints, and thus larger FO Vs, with a re- duced element count. Extensiv e work has been done to explore sparse arrays as an option for reducing the ele- ment count of 2D ultrasound arrays [14]. W ith sparsity , element spacing is allowed to increase far beyond the usual 1 / 2 - 1 wa velength in some kind of optimized pat- tern. The most serious dra wback of sparse arrays is that they ha ve reduced transmit power and SNR because of the low element count and small elements. Larger ele- ments have been used in sparse, periodic arrays used for ultrasound localization microscopy (ULM) [15, 16, 17]. W ith very large elements, the directivity of elements raises the minimum F-number of the array , reducing res- olution. For ULM, researchers overcame this issue by placing a diver ging acoustic lens ov er each element to widen the directivity [17]. A wider directivity will intro- duce grating lobes into the point-spread function (PSF) of the imaging system when the array pitch is greater than one wav elength. The lens approach is successful with ULM because ULM displays the tracks of mov- ing microb ubbles over many thousands of frames, and the microbubble localization and tracking algorithms are able to ignore grating lobes. Howe ver , for B-mode imaging, grating lobes create artifacts that can obscure other image details [18, 19], meaning a lens which al- lows grating lobes to appear w ould not be optimal. Another example of a reduced element count ar- ray design is a row-column addressed (RCA) array [20]. These are fully populated arrays where instead of addressing individual elements, entire rows or entire columns are accessed at once. This e ff ectiv ely creates two orthogonal arrays of long, line elements, reducing the element count from N x N to N + N. Howe ver , the main limitation of these arrays are that they cannot fo- cus or steer along diagonals, they can only focus or steer in a cross (rows, then columns, or vice-versa). This means that resolution will be lo wer with RCA arrays compared to fully addressed arrays which can steer and focus along both directions at once, in both transmit and receiv e. The design we propose is a periodic matrix array with large, square elements. Rather than using a lens, we also propose the use of adapti ve or non-linear beamformers for the task of re gaining resolution lost by increased el- ement directi vity . Sev eral beamformers ha ve been stud- ied for their e ff ectiv eness at improving resolution and contrast. The Null Subtraction Imaging (NSI) beam- former has been shown to greatly improv e resolution in 2D B-mode images [21], and reduce grating lobe arti- facts [22, 19]. The NSI beamformer was also recently implemented on a 2D array for 3D imaging [23]. Co- herence factor beamformers hav e often been used for improving contrast, and the y can hav e benefits to reso- lution as well. Recently , a Directional Coherence F actor (DCF) beamformer was proposed specifically for use in matrix arrays [24], which calculates coherence on direc- tional projections of matrix data, then combines direc- tions for improved resolution. The Minimum V ariance (MV) beamformer has also been shown to greatly im- prov e resolution for 1D arrays [25], but has not been extensi vely used on 2D matrix arrays due to computa- tional complexity . Therefore, we implemented MV on the same directional projections as DCF to reduce the data size [24]. The resolution improvement makes these beamformers attracti ve for solving the issues one would encounter when creating an array with a lar ger pitch and element size. The goals of this paper are to demonstrate how larger elements can lead to larger apertures without increase to element count and to determine which beamformer is best for maintaining resolution with large elements. T o test our approach of using large, square elements, we used a commercial matrix array and electronically coupled adjacent elements to act as if the y were one el- ement. T o demonstrate a lar ger aperture, we also used a positioning system to acquire data from a virtual aper- ture where data was collected in quadrants. Images were beamformed with D AS, NSI, DCF , and MV beamform- ers for quality comparison. 2 2. Background Theory This section will give a brief analysis on the e ff ects of increased element size on the beam pattern of a 2D ma- trix array . The array beam pattern is the product of tw o factors: an array factor arising from element spacing, and an element factor (often referred to as the directiv- ity) arising from element size. Mathematically , this can be expressed as B ( θ ) = H ( θ ) G ( θ ) (1) where B ( θ ) is the beam pattern, while H ( θ ) and G ( θ ) are the array f actor and element directi vity respectiv ely . The far -field directivity of an element is giv en by its spa- tial Fourier transform. A square element can be repre- sented by two orthogonal rectangle functions, as in A ( x , y ) = rect  x W  rect  y H  (2) where A ( x , y ) is the aperture function, W is the width and H is the height of the element. Because these two rectangle functions are independent, we can separate the 2D spatial Fourier transform of the aperture into two 1D spatial F ourier transforms along the x- and y- axes, leading to F { A ( x , y ) } = F x  rect  x W  × F y  rect  y H  = W sinc W k x 2 ! × H sinc H k y 2 ! (3) where F {·} denotes the spatial F ourier transform, and k x and k y are the wave numbers in the x and y directions. From Equation (3), increased W and H values create more narro w sinc functions, representing ho w larger el- ement are less sensitiv e to o ff -axis echoes. This direc- tivity raises the achie vable F-number of the array , which is giv en by [24]: F # = z a per tur e = 1 2 tan( α ) (4) where z is imaging depth and α is the -3 dB angle of the directivity . This raised F-number will reduce the lateral resolution that can be achiev ed with the array when us- ing D AS beamforming. Lastly , we are using coupled elements as an approxi- mation of large elements. A block of coupled elements will hav e gaps, called the kerf, from the sa w used to cut the piezo-electric elements. Figures 1a and 1c display representations of a coupled element (with kerf gaps) and a large element (without kerf gaps), while Figures 1b and 1d display the corresponding 2D spatial Fourier Figure 1: The left column displays representations of (a) a block of coupled elements and (c) and a large element, with corresponding Fourier transforms representing their directivities in (b) and (d). (e) is the cross section of the 2D Fourier transforms at an elev ation of 0. Note there are only small di ff erences in the side lobes, meaning cou- pled elements are a good approximation of a large element, despite kerf gaps. transforms. Figure 1c displays an azimuthal cross sec- tion of both directi vities at an elev ation angle of zero. From these figures, it can be observ ed that the kerf gaps only create small di ff erences in the side lobes of the di- rectivities, as well as nearly identical main lobes, sug- gesting that a coupled element is a close approximation of a large element despite the k erf gaps. 3. Methods 3.1. Simulations Simulations were conducted in K-wa ve [26, 27] to examine the point-spread functions of the di ff erent beamformers on a large-element array . The simulated array had 64 elements in an 8 x 8 grid, with a pitch of 1.25 mm in both lateral and elev ational directions and an element size of 1.20 mm in both directions. The trans- mitted frequency was 7.81 MHz with a medium sound speed of 1540 m / s, making the w avelength 197 µ m. The array parameters then translate to a pitch of 6.7 λ , and an element size of 6.5 λ , where λ is the wa velength. This element size yields a minimum F-number of 7.07. The plane wav e transmission sequence consisted of 13 angles, arranged in a star pattern over azimuth and ele- vation as displayed in Figure 2. 3.2. V irtual Large Apertur e A V ermon 1024-element 8 MHz matrix probe (V er - mon S.A., T ours, France) was used to create the virtual aperture. The V ermon probe is made up of four pan- els of 8x32 elements each. Element coupling on this 3 Figure 2: Diagram of transmission angles used for simulations and virtual large aperture experiment. Fewer diagonal angles (i.e. those with both azimuth and ele vation components) were included due to increased element width, and thus a narrowed directivity , on the diag- onals. probe was performed by grouping blocks of adjacent el- ements, av eraging their transmission delays, then sum- ming their received channel data. Figure 3 illustrates which elements would be coupled for a coupling factor of 4. Because of the panel separation and multiplex- ing on the V ermon probe, coupling factors were lim- ited to 1 (i.e. no coupling), 2, and 4 to avoid coupling across panels. Within each block of coupled elements, for a given transmission (e.g. plane-wave angle), the de- lays would be a veraged so that all elements in that block fire at the same time. Then on the recei ve side, the RF traces from each element in a block would be summed together without applying any time delay . That way re- ceiv ed signals superimpose like they would for a single element. By coupling on transmit and receive, each cou- pled block acted as if it was a single, lar ge element. T o create the virtual large aperture, the probe was attached to a positioning system, and programmed to mov e into 4 quadrants, as displayed in Figure 4. T rans- mission delays for 13 plane-wave angles, arranged in a star as displayed in Figure 2, were calculated for each element in the virtual aperture. These angles were cho- sen based on the directivity of elements with a cou- pling factor of 4 because the largest elements had the most limited steering capabilities. The width of such Figure 3: Example of element coupling for a coupling factor of four . Black dots represent the original element positions for the 32 x 32 array . Red boxes represent the larger elements made from coupling. Only coupling factors 1, 2 and 4 were used to divide 32x32 elements ev enly and avoid coupling across panels. elements was 6.32 λ in the lateral and ele vational di- rections, leading to a max steering angle of about 4°. Also, these elements were 6 . 32 λ × √ 2 = 8 . 98 λ along their diagonals, leading to a max steering angle of about 3° diagonally . Therefore, we chose the star pattern of transmission angles, as opposed to a full square, so that fewer diagonal angles (i.e. those with both azimuth and elev ation components) were included. W e opted to use the same angle set for all coupling factors for a more direct comparison of quality based only on element size and recei ve beamforming, ev en though a larger steering range was a vailable with less coupling. Once the transmission profiles were calculated, the probe was set in one position, fired the transmission pro- file associated with the VLA elements in that position, and collected data. Then, the probe was mov ed to the next position, fired with the corresponding transmission delays, collected data, and so on until the full aperture data was acquired. 4 Figure 4: (a) Photograph of the positioning system setup with the gen- eral purpose phantom. (b) Diagram of the 4 quadrants making up the virtual large aperture. Each colored square represents one position of the V ermon probe, and its corresponding section of the virtual aper- ture. 3.3. Beamforming 3.3.1. Delay-and-Sum Delay-and-sum (DAS) beamforming is the con- ventional method of beamforming which simply de- lays channel data based on transmit and receive dis- tances and sums the delayed data [28]. Receiv ed radio-frequency (RF) data were demodulated into In-phase / Quadrature (IQ) data before beamforming. Beamforming was done on a pixel grid with spacing 200 µ m in lateral, ele vational, and axial directions. Once time delays were applied, transmissions were com- pounded coherently and elements were summed y (  r ) = N X i = 1 M X k = 1 I Q i , k  τ (  r ,  x i )  e − j 2 π f 0 τ (  r ,  x i ) . (5) In the abov e equation, y   r  is the beamformed pixel value for pixel  r , M and N are the number of transmis- sions and elements respectively , I Q i , k is the delayed IQ data from element i and transmission k , τ   r ,  x i  is the time delay applied to the IQ data, and f 0 is the center frequency of transmission. Receiv e sub-apertures were determined using the minimum constant F-number approach described in [28], where the -3 dB point of the element directivity determined which pix els a giv en element should con- tribute to. Elements in the V ermon probe are square, and the lateral width of the elements was used to esti- mate the directivity and minimum F-number . 3.3.2. Null Subtraction Imaging Null subtraction imaging is implemented by beam- forming with three di ff erent apodizations in parallel and incoherently summing the results [21]. These three apodizations are a zero-mean (ZM) apodization, and two direct current (DC) apodizations which are o ff set versions of the ZM. The ZM apodization creates a null at 0°, or broadside to the receive sub-aperture. The DC apodizations bridge the null but create similar side lobes in the beam pattern. Once beamforming is done with the three apodizations, the null is subtracted leav- ing a narrow main lobe and low side lobes [21]. The amount of o ff set, referred to as the DC o ff set, is a tun- able parameter of NSI typically in the range 0.1-1. In our simulations and experiments, we set the DC o ff - set to 0.5. Lower DC o ff sets result in more narrow main lobes, more suppressed side lobes, reduced grat- ing lobes, and speckle that is no longer fully dev eloped [21, 22, 19, 29]. T o implement NSI on a 2D array , a directional ap- proach was taken as described in [23]. In this ap- proach, tw o zero-mean apodizations are made by taking an equal number of + 1 and -1 along row and column di- rections of the receiv e sub-aperture. Once beamforming and en velope detection are done for each ZM apodiza- tion, the maximum value over each result is taken as a directional zero-mean, as in E Z M = max  E Z M row , E Z M col  (6) where E represents the en velope of a beamformed sig- nal. Beamforming is also performed with two DC o ff set apodizations for each ZM apodization. For a gi ven ZM apodization, the DC o ff set apodizations are giv en by DC 1 X = Z M X + d c DC 2 X = Z M X − d c (7) where Z M represents the ZM apodization, dc represents the DC o ff set, DC 1 and D C 2 represent the DC apodiza- tions, and X is a placeholder for direction (i.e. ro w or column). Once beamforming is performed for each DC o ff set apodization, the en velopes of the beamformed signal in each direction are av eraged, as in E DC X = E DC 1 X + E DC 2 X 2 . (8) Then, a directional DC env elope is taken as the pixel- by-pixel maximum over both directions, same as with the ZM E DC = max  E DC row , E DC col  (9) Finally , the NSI image is obtained by subtracting the ZM en velope from the DC en velope, as E N S I = | E DC − E Z M | (10) where the absolute v alue is taken to av oid negati ve en- velope v alues. 5 3.3.3. Dir ectional Coher ence F actor The Directional Coherence Factor (DCF) in volv es projecting matrix array data along azimuth and elev a- tional dimensions, calculating a coherence factor for each direction, then multiplying those factors together to get the directional coherence factor [24]. Mathemat- ically , the projected data can be written as P i = Q X j = 1 S i , j = h P AZ 1 , ..., P AZ Q i P j = Q X i = 1 S i , j = h P E L 1 , ..., P E L Q i T (11) where P i is the projection onto azimuth, P j is the projec- tion onto elev ation, Q = √ N is the number of elements in ele vation and azimuth, and S i , j is time-delayed chan- nel data from the element at row i and column j . From these projected vectors, a coherence factor is calculated for each of them as C F AZ =    P Q i = 1 P i    2 Q P Q i = 1 | P i | 2 C F E L =     P Q j = 1 P j     2 Q P Q j = 1    P j    2 . (12) Finally , the DCF is calculated by multiplying C F AZ and C F E L : DC F = C F AZ × C F E L =     P Q i = 1 P Q j = 1 P i P j     2 Q 2 P Q i = 1 P Q j = 1    P i P j    2 . (13) After the DCF value is calculated for each pixel, the pixel is multiplied by the DCF . W e note that in the orig- inal published method, the authors included diagonal projections [24]. Howe ver , for reduced computation, we opted to only include the row and column projections. 3.3.4. Minimum V ariance Inspired by the directional projections from DCF , we also performed Minimum V ariance (MV) weighting on the projected vectors in Equation 11 [25]. The optimal MV weights are giv en by w = ˆ R − 1 a a H ˆ R − 1 a (14) where w is the desired weighting vector , ˆ R is an esti- mate of the spatial co variance matrix, and a is a steering vector set to all ones. The spatial cov ariance matrix w as estimated on either direction using spatial smoothing as in R AZ = 1 Q − L + 1 Q − L X l = 0 P i , l P H i , l R E L = 1 Q − L + 1 Q − L X l = 0 P j , l P H j , l (15) where P i , l and P j , l were sub-vectors of the directional projections in Equation 11, defined as the following: P i , l = h P AZ l , ..., P AZ l + L − 1 i P j , l = h P E L l , ..., P E L l + L − 1 i . (16) The sub-array size L was set to Q / 2 for both azimuth and elev ation projections. In both cases, diagonal load- ing was also performed on the covariance matrix esti- mates: ˆ R = R + ϵ I (17) where I is the identity matrix and ϵ = 1 10 L · tr ace { R } is the diagonal loading factor . After estimating the cov ari- ance matrices for either direction, weights were calcu- lated with Equation 14, resulting in a w AZ and a w E L . Then, projected MV outputs were calculated as y AZ = 1 Q − L + 1 Q − L X l = 0 w H AZ P i , l y E L = 1 Q − L + 1 Q − L X k = 0 w H E L P j , l (18) and finally , directions were combined with y = q y AZ y ∗ E L (19) where y is the final directional MV output, and ∗ denotes complex conjugate. 3.3.5. Quality Metrics T o assess the spatial resolution, the full width at half- maximum (FWHM) of the lateral profiles was estimated for images with wire targets. Lateral profiles were taken across the X / Z (lateral / depth slice) of the 3D volumes, av eraged axially with depth 0.5 mm and elev ationally with width 1 mm. Profiles were then interpolated by a factor of 64 using spline interpolation for FWHM esti- mates. The contrast of anechoic regions was estimated using the contrast ratio (CR), the contrast-to-noise ra- tio (CNR) and the generalized contrast-to-noise ratio 6 (gCNR). These metrics are given by the following equa- tions [30]: C R = 10 log 10 µ i µ o ! (20) C N R = µ i − µ o q σ 2 i + σ 2 o (21) gC N R = 1 − Z ∞ −∞ min x { p i ( x ) , p o ( x ) } d x (22) In the abov e equations, µ i is the a verage en velope level inside some region of interest (ROI), and µ o is av erage en velope level of the background. Likewise, σ 2 i and σ 2 o are the variances of the en velopes in those regions. Fi- nally , p i ( x ) and p o ( x ) are the probability distributions of the en velopes in the same regions. Histograms with 256 bins were used to estimate the distributions. Finally , to assess speckle quality , the speckle signal- to-noise ratio (sSNR) was also estimated. The sSNR is giv en by sS N R =    ¯ A    / √ var ( A ) , where A is the en- velope lev el in some region of interest, the v alue ¯ A is the mean and var ( A ) is the variance of the en velope in that region [31]. A higher sSNR indicates better speckle quality , with fully dev eloped speckle approaching an sSNR value of 1.91 [31]. 4. Results 4.1. Simulations Resulting max projections over depth of the simu- lated point-spread functions are displayed in Figure 5, along with lateral and elev ational profiles. Quality met- rics are giv en in T able 1. W e observed the highest side lobes with D AS, with main lobe to side lobe contrast of 26.3 dB. The NSI beamformer had the lowest side lobes, with their main lobe to side lobe contrast going down to 41.3 dB. The next lowest side lobes were DCF followed by MV , with contrasts of 30.2 dB and 28.8 dB respectiv ely . In spatial resolution, NSI and DCF achiev ed similar resolutions in both directions. The NSI and DCF beamformers achie ved 1.22 mm and 1.26 mm FWHM estimates respecti vely in the lateral direction. Then, they achieved 1.25 mm and 1.23 mm in the ele- vational direction. The MV beamformers demonstrated lowest resolution performance, with FWHM estimates of 1.49 mm and 1.33 mm in lateral and ele vational di- rections respectiv ely . Lastly , DAS displayed the lo west resolution, with FWHM estimates around 2 mm in ei- ther direction. Figure 5: Max depth projection of simulated PSFs of data from an 8 x 8 matrix array with a pitch of 6.7 wav elengths and element size 6.5 wav elengths. (a) D AS, (b) NSI with DC = 0.5, (c) DCF , (d) MV with L = Q / 2. (e) Lateral profile across elev ation of 0 mm of each max projection. (f) Elev ation profile across lateral position of 0 mm of each max projection. The red and yellow boxes in (a) illustrate the R OIs used for contrast estimates in T able 1. The lateral profiles were used for FWHM estimates. Beamformer D AS NSI DCF MV Lateral 2.29 mm 1.22 mm 1.26 mm 1.49 mm Elev ation 1.81 mm 1.25 mm 1.23 mm 1.33 mm Contrast 26.3 dB 41.3 dB 30.2 dB 28.8 dB T able 1: The first two rows are FWHM values along elev ation and lateral profiles of simulated PSFs with di ff erent beamformers. The last row are contrast estimates of main lobe versus side lobe lev els, with R OIs defined in Figure 5a. 4.2. V irtual Large Apertur e T able 2 lists coupling factors with corresponding el- ement counts, widths, and F-numbers for the virtual large aperture. W ith a coupling factor of two, the to- tal element count of the virtual aperture was the same as the V ermon matrix array , meaning the aperture size and FO V have increased by a factor of two in both dimen- sions with no increase to the element count. Figs. 6-9 display B-mode images acquired with the virtual large aperture. First, Figs. 6-7 display the constant depth slices of wire and cyst phantoms respectively , with each coupling factor and beamformer . Then, Figs. 8-9 dis- play the lateral and ele vational slices of the wire and cyst phantoms, ag ain with each coupling f actor . Fig. 10 displays lateral profiles for the wire targets and cyst re- gions. Finally , Fig. 11 displays the charts of quality metrics for each beamformer and coupling factor . W e observed a decrease to resolution with increas- ing coupling factors, with FWHM estimates going from 0.87 mm with no coupling to 1.88 mm with coupling by four, both with DAS. W ith a coupling factor of two, FWHM estimates improved from 1.17 mm with D AS to 0.81 mm with NSI beamforming. Overall, the NSI 7 Coupling Factor Element Count W idth ( λ ) F-number 1 4096 1.48 1.70 2 1024 3.10 3.51 4 256 6.32 7.16 T able 2: Element counts, widths, and minimum F-numbers for each coupling factor on the virtual large aperture. Figure 6: Lateral / elev ational slices of the wire phantom from the vir- tual large aperture. All images are displayed with a dynamic range of 60 dB. (a)-(d) Coupling factor 1 (i.e. no coupling). (e)-(h) Coupling factor 2. (i)-(l) Coupling factor 4. The green line in (a) illustrates the region for the lateral profiles and FWHM estimates. beamformer performed best for resolution. W ith con- trast metrics, the coupling factors had little e ff ect, while the beamformers influenced the estimates immensely . The D AS beamformer had contrast estimates around - 10 dB, CNR estimates near -1, and gCNR estimates near 0.6 for all coupling factors. Overall, D AS performed the best for all contrast metrics and sSNR across all cou- pling factors. Every other beamformer reduced one or more of the contrast metrics. The DCF beamformer was the worst for CNR (all abov e -0.6) and sSNR (all below 1), but the best for basic contrast for all coupling factors. The NSI beamformer had the lowest gCNR estimates, consistently around 0.4. The best CNR and gCNR met- rics behind D AS came from MV , with those metrics be- ing within 13% of corresponding estimates from DAS across all coupling factors. 5. Discussion Our results demonstrate how larger elements can lead to larger apertures, as well as highlight trade-o ff s be- tween di ff erent beamformers used to maintain resolu- Figure 7: Lateral / elev ational slices of the cyst phantom from the vir- tual large aperture. All images are displayed with a dynamic range of 60 dB. (a)-(d) Coupling factor 1 (i.e. no coupling). (e)-(h) Coupling factor 2. (i)-(l) Coupling factor 4. The red and yellow boxes in (a) illustrate the R OIs for contrast metrics. tion. W e observed that a coupling factor of two on the V ermon array still produced high quality images with advanced beamformers, where the resolution estimates were comparable to resolution estimates of uncoupled D AS. These results provide evidence that matrix ar- rays with doubled FO V could be produced by doubling the element size, and advanced beamformers can main- tain resolution. Howe ver , the increased resolution with advanced beamformers also came with a decrease to speckle or contrast metrics. Each of these beamform- ers has methods of tuning this trade-o ff , which can help maximize quality for various applications. NSI is tuned by the DC o ff set [21], DCF could be tuned if performing a GCF calculation on the projections instead of merely CF [24], and MV can be tuned by the diagonal loading factor and sub-array size [25]. Also, a coupling factor of four produced images that were grainy , blurry , and smeared, with a bigger dead zone near the transducer surface, compared to the lower coupling factors. W ith such a high coupling factor , too few elements are able to contribute to pixel values for any beamforming method to be e ff ecti ve. A coupling factor of four represented an element width of 6.3 λ , while a coupling factor of two represented a width of 3.1 λ (T able 2). Therefore, the element size cannot be pushed too far abo ve 3 λ if using large, square elements in a periodic array . Comparing the beamformers, we observed that NSI was the best for maintaining resolution compared to D AS, DCF , and MV , based on estimated FWHM val- 8 Figure 8: W ire phantom image slices from the virtual large aperture. (a)-(c) Lateral / depth slices. (d)-(f) Elevational / depth slices. (a), (d) Coupling factor 1. (b), (e) Coupling factor 2. (c), (f) Coupling factor 4. T op to bottom: DAS, NSI, DCF , MV beamforming. All images displayed with a 60 dB dynamic range. Figure 9: Cyst phantom image slices from the virtual large aperture. (a)-(c) Lateral / depth slices. (d)-(f) Elevational / depth slices. (a), (d) Coupling factor 1. (b), (e) Coupling factor 2. (c), (f) Coupling factor 4. T op to bottom: DAS, NSI, DCF , MV beamforming. All images displayed with a 60 dB dynamic range. 9 Figure 10: Lateral profiles from the virtual large aperture. (a)-(c) Lat- eral profiles of a wire target. (d)-(f) Lateral profiles of anechoic re- gions. (a), (d) Coupling factor 1. (b), (e) Coupling factor 2. (c), (f) Coupling factor 4. Figure 11: Quality metrics for the virtual aperture. (a) FWHM. (b) Contrast. (c) CNR. (d) gCNR. (e) sSNR. ues. Howe ver , it also reduced the speckle signal sig- nificantly , producing the second-lo west sSNR estimates behind DCF . The MV beamformer had the best CNR, gCNR, and sSNR metrics behind D AS for each cou- pling factor , but it did not improve resolution for cou- pling factors of two or four , making it the least suit- able for large-element array designs. This di ff erence in performance can be explained by considering the fact these beamformers were used in a low-element count scenario. W ith NSI, the fundamental idea is to form and in vert a null [21]. A null can be created using as little as two elements, one positi ve and one negati ve, mean- ing this method should still be e ff ective with fewer ele- ments. Ho wever , with MV beamforming, the idea is to calculate apodization weights that minimize the po wer of the apodized signal while maintaining unit gain at the focal point [25]. W ith fewer elements, and thus fewer apodization weights, the rejection of signals from out- side the focal point is reduced, making MV less e ff ective for lo w-element count scenarios. The DCF beamformer estimates the coherence of signals across multiple el- ements, in this case directional projections of the ma- trix elements. With low element counts, the coherence estimate might be less accurate, making this method less e ff ectiv e with fewer elements as well. It also re- duces speckle the most out of all beamformers because speckle signal is not coherent. Another important consideration for each of these beamformers is the computation that each one re- quires with respect to the element count. D AS is the most e ffi cient, as it only requires a sum ov er the elements. In this case, the complexity is O(N) where N is the total number of elements. NSI on a 2D array essentially consists of apodizing de- layed channel data six times, one for each apodization ( Z M AZ , Z M E L , DC 1 AZ , DC 1 E L , DC 2 AZ , DC 2 E L ). The remaining operations, such as the extra en velope detec- tions, the DC averaging (Eq. 8), the max over directions (Eqs. 6 and 9), and subtracting the null (Eq. 10) all run in constant time with respect to the elements. There- fore, NSI is also of complexity class O(N), with some added constant factors. The DCF beamformer first re- quires directional projections (Eq. 11), which are sums ov er ro ws or columns of channel data. In other words, it sums √ N values √ N times, making the projection calculation O(N) comple xity . The actual coherence fac- tor calculation then only in volv es summing √ N values from these projected vectors. Therefore, DCF is ov erall in complexity class O(N). Finally , our MV implemen- tation starts with the same projections, then estimates cov ariance matrices using spatial smoothing on the pro- jections (Eqs. 15 and 16). If the subaperture length for 10 spatial smoothing is denoted by L, then L ≤ √ N / 2, and the size of the covariance matrix is L x L. The most computationally expensi ve operation in MV is then in- verting that cov ariance matrix for the weight calculation (Eq. 14). Inv erting an L x L matrix is of complexity O( L 3 ). Therefore, our MV implementation has com- plexity O( L 3 ), or O( N 3 / 2 ). Based on the complexities, MV will require the most computation, e ven for low el- ement counts. T o examine the di ff erences between the other three beamformers, we include the run times for coupling factors of two and four with the virtual large aperture in T able 3. From the T able, we observe NSI be- ing the fastest of the adv anced beamformers, then DCF , then MV . Coupling factor D AS NSI DCF MV 2 86 µ s 3.1 ms 41 ms 1.42 s 4 91 µ s 1.2 ms 20 ms 490 ms T able 3: Run times for di ff erent beamformers in low-element count scenarios. 6. Conclusion W e have demonstrated in a phantom ho w a doubled element size (i.e. around 3 λ ) can allow for twice the FO V in 2D matrix arrays with no increase to the element count. The resulting loss in resolution from increased element size can be mitigated by advanced beamform- ers, making large elements more viable. The best beam- former we observed for maintaining resolution was NSI for both simulations and phantom experiments, and it was the most computationally e ffi cient of the advanced beamformers. Howe ver , NSI also reduced speckle and contrast metrics. Our approach to array design can al- low for larger 2D matrix arrays with an increased FO V with no added electronics. Such an improvement can make 3D ultrasound more viable and useful for clinics without increasing the cost or complexity of 3D ultra- sound systems. Declaration of Competing Interest The Authors declare that they hav e no competing fi- nancial interests or personal relationships that appear to hav e influenced the work reported in this paper . Data A vailability Data will be made av ailable upon request. Acknowledgments The authors would like to acknowledge and thank Matt Lowerison and Y iRang Shin for performing ele- ment calibration on the V ermon matrix array . This work was supported by the National Institutes of Health [grant numbers R21EB035714, R01CA273700, R01EB036800]. CRediT A uthorship Contribution Statement Mick Gardner: Writing - original draft; Writing - re- view and editing; In vestigation; Formal Analysis; Soft- ware. Michael L. Oelze: Writing - re view and editing; Supervision; Project Administration. References [1] Q. Huang, Z. Zeng, A Revie w on Real-T ime 3D Ultrasound Imaging T echnology , BioMed Research International 2017 (2017) 6027029. doi:doi:10.1155 / 2017 / 6027029. [2] J. Y u, H. Y oon, Y . M. Khalifa, S. Y . Emelianov , Design of a V olumetric Imaging Sequence Using a V antage-256 Ultrasound Research Platform Mul- tiplex ed With a 1024-Element Fully Sampled Ma- trix Array , IEEE T ransactions on Ultrasonics, Fer - roelectrics, and Frequency Control 67 (2020) 248– 257. doi:doi:10.1109 / TUFFC.2019.2942557. [3] A. Chavignon, B. Heiles, V . Hingot, C. Orset, D. V ivien, O. Couture, 3D T ranscranial Ultra- sound Localization Microscop y in the Rat Brain W ith a Multiplexed Matrix Probe, IEEE Transac- tions on Biomedical Engineering 69 (2022) 2132– 2142. doi:doi:10.1109 / TBME.2021.3137265. [4] B. Savord, R. Solomon, Fully sampled ma- trix transducer for real time 3D ultrasonic imag- ing, in: IEEE Symposium on Ultrasonics, 2003, v olume 1, 2003, pp. 945–953 V ol.1. doi:doi:10.1109 / UL TSYM.2003.1293556. [5] J. M. Rothberg, T . S. Ralston, A. G. Rothberg, J. Martin, J. S. Zahorian, S. A. Alie, N. J. Sanchez, K. Chen, C. Chen, K. Thiele, D. Gros- jean, J. Y ang, L. Bao, R. Schneider, S. Schaetz, C. Meyer , A. Neben, B. Ryan, J. R. Petrus, J. Lut- sky , D. McMahill, G. Corteville, M. R. Hage- man, L. Miller , K. G. Fife, Ultrasound-on-chip platform for medical imaging, analysis, and col- lectiv e intelligence, Proceedings of the National 11 Academy of Sciences 118 (2021) e2019339118. doi:doi:10.1073 / pnas.2019339118. [6] P . Acar , L. Battle, Y . Dulac, M. Pe yre, H. Dubourdieu, S. Hascoet, M. Groussolles, C. V ayssière, Real-time three-dimensional foetal echocardiography using a new transab- dominal xMA TRIX array transducer, Archi ves of Cardiov ascular Diseases 107 (2014) 4–9. doi:doi:10.1016 / j.acvd.2013.10.003. [7] K.-Q. Zhao, T . G. Bjastad, K. Kristo ff ersen, Er- ror analysis of subaperture processing in 1-D ultra- sound arrays, IEEE T ransactions on Ultrasonics, Ferroelectrics, and Frequency Control 62 (2015) 663–672. doi:doi:10.1109 / TUFFC.2014.006822. [8] Z. Y u, M. A. Pertijs, G. C. M. Meijer , A programmable analog delay line for Micro- beamforming in a transesophageal ultrasound probe, in: 2010 10th IEEE International Conference on Solid-State and Integrated Circuit T echnology , 2010, pp. 299–301. doi:doi:10.1109 / ICSICT .2010.5667749. [9] L. Castrignano, P . T ortoli, G. Matrone, M. Crocco, A. S. Sav oia, A. Ramalli, On the Impact of Microbeamformers in 3-D High Frame Rate Ultrasound Imaging: A Simulation Study, IEEE Transactions on Biomedical Engineering 72 (2025) 1941–1950. doi:doi:10.1109 / TBME.2025.3529198. [10] J. Kang, D. Go, I. Song, Y . Y oo, W ide Field- of-V iew Ultrafast Curved Array Imaging Us- ing Diver ging W aves, IEEE Transactions on Biomedical Engineering 67 (2020) 1638–1649. doi:doi:10.1109 / TBME.2019.2942164. [11] S. H. Kim, B. I. Choi, K. W . Kim, K. H. Lee, J. K. Han, Extended Field-of-V iew Sonography , Jour- nal of Ultrasound in Medicine 22 (2003) 385–394. doi:doi:10.7863 / jum.2003.22.4.385. [12] T . C. Poon, R. N. Rohling, Three-dimensional extended field-of-view ultrasound, Ultrasound in Medicine & Biology 32 (2006) 357–369. doi:doi:10.1016 / j.ultrasmedbio.2005.11.003. [13] C. W achinger , W . W ein, N. Nav ab, Three- Dimensional Ultrasound Mosaicing, in: N. A y- ache, S. Ourselin, A. Maeder (Eds.), Medical Im- age Computing and Computer-Assisted Interven- tion – MICCAI 2007, Springer , Berlin, Heidel- berg, 2007, pp. 327–335. doi:doi:10.1007 / 978-3- 540-75759-7_40. [14] A. Ramalli, E. Boni, E. Roux, H. Liebgott, P . T ortoli, Design, Implementation, and Med- ical Applications of 2-D Ultrasound Sparse Ar- rays, IEEE T ransactions on Ultrasonics, Ferro- electrics, and Frequency Control 69 (2022) 2739– 2755. doi:doi:10.1109 / TUFFC.2022.3162419. [15] H. Fa vre, M. Pernot, M. T anter, C. Papadacci, Boosting transducer matrix sensiti vity for 3D large field ultrasound localization microscopy using a multi-lens di ff racting layer: a simulation study , Physics in Medicine & Biology 67 (2022) 085009. doi:doi:10.1088 / 1361-6560 / ac5f72. [16] H. Favre, M. Pernot, M. T anter, C. Pa- padacci, T ranscranial 3D ultrasound localization microscopy using a large element matrix array with a multi-lens di ff racting layer: an in vitro study , Physics in Medicine & Biology 68 (2023) 075003. doi:doi:10.1088 / 1361-6560 / acbde3. [17] N. Haidour , H. Favre, P . Mateo, J. Reydet, A. Bizé, L. Sambin, J. Dai, P .-M. Chiaroni, B. Ghaleh, M. Pernot, M. T anter , C. Papadacci, Multi- lens ultrasound arrays enable large scale three- dimensional micro-vascularization characteriza- tion ov er whole organs, Nature Communica- tions 16 (2025) 9317. doi:doi:10.1038 / s41467- 025-64911-z. [18] Y . Paul, D. Barthez, R. Léveillé, V . Peter, D. Scriv ani, Side Lobes and Grating Lobes Artifacts in Ultrasound Imaging, V eterinary Radiology & Ultrasound 38 (1997) 387–393. doi:doi:10.1111 / j.1740-8261.1997.tb02104.x. [19] M. Gardner, R. J. Miller, M. L. Oelze, Grat- ing lobe mitigation on large-pitch arrays using null subtraction imaging, Ultrasonics 140 (2024) 107302. doi:doi:10.1016 / j.ultras.2024.107302. [20] J. A. Jensen, M. Schou, L. T . Jørgensen, B. G. T omov , M. B. Stuart, M. S. T raberg, I. T aghavi, S. H. Øygaard, M. L. Ommen, K. Steenberg, E. V . Thomsen, N. S. Panduro, M. B. Nielsen, C. M. Sørensen, Anatomic and Functional Imaging Using Row–Column Arrays, IEEE T ransactions on Ultrasonics, Ferroelectrics, and Frequency Control 69 (2022) 2722–2738. doi:doi:10.1109 / TUFFC.2022.3191391. [21] A. Agarw al, J. Reeg, A. S. Podkow a, M. L. Oelze, Improving Spatial Resolution Using Incoherent Subtraction of Receive Beams Having Di ff erent 12 Apodizations, IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control 66 (2019) 5–17. doi:doi:10.1109 / TUFFC.2018.2876285. [22] Z. Kou, R. J. Miller, M. L. Oelze, Grating Lobe Reduction in Plane-W ave Imaging W ith Angular Compounding Using Subtraction of Coherent Sig- nals, IEEE T ransactions on Ultrasonics, Ferro- electrics, and Frequency Control 69 (2022) 3308– 3316. doi:doi:10.1109 / TUFFC.2022.3217993. [23] M. Y ociss, K. Brown, K. Hoyt, Null Sub- traction Beamforming for Improved V essel Res- olution in V olumetric Contrast-Enhanced Ultra- sound, in: 2021 IEEE International Ul- trasonics Symposium (IUS), 2021, pp. 1–4. doi:doi:10.1109 / IUS52206.2021.9593368. [24] X. W u, W .-N. Lee, Directional Coherence Fac- tor for V olumetric Ultrasound Imaging with Ma- trix Arrays, IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control (2025) 1–1. doi:doi:10.1109 / TUFFC.2025.3557519. [25] J. F . Synne vag, A. Austeng, S. Holm, Adap- tiv e Beamforming Applied to Medical Ultrasound Imaging, IEEE T ransactions on Ultrasonics, Ferroelectrics, and Frequency Control 54 (2007) 1606–1613. doi:doi:10.1109 / TUFFC.2007.431. [26] E. Martin, Y . T . Ling, B. E. T reeby , Sim- ulating Focused Ultrasound Transducers Using Discrete Sources on Regular Cartesian Grids, IEEE T ransactions on Ultrasonics, Ferroelectrics, and Frequency Control 63 (2016) 1535–1542. doi:doi:10.1109 / TUFFC.2016.2600862. [27] B. E. T reeby , J. Budisk y , E. S. W ise, J. Jaros, B. T . Cox, Rapid calculation of acoustic fields from arbitrary continuous-w av e sources, The Journal of the Acoustical Society of America 143 (2018) 529–537. doi:doi:10.1121 / 1.5021245. [28] V . Perrot, M. Polichetti, F . V array , D. Garcia, So you think you can DAS? A vie wpoint on delay- and-sum beamforming, Ultrasonics 111 (2021) 106309. doi:doi:10.1016 / j.ultras.2020.106309. [29] Z. Kou, M. R. Lowerison, Q. Y ou, Y . W ang, P . Song, M. L. Oelze, High-Resolution Power Doppler Using Null Subtraction Imaging, IEEE T ransactions on Medical Imaging 43 (2024) 3060– 3071. doi:doi:10.1109 / TMI.2024.3383768. [30] A. Rodriguez-Molares, O. M. H. Rindal, J. D’hooge, S.-E. Måsøy , A. Austeng, M. A. Lediju Bell, H. T orp, The Generalized Contrast- to-Noise Ratio: A F ormal Definition for Lesion Detectability , IEEE T ransactions on Ultrasonics, Ferroelectrics, and Frequency Control 67 (2020) 745–759. doi:doi:10.1109 / TUFFC.2019.2956855. [31] R. W agner , S. Smith, J. Sandrik, H. Lopez, Statis- tics of Speckle in Ultrasound B-Scans, IEEE T ransactions on Sonics and Ultrasonics 30 (1983) 156–163. doi:doi:10.1109 / T -SU.1983.31404. 13