Eigenspace generalized sidelobe canceller combined with SNR dependent coherence factor for plane wave imaging

Background The eigenspace generalized sidelobe canceller (EGSC) beamformer combined with a signal-to-noise ratio (SNR) dependent coherence factor (CF) is suggested for coherent plane wave compounding (PW) imaging. Conventional CF based methods such as generalized CF and subarray CF can improve the image quality, however, they are not suitable for low SNR. On the other hand, the EGSC CF based approach can introduce improvements in image quality, however, in PW imaging is susceptible to suffer from degradation due to low SNR which leads to a poor image quality. To overcome this limitation, the SNR dependent CF method is suggested for application in such situations due to its ability to control the SNR levels. Methods The Field II and the Verasonics ultrasound imaging system with a L11-4v array transducer with a contrast resolution phantom were used to capture the plane wave sequences of simulation and experimental data, respectively. The performance evaluation using full width at half maximum (FWHM), contrast (CR and CNR) and the speckle statistics by using the signal to noise ratio (SNR) complemented by the Rayleigh distribution analysis was performed. In order to evaluate the performance of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_{3}$$\end{document}EGSC3 (the SNR CF) beamformer, the comparison is done with particular importance to other CF-based approaches such as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_{1}$$\end{document}EGSC1 (the generalized CF) and, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_{2}$$\end{document}EGSC2 (the subarray CF) respectively. Results Taking DAS as reference, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_3$$\end{document}EGSC3 showed 30.3 and 39.5% of improvement for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {CR(dB)}$$\end{document}CR(dB) and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {CNR}$$\end{document}CNR, respectively, when using experimental data. The proposed method also slightly outperforms the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_1$$\end{document}EGSC1 and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_2$$\end{document}EGSC2 methods for \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {CR(dB)}$$\end{document}CR(dB), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {CNR}$$\end{document}CNR, and speckle statistics assessment. Conclusion The \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {EGSC}_3$$\end{document}EGSC3 is, therefore, suitable for CPWC by improving the spatial resolution and contrast while preserving the speckle pattern.


Background
Medical ultrasonic imaging is a noninvasive and low-cost technology widely used for diagnosis. Several techniques have been recently introduced for medical ultrasonic imaging in addition to the traditional Delay and Sum (DAS) beamforming method. DAS is a non adaptive beamforming method since it applies a fixed weight function (e.g., boxcar and humming) for array data summation. Adaptive weighting is, therefore, expected to result in increased image quality as introduced later. In addition, image reconstruction techniques have been proposed to increase frame rate as in plane wave imaging (PWI) [1,2]. In PWI the medium is illuminated using all the aperture elements simultaneously to create a plane ultrasound wave. The signal to noise ratio (SNR) is commonly smaller compared to DAS method because PWI uses a single or a few ultrasound plane wavefronts to interrogate the medium in order to create a frame, therefore reducing the ultrasound power per generated frame opposed to the focused firing in DAS. On the other hand, since frame rate is generally limited by the time of flight of the ultrasound signal, PWI can significantly contribute to increase the frame rate of the system [1,2].
Besides the lower interrogating power, contrary to DAS, in PWI the received signal in all transducer elements is composed by backscattered echoes from many different points from the interrogated medium resulting in lower SNR [1][2][3][4] being expected a tradeoff between the frame rate and the image quality in terms of spatial resolution and contrast [1,2]. The possibility of compounding a set of low-resolution images to get an image with improved spatial resolution and contrast is the basis of spatial compounding (SC) as described by Berson et al. [5] to reduce speckle noise. Additionally, Montaldo [1] proposed a SC method called coherent plane wave compounding (CPWC) for further improvements.
As mentioned earlier, data-dependent weights can be created to also improve image quality in terms of spatial resolution and contrast by maintaining the main lobe while reducing the side lobes [6][7][8][9]. One of the most representative adaptive method is the minimum variance (MV) beamformer proposed by Capon [10]. The aim of the MV beamformer is minimizing the output power of a beamformer while keeping the signal undistorted [7,8] and was previously evaluated into ultrasound imaging [6][7][8]. In MV-based methods, the covariance matrix (CM) estimate is one of the key steps which determine the performance of the algorithm [7]. Additionally, the adaptive weights are obtained at the cost of additional computation which comes from the need of inverting the CM [6,7].
In order to obtain a robust CM estimation, Sasso and Cohen-Bacrie [11] proposed the spatial smoothing method, Synnevag et al. [7] introduced the diagonal loading (DL) technique and, Asl and Mohammadzadeh [12] introduced forward-backforward method. Likewise, Mehdizadeh [13] proposed the Eigenspace-based MV (EMV) beamformer to improve the image quality. The EMV consists of projecting the MV weights onto the signal subspace where a considerable part of the noise is removed from the received data. The signal subspace is constructed from the eigen decomposition of the CM. This procedure consists of sorting the eigenvalues of CM such that the eigenvectors associated to the larger eigenvalues determine the basis for the signal subspace [13,14].
Another robust representation of MV beamformer is the generalized sidelobe canceller (GSC) beamformer proposed by Applebaum and Chapman [15] and later popularized by Jim and Griffiths [16] in the area of array signal processing. The GSC has been implemented in ultrasound imaging and compared to DAS and MV beamformers [17][18][19][20]. Similarly to EMV compared to MV, projecting the GSC weight onto the signal subspace lead to the eigenspace-based GSC (EGSC) beamformer implemented by Aliabadi et al. [19] resulting in improved image quality at increased computational cost when compared to GSC beamformer.
A coherence factor (CF) may also be combined into adaptive beamformers such as MV and was previously used for aberration correction and sidelobe suppression in ultrasound imaging [21].
Li and Li [22] have proposed the implementation of generalized CF (GCF) approach in MV beamformer to improve imaging resolution by protecting the main lobe. To improve the GCF performance, Zhau et al. [23] have suggested the subarray CF (SCF) approach. The SCF has optimized the performance of the EMV beamformer [23] at the cost of increased computational complexity. As in GCF, the SCF results from a combination of the incoherent sum and the coherent sum of signal energy [24,25]. To address the limitations of GCF and SCF, Wang et al. proposed the SNR dependent CF (SNR-CF). In the SNR-CF method, a weighting scaled factor determined by the sigmoid function is combined with the coherent and incoherent summations of signal energy in order to control local levels of SNR [23]. Due to the SNR-CF features, the effects of signal cancellation at low SNR were minimized while the effects of self-cancellation were controlled at high SNR, which may introduce improvements in spatial resolution and contrast [24].
Considering SNR-CF based method is suitable for applications in signals with low SNR and recognizing that PWI typically presents low SNR [23], the SNR-CF method could aggregate some benefits [24,25] into PWI-based methods. In this work, we suggest a combination between the EGSC with SNR dependent CF ( CF 3 ) method to form the ( EGSC 3 ) beamformer. For comparison, the GCF ( CF 1 ) is also combined with EGSC to obtain the ( EGSC 1 ) beamformer while the SCF ( CF 2 ) combines with EGSC to form ( EGSC 2 ) beamformer. The tests were performed using simulated and phantom data.

Methods
The coherent plane wave compounding principle is firstly presented followed by evolved beamforming methods including different adaptive techniques such as MV, EMV, GSC, EGSC, EGSC 1 , EGSC 2 and EGSC 3 . These techniques are the basis for the proposed SNRCF-based EGSC method and are later used for beamforming performance assessment.

Coherent plane wave compounding
In Eq. (1), we assume that a sampled time series X(k) is represented in a measurement grid as X CPWC (x, z) whose geometrical entities are depicted in Fig. 1, is measured for each element in a linear array where, M ( r = 1, . . . , M ) is the element number, N ( i = 1, . . . , N ) is the number of wave emissions, w represents the receive apodization window, u is the angular apodization, h ir is the impulse response of element r to plane wave i. For each time sample, τ i = (di + dr)/c represents the total travel time for ith wave where d i = zcosα i + xsinα i , is the distance from wave emission to point (x, z), dr is the distance from point (x, z) to the receive rth array element and, c is the sound speed set to be 1540 m/s [26].
The angle sequence α can be formulated as in Eq. (2) where, is the wavelength of the system and, A represents the probe aperture [26].
Individual signals collected in Eq. (1) follow the array geometry representation in Fig. 1 and can be represented in the form of a 2-D data matrix X(k) as in Eq. (3), where k = 1, . . . , K . In this context, K represents the total amount of time samples recorded in each array element, where, x ir (k) , is the received signal in element r for the ith wave emission. In CPWC, the beamformer output Eq. (4) is obtained by averaging the data in Eq. (3) as follows: x ir (k) In CPWC with a sequence of tilted plane waves, the angular and the receive apodization effects are taken into account in order to form Eq. (4). Additionally, in the conventional adaptive beamforming, the receive apodization window is placed by the data-dependent weight [4].

Minimum variance
In MV beamformer (BF) the weight is found by minimizing the Power ( w H Rw ) of the BF output (s.t) subject to the directional constraints as follows: is the expectation operator and, [·] H is the Hermitian operator, respectively.
The steering vector a (i.e., a vector of ones) follows a structure designed for the cases when time-delay steering is used such that the desired signal is approximately in phase at the steered outputs [7,16]. The solution of the minimization problem in Eq. (6) using the Lagrangian multipliers [7,9] is: In accordance with the presented in Eq. (5), the MV BF aims to minimize the output power of the output signal while keeping it undistorted [7,9]. In practice, the covariance matrix (CM) R Eq. (7) is estimated from the data by averaging in the temporal domain, averaging in the spatial domain or, the combining temporal and spatial domains [7].
While the former is performed taking into account F = 2T + 1 , where T is the number of temporal samples used in covariance matrix (CM) estimation, the latter considers the G = M − L + 1 overlapping subarrays X l (k) whose length is limited to ( L ≤ M/2 ) [7].
Subarray averaging aims to decorrelates the coherence between ultrasound signals and consists of dividing the dataset into overlapped submatrices as follows: X l (k) = [x l (k), x l+1 (k), . . . , x l+L−1 (k)] T and, (·) T denotes the transpose. After performing this procedure, the CM estimate R, can be subject to diagonal loading (DL) method in order to get a robust estimate [7]. The formulation of CM estimate in Eq. (7) combines the previous aspects as follows [7]: where µ = 1/(�L)tr R , R corresponds to the first term of Eq. (7), tr is the trace and, stands for a DL factor [7,27]. The value of the DL factor adopted for data processing is included in Table 2.
In minimum variance (MV) BF, the output in Eq. (8) is obtained by applying a set of weights to the received signals as follow.

Eigenspace minimum variance
In the EMV beamformer, the eigen analysis of the CM is performed in order to construct the signal and the noise subspaces as in Eq. (9) [13]: where R s and R n are the signal and the noise the CM , � = diag 1 , 2 , . . . , L are the eigenvalues in decreasing order such that 1 ≥ 2 ≥, . . . , ≥ L and U = v 1 , v 2 , . . . , v L are the eigenvectors v i ( i = 1, 2, . . . , L ) of i . The subspace estimation is obtained based on a threshold ( th ) such that U s = v 1 , v 2 , . . . , v th is the basis for the signal subspace (i.e. for the largest eigenvalues) and U n = v th+1 , v th+2 , . . . , v L makes the basis for the noise subspace for the remaining (smaller) eigenvalues [13]. Moreover, the number of eigenvalues used to construct the signal subspace is data dependent and varies in accordance with the amount of signal energy attained at a specific point. The largest eigenvalues used to identify the eigenvectors for constructing the signal subspace are obtained as follows: i ≥ α th 1 , for ( α th > 1 ) or i ≤ β th L for ( β th < 1 ), where the i is the eigenvalue corresponding to eigenvector used to construct the signal subspace [13,19,28].
The EMV weight w EMV is found by projecting the MV weight w MV onto the signal subspace ( P s = U s U H s ) as in Eq. (10) [13].

The generalized side lobe canceller
The weight of GSC BF is formulated by minimizing the output power of adaptive beamformer [15,16,19]. The GSC weight in Eq. (13) is found by partitioning the optimal weight of Eq. (6) in terms of orthogonal (adaptive) branch as in Eq. (11) and the nonadaptive branch as in Eq. (12), respectively [16]: where H = 2B H RB is the Hessian matrix (HM) [17], B the blocking matrix (i.e., Eq. 14) of dimensions L × ( L − 1) and the steering vector a with 1 × L such that: a H B = 0 [19,27].

The eigenspace GSC
The EGSC BF is found by projecting w GSC onto the signal subspace [19,27].

Generalized coherence factor (GCF) and subarray-based coherent factor (SCF)
The ratio between the coherent sum (CS) and incoherent sum (IS) is termed as coherence factor (CF) [22]. The GCF or CF 1 is formulated by combining the CS and IS as follow [24]: Accordingly, the SCF of CF 2 is formulated as in Eq. (17) [23,24].

The SNR dependent CF
One problem of using coherent factors is that the imaging system can give suboptimal results because in such applications the noise lowers the signal coherence and thus the noise suppression capabilities of other CF-based methods are compromised. Different to the CF 1 and CF 2 methods, the CF 3 method takes into account the local SNR in the CF formulation so that contrast of the imaging system can be restored despite the low SNR [23,24]. In practice, such local SNR values are averaged for half or one wavelength to improve the robustness [24]. In order to compute the weighting value η(SNR) = η(P s , P n ) the reference ratio is formulated as follows [24]: where L represents the subarray length.
The values defined on Eq. (18), range from 1/L to 1 allowing control the SNR levels while performing adaptive processing avoiding signal cancelation at low SNR while controlling self-cancelation at high SNR. A sigmoid structure for different values of α snr (see Fig. 2a) and β snr (see Fig. 2b) is depicted. For data processing, values of α snr and β snr are presented in Table 2.

Proposed SNR CF-based EGSC
In this work, the core of the implementation of the CF-based methods is the EGSC beamformer. More emphasis is given to the combination with the SNR CF ( CF 3 ) since it appears to offer better performance over low SNR data [23,24]. Moreover, the EGSC ideally provides the better performance than the DAS, MV, EMV and GSC [27] beamformers, so that the CF-based methods are suggested to be combined with the EGSC beamformer to obtain the EGSC-CF based beamformers. In addition, the (SNR) coherence factor (CF) method ( EGSC 3 ) is then compared with the following categories of CF based methods: the GCF ( CF 1 ) and the SCF ( CF 2 ) [23,24], respectively.

Implementation procedure
The key steps for implementation procedures are summarized in Table 1. This process is repeated for each imaging point to get the final image which is finally displayed in dynamic range of 60 dB.

Evaluation metrics
The evaluation metrics for resolution analysis were FWHM and PSL while CR, CNR, SNR, and the theoretical Rayleigh distribution (TRD), were used for image contrast and speckle assessment as follow.

Full-width half maximum (FWHM) and peak side lobe levels (PSL)
The FWHM [29] and the PSL [27] values are computed for evaluation of spatial resolution and interference. The former defines the beam width of the main lobe at − 6 dB while the latter defines the level of the first side lobe level for evaluating the interference and noise suppression abilities of a beamformer. Normally, the evaluation of FWHM and PSL is complemented by plots of normalized profiles and, improvements are presented by the narrower beam width and lower PSL [29]. (20)

Contrast ratio (CR), contrast to noise ratio (CNR) and, signal to noise ratio
For anechoic cyst and hyperechoic target to background evaluation, the contrast ratio (CR) and contrast to noise ratio (CNR) was assessed using Eqs. (24) and (25) [23]. In this context, an increase of CR and CNR indicates improvements in image quality.
where cyst , is the mean amplitude (before log-compression) in the cyst and bck is the mean amplitude in the speckle, σ 2 cyst variances of intensities inside cyst and σ 2 bck in the background, respectively. In Eq. (25), and σ respectively the mean value and standard deviation. The fully developed speckle is approximated by the RD which is expressed by the probability density function Eq. (26) for the envelope amplitude in an image as follows [29]: where 2 is the mean of the squared amplitude.
Since the speckle pattern of DAS beamformer follows closely the RD, the TRD with SNR = 1.91 [7,29] was included for illustration in Fig. 3 where the TRD together with the histogram of echo brightness (amplitudes) for simulated data are presented.
One can examine if the speckle pattern of the beamformed data follows the RD by performing the hypothesis test (K-test 5% significance). The test is decided when the computed p-value is compared with the maximal extreme of the significance interval [23]. Furthermore, a decrease of the amplitude level of the histogram indicates the speckle degradation effect [29].

Simulations and experiments
In order to test the applicability of the proposed method in CPWC imaging, we used the simulated and phantom dataset. The simulated data were acquired using simulated phantoms in Field II program for individual point scatters and the circular anechoic cysts [30]. A 128-element linear array transducer model L11-4v with center frequency of 6.25 MHz and fractional bandwidth of 60% was simulated. A two-cycle sine wave excitation pulse at the central frequency and sampling frequency of 40 MHz were used.
A set of point targets were simulated in a homogeneous medium. The point target simulation phantom with fifty scattering points was created. The points were laterally distributed (x = 0 mm, x = ±5 mm and x = ±10 mm) from an axial depth of (z = 20 mm) to (z = 45 mm) spaced 5 mm apart from each other.
Two circular anechoic cysts of 5 mm radius were simulated in a speckle pattern. The anechoic cysts with a diameter of 10 mm centered at lateral positions ( x = +10 ) and ( x = −10 mm) at the axial position of (z = 42.5 mm) were placed in a simulated homogenous media with 80,000 scattering points.
For experiments, the contrast detail resolution ultrasound phantom ATS model 532A containing several circular anechoic cysts and hyperechoic targets with different values of contrast was used. In accordance with Fig. 4, we decided recording data positioned at 8 mm diameter along the phantom.
Data corresponding to a pair of closely set anechoic cysts (label A) with − 12 and − 9 dB and another pair of hyperechoic cysts (label B) with + 9 and + 12 dB (see Fig. 4) were acquired. The background (speckle) and the regions inside the cysts have been used as the reference to compare the performance of the different beamformers.
A 128 element linear array transducer working at 6.25 MHz central frequency and pitch of 0.3 mm was used. For data processing, the sampling frequency was set to be 25 MHz. Verasonics imaging system using an L11-4v (Verasonics Ltd, Kirkland-WA, USA) transducer was used to acquire data with a sampling frequency of 25 MHz and transducer center frequency of 6.25 MHz. For both simulation and the experimental studies, all the array elements were used to emit the plane waves and record the RF data. Furthermore, either in emission or the reception, the acquisition system was not apodized.
Additionally, for simulation and the experiments, the plane waves were emitted with different steering angles for multiple scanning. All the transmitting steering angles were set to range from −5 • to 5 • (21 angles with an interval of 0.5 • ) and, the f-number of 1.75 was used.

Results and discussion
The simulated and phantom dataset was used to evaluate the performance of the proposed beamformers.

Simulation: point targets
Fifteen point targets were simulated in order to test the performance of different beamformers as shown in Fig. 5 with a dynamic range of 60 dB. The point targets are respectively located at depths of 20 mm, 25 mm, 30 mm, 35 mm and 40 mm. The lateral profiles for FWHM and PSL assessment are shown in Fig. 6 using the point target highlighted in green box at (x, z) = (0, 30) mm in Fig. 5a as reference. Figure 6 presents the improvements introduced by adaptive techniques over the traditional DAS method. Figure 7 presents the 0 to − 6 dB region in order to better the visualize the lateral resolution improvement.
From displayed images, it can be seen that DAS has a wide main lobe and higher side lobes compared to the MV beamformer. The EMV provides a slightly narrower main lobe compared to the MV beamformer.
Additionally, the EMV presents the lower sidelobe levels compared to MV beamformer meaning that the EMV outperforms the MV beamformer in terms of lateral resolution by providing a lower FWHM and reduced PSL, respectively.
The GSC technique presented slightly narrower main lobe and lower side lobe level compared to EMV beamformer. The EGSC presented a narrower main lobe while lowering down the side lobe energy compared to GSC beamformer, which justifies the welldefined scatters on the displayed images. Moreover, the CF-based methods came with remarkable improvements in image quality so that by combining the CF 1 with EGSC resulted in the EGSC 1 beamformer which outperforms the EGSC beamformer in terms of FWHM and PSL, respectively. Analogously, by combining the CF 2 with EGSC we formulated the EGSC 2 beamformer which appears to outperform EGSC 1 . Finally, on combining the CF 3 with EGSC we formulated the EGSC 3 which in turn, outperforms EGSC 2 , respectively.
Beyond the analysis of the displayed images, all the FWHM and PSL values are presented in Table 3. Remarkably, the EGSC 3 presented the improved lateral resolution by providing a narrower FWHM and lower PSL compared to the EGSC 1 and EGSC 2 so that the margins of difference in PSL (dB) were of 19.4 and 3, respectively. However, the FWHM was not as highly improved as the PSL meaning that the CF 3 appears to suppress more side lobes rather than narrowing the main lobe compared to the CF 1 .
Beyond the analysis of the displayed images, all the FWHM and PSL values are presented in Table 2. Remarkably, the EGSC 3 presented the improved lateral resolution by providing a narrower FWHM and lower PSL compared to the EGSC 1 and EGSC 2 so that the margins of difference in PSL were of 19.4 and 3 dB, respectively. However, the FWHM was not as highly improved as the PSL meaning that the CF 3 appears to suppress more side lobes rather than narrowing the main lobe compared to the CF 1 .

Simulation: circular anechoic cysts
The simulated speckle pattern containing two circular anechoic cysts closely spaced was used for contrast and speckle statistics evaluation. For CR, CNR and SNR calculation, the background region marked in a white box and the interior of the cysts marked with green as shown in Fig. 8a were used for calculating the mean intensity in the background and inside the cyst, respectively. The responses for different techniques are displayed in Fig. 8.
In Fig. 8, it is shown that DAS, MV and GSC have a poor contrast due to higher side lobes level while the EMV and the EGSC can provide a better contrast. The EGSC with all versions of coherent factors has a relatively higher contrast and also exhibits clearly the cyst edge. Additionally, we have measured the speckle statistics for the different beamformers as presented in Fig. 9 complemented by different value of SNR presented in Table 4. For this purpose, the normalized pdf of the speckle region (see Fig. 8a white box) was computed for different beamformers.
Since the speckle pattern of DAS beamformer follows closely the TRD [7,29], was included for illustration. All the beamformers were subject to the hypothesis test separately and the results showed that they followed the RD [7,23,29].
To complement the displayed images, the results are presented in Table 4. From Table 4, we can see that all the adaptive beamformers outperformed the DAS beamformer in terms of CR (dB)/CNR with values of 1.5/0.08, 2.8/0.11, 3.0/0.16, 4.2/0.93, 6.1/1.07, 7.6/1.10 and, 9.1/1.12, were respectively obtained. These values represent  . 9 The speckle statistics of simulated data of Fig. 8 Fig. 9, the MV, EMV, GSC and EGSC speckle pattern appears to be somewhat damaged compared to that presented by DAS so that the plots of the speckle patterns are slightly disparate. In GSC method, a slightly dark image can be seen in Fig. 8d which justifies the reduced SNR value, however, the EGSC 1 , EGSC 2 and EGSC 3 , respectively presented an increased SNR and hence, their corresponding plots of speckle pattern in Fig. 9 exhibited that their intensity levels in dB increased. This effect was in agreement with the SNR values presented in Table 4 in the sense that, among other beamformers, they were more similar compared to DAS beamformer. Table 7 presents the speckle magnitude of simulated data whose plots are shown in Fig. 9. The large values of speckle magnitude degradation for Fig. 9 in dB were − 20.0, − 16.2, − 17.3, − 17.3 for MV, EMV, GSC and EGSC respectively, contrary to − 14.8, − 15.7, − 13.8 for EGSC 1 , EGSC 2 and EGSC 3 , respectively. We notice that the speckle pattern of the EGSC-CF based methods is colse to that provided by DAS beamformer while providing improved lateral resolution and contrast.

Experiments: phantom circular anechoic cyst
The cysts illustrated in Fig. 4 are ordered such that the rightmost cyst has − 9 dB contrast whereas, the cyst on the left has − 12 dB contrast, respectively. The images for the anechoic cyst using the various techniques are shown in Fig. 10.
For CR, CNR and SNR calculation, the background has been marked in a white box while the interior of the cysts was marked with the green circles as shown in Fig. 10a. These positions were used as the references for calculating the mean intensity in the background and inside the cyst, respectively. The DAS, MV, and GSC exhibited a poor contrast due to higher side lobes level while the EMV and the EGSC can provide a better contrast than DAS, MV, and GSC, respectively, but the visibility of the edge of the cyst is still difficult. However, the EGSC with all versions of CF outperforms in terms of contrast and exhibits clearly the improvements of the cyst definition. Taking DAS as reference, the CR(dB)/CNR presented the following improvements 1. Additionally, we have measured the speckle statistics for the different beamformers as presented in Fig. 11 complemented by different value of SNR presented in Table 5. For this purpose, the normalized pdf of the speckle region (see Fig. 11a white box) was computed for different beamformers. Similarly to the simulated data, all the beamformers passed the hypothesis tests.
On observing Fig. 11, the MV, EMV, GSC and EGSC speckle pattern appears to be somewhat damaged compared to that presented by DAS so that the plots of the speckle patterns are slightly different. In GSC method, a slightly dark image can be seen in Fig. 10d which justifies the reduced SNR value, however, the EGSC 1 , EGSC 2 and EGSC 3 , respectively presented an increased SNR and hence, their corresponding plots of speckle pattern in Fig. 11 exhibited that their intensity levels in dB increased. This effect was in agreement with the SNR values presented in Table 5 in the sense that, among other beamformers, they were similar to DAS beamformer. Table 7 presents values of plots depicted in Fig. 11. The large values of speckle magnitude degradation for Fig. 11 in dB were − 21.3, − 20.0, − 19.0, − 16.6 for MV, EMV, GSC and EGSC respectively, contrary to − 15.2, − 14.4, − 13.6 for EGSC 1 , EGSC 2 and EGSC 3 , respectively. We remark that the speckle pattern of the EGSC-CF based methods is close to that provided by DAS while providing improved lateral resolution and contrast similarly to the simulated data. Fig. 11 The speckle statistics of real data of Fig. 10 for different beamformers. The TRD fit to DAS is included as the reference Page 18 of 23 Zimbico et al. BioMed Eng OnLine (2018) 17:109 Experiments: phantom circular hyperechoic targets Circular hyperechoic targets were also assessed and the images are shown in Fig. 12. For CR, CNR and SNR calculation, the background and the interior of the cysts was marked as shown in Fig. 12a.
The beamformed images reveal that the DAS, MV and GSC shows a poor contrast due to the higher side lobes level while the EMV and the EGSC provide an improved contrast compared to MV and the GSC.
The imaging contrast performance is quantified by using the contrast ratio (CR) and the analysis of the beamformed image data is done in terms of the quantitative values of intensities inside the target, intensities outside target, the CR, and the SNR as presented in Table 6.
Taking DAS as reference, the CR(dB)/CNR presented the following improvements 2.  Additionally, the speckle statistics was evaluated as presented in Fig. 13 and Table 6. Similarly to the simulated data, all the beamformers passed the hypothesis tests. On observing Fig. 13, the MV, EMV, GSC and EGSC speckle pattern appears to be somewhat damaged compared to that presented by DAS so that the plots of the speckle patterns are slightly distinct.
In GSC method, a slightly dark image can be seen in Fig. 11d which justifies the reduced SNR value, however, the EGSC 1 , EGSC 2 and EGSC 3 , respectively presented an increased SNR and hence, their corresponding plots of speckle pattern in Fig. 13 exhibited that their intensity levels in dB increased. This effect was in agreement with the SNR values presented in Table 5 in the sense that, among other beamformers, they were relatively close compared to DAS beamformer. Table 7 present values of plots depicted in Fig. 13. The large values of speckle magnitude degradation for Fig. 13 in dB were − 18.3, − 18.1, − 19.8, − 17.0 for MV, EMV, GSC and EGSC respectively, contrary to − 16.4, − 15.8, − 14.8 for EGSC 1 , EGSC 2 and EGSC 3 , respectively. We can see that the speckle pattern of the EGSC-CF based methods is close to that provided by DAS beamformer while providing improved lateral resolution and contrast.  . 13 The speckle statistics of real data of Fig. 11 for different beamformers. The TRD fit to DAS is included as the reference

Computational complexity analysis
The main purpose of CPWC is to reduce the computational complexity (CC) aiming to increase the frame rate of the imaging system hence, the analysis of the CC is an important task. The indicative CC required to perform the DAS BF goes under O(M) where M is the array length set to be 128 and hence, it will need O(M) floating operations.
In adaptive methods the CC increases because the inversion and eigen decomposition operations of the CM. The CM inversion for MV require 2/3 × L 3 floating operations when applying Gauss-Jordan eliminations [4,31] while the Hessian matrix inversion in Eq. (11) for GSC [17] needs L 3 floating operations. Note that the HM is more demanding than the CM. However, EMV and EGSC undergo the eigen decomposition of CM which in accordance with the Golub-Reinsch algorithm [32], need 21 × L 3 floating operations. Hence, the EMV will need 2/3 × L 3 + 21 × L 3 while the EGSC L 3 + 21 × L 3 = 22 × L 3 [4,17,31].
For GCF computation in (16) an additional computational burden lead to the computation of the coherent and the incoherent energy [23] which is L 2 . In the formulation of SCF ( CF 2 ) Eq. (17) and the SNR-CF ( CF 3 ) Eq. (19), the difference of the IS and CS is weighted using the scaling factors of 1/L and Eq. (18), respectively. In SNR-CF ( CF 3 ) formulation, the weight η(SNR) Eq. (18) is more complex compared to 1/L in SCF ( CF 2 ) formulation so that the computation of CF 3 will be more complex than CF 2 .
In order to examine the processing flow, we have measured the time cost (in seconds) for each technique without any special implementation optimization in MatLab. The results are presented in Table 8 in terms of (mean ± standard deviation) reflecting the averaged time over 10 tests for simulated and the real data respectively, including the sizes of the reconstructed images. Computational time was measured for total image formation (i.e. for 21 steering angles for one image). For example, in accordance with Table 8, the proposed method increased the time cost compared to DAS , MV , EMV , GSC , EGSC , EGSC 1 and EGSC 2 , respectively. For EGSC 3 compared to EGSC 1 and EGSC 2 the time in percentage (%) increased with the following values 7.0, 5.0 and 11.6, 6.5 for point targets and cyst simulation and, 13.5, 11.0 and 11.8, 8.2 for real data involving the anechoic cysts and the hyperchoic targets, respectively.
Among the CF based methods, the processing flow indicated values of time somewhat close one another however, there was an agreement with the indicative CC analysis. In general, the CF-based methods increase the computational complexity so that the EGSC 3 method will be more time demanding than EGSC 2 , and EGSC 2 more time consuming than the EGSC 1 technique. Therefore, the improvements in image quality introduced by the EGSC-CF beamformers were at the cost of an extra computational load compared to the different adaptive methods. For real-time applications, the most significant computational complexity added (i.e., the matrix inversion computation of the CM or HM) can be executed applying appropriate algorithms and architectures such as recursive updating [14], and graphics processing units [4], respectively.

Conclusion
The eigenspace generalized sidelobe canceller EGSC beamformer combined with the SNR dependent coherence factor SNRCF method has been proposed for coherent plane wave compounding imaging. The technique, here called EGSC 3 beamformer was compared with various adaptive beamformers (i.e., MV, EMV, GSC and EGSC) and other coherence factor based adaptive beamformers such as generalized coherence factor (GCF) and subarray coherence factor (SCF) using simulated and experimental data.
Taking DAS as reference, EGSC 3 showed improvements of 30.3 and 39.5% for CR (dB) and CNR , respectively, using experimental data. As compared to EGSC 1 and EGSC 2 , improvements of 2.9 and 1.1% for CR (dB) and 0.27 and 0.16 in CNR were found, respectively. Even thought the image quality improvement compared to other EGSC methods is small, the values of speckle statistics of EGSC 3 outperforms the EGSC 1 and EGSC 2 methods and remained close to DAS contrary to the remaining adaptive beamformers. The EGSC 3 is, therefore, suitable for CPWC by improving the spatial resolution and contrast while preserving the speckle pattern different from those beamforming methods that do not use coherence factors.