- Email: [email protected]
A decimated minimum variance beamformer applied to ultrasound imaging
Accepted Manuscript A decimated minimum variance beamformer applied to ultrasound imaging Sayed Mahmoud Sakhaei PII: DOI: Reference:
S0041-624X(15)00038-4 http://dx.doi.org/10.1016/j.ultras.2015.02.005 ULTRAS 5000
To appear in:
Ultrasonics
Received Date: Revised Date: Accepted Date:
21 June 2013 26 January 2015 5 February 2015
Please cite this article as: S.M. Sakhaei, A decimated minimum variance beamformer applied to ultrasound imaging, Ultrasonics (2015), doi: http://dx.doi.org/10.1016/j.ultras.2015.02.005
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
A decimated minimum variance beamformer applied to ultrasound imaging Sayed Mahmoud Sakhaei1∗ ABSTRACT Minimum variance beamforming has performed significant improvement in the resolution of the ultrasound images. However, its computational complexity is a serious problem. This paper introduces a new implementation of the minimum variance beamformer for ultrasound imaging with a focused transmit beam. In this method, a decimated aperture data instead of full of it, is used as the beamformer input, on which the minimum variance beamforming is applied, with the covariance matrix estimated using the full aperture data. In this way, the method can give a linear complexity while it can show a performance comparable to that of the full array implementation of the minimum variance beamforming, as the simulation and experimental results confirm this. Therefore, this adaptive beamforming method can be viewed as an approximate implementation of the minimum variance beamforming with a linear computational complexity.
1. Introduction Medical ultrasound imaging is conventionally done through a delay and sum (DAS) beamformer. This data independent beamformer applies a preset apodization on the signals transmitted or received by the array elements to control the beampattern properties. In contrast, data dependent adaptive beamformers
∗
S. M. Sakhaei is with the department of electrical and computer engineering at Babol University of Technology, Babol, Iran, (e-mail: [email protected])
1
update the apodization vector for each point of image, so that the signals received from the points far from the point of interest are lowered. The minimum variance (MV) beamforming is one of the most widely studied adaptive beamforming methods for ultrasound imaging [1]- [6], which has shown significant improvement in the resolution of the ultrasound images. In spite of its capability in the enhancement of the image, the computational complexity is a major drawback of it, which makes its application in real time ultrasound imaging a serious problem. While the conventional delay and sum (DAS) beamformer requires the computation load at the order of N ( i.e. N) at each time, that of the MV is N , where N is the number of array elements. This complexity is a cost paid for the enhancement achieved by the MV method. One way proposed to reduce the computational complexity of the MV is the beamspace domain MV beamformer [7]-[8]. This method employs the fact that the signals received from the points far from the focal points are very small and negligible, as a result of the focusing applied on the transmission. In this way, a good approximate implementation of the MV beamforming is reached. In [9], the authors proposed a method to obtain a covariance matrix in a Toeplitz structure, based on which the inverse of the matrix can be computed using a fast algorithm. Both of these methods reduce the complexity to N , while performing a performance comparable to that of MV. Synnevag et. al. [10] have proposed a low complexity adaptive (LCA) beamforming method based on the MV beamforming idea with a complexity of N. The method is based on the evaluation of the beamformer output for a set of predefined apodization and applying the best one. The evaluation criterion was an approximate value of the output power acquired through a temporal averaging. The method gives significant improvement in image resolution compared with the DAS. To retain the complexity to N, they have ignored the spatial averaging and this restricts its capability to enhance the resolution.
2
Moreover, the window design should be carefully done and may be application-dependent. The LCA also shows some artifacts on the image because of the discrete solution space used. This paper describes a new method to reduce the complexity of the MV. This method named as decimated minimum variance (DMV) beamformer uses the same fact used in the beamspace MV [7]; but applying it in array domain, instead of transform domain (beamspace). The key feature of the method is the decimation of the array without significant degradation of the obtained image. This feature reduces the number of weights which are adaptively calculated; however the covariance matrix is estimated through using all data received by the array. This property distinguishes the new method from applying the MV on a sparser array. In this way, the DMV beamformer is an approximate implementation of the MV and is capable to reduce the complexity to N, even with spatial and temporal averaging. The paper is organized as follows. The next section describes the new method through mathematical expression. Section 3 is devoted to the way of estimating the covariance matrix and then some designing considerations are discussed in section 4. The simulation results are given in sections 5, including some discussions on resolution, contrast and robustness properties and the comparison with the DAS and MV methods. The performance of the method on the experimental data is presented in section 6 and concluding remarks are given in section 7. 2. The method In array imaging systems, echo signals from the focal point arrived at individual transducer array elements are coherent and those emanating from the points far from the focal point are incoherent. It means that in the spatial frequency domain, the low frequency components are related to the signals received from the focal point [11]. Assuming dynamic receive focusing in the ultrasound imaging; this fact is valid for each imaging point. Therefore, high spatial frequency components are thrown away to obtain a good viewing image. This filtering process can simply be performed by smoothing the array
3
signals through summation or a smooth weighting and then summation, as it is done in the conventional beamforming methods. The aforesaid fact can be utilized in the MV beamforming to reduce its computational complexity. Assume that the array signals are filtered by a spatial frequency low pass filter to obtain a low pass array signal. Considering this filter as a decimation filter, the filtered array signal can be down-sampled or decimated in space with no significant aliasing distortion which means the filtered data can be recovered from its decimated version. Therefore, the minimum variance beamforming can be applied on the decimated- filtered data instead of the full data of the primary array, while the performance is almost preserved. It is worth to note that the transmit focusing makes the desired signal be a narrow band in the spatial frequency domain. This feature permits a high value for the down sampling rate, even with no decimation filter. This is in agreement with the results reported by Nilsen et.al [7] where the beamspace adaptive beamforming through using 3 beams have shown an acceptable performance. Consider an array of N equally spaced elements, and assume that [ ] is the primary array signal at time k after applying delays for dynamic receive focusing and is the tap weight vector of a J-length spatial FIR filter used for smoothing the signal (J
[ ] = [ ]
(1)
where, F is the convolution matrix for vector and [ ] a = − + 1 length vector, is the filtered array signal. Now, assume that [ ] is decimated to a N -length vector called [ ], such that [ ] and [ ] are related through the following equation:
, [ ] = , [ ]
(2)
4
where, the second subscript represents the element number and = !
"
is the decimation factor, which
is assumed to be an integer. The element-wise relation in (2) can be expressed in a matrix form as:
[ ] = # [ ]
(3)
where, Id is a × matrix obtained through decimation on column of a × identity matrix:
[# ]%, = &
1 , ' = ( − 1 + 1 1 0 , *+ℎ- ./0-
(4)
The beamforming process is applied on the decimated data:
2[ ] = . 334 [ ] [ ]
(5)
where, . 33 is the apodization vector and y is the beamformer output. In the minimum variance beamforming approach, the apodization vector is calculated to minimize the mean power of output which is obtained through the following optimization problem:
/( 56[|2[ ]| ] = . 334 8 . 33 9 4 3 0:;<-=+ +*: . 33 1 = 1
(6)
3 is a N × 1 vector of one, assuming dynamic receive focusing, and 8 is the where, the steering vector 1 covariance matrix of the decimated data vector:
8 = 6[ 4 ]
(7)
The optimal vector . 33 can be expressed as follows:
. 33 =
3 8 1 3? 8 1 3 1
(8)
It is seen that the decimated MV is an algorithm with computational complexity for the matrix inversion, instead of in the standard MV. It is worth to note that if the main sources of the interference signals are assumed those near to the point of interest, the weights on the decimated array can be calculated such that the output power of the primary array is a scaled version of that of the decimated array (see Appendix A). As a result, the
5
minimum power of MV beamformer on the primary array can be reached through applying optimal weights in (8) on the decimated one. Hence, using MV on the decimated array will obtain an image similar to that obtained by the primary full array. To get a robust response, the MV beamformer usually uses a diagonally load covariance matrix. This technique can be utilized in the DMV and the covariance matrix 8 can diagonally be loaded by a factor @ as follows:
8 = 8 +
@ + A=- 8 #
(9)
3. The spatial and temporal averaging In MV beamformer, the usual way to estimate the covariance matrix through using the single snapshot data is the subarray averaging method. In this method, estimation is done by dividing the aperture into overlapping subarrays and averaging the covariance matrices of each subarray. For the decimated MV, the subarrays are built on the filtered array and then each subarray is decimated into elements, where averaging is done over these decimated subarrays. Assuming each subarray has length L and M is the length of the filtered array, P=M-L+1 subarrays can be used in averaging. Also, consider B as a × 1 vector containing the decimated data of p’th subarray. Then the decimated covariance matrix is estimated as: F
1 4 8C [ ] = E B [ ] B [ ] D
(10)
BG
This technique for estimating the covariance matrix can be considered as the decimation in columns and rows of the covariance matrix of the full array data (see Appendix A). In a similar way, the covariance matrix can be estimated through the temporal and spatial averaging:
6
J
F
1 4 8C [ ] = E E B [ − / ] B [ − /] D 2I + 1
(11)
KGJ BG
It is apparent that spatial (subarray) averaging is identical to the temporal- spatial averaging with K=0. The beamformer output, in both cases, is calculated as follows: F
1 2[ ] = E . 334 [ ]B [ ] D
(12)
BG
where, . 33 is given by (8). The computational complexity for the DMV, assuming spatial averaging is only applied, can be determined as follows: • Calculating the output of the filters requires MJ multiplications. • The estimation of the covariance matrix is done through P times multiplying an × 1 vector by its transpose conjugate (according to (10)), so it needs the D multiplications. • The computation of the inverse of the covariance matrix is . • The computation of . 33 and then calculating the beamformer output is done through + and D multiplications, respectively. Hence, assuming the inverse of the covariance matrix is available, calculating the beamformer output can be done with a computational complexity of D + . Therefore, the complexity for the DMV is + + D . As explained before, in the transmit focusing mode, the filters can be ignored without appearing serious aliasing distortion. Moreover, it is possible to use filters with the uniform weights. These show that the complexity for the DMV can be reduced to + D . By contrast, that for the MV is L + DL , and for the beamspace MV is M + DM + DLM , while the LCA needs NN multiplications, where M and N are the number of
7
beams and the number of windows used. This comparison shows that if is set to a small value (i.e. ≪ L), then the DMV outperforms the three others in the computational complexity aspect. 4. Designing considerations This section is devoted to describe some points useful for designing the decimation filter and determining the proper value for . As stated above, the decimation filter is mainly used to cancel out high spatial frequency portion of the received data, as well as it permits the signal be down-sampled without aliasing distortion. It also suppresses the noise in array signal and prevents the noise cause the aliasing distortion. Hence, if the noise on the received array signal is high, applying a proper decimation filter would be needed before down-sampling. The decimation filter may cause some degradation in the resolving capability of the MV beamformer. This degradation is mainly affected by the length of the decimation filter used, as its increased length is resulted in the shortened spatial length of the filtered array. For example, if a filter of length J=32 is applied on a 64-elements array, the length of the filtered array is approximately equal to the half of the primary array. This shortened length leads to the reduced resolving capability. Therefore, the decimation filter should carefully be designed to have a beamformer with a good resolution property and no significant aliasing effects. Nevertheless, the decimation filter may result in a better contrast property assuming a decreased sidelobe level of the decimation filter at the price of its mainlobe width.
8
The most important parameter which should carefully be determined is , as the main advantage of the DMV over the MV method is the reduced computational complexity, while the performance is comparable. It is clear that a high value of can produce an enhanced image at the price of complexity,
1
Full array Decimated array 0.8
Normalized Energy
Normalized Beampattern
1
0.6
0.4
0.2
0 0
0.2
0.4
0.6
0.8
0.9 0.85
0.8
0.75 0
1
M=32 M=64 M=128
0.95
5
10
15
20
25
30
Nd u Fig. 1. (a) One way beampattern of a 9-elements full array (solid line) and that of decimated to 3 elements (dash line), (b) the energy of 2-way beampattern versus for three values of number of array elements. The energy is normalized to the energy of full aperture array.
which is an undesired result. On the other hand, using a very small value of may result in a highly degraded image. To investigate the effects of on the performance of the proposed beamformer, we assume that the focal point in the transmission and the reception is the same. The filed pattern at the focal point is well described by the far filed beampattern [12]. Moreover, to simplify our discussions, the analysis is restricted to the continuous wave (CW) scenario because of its mathematical tractability. In fact, as the analysis is used to compare the performance of the DAS beamformer with the DMV, the conclusions would be valid in the practical situations where the pulsed wave (PW) excitation is applied. Furthermore, the array elements are assumed to be spaced by P/2 and apodized by a rectangular window in the transmission.
9
The steered response (i.e. the received signal by the array at various angles) is described by the twoway beampattern, which is the production of the transmit and receive beampattern:
RD: = RD? :. RDT :
(13)
Here, RD, RD? and RDT are 2-way, transmit and receive beampattern, respectively, and : = 0/(U where U is the observation angle. Since the rectangular apodization is assumed in transmission, the transmit beampattern of an array of M elements can be represented as follows:
RD? : =
sin Y:/2 sin Y:/2
(14)
Assuming the receive aperture be decimated to elements and with uniform weighting, the receive beampattern is:
RDT :; = where = !
sin Y :/2 sin Y :/2
(15)
. An exemplary transmit and receive beampattern is plotted in Fig. 1(a), for M=9 and
"
=3. It can be seen that the 2-way beampattern would have small magnitude at angles far from the main axis even the receive beampattern has high magnitude. To investigate the energy lost due to decimating the receive aperture, the energy of 2-way beampattern is calculated as:
6 = [ |RD? :. RDT :; | \:
(16)
and, after normalizing to E(M), is depicted in Fig. 1(b) versus for different values of M (M=32, 64 and 128). The figure shows that the energy received by the decimated aperture with =3 can be almost 84% of the full aperture and that for =5 is higher than 90%. Therefore, it can be concluded that decimating the receive aperture to 3 or 5 elements would not significantly lose the energy. Moreover, the difference between the performance of the decimated and full aperture may be eliminated by applying an adaptive beamforming method on the decimated aperture. It is because,
10
according to (14), the transmit beampattern has nulls at directions specified by sin Y:/2 = 0, except for u=0. This means that the first right and left nulls occur at : = ±2/ and the second at : = ±4/ and so on. Moreover, the first sidelobes are centered at : = ±3/ with a level of approximately M/4.7 and the other sidelobes would have levels that are monotonically decreased approximately by 1/u. These properties of the transmit beampattern can be observed in Fig. 1. This discussion emphasizes that the main source of interference is the first sidelobes and the energy interfered form the other sidelobes is usually negligible. Hence, to catch a qualified image, the receive beampattern should have small magnitude especially at the first sidelobe regions of the transmit beampattern. Consequently, the adaptive receive beamformer would show a performance comparable or better than that of DAS if the receive beampattern can apply at least two nulls at the first sidelobes. This implies that the receive beamformer should have three degrees of freedom to be capable of imposing two nulls while satisfying the condition of the distortionless response (the constraint in the optimization problem (6)). In the other words, a beamformer with three weighting parameters and an aperture of the same length as the transmit aperture can afford the mentioned conditions. Therefore, by setting = 3 and assuming the signal received from the sidelobes other than the nearest ones is negligible, the DMV beamformer would outperform the DAS. Correspondingly, setting = 5 enables the beamformer to cancel the interfering signals from the four nearest sidelobes and hence a better performance is reached. It is worth to emphasize that selected according to the above discussion is such that the DMV method outperforms the DAS in terms of resolution assuming that the main sources of interferences are the points at the vicinity of the point of interest; an assumption that is usually valid in the focused transmit mode. In this way, the above idea for determining can be used for every uniformly spaced array transducer such as a linear array with an element space of P.
11
5. Simulation results The proposed beamformer was tested on the simulated ultrasound data, obtained using Field II [13]. Two phantoms were simulated: one consisting of pairs of wire targets and the other was a circular cyst in a speckle medium. For the simulations, a 3.5MHz, 64-element linear array transducer with λ/2 spacing was used. Excitation pulse had a Gaussian envelop with -6 dB relative bandwidth equal to 0.5. Moreover, the transmit focus was fixed at distance 45 mm and dynamic focusing was applied in the reception. To fully evaluate the performance of the DMV method, ten different beamforming methods were simulated according to the following list: BF1: a DAS beamformer with rectangular weighting. BF2: a MV beamformer with subarray length of L=N/2. BF3: a beamspace beamformer with three central beams and subarray length of L=N/2. BF4: an LCA beamformer with 12 windows defined according to TABLE I of [10]. BF5: a DMV beamformer with = 3 without spatial filtering, i.e. J=1. BF6: a DMV beamformer with = 5 without spatial filtering, i.e. J=1.
12
BF2
BF3
BF4
BF5
30
30
35
35
35
35
35
40
40
40
40
40
45
45
45
45
45
50
50
50
50
50
55 60
55 60
55 60
Depth [mm]
30
Depth [mm]
30
Depth [mm]
30
Depth [mm]
Depth [mm]
BF1
55 60
55 60
65
65
65
65
65
70
70
70
70
70
75
75
75
75
75
80
80
80
80
80
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
Fig. 2. Simulated point targets; BF1: DAS with uniform weighting; BF2: the MV beamformer; BF3: Beamspace beamformer; BF4: LCA beamformer; BF5: the DMV with J=1 and =3. The images are displayed using 50 dB dynamic range.
BF7: a DMV beamformer with = 3 and a uniform weighting decimation filter with length of J=16. BF8: a beamformer similar to BF7 but with J=32. BF9: a beamformer similar to BF8 but without decimation (J=32 and = 33. BF10: a MV beamformer applied on the 32 central elements of the array with subarray length of 16. The subarrays length L were set to [N/2] for the MV beamformers and to [M/2] for the DMVs, where [ ] is the largest integer value not greater than . Moreover, the covariance matrices were diagonally loaded by a factor of 0.01 for all MV- based methods. A. Simulated Point Targets
13
0
0 BF1 BF2 BF3 BF4 BF5
-10
BF1 BF2 BF3 BF4 BF5
-10
-20
Magnitude [dB]
Magnitude [dB]
-20
-30
-30
-40
-40
-50 -50
-60 -60
-4
-2
0
2
-8
4
-6
-4
-2
0
2
4
6
8
Lateral distance [mm]
Lateral distance [mm]
Fig. 3. Steered responses of DAS, LCA and other MV-based beamformers at depths 40 mm (the left figure) and 70 mm (the right one). Transmit focus was at 45 mm and dynamic receive focusing was applied.
The point targets were 6 paired wires separated by 3 mm laterally at 6 different depths. Fig. 2 displays images of the point targets using the first five beamformers (BF1 to BF5). We see that all the MV-based beamformers (BF2, BF3 and BF5) as well as LCA beamformer (BF4) give significantly better resolution than DAS (BF1). In comparison to the LCA, the figure shows that the image produced by the DMV beamformer is more similar to that of MV. These properties are verified by looking at Fig. 3, where the steered responses at two different depths are depicted. Fig. 3(a) shows the responses at depth 50 mm and Fig. 3(b) shows those at depth 80 mm. It is obvious that the DMV beamformer is capable of resolving the two points at depth 80mm, where the DAS beamformer is not. Furthermore, the figures illustrate that the performance reconstruction of the DMV methods is depth dependent such that at depth 50 mm the accuracy of reconstruction of DMV is similar to that of MV, but at depth 80 mm, the advantage of MV over DMV is obvious. This property may be considered as a result of aliasing effect which causes the coherent portion of the array signal be increased. Although the amount of the aliased energy is generally low, it may be considerable at the vicinity of a high echoic reflector. This leads to a reduced resolving capability of two point reflectors at deeper distances, where the two near points are seen at a narrower angle and the aliasing distortion is appeared more. It seems that one way to maintain the resolving
14
capability for different depth's is to use a small value of (e.g. 3) at small distances and increasing it (e.g. 5 or 7) at deeper distances. More results of the DMV beamformers are displayed in Fig. 4, where the images obtained through using different parameters of and J are shown. The image of BF6 beamformer emphasizes that increasing from 3 to 5 can slightly enhance the DMV performance. The figure also indicates the effect of the length of the spatial filter on the obtained image: by increasing J, the capability of resolving two near point reflectors is decreased, while the amplitude of the image at points laterally far from the reflectors are diminished. As stated in the previous section, these properties are the results of the reduced spatial length of the filtered array and attenuated high spatial frequency contents of the signal. Furthermore, the figure shows that BF8 (the DMV with a spatial filter of length J=32) suffers from an artifact problem for deeper distances. This artifact is a peak response between two point reflectors and appeared as a non-real point in the image. It could not be considered as a result of decimation, as it is seen in BF9 (which is the same as BF8 but without decimation). This artifact is also appeared in the standard MV beamformer of BF10 and could be appeared everywhere the two point reflectors are such close that the MV cannot resolve them. As mentioned in the previous section, applying the low pass filter restricts the resolving capability of the beamformer and this capability is more reduced by increasing the length of the filter. This restricted resolution is responsible for the aforementioned artifact. This fact can be verified by Fig. 5, showing that the steered responses at depth 80 mm for beamformers BF8, BF9 and BF10 are very similar. In other words, the performance of a DMV with J=32 is similar to that of a MV on a 32elements array.
15
BF7
BF8
BF9
BF10
30
30
35
35
35
35
35
40
40
40
40
40
45
45
45
45
45
50
50
50
50
50
55 60
55 60
55 60
Depth [mm]
30
Depth [mm]
30
Depth [mm]
30
Depth [mm]
Depth [mm]
BF6
55 60
55 60
65
65
65
65
65
70
70
70
70
70
75
75
75
75
75
80 -5
80 0
5
-5
Lateral Distance [mm]
80 0
5
Lateral Distance [mm]
-5
80 0
5
Lateral Distance [mm]
-5
80 0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
Fig. 4. Simulated point targets; BF6: the DMV with J=1 and =5; BF7: the DMV with J=16 and =3; BF8: the DMV with J=32 and =3; BF9: the DMV with J=32 and =33; BF10: the MV on a 32-element array. The images are displayed using 50 dB dynamic range.
0 BF6 BF7 BF8 BF9 BF10
Magnitude [dB]
-10 -20 -30 -40 -50 -60 -10
-5
0 5 Lateral distance [mm]
10
Fig. 5. Steered responses of DMV beamformers (BF6 to BF9) and MV on a 32-element array (BF10) at depth 80 mm. Transmit focus was at 45 mm and dynamic receive focusing was applied.
16
BF2
50 -5
5
0
BF6
5
45
50 -5
Lateral Distance [mm]
0
5
Lateral Distance [mm]
Depth [mm]
Depth [mm]
0
50 -5
5
BF9
45
0
5
Lateral Distance [mm]
0
5
Lateral Distance [mm]
BF10
40
50 -5
45
Lateral Distance [mm]
40 Depth [mm]
Depth [mm]
45
50 -5
5
BF8
40
0
0
40
45
Lateral Distance [mm]
BF7
40
50 -5
45
50 -5
5
Lateral Distance [mm]
BF5
40
Depth [mm]
0
45
Lateral Distance [mm]
BF4
40 Depth [mm]
Depth [mm]
Depth [mm]
45
50 -5
Depth [mm]
BF3
40
40 Depth [mm]
BF1 40
45
50 -5
0
5
Lateral Distance [mm]
45
50 -5
0
5
Lateral Distance [mm]
Fig. 6. Simulated cyst phantom; BF1: DAS with uniform weighting; BF2: the MV beamformer; BF3: Beamspace beamformer; BF4: LCA beamformer; BF5: the DMV with J=1 and =3; BF6: the DMV with J=1 and =5; BF7: the DMV with J=16 and =3; BF8: the DMV with J=32 and =3; BF9: the DMV with J=32 and =33; BF10: the MV on a 32-element array. The images are displayed using 50 dB dynamic range.
B. Simulated Cyst Phantom The cyst phantom consisted of a 4 mm cyst located at depth of 45 mm in a speckle medium. The speckle pattern is simulated with 10 randomly placed scatterers within a resolution cell of λ3, where λ is the central wavelength of the propagating waveform in the medium, to ensure fully developed speckle. The scattering amplitudes are Gaussian distributed. For all MV-based beamformers, the estimation of the covariance matrix is done through spatial and temporal averaging, where for temporal averaging K is set to 10 at the sampling rate of 20 MHz. The images of the cyst target are shown in Fig. 6. The contrast parameters (CNR and CR) are listed in Table I. The CR is Contrast Ratio and defined as the ratio of the mean value in the background to the mean value in the cyst region [14] and the CNR is Contrast to Noise Ratio and defined as the CR divided by the standard deviation of image intensity in the background region [15]. It illustrates that except for the LCA beamformer (BF4) which gives the best contrast
17
properties, the DMV beamformers (BF5 to BF9) show better contrast than DAS and other MV-based beamformers. The Table also indicates that the contrast properties of the DMV are enhanced by increasing from 3 to 5, which is similar to that seen for the resolution properties mentioned before. Moreover, the Table shows the effects of spatial filtering on the contrast enhancement, which are the results of the reduced high spatial frequency of the signal, i.e. more attenuating the signals received from far sidelobes. TABLE I. CONTRAST PARAMETERS OF THE SIMULATED CYST PHANTOM FOR DIFFERENT
BF1 BF2 BF3 BF4 BF5 BF6 BF7 BF8 BF9 BF10
CR(dB)
CNR
34.1 34.9 34.3 41.7 36.5 38.3 37.9 38.6 38.7 32.1
4.9 4.9 4.4 5 4.6 4.8 4.7 4.6 4.8 4
C. Sensitivity to the sound speed inhomogeneities A potential weakness of an adaptive beamformer is the lack of robustness against errors in the assumed wavefield parameters [2]. To investigate the robustness of the DMV, we assumed that the acoustic velocity is known with an overestimation by 5%. The steered responses of DMV beamformers for the points at depth 40 mm of the point phantom are depicted in Fig. 7 and compared with those of DAS and other MV-based methods. It is seen that the LCA is the best adaptive beamformer in robustness aspect and the MV is the worst one. The figure also shows that the DMV methods estimate the amplitude of points near 6 dB better than the MV beamformer does. This looking strange result can be explained by examination of the MV properties. It is known that the MV beamformer cancels out the interference
18
0 BF1 BF2 BF3 BF4 BF5 BF7
Magnitude [dB]
-10 -20 -30 -40 -50 -60 -4
-2
0
2
4
Lateral distance [mm]
Fig. 7. Steered responses of DAS, LCA and other MV-based beamformers at depth 40mm with 5% overestimation of the acoustic velocity. Transmit focus was at 45 mm and dynamic receiver focusing was applied.
signals which are incoherent, while the coherent portion of the received signal is considered as the desired one and retained by the beamformer. But at the presence of the sound speed error, the signals received from the point of interest would not be at the same phase and the property of coherency would slightly be lost. As a result, the MV beamformer attenuates that signal. This property would be appeared more obviously when the resolving capability of the beamformer is increased, where a received signal with no complete coherency may be considered as the interference and be completely cancelled out. This means that the robustness of the MV beamformer is inversely related to the number of weights adaptively determined. This fact is verified by [2], where the authors showed that increasing the subarray length in applying MV leads to the more sensitivity against errors. A similar situation occurs in DMV when compared to the MV and hence, the DMV is more robust than MV. Fig. 7 also shows that increasing J from 1 to 16 in DMV beamformer results in decreased sensitivity to the velocity error. This decreased sensitivity is similar to that appeared in DMV compared to MV and can be explained by the same facts.
19
6. Experimantal Results The proposed beamformer was applied on two experimental ultrasound data acquired at the Biomedical Ultrasonics Laboratory, University of Michigan (http://www.bme.umich.edu/labs/bul) using Synthetic Aperture Focusing Technique (SAFT) [16]. The first data was for a phantom of 6 wires at different depths. It was a complete data set from 128 transmit-receive measurements using a 128element linear array with the center frequency of 3.5 MHz sampled at rate 13.8889MHz. To compare the performance of different beamformers and for the sake of simplicity, we utilized the data only on the 64 central elements of the array such that the number of elements in transmit and receive modes was 64. Moreover, the transmit array was apodized by a rectangular window and focused at depth 60mm, as well as dynamic focusing was applied in the receiving mode. To do focusing, the data were firstly upsampled by a factor of 10 and then focusing delays were applied. Because of high resolution of MV-based beamformers, aliasing distortion could be occurred if sufficient number of lines per sector was not existed. Hence, we used 400 lines for a sector of 90 degrees for generating the images. The images obtained through three different beamformers (DAS, MV and DMV) are displayed in Fig. 8. It illustrates that the images obtained through DMV and MV beamformers have better resolution properties than DAS. The second experimental data used to evaluate the proposed beamformer was for a heart phantom acquired by a 64-element linear array with center frequency of 3.333 MHz and 17.76 MHz sampling rate. The focusing operation was completely the same as that of the first data. Figure 9 displays the images obtained by different beamformers and shows that the beamformers have similar contrast properties, while the walls are better defined in MV and DMV compared to DAS.
20
30 40
20
50
40
60
60
70 80
80 90
100
100
120
110
-50
0
120 -50
50
0
50
0
50
0
50
30 40
20
50
40
60
60
70 80
80
90
100
100 110
120 -50
0
120 -50
50
30 40
20
50
40
60
60
70 80
80 90
100
100
120
110
-50
0
120 -50
50
Fig. 9. Images of the second experimental data using DAS (top figure), MV (the middle one) and DMV beamformer (the bottom figure). The images are displayed using 50 dB dynamic range. Horizontal and vertical axes are lateral and axial distances, respectively, both in millimeters.
Fig. 8. Images of the first experimental data using DAS (top figure), MV (the middle one) and DMV beamformer (the bottom figure). The images are displayed using 50 dB dynamic range. Horizontal and vertical axes are lateral and axial distances, respectively, both in millimeters. 21
7. Conclusions This paper applied a new method for implementing the MV beamformer in ultrasound imaging systems with focused transmit beam. This method referred to as the decimated minimum variance (DMV) beamformer, uses a decimated aperture data, instead of full aperture data to apply the apodization vector, while utilizing full data to estimate the covariance matrix. To reduce the aliasing effects caused by the decimation, some antialiasing filters can be applied before decimation. The most important feature of the new method is its reduced computational complexity compared to the MV beamformer. The simulations have indicated that the elements of the subarrays used in the DMV can be lowered to a smaller value such as 3, while the performance remains almost comparable to the full MV beamformer. While the MV shows better resolution and contrast properties in no modeling errors conditions, the DMV behaves more robust at the presence of errors. This fact denotes that further works are needed to compare the DMV with MV. The most advantage of DMV over MV is its complexity. Assuming focused transmit beam and with no decimation filter or with a filter with uniform weighting, the computational complexity of the DMV would be at N, compared to N for the MV method. Thus, the DMV is an adaptive beamforming whose complexity is in order of the complexity of the DAS, while outperforming the DAS in contrast and resolution properties. As the complexity is most restrictive feature for applying MV on a 2D array, using the DMV on a 2D array can obviate this problem and hence, the DMV is a good candidate for applying an MV-based adaptive beamforming on a 2D array.
Appendix A In this appendix, it is shown that under some reasonable assumptions, the beamformer output of a decimated array is approximately equal to the scaled output of the primary array. Moreover, the
22
covariance matrix of the decimated array is related to that of the primary array by the decimation. The first assumption is that the main sources of the interference signals are those near to the point of interest. The second one is that the signals arrived to the array are spatially stationary. This assumption results in a Toeplitz structure for the covariance matrix [9] and also permits analyzing it in spatial frequency domain. The mean power of the beamformer output can be written as:
6 [|2[
]| ]
1 e 1 e [ Φ d\d = [ Φf d |g d| \d = 2Y e c 2Y e
(A-1)
where Φc d and Φf d are respectively the Fourier transforms of the correlation functions of the y[k] and x[k], which called the power spectral density and gd is the Fourier transform of the applied weighting on the array data x[k]. By the first assumption, Φf d has negligible values outside the interval (−Y/h, Y/h for some values of h > 1 and hence: e
1 k [ Φ d |gd| \d 6 [|2[ ]| ] ≅ 2Y e f k
=
e 1 m m [ Φf l n og l no \m 2Y h e h h
(A-2)
If we set Φf d = hΦf d/h) and g d = hgd/h, (A-2) can be written as:
6 [|2[ ]| ] ≅
h e [ Φ d |g d| \d 2Y e f
(A-3)
Equation (A-3) states that the power at the beamformer output of a full array would approximately be equal to h times of a decimated array if two conditions are satisfied. The first one is that the weighting applied on the decimated array is such that its Fourier transform is an expansion of that of the full array: g d = hgd/h. Note that this condition does not necessarily mean that . [ ] is a decimated version of .[ ], because of possible aliasing which can be occurred in decimating .[ ].
23
The second condition is that the correlation function of input is the decimated correlation function of the primary input x[k] by a factor of L. It means that, the decimated array has an input of [ ] which is spatially decimated of x[n] and its covariance matrix is obtained by applying the decimation on the columns and rows of the covariance matrix of x[k]. In other words, the covariance matrix of [ ] can be calculated through spatial averaging of the decimated subarrays.
REFERENCES
[1] F. Vignon, and M. R. Burcher “Capon beamforming in medical ultrasound imaging with focused beams,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 55, no. 3, pp. 619-628, March. 2008. [2]
F. Synnevag, A. Austeng, and S. Holm, “Adaptive beamforming applied to medical ultrasound imaging,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 54, no.8, pp. 1606-1613, Aug. 2008.
[3]
I. Holfort, F. Gran, and J. Jensen, “Broadband minimum variance beamforming for ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 2, pp. 314-325, Feb. 2009.
[4]
B. Mohammadzadeh Asl, and A. Mahloojifar, “Minimum variance beamforming combined with adaptive coherence weighting applied to medical ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 9, pp. 1923-1931, Sep. 2009.
[5] B. Mohammadzadeh Asl, and A. Mahloojifar, “Contrast enhancement and robustness improvement of adaptive ultrasound imaging using forward-backward minimum variance beamforming,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 58, no.4, pp. 858-867, Appr. 2011. [6] X. Zeng, C. Chen, and Y. Wang, “Eigenspace-based minimum variance beamformer combined with Wiener postfilter for medical ultrasound imaging,” Ultrasonic., vol. 52, pp. 996-1004, 2012. [7]
C. C. Nilsen, and I. Hafizovic, “Beamspace adaptive beamforming for ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 10, pp. 2187–2197, Oct. 2009.
[8]
X. Zeng, Y. Wang, J. Yu, and Y. Guo, “Beam-domain eigenspace-based minimum variance beamformer for medical ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 60, no. 12, pp. 2670–2676, 2013.
24
[9] C B. Mohammadzadeh Asl, and A. Mahloojifar, “A low-complexity adaptive beamformer for ultrasound imaging using structured covariance matrix,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 59, no.4, pp. 660-667, Appr. 2012. [10]
J.-F. Synnevåg, A. Austeng, and S. Holm, “A low-complexity data-dependent beamformer,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 57, no. 2, pp. 281–289, Feb. 2010.
[11] P. S. Naidu, Sensor array signal processing, CRC Press, 2001. [12] J. L. Prince, and J. M. Links, Medical imaging signals and systems, Pearson Prentice Hall, 2006. [13]
J.A. Jensen, “Field: a program for simulating ultrasound systems,” Med. Biol. Eng. Comp. 10th Nordic-Baltic Conference on Biomedical Imaging 4 (Suppl. 1), 1996, pp. 351–353, Part 1.
[14] M. O’donnell and S. W. Flax, “ Phase aberration correction using signals from point reflectors and diffuse scatterers: Measurements,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 35, no. 6, pp. 768–774, 1988. [15] S . Krishnan, K. W. Rigby, and M. O’donnell, “Improved estimation of phase aberration profile,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 44, no. 3, pp. 701–713, 1997. [16] M. Karaman, P. C. Li, and M. O'Donnell, “Synthetic aperture imaging for small scale systems,” IEEE Trans. Ultrason., Ferroelect., Freq. Contr. , vol. 42, pp. 429-442, May 1995.
25
S0041-624X(15)00038-4 http://dx.doi.org/10.1016/j.ultras.2015.02.005 ULTRAS 5000
To appear in:
Ultrasonics
Received Date: Revised Date: Accepted Date:
21 June 2013 26 January 2015 5 February 2015
Please cite this article as: S.M. Sakhaei, A decimated minimum variance beamformer applied to ultrasound imaging, Ultrasonics (2015), doi: http://dx.doi.org/10.1016/j.ultras.2015.02.005
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
A decimated minimum variance beamformer applied to ultrasound imaging Sayed Mahmoud Sakhaei1∗ ABSTRACT Minimum variance beamforming has performed significant improvement in the resolution of the ultrasound images. However, its computational complexity is a serious problem. This paper introduces a new implementation of the minimum variance beamformer for ultrasound imaging with a focused transmit beam. In this method, a decimated aperture data instead of full of it, is used as the beamformer input, on which the minimum variance beamforming is applied, with the covariance matrix estimated using the full aperture data. In this way, the method can give a linear complexity while it can show a performance comparable to that of the full array implementation of the minimum variance beamforming, as the simulation and experimental results confirm this. Therefore, this adaptive beamforming method can be viewed as an approximate implementation of the minimum variance beamforming with a linear computational complexity.
1. Introduction Medical ultrasound imaging is conventionally done through a delay and sum (DAS) beamformer. This data independent beamformer applies a preset apodization on the signals transmitted or received by the array elements to control the beampattern properties. In contrast, data dependent adaptive beamformers
∗
S. M. Sakhaei is with the department of electrical and computer engineering at Babol University of Technology, Babol, Iran, (e-mail: [email protected])
1
update the apodization vector for each point of image, so that the signals received from the points far from the point of interest are lowered. The minimum variance (MV) beamforming is one of the most widely studied adaptive beamforming methods for ultrasound imaging [1]- [6], which has shown significant improvement in the resolution of the ultrasound images. In spite of its capability in the enhancement of the image, the computational complexity is a major drawback of it, which makes its application in real time ultrasound imaging a serious problem. While the conventional delay and sum (DAS) beamformer requires the computation load at the order of N ( i.e. N) at each time, that of the MV is N , where N is the number of array elements. This complexity is a cost paid for the enhancement achieved by the MV method. One way proposed to reduce the computational complexity of the MV is the beamspace domain MV beamformer [7]-[8]. This method employs the fact that the signals received from the points far from the focal points are very small and negligible, as a result of the focusing applied on the transmission. In this way, a good approximate implementation of the MV beamforming is reached. In [9], the authors proposed a method to obtain a covariance matrix in a Toeplitz structure, based on which the inverse of the matrix can be computed using a fast algorithm. Both of these methods reduce the complexity to N , while performing a performance comparable to that of MV. Synnevag et. al. [10] have proposed a low complexity adaptive (LCA) beamforming method based on the MV beamforming idea with a complexity of N. The method is based on the evaluation of the beamformer output for a set of predefined apodization and applying the best one. The evaluation criterion was an approximate value of the output power acquired through a temporal averaging. The method gives significant improvement in image resolution compared with the DAS. To retain the complexity to N, they have ignored the spatial averaging and this restricts its capability to enhance the resolution.
2
Moreover, the window design should be carefully done and may be application-dependent. The LCA also shows some artifacts on the image because of the discrete solution space used. This paper describes a new method to reduce the complexity of the MV. This method named as decimated minimum variance (DMV) beamformer uses the same fact used in the beamspace MV [7]; but applying it in array domain, instead of transform domain (beamspace). The key feature of the method is the decimation of the array without significant degradation of the obtained image. This feature reduces the number of weights which are adaptively calculated; however the covariance matrix is estimated through using all data received by the array. This property distinguishes the new method from applying the MV on a sparser array. In this way, the DMV beamformer is an approximate implementation of the MV and is capable to reduce the complexity to N, even with spatial and temporal averaging. The paper is organized as follows. The next section describes the new method through mathematical expression. Section 3 is devoted to the way of estimating the covariance matrix and then some designing considerations are discussed in section 4. The simulation results are given in sections 5, including some discussions on resolution, contrast and robustness properties and the comparison with the DAS and MV methods. The performance of the method on the experimental data is presented in section 6 and concluding remarks are given in section 7. 2. The method In array imaging systems, echo signals from the focal point arrived at individual transducer array elements are coherent and those emanating from the points far from the focal point are incoherent. It means that in the spatial frequency domain, the low frequency components are related to the signals received from the focal point [11]. Assuming dynamic receive focusing in the ultrasound imaging; this fact is valid for each imaging point. Therefore, high spatial frequency components are thrown away to obtain a good viewing image. This filtering process can simply be performed by smoothing the array
3
signals through summation or a smooth weighting and then summation, as it is done in the conventional beamforming methods. The aforesaid fact can be utilized in the MV beamforming to reduce its computational complexity. Assume that the array signals are filtered by a spatial frequency low pass filter to obtain a low pass array signal. Considering this filter as a decimation filter, the filtered array signal can be down-sampled or decimated in space with no significant aliasing distortion which means the filtered data can be recovered from its decimated version. Therefore, the minimum variance beamforming can be applied on the decimated- filtered data instead of the full data of the primary array, while the performance is almost preserved. It is worth to note that the transmit focusing makes the desired signal be a narrow band in the spatial frequency domain. This feature permits a high value for the down sampling rate, even with no decimation filter. This is in agreement with the results reported by Nilsen et.al [7] where the beamspace adaptive beamforming through using 3 beams have shown an acceptable performance. Consider an array of N equally spaced elements, and assume that [ ] is the primary array signal at time k after applying delays for dynamic receive focusing and is the tap weight vector of a J-length spatial FIR filter used for smoothing the signal (J
[ ] = [ ]
(1)
where, F is the convolution matrix for vector and [ ] a = − + 1 length vector, is the filtered array signal. Now, assume that [ ] is decimated to a N -length vector called [ ], such that [ ] and [ ] are related through the following equation:
, [ ] = , [ ]
(2)
4
where, the second subscript represents the element number and = !
"
is the decimation factor, which
is assumed to be an integer. The element-wise relation in (2) can be expressed in a matrix form as:
[ ] = # [ ]
(3)
where, Id is a × matrix obtained through decimation on column of a × identity matrix:
[# ]%, = &
1 , ' = ( − 1 + 1 1 0 , *+ℎ- ./0-
(4)
The beamforming process is applied on the decimated data:
2[ ] = . 334 [ ] [ ]
(5)
where, . 33 is the apodization vector and y is the beamformer output. In the minimum variance beamforming approach, the apodization vector is calculated to minimize the mean power of output which is obtained through the following optimization problem:
/( 56[|2[ ]| ] = . 334 8 . 33 9 4 3 0:;<-=+ +*: . 33 1 = 1
(6)
3 is a N × 1 vector of one, assuming dynamic receive focusing, and 8 is the where, the steering vector 1 covariance matrix of the decimated data vector:
8 = 6[ 4 ]
(7)
The optimal vector . 33 can be expressed as follows:
. 33 =
3 8 1 3? 8 1 3 1
(8)
It is seen that the decimated MV is an algorithm with computational complexity for the matrix inversion, instead of in the standard MV. It is worth to note that if the main sources of the interference signals are assumed those near to the point of interest, the weights on the decimated array can be calculated such that the output power of the primary array is a scaled version of that of the decimated array (see Appendix A). As a result, the
5
minimum power of MV beamformer on the primary array can be reached through applying optimal weights in (8) on the decimated one. Hence, using MV on the decimated array will obtain an image similar to that obtained by the primary full array. To get a robust response, the MV beamformer usually uses a diagonally load covariance matrix. This technique can be utilized in the DMV and the covariance matrix 8 can diagonally be loaded by a factor @ as follows:
8 = 8 +
@ + A=- 8 #
(9)
3. The spatial and temporal averaging In MV beamformer, the usual way to estimate the covariance matrix through using the single snapshot data is the subarray averaging method. In this method, estimation is done by dividing the aperture into overlapping subarrays and averaging the covariance matrices of each subarray. For the decimated MV, the subarrays are built on the filtered array and then each subarray is decimated into elements, where averaging is done over these decimated subarrays. Assuming each subarray has length L and M is the length of the filtered array, P=M-L+1 subarrays can be used in averaging. Also, consider B as a × 1 vector containing the decimated data of p’th subarray. Then the decimated covariance matrix is estimated as: F
1 4 8C [ ] = E B [ ] B [ ] D
(10)
BG
This technique for estimating the covariance matrix can be considered as the decimation in columns and rows of the covariance matrix of the full array data (see Appendix A). In a similar way, the covariance matrix can be estimated through the temporal and spatial averaging:
6
J
F
1 4 8C [ ] = E E B [ − / ] B [ − /] D 2I + 1
(11)
KGJ BG
It is apparent that spatial (subarray) averaging is identical to the temporal- spatial averaging with K=0. The beamformer output, in both cases, is calculated as follows: F
1 2[ ] = E . 334 [ ]B [ ] D
(12)
BG
where, . 33 is given by (8). The computational complexity for the DMV, assuming spatial averaging is only applied, can be determined as follows: • Calculating the output of the filters requires MJ multiplications. • The estimation of the covariance matrix is done through P times multiplying an × 1 vector by its transpose conjugate (according to (10)), so it needs the D multiplications. • The computation of the inverse of the covariance matrix is . • The computation of . 33 and then calculating the beamformer output is done through + and D multiplications, respectively. Hence, assuming the inverse of the covariance matrix is available, calculating the beamformer output can be done with a computational complexity of D + . Therefore, the complexity for the DMV is + + D . As explained before, in the transmit focusing mode, the filters can be ignored without appearing serious aliasing distortion. Moreover, it is possible to use filters with the uniform weights. These show that the complexity for the DMV can be reduced to + D . By contrast, that for the MV is L + DL , and for the beamspace MV is M + DM + DLM , while the LCA needs NN multiplications, where M and N are the number of
7
beams and the number of windows used. This comparison shows that if is set to a small value (i.e. ≪ L), then the DMV outperforms the three others in the computational complexity aspect. 4. Designing considerations This section is devoted to describe some points useful for designing the decimation filter and determining the proper value for . As stated above, the decimation filter is mainly used to cancel out high spatial frequency portion of the received data, as well as it permits the signal be down-sampled without aliasing distortion. It also suppresses the noise in array signal and prevents the noise cause the aliasing distortion. Hence, if the noise on the received array signal is high, applying a proper decimation filter would be needed before down-sampling. The decimation filter may cause some degradation in the resolving capability of the MV beamformer. This degradation is mainly affected by the length of the decimation filter used, as its increased length is resulted in the shortened spatial length of the filtered array. For example, if a filter of length J=32 is applied on a 64-elements array, the length of the filtered array is approximately equal to the half of the primary array. This shortened length leads to the reduced resolving capability. Therefore, the decimation filter should carefully be designed to have a beamformer with a good resolution property and no significant aliasing effects. Nevertheless, the decimation filter may result in a better contrast property assuming a decreased sidelobe level of the decimation filter at the price of its mainlobe width.
8
The most important parameter which should carefully be determined is , as the main advantage of the DMV over the MV method is the reduced computational complexity, while the performance is comparable. It is clear that a high value of can produce an enhanced image at the price of complexity,
1
Full array Decimated array 0.8
Normalized Energy
Normalized Beampattern
1
0.6
0.4
0.2
0 0
0.2
0.4
0.6
0.8
0.9 0.85
0.8
0.75 0
1
M=32 M=64 M=128
0.95
5
10
15
20
25
30
Nd u Fig. 1. (a) One way beampattern of a 9-elements full array (solid line) and that of decimated to 3 elements (dash line), (b) the energy of 2-way beampattern versus for three values of number of array elements. The energy is normalized to the energy of full aperture array.
which is an undesired result. On the other hand, using a very small value of may result in a highly degraded image. To investigate the effects of on the performance of the proposed beamformer, we assume that the focal point in the transmission and the reception is the same. The filed pattern at the focal point is well described by the far filed beampattern [12]. Moreover, to simplify our discussions, the analysis is restricted to the continuous wave (CW) scenario because of its mathematical tractability. In fact, as the analysis is used to compare the performance of the DAS beamformer with the DMV, the conclusions would be valid in the practical situations where the pulsed wave (PW) excitation is applied. Furthermore, the array elements are assumed to be spaced by P/2 and apodized by a rectangular window in the transmission.
9
The steered response (i.e. the received signal by the array at various angles) is described by the twoway beampattern, which is the production of the transmit and receive beampattern:
RD: = RD? :. RDT :
(13)
Here, RD, RD? and RDT are 2-way, transmit and receive beampattern, respectively, and : = 0/(U where U is the observation angle. Since the rectangular apodization is assumed in transmission, the transmit beampattern of an array of M elements can be represented as follows:
RD? : =
sin Y:/2 sin Y:/2
(14)
Assuming the receive aperture be decimated to elements and with uniform weighting, the receive beampattern is:
RDT :; = where = !
sin Y :/2 sin Y :/2
(15)
. An exemplary transmit and receive beampattern is plotted in Fig. 1(a), for M=9 and
"
=3. It can be seen that the 2-way beampattern would have small magnitude at angles far from the main axis even the receive beampattern has high magnitude. To investigate the energy lost due to decimating the receive aperture, the energy of 2-way beampattern is calculated as:
6 = [ |RD? :. RDT :; | \:
(16)
and, after normalizing to E(M), is depicted in Fig. 1(b) versus for different values of M (M=32, 64 and 128). The figure shows that the energy received by the decimated aperture with =3 can be almost 84% of the full aperture and that for =5 is higher than 90%. Therefore, it can be concluded that decimating the receive aperture to 3 or 5 elements would not significantly lose the energy. Moreover, the difference between the performance of the decimated and full aperture may be eliminated by applying an adaptive beamforming method on the decimated aperture. It is because,
10
according to (14), the transmit beampattern has nulls at directions specified by sin Y:/2 = 0, except for u=0. This means that the first right and left nulls occur at : = ±2/ and the second at : = ±4/ and so on. Moreover, the first sidelobes are centered at : = ±3/ with a level of approximately M/4.7 and the other sidelobes would have levels that are monotonically decreased approximately by 1/u. These properties of the transmit beampattern can be observed in Fig. 1. This discussion emphasizes that the main source of interference is the first sidelobes and the energy interfered form the other sidelobes is usually negligible. Hence, to catch a qualified image, the receive beampattern should have small magnitude especially at the first sidelobe regions of the transmit beampattern. Consequently, the adaptive receive beamformer would show a performance comparable or better than that of DAS if the receive beampattern can apply at least two nulls at the first sidelobes. This implies that the receive beamformer should have three degrees of freedom to be capable of imposing two nulls while satisfying the condition of the distortionless response (the constraint in the optimization problem (6)). In the other words, a beamformer with three weighting parameters and an aperture of the same length as the transmit aperture can afford the mentioned conditions. Therefore, by setting = 3 and assuming the signal received from the sidelobes other than the nearest ones is negligible, the DMV beamformer would outperform the DAS. Correspondingly, setting = 5 enables the beamformer to cancel the interfering signals from the four nearest sidelobes and hence a better performance is reached. It is worth to emphasize that selected according to the above discussion is such that the DMV method outperforms the DAS in terms of resolution assuming that the main sources of interferences are the points at the vicinity of the point of interest; an assumption that is usually valid in the focused transmit mode. In this way, the above idea for determining can be used for every uniformly spaced array transducer such as a linear array with an element space of P.
11
5. Simulation results The proposed beamformer was tested on the simulated ultrasound data, obtained using Field II [13]. Two phantoms were simulated: one consisting of pairs of wire targets and the other was a circular cyst in a speckle medium. For the simulations, a 3.5MHz, 64-element linear array transducer with λ/2 spacing was used. Excitation pulse had a Gaussian envelop with -6 dB relative bandwidth equal to 0.5. Moreover, the transmit focus was fixed at distance 45 mm and dynamic focusing was applied in the reception. To fully evaluate the performance of the DMV method, ten different beamforming methods were simulated according to the following list: BF1: a DAS beamformer with rectangular weighting. BF2: a MV beamformer with subarray length of L=N/2. BF3: a beamspace beamformer with three central beams and subarray length of L=N/2. BF4: an LCA beamformer with 12 windows defined according to TABLE I of [10]. BF5: a DMV beamformer with = 3 without spatial filtering, i.e. J=1. BF6: a DMV beamformer with = 5 without spatial filtering, i.e. J=1.
12
BF2
BF3
BF4
BF5
30
30
35
35
35
35
35
40
40
40
40
40
45
45
45
45
45
50
50
50
50
50
55 60
55 60
55 60
Depth [mm]
30
Depth [mm]
30
Depth [mm]
30
Depth [mm]
Depth [mm]
BF1
55 60
55 60
65
65
65
65
65
70
70
70
70
70
75
75
75
75
75
80
80
80
80
80
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
Fig. 2. Simulated point targets; BF1: DAS with uniform weighting; BF2: the MV beamformer; BF3: Beamspace beamformer; BF4: LCA beamformer; BF5: the DMV with J=1 and =3. The images are displayed using 50 dB dynamic range.
BF7: a DMV beamformer with = 3 and a uniform weighting decimation filter with length of J=16. BF8: a beamformer similar to BF7 but with J=32. BF9: a beamformer similar to BF8 but without decimation (J=32 and = 33. BF10: a MV beamformer applied on the 32 central elements of the array with subarray length of 16. The subarrays length L were set to [N/2] for the MV beamformers and to [M/2] for the DMVs, where [ ] is the largest integer value not greater than . Moreover, the covariance matrices were diagonally loaded by a factor of 0.01 for all MV- based methods. A. Simulated Point Targets
13
0
0 BF1 BF2 BF3 BF4 BF5
-10
BF1 BF2 BF3 BF4 BF5
-10
-20
Magnitude [dB]
Magnitude [dB]
-20
-30
-30
-40
-40
-50 -50
-60 -60
-4
-2
0
2
-8
4
-6
-4
-2
0
2
4
6
8
Lateral distance [mm]
Lateral distance [mm]
Fig. 3. Steered responses of DAS, LCA and other MV-based beamformers at depths 40 mm (the left figure) and 70 mm (the right one). Transmit focus was at 45 mm and dynamic receive focusing was applied.
The point targets were 6 paired wires separated by 3 mm laterally at 6 different depths. Fig. 2 displays images of the point targets using the first five beamformers (BF1 to BF5). We see that all the MV-based beamformers (BF2, BF3 and BF5) as well as LCA beamformer (BF4) give significantly better resolution than DAS (BF1). In comparison to the LCA, the figure shows that the image produced by the DMV beamformer is more similar to that of MV. These properties are verified by looking at Fig. 3, where the steered responses at two different depths are depicted. Fig. 3(a) shows the responses at depth 50 mm and Fig. 3(b) shows those at depth 80 mm. It is obvious that the DMV beamformer is capable of resolving the two points at depth 80mm, where the DAS beamformer is not. Furthermore, the figures illustrate that the performance reconstruction of the DMV methods is depth dependent such that at depth 50 mm the accuracy of reconstruction of DMV is similar to that of MV, but at depth 80 mm, the advantage of MV over DMV is obvious. This property may be considered as a result of aliasing effect which causes the coherent portion of the array signal be increased. Although the amount of the aliased energy is generally low, it may be considerable at the vicinity of a high echoic reflector. This leads to a reduced resolving capability of two point reflectors at deeper distances, where the two near points are seen at a narrower angle and the aliasing distortion is appeared more. It seems that one way to maintain the resolving
14
capability for different depth's is to use a small value of (e.g. 3) at small distances and increasing it (e.g. 5 or 7) at deeper distances. More results of the DMV beamformers are displayed in Fig. 4, where the images obtained through using different parameters of and J are shown. The image of BF6 beamformer emphasizes that increasing from 3 to 5 can slightly enhance the DMV performance. The figure also indicates the effect of the length of the spatial filter on the obtained image: by increasing J, the capability of resolving two near point reflectors is decreased, while the amplitude of the image at points laterally far from the reflectors are diminished. As stated in the previous section, these properties are the results of the reduced spatial length of the filtered array and attenuated high spatial frequency contents of the signal. Furthermore, the figure shows that BF8 (the DMV with a spatial filter of length J=32) suffers from an artifact problem for deeper distances. This artifact is a peak response between two point reflectors and appeared as a non-real point in the image. It could not be considered as a result of decimation, as it is seen in BF9 (which is the same as BF8 but without decimation). This artifact is also appeared in the standard MV beamformer of BF10 and could be appeared everywhere the two point reflectors are such close that the MV cannot resolve them. As mentioned in the previous section, applying the low pass filter restricts the resolving capability of the beamformer and this capability is more reduced by increasing the length of the filter. This restricted resolution is responsible for the aforementioned artifact. This fact can be verified by Fig. 5, showing that the steered responses at depth 80 mm for beamformers BF8, BF9 and BF10 are very similar. In other words, the performance of a DMV with J=32 is similar to that of a MV on a 32elements array.
15
BF7
BF8
BF9
BF10
30
30
35
35
35
35
35
40
40
40
40
40
45
45
45
45
45
50
50
50
50
50
55 60
55 60
55 60
Depth [mm]
30
Depth [mm]
30
Depth [mm]
30
Depth [mm]
Depth [mm]
BF6
55 60
55 60
65
65
65
65
65
70
70
70
70
70
75
75
75
75
75
80 -5
80 0
5
-5
Lateral Distance [mm]
80 0
5
Lateral Distance [mm]
-5
80 0
5
Lateral Distance [mm]
-5
80 0
5
Lateral Distance [mm]
-5
0
5
Lateral Distance [mm]
Fig. 4. Simulated point targets; BF6: the DMV with J=1 and =5; BF7: the DMV with J=16 and =3; BF8: the DMV with J=32 and =3; BF9: the DMV with J=32 and =33; BF10: the MV on a 32-element array. The images are displayed using 50 dB dynamic range.
0 BF6 BF7 BF8 BF9 BF10
Magnitude [dB]
-10 -20 -30 -40 -50 -60 -10
-5
0 5 Lateral distance [mm]
10
Fig. 5. Steered responses of DMV beamformers (BF6 to BF9) and MV on a 32-element array (BF10) at depth 80 mm. Transmit focus was at 45 mm and dynamic receive focusing was applied.
16
BF2
50 -5
5
0
BF6
5
45
50 -5
Lateral Distance [mm]
0
5
Lateral Distance [mm]
Depth [mm]
Depth [mm]
0
50 -5
5
BF9
45
0
5
Lateral Distance [mm]
0
5
Lateral Distance [mm]
BF10
40
50 -5
45
Lateral Distance [mm]
40 Depth [mm]
Depth [mm]
45
50 -5
5
BF8
40
0
0
40
45
Lateral Distance [mm]
BF7
40
50 -5
45
50 -5
5
Lateral Distance [mm]
BF5
40
Depth [mm]
0
45
Lateral Distance [mm]
BF4
40 Depth [mm]
Depth [mm]
Depth [mm]
45
50 -5
Depth [mm]
BF3
40
40 Depth [mm]
BF1 40
45
50 -5
0
5
Lateral Distance [mm]
45
50 -5
0
5
Lateral Distance [mm]
Fig. 6. Simulated cyst phantom; BF1: DAS with uniform weighting; BF2: the MV beamformer; BF3: Beamspace beamformer; BF4: LCA beamformer; BF5: the DMV with J=1 and =3; BF6: the DMV with J=1 and =5; BF7: the DMV with J=16 and =3; BF8: the DMV with J=32 and =3; BF9: the DMV with J=32 and =33; BF10: the MV on a 32-element array. The images are displayed using 50 dB dynamic range.
B. Simulated Cyst Phantom The cyst phantom consisted of a 4 mm cyst located at depth of 45 mm in a speckle medium. The speckle pattern is simulated with 10 randomly placed scatterers within a resolution cell of λ3, where λ is the central wavelength of the propagating waveform in the medium, to ensure fully developed speckle. The scattering amplitudes are Gaussian distributed. For all MV-based beamformers, the estimation of the covariance matrix is done through spatial and temporal averaging, where for temporal averaging K is set to 10 at the sampling rate of 20 MHz. The images of the cyst target are shown in Fig. 6. The contrast parameters (CNR and CR) are listed in Table I. The CR is Contrast Ratio and defined as the ratio of the mean value in the background to the mean value in the cyst region [14] and the CNR is Contrast to Noise Ratio and defined as the CR divided by the standard deviation of image intensity in the background region [15]. It illustrates that except for the LCA beamformer (BF4) which gives the best contrast
17
properties, the DMV beamformers (BF5 to BF9) show better contrast than DAS and other MV-based beamformers. The Table also indicates that the contrast properties of the DMV are enhanced by increasing from 3 to 5, which is similar to that seen for the resolution properties mentioned before. Moreover, the Table shows the effects of spatial filtering on the contrast enhancement, which are the results of the reduced high spatial frequency of the signal, i.e. more attenuating the signals received from far sidelobes. TABLE I. CONTRAST PARAMETERS OF THE SIMULATED CYST PHANTOM FOR DIFFERENT
BF1 BF2 BF3 BF4 BF5 BF6 BF7 BF8 BF9 BF10
CR(dB)
CNR
34.1 34.9 34.3 41.7 36.5 38.3 37.9 38.6 38.7 32.1
4.9 4.9 4.4 5 4.6 4.8 4.7 4.6 4.8 4
C. Sensitivity to the sound speed inhomogeneities A potential weakness of an adaptive beamformer is the lack of robustness against errors in the assumed wavefield parameters [2]. To investigate the robustness of the DMV, we assumed that the acoustic velocity is known with an overestimation by 5%. The steered responses of DMV beamformers for the points at depth 40 mm of the point phantom are depicted in Fig. 7 and compared with those of DAS and other MV-based methods. It is seen that the LCA is the best adaptive beamformer in robustness aspect and the MV is the worst one. The figure also shows that the DMV methods estimate the amplitude of points near 6 dB better than the MV beamformer does. This looking strange result can be explained by examination of the MV properties. It is known that the MV beamformer cancels out the interference
18
0 BF1 BF2 BF3 BF4 BF5 BF7
Magnitude [dB]
-10 -20 -30 -40 -50 -60 -4
-2
0
2
4
Lateral distance [mm]
Fig. 7. Steered responses of DAS, LCA and other MV-based beamformers at depth 40mm with 5% overestimation of the acoustic velocity. Transmit focus was at 45 mm and dynamic receiver focusing was applied.
signals which are incoherent, while the coherent portion of the received signal is considered as the desired one and retained by the beamformer. But at the presence of the sound speed error, the signals received from the point of interest would not be at the same phase and the property of coherency would slightly be lost. As a result, the MV beamformer attenuates that signal. This property would be appeared more obviously when the resolving capability of the beamformer is increased, where a received signal with no complete coherency may be considered as the interference and be completely cancelled out. This means that the robustness of the MV beamformer is inversely related to the number of weights adaptively determined. This fact is verified by [2], where the authors showed that increasing the subarray length in applying MV leads to the more sensitivity against errors. A similar situation occurs in DMV when compared to the MV and hence, the DMV is more robust than MV. Fig. 7 also shows that increasing J from 1 to 16 in DMV beamformer results in decreased sensitivity to the velocity error. This decreased sensitivity is similar to that appeared in DMV compared to MV and can be explained by the same facts.
19
6. Experimantal Results The proposed beamformer was applied on two experimental ultrasound data acquired at the Biomedical Ultrasonics Laboratory, University of Michigan (http://www.bme.umich.edu/labs/bul) using Synthetic Aperture Focusing Technique (SAFT) [16]. The first data was for a phantom of 6 wires at different depths. It was a complete data set from 128 transmit-receive measurements using a 128element linear array with the center frequency of 3.5 MHz sampled at rate 13.8889MHz. To compare the performance of different beamformers and for the sake of simplicity, we utilized the data only on the 64 central elements of the array such that the number of elements in transmit and receive modes was 64. Moreover, the transmit array was apodized by a rectangular window and focused at depth 60mm, as well as dynamic focusing was applied in the receiving mode. To do focusing, the data were firstly upsampled by a factor of 10 and then focusing delays were applied. Because of high resolution of MV-based beamformers, aliasing distortion could be occurred if sufficient number of lines per sector was not existed. Hence, we used 400 lines for a sector of 90 degrees for generating the images. The images obtained through three different beamformers (DAS, MV and DMV) are displayed in Fig. 8. It illustrates that the images obtained through DMV and MV beamformers have better resolution properties than DAS. The second experimental data used to evaluate the proposed beamformer was for a heart phantom acquired by a 64-element linear array with center frequency of 3.333 MHz and 17.76 MHz sampling rate. The focusing operation was completely the same as that of the first data. Figure 9 displays the images obtained by different beamformers and shows that the beamformers have similar contrast properties, while the walls are better defined in MV and DMV compared to DAS.
20
30 40
20
50
40
60
60
70 80
80 90
100
100
120
110
-50
0
120 -50
50
0
50
0
50
0
50
30 40
20
50
40
60
60
70 80
80
90
100
100 110
120 -50
0
120 -50
50
30 40
20
50
40
60
60
70 80
80 90
100
100
120
110
-50
0
120 -50
50
Fig. 9. Images of the second experimental data using DAS (top figure), MV (the middle one) and DMV beamformer (the bottom figure). The images are displayed using 50 dB dynamic range. Horizontal and vertical axes are lateral and axial distances, respectively, both in millimeters.
Fig. 8. Images of the first experimental data using DAS (top figure), MV (the middle one) and DMV beamformer (the bottom figure). The images are displayed using 50 dB dynamic range. Horizontal and vertical axes are lateral and axial distances, respectively, both in millimeters. 21
7. Conclusions This paper applied a new method for implementing the MV beamformer in ultrasound imaging systems with focused transmit beam. This method referred to as the decimated minimum variance (DMV) beamformer, uses a decimated aperture data, instead of full aperture data to apply the apodization vector, while utilizing full data to estimate the covariance matrix. To reduce the aliasing effects caused by the decimation, some antialiasing filters can be applied before decimation. The most important feature of the new method is its reduced computational complexity compared to the MV beamformer. The simulations have indicated that the elements of the subarrays used in the DMV can be lowered to a smaller value such as 3, while the performance remains almost comparable to the full MV beamformer. While the MV shows better resolution and contrast properties in no modeling errors conditions, the DMV behaves more robust at the presence of errors. This fact denotes that further works are needed to compare the DMV with MV. The most advantage of DMV over MV is its complexity. Assuming focused transmit beam and with no decimation filter or with a filter with uniform weighting, the computational complexity of the DMV would be at N, compared to N for the MV method. Thus, the DMV is an adaptive beamforming whose complexity is in order of the complexity of the DAS, while outperforming the DAS in contrast and resolution properties. As the complexity is most restrictive feature for applying MV on a 2D array, using the DMV on a 2D array can obviate this problem and hence, the DMV is a good candidate for applying an MV-based adaptive beamforming on a 2D array.
Appendix A In this appendix, it is shown that under some reasonable assumptions, the beamformer output of a decimated array is approximately equal to the scaled output of the primary array. Moreover, the
22
covariance matrix of the decimated array is related to that of the primary array by the decimation. The first assumption is that the main sources of the interference signals are those near to the point of interest. The second one is that the signals arrived to the array are spatially stationary. This assumption results in a Toeplitz structure for the covariance matrix [9] and also permits analyzing it in spatial frequency domain. The mean power of the beamformer output can be written as:
6 [|2[
]| ]
1 e 1 e [ Φ d\d = [ Φf d |g d| \d = 2Y e c 2Y e
(A-1)
where Φc d and Φf d are respectively the Fourier transforms of the correlation functions of the y[k] and x[k], which called the power spectral density and gd is the Fourier transform of the applied weighting on the array data x[k]. By the first assumption, Φf d has negligible values outside the interval (−Y/h, Y/h for some values of h > 1 and hence: e
1 k [ Φ d |gd| \d 6 [|2[ ]| ] ≅ 2Y e f k
=
e 1 m m [ Φf l n og l no \m 2Y h e h h
(A-2)
If we set Φf d = hΦf d/h) and g d = hgd/h, (A-2) can be written as:
6 [|2[ ]| ] ≅
h e [ Φ d |g d| \d 2Y e f
(A-3)
Equation (A-3) states that the power at the beamformer output of a full array would approximately be equal to h times of a decimated array if two conditions are satisfied. The first one is that the weighting applied on the decimated array is such that its Fourier transform is an expansion of that of the full array: g d = hgd/h. Note that this condition does not necessarily mean that . [ ] is a decimated version of .[ ], because of possible aliasing which can be occurred in decimating .[ ].
23
The second condition is that the correlation function of input is the decimated correlation function of the primary input x[k] by a factor of L. It means that, the decimated array has an input of [ ] which is spatially decimated of x[n] and its covariance matrix is obtained by applying the decimation on the columns and rows of the covariance matrix of x[k]. In other words, the covariance matrix of [ ] can be calculated through spatial averaging of the decimated subarrays.
REFERENCES
[1] F. Vignon, and M. R. Burcher “Capon beamforming in medical ultrasound imaging with focused beams,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 55, no. 3, pp. 619-628, March. 2008. [2]
F. Synnevag, A. Austeng, and S. Holm, “Adaptive beamforming applied to medical ultrasound imaging,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 54, no.8, pp. 1606-1613, Aug. 2008.
[3]
I. Holfort, F. Gran, and J. Jensen, “Broadband minimum variance beamforming for ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 2, pp. 314-325, Feb. 2009.
[4]
B. Mohammadzadeh Asl, and A. Mahloojifar, “Minimum variance beamforming combined with adaptive coherence weighting applied to medical ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 9, pp. 1923-1931, Sep. 2009.
[5] B. Mohammadzadeh Asl, and A. Mahloojifar, “Contrast enhancement and robustness improvement of adaptive ultrasound imaging using forward-backward minimum variance beamforming,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 58, no.4, pp. 858-867, Appr. 2011. [6] X. Zeng, C. Chen, and Y. Wang, “Eigenspace-based minimum variance beamformer combined with Wiener postfilter for medical ultrasound imaging,” Ultrasonic., vol. 52, pp. 996-1004, 2012. [7]
C. C. Nilsen, and I. Hafizovic, “Beamspace adaptive beamforming for ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 56, no. 10, pp. 2187–2197, Oct. 2009.
[8]
X. Zeng, Y. Wang, J. Yu, and Y. Guo, “Beam-domain eigenspace-based minimum variance beamformer for medical ultrasound imaging,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 60, no. 12, pp. 2670–2676, 2013.
24
[9] C B. Mohammadzadeh Asl, and A. Mahloojifar, “A low-complexity adaptive beamformer for ultrasound imaging using structured covariance matrix,” IEEE Trans. Ultras., Ferroelec., Freq. Contr., vol. 59, no.4, pp. 660-667, Appr. 2012. [10]
J.-F. Synnevåg, A. Austeng, and S. Holm, “A low-complexity data-dependent beamformer,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 57, no. 2, pp. 281–289, Feb. 2010.
[11] P. S. Naidu, Sensor array signal processing, CRC Press, 2001. [12] J. L. Prince, and J. M. Links, Medical imaging signals and systems, Pearson Prentice Hall, 2006. [13]
J.A. Jensen, “Field: a program for simulating ultrasound systems,” Med. Biol. Eng. Comp. 10th Nordic-Baltic Conference on Biomedical Imaging 4 (Suppl. 1), 1996, pp. 351–353, Part 1.
[14] M. O’donnell and S. W. Flax, “ Phase aberration correction using signals from point reflectors and diffuse scatterers: Measurements,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 35, no. 6, pp. 768–774, 1988. [15] S . Krishnan, K. W. Rigby, and M. O’donnell, “Improved estimation of phase aberration profile,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control, vol. 44, no. 3, pp. 701–713, 1997. [16] M. Karaman, P. C. Li, and M. O'Donnell, “Synthetic aperture imaging for small scale systems,” IEEE Trans. Ultrason., Ferroelect., Freq. Contr. , vol. 42, pp. 429-442, May 1995.
25