Journal of Medical Ultrasonics (2020) 47:179–192
https://doi.org/10.1007/s10396-019-01004-7

ORIGINAL ARTICLE–PHYSICS & ENGINEERING

Generalized coherence factor estimated from real signals
in ultrasound beamforming
Masanori Hisatsu1,2· Shohei Mori3 · Mototaka Arakawa2,3 · Hiroshi Kanai2,3Received: 18 October 2019 / Accepted: 5 December 2019 / Published online: 31 January 2020
© The Japan Society of Ultrasonics in Medicine 2020

Abstract
Purpose The generalized coherence factor (GCF), an adaptive beamforming technique, can reduce unnecessary signals froman unfocused position without reducing the contrast-to-noise ratio. However, the computational complexity of this method
is large compared to the conventional delay-and-sum (DAS) beamformer. In the present paper, we propose a novel method
to achieve the same reduction effect of unnecessary signals with a smaller computational load than that of the conventional
GCF approach.
Methods One of the factors increasing the computational complexity of the GCF-based beamformer compared with DAS
is the generation of analytic signals at receiving elements. We clarified the mechanism of generating unnecessary signal
components to enable the calculation of the GCF value directly from real signals without generating analytic signals. Furthermore, we proposed a method to filter out these components without generating analytic signals.
Results The GCF values obtained using the proposed and conventional methods were compared and verified using theactual data acquired from a phantom with an ultrasound diagnostic system. We also compared the B-mode images. As aresult, equivalent GCF values and similar B-mode image quality were achieved with the proposed method with reduced
computational complexity.
Conclusion With the proposed method, generation of analytic signals at receiving elements can be omitted, and as a result,the computational load of the GCF method can be greatly reduced, while preserving the effect of reducing unnecessary
signals like with the conventional method.
Keywords Ultrasound imaging · Adaptive beamforming · Generalized coherence factor

Introduction
Ultrasound imaging is widely used for medical diagnosis
purposes as it provides real-time cross-sectional images of
the human body noninvasively. In conventional ultrasounddiagnostic systems, the received signals are focused at any
position in the B-mode image using delay-and-sum (DAS)
beamforming. However, even with the use of DAS beamforming, reflected signals from positions other than the* Masanori Hisatsu
masanori.hisatsu.uc@hitachi.com
1

Healthcare Business Unit, Hitachi, Ltd., 3‑1‑1,
Higashikoigakubo, Kokubunji, Tokyo 185‑0014, Japan

2

Graduate School of Biomedical Engineering, Tohoku
University, Sendai, Miyagi 980‑8579, Japan

3

Graduate School of Engineering, Tohoku University, Sendai,
Miyagi 980‑8579, Japandesired one, such as, for example, a signal generated by thesidelobe, cannot be completely removed. Due to such unnecessary signals, artifacts are generated on a B-mode image,and the image quality deteriorates. In ultrasound image diagnosis, the examiner often diagnoses the lesion based on aslight change in brightness. Consequently, the deteriorationin the image quality affects the diagnosis and might lead to
misdiagnosis and stress on the examiner.
Many adaptive beamformers have been proposed to
improve image quality by reducing unnecessary signals
contained in received signals. For example, minimum variance beamforming (MVB), which adaptively calculates theweight values of the signals in receiving elements, has been
actively studied [1-3]. However, the computational complexity of MVB is large, and it is difficult to implement ingeneral ultrasound systems. Coherence-based beamforming(CBB) using the coherence factor (CF) [4, 5] has been proposed as an adaptive beamforming technique with relatively

13Vol.:(0123456789)

180

small computational complexity. CF is a factor representing
the coherence of the received signals after the delay compensation of individual elements. In CBB, by weighting the
CF value to the received signal after DAS, the brightness
value of the pixel where the unnecessary signal is dominant
is reduced. Moreover, various factors other than CF havebeen proposed in the literature [6-9].
Speckle noises [10, 11] are generated in B-mode images
as scatterers are uniformly distributed in most parts of thehuman body. Speckle noise is caused by the interferenceof sound waves from many scatterers and is a variation ofa brightness value that is not directly related to the structure in the human body. In such a region, the phases of
the received signals are not completely aligned with eachother, and minor phase fluctuations occur. If CF is applied
to such a region, it can result in a problem where the average
brightness in the region of the diffuse scattering medium islowered, and speckle noise is emphasized [6, 12]. Although
there have been studies seeking to clinically obtain useful
information from the characteristics of the speckle noise[13] on B-mode images, enhancement of the speckle noise
generally causes degradation of the contrast-to-noise ratio(CNR), which hinders making a correct diagnosis. Therefore, speckle reduction methods, such as spatial compound[14], frequency compound [15], or image processing [16,17], have been studied and developed. Therefore, when using
CBB for imaging of a diffuse scattering medium, it is important not to emphasize speckle noises. The generalized coherence factor (GCF) [6], which is one of the factors used for
CBB, can reduce unnecessary signals without emphasizing
speckle noise and has CNR higher than that of other factors[12]. As another approach that focuses on the imaging ability
in diffuse scattering media, a phase coherence factor (PCF),
in which phase variation is suppressed by applying DAS ineach divided sub-aperture, has been proposed [18].
Recently, with the improvement in graphical processing units (GPU), high-performance beamforming including
CBB can be realized in real time [19, 20]. On the other hand,miniaturization and portability of ultrasonic diagnostic apparatuses have been progressing as well. In such small models
and low-end models, it is important to achieve high image
quality with limited computing power. Therefore, there is aneed for a method with low computational complexity when
applying CBB, so that CBB can be expanded to a wider
range of models using this technique.
In the present study, we aimed to reduce the computational complexity of the GCF [6], which has CNR higher
than that of other factors. The computational complexityof the GCF is higher than that of the DAS beamformer
as it requires generation of analytic signals for individual
elements and calculation of the discrete Fourier transform(DFT) of the element direction for each pixel. The generation of analytic signals for individual elements is required in

13

Journal of Medical Ultrasonics (2020) 47:179–192

many adaptive beamformers including the PCF. To omit thecalculation of the DFT of the element direction, a method
has been proposed in which the GCF value is estimatedfrom autocorrelation [21]. In the present paper, we propose
a method to calculate the GCF values from real signals,
omitting generation of analytic signals. It has been reported
that the correct GCF value cannot be calculated using real
signals due to the presence of the carrier frequency [7].
However, the specific problems related to the presence ofthe carrier frequency have not been reported in detail. In
the present paper, this problem is clarified by theoreticallycomparing cases where the GCF values are calculated fromreal signals and where the GCF values are calculated fromanalytic signals following the conventional method. Then,
we propose a method to calculate the value equivalent tothe conventional GCF using real signals. To conclude on
the applicability of the proposed method, we evaluated the
GCF values obtained using the channel data acquired froma phantom with the ultrasonic diagnostic equipment and the
B-mode images calculated by the conventional and the proposed methods.

Conventional methods
Coherence‑based beamforming [4, 6–9]
Figure 1a shows the system block diagram for CBB using
CF. The focused beam is generated from received signals in
the number of channels connected to the ultrasound diagnostic apparatus. The received signals in channels are converted
to digital signals discretized in the time direction by analogto-digital converters (ADC) and delayed to be in phase forwaves coming from the focus point. The focused signal
xin (n) is obtained by the sum of delayed signals s(m, n) ,where m(m = 0, 1, … , M − 1) represents the channel number
corresponding to the element, M represents the number ofreceived channels, and n is the sampling number in the time
direction. CF is calculated from s(m, n) and then used as aweighting factor for xin (n), as described below. As shownin Fig. 1, the calculated CF value can be adjusted using alook-up table (LUT).

Coherence factor (CF) [4]
The CF value CF(n) [4, 5] is obtained as follows:

�∑M−1
�2
� m=0 s(m, n)�
� .
CF(n) = � ∑M−1
M ⋅ m=0 �s(m, n)�2

(1)

The denominator is the total power value (incoherent sum) of the received signals in the channel direction,

Journal of Medical Ultrasonics (2020) 47:179–192

181

Fig. 1  a System block diagram for coherence-based imaging. b GCF estimator. c The proposed method ­(GCFreal)

and the numerator represents the power value of the DC
component (coherent sum). As the difference between
phases of the received signals among channels is small,the CF value increases because the numerator, which corresponds to the power of the DC component, becomeslarge. If the received signals are completely in phase andtheir amplitudes are equal for all channels, the CF value
gets the maximum value of 1. However, when the difference between phases of the received signals among channels is large, the CF value decreases, because the powerof the DC component becomes small. CF is large whena reflected or scattered wave is received only from thefocus point. Therefore, by weighting CF to the signal xin (n)after applying DAS, xin (n) including unnecessary signals
is suppressed. However, as CF becomes too small aroundthe strong scatterer, the signals from the surrounding diffuse scattering medium are excessively reduced, and a dark
region artifact is generated [22]. Therefore, to apply CF
to a B-mode image, it is necessary to adjust the reduction
effect. In the present paper, similarly to the sign coherence factor (SCF) [8], the index p in Eq. (2) was used asa weight:

13

182

Journal of Medical Ultrasonics (2020) 47:179–192

xout (n) = [CF(n)]p xin (n).

(2)

The adjustment by the power p is implemented using
LUT, as outlined in Fig. 1a.

Generalized coherence factor (GCF) [6]k=−K

GCF(n; K0 ) = ∑K−1 0�c(k, n)�2

k=−K �c(k, n)�

2

,

(3)

where k is the frequency index in the channel direction, andk = −K, −K + 1, … , 0, … , K − 1, (K = M∕2) for an evennumber M . Here, the CF value is calculated from the GCF
value as CF(n) = GCF(n;0). While CF is calculated fromthe power value of the DC component c(0, n), GCF is calculated from the[ power valueof the frequency components]
in the range of −K0 , K0 . Therefore, the decrease in GCF
caused by the phase fluctuation due to the diffuse scattering
media is suppressed. The effect on GCF can be adjusted by
K0 . The GCF value approaches the maximum value of 1
if the frequency component around DC determined by K0
is dominant over the entire power spectrum in the channel
direction. As GCF aims to reduce unnecessary signals, it isnecessary to set an appropriate value K0 to avoid inclusion of
frequency components of signals outside of the focal point.
As K0 becomes larger, unnecessary signal components far
from DC are included in the numerator of Eq. (3). As aresult, the GCF value increases, and unnecessary signals
cannot be suppressed.
For example, the GCF value can be obtained from
Eq. (3) after calculating c(k, n) for analytic signals in thechannel direction using the fast Fourier transform (FFT).
However, as the value of K0 is generally a small number,

13

M−1

K−1

In a signal received from the diffuse scattering medium,the incoherent components are included due to the interference of ultrasonic waves from many scatterers. Therefore, the CF value decreases, and weighting of CF reduces
the brightness value in the region of the diffuse scattering medium accordingly. To deal with this problem, GCF,
which is an improvement of CF, has been proposed. Figure 1b shows the configuration, where the CF estimator
unit in Fig. 1a is changed to that of GCF. GCF is calculated from the following equation using the power valueof the Fourier coefficient c(k, n) obtained using DFT of theanalytic signal I(m, n) + jQ(m, n) in the channel direction:∑ K0

the computational complexity can be reduced by DFT only
for necessary components. Also, there is no need to calculate DFT for the denominator, because the value can becalculated from the squared sum of the received analytic
signals in the channel direction based on Parseval’s theorem as follows:

∑

k=−K

|c(k, n)|2 = M ⋅

∑

m=0

(4)

|I(m, n) + jQ(m, n)|2 ,where I(m, n) and Q(m, n) are the real and imaginary partsof the analytic signal, respectively. Substituting Eq. (4) into
Eq. (3) provides the follows:

GCF(n; K0 ) =

∑ K0

k=−K0

M⋅

∑M−1

�c(k, n)�2

m=0 �I(m, n) + jQ(m, n)�

2

.

(5)

As in the case of CF, the signal xin (n) after applying DAS
is weighted using the GCF value.

Proposed method
Generalized coherence factor estimated from real
signal ­(GCFreal)
According to the conventional formula, GCF is calculatedusing analytic signals generated from the received signals
in individual channels. However, as multipliers such as
mixers and low pass filters are required at a high sampling frequency, the computational complexity becomessignificantly large when generation of the analytic signal
is applied to the signal in each element. Therefore, we
propose a method to calculate GCF based on real signals
to omit the generation of analytic signals.
First, let us consider the conventional GCF calculation.
Although we derive analytic signals by Hilbert transform
here, the same approach can also be applied to the case of
baseband demodulation. For the delayed signal in channel
m , the analytic signal calculated by the Hilbert transformin the n th direction is expressed as follows:[]I(m, n) + jQ(m, n) = A(m, n)exp j2𝜋fnT ,(6)
where A(m, n) is the complex amplitude including the phasechange in the channel direction, and its absolute value represents the signal envelope. Here, f represents the carrier

Journal of Medical Ultrasonics (2020) 47:179–192

183

frequency of the received signal, and T represents the reception sampling period. When DFT is applied to the analytic
signal in the channel direction, the Fourier coefficient c(k, n)
is expressed as follows:

[

=

[
]
1 ∑
2mk𝜋
A(m, n) exp j2𝜋fnT exp −j
M m=0
M

]

(7)

∑ K0

GCFI (n; K0 ) = ∑K−1

�d(k, n)�22
k=−K �d(k, n)�

Consequently, the power value of c(k, n) is calculated by
2
|M−1
]|
[
2mkπ ||
1 || ∑,A(m, n)exp −j|c(k, n)| = 2 |M ||M | m=0||

Substituting d(k, n) instead of c(k, n) in Eq. (3) provides
the following:
k=−K0

[
]
∑
[
] M−1
2mk𝜋
1
A(m, n) exp −j.
exp j2𝜋fnT
=
M
M
m=0

2

(11)

}

+c∗ (k, n)c∗ (−k, n) + |c∗ (−k, n)|2 .

M−1
[
]
2mk𝜋
1 ∑
c(k, n) ={I(m, n) + jQ(m, n)} exp −j
M m=0
M
M−1

1
|c(k, n) + c∗ (−k, n)|2
4{
1=
|c(k, n)|2 + c(k, n)c(−k, n)4|d(k, n)|2 =

(8)

=

∑ K0

k=−K0

M⋅

�d(k, n)�2

,
2
n)�
�s(m,
m=0
(12)

∑M−1

where notation GCFI (n;K0 ) is used to distinguish from theconventional GCF calculated from the analytic signals.
The right hand side of Eq. (12) shows the case when the
denominator is calculated from the real signal s(m, n) similarly as in Eq. (5). Substituting Eq. (11) into Eq. (12) gives
the following:���c(k, n)�2 + c(k, n)c(−k, n) + c∗ (k, n)c∗ (−k, n) + �c∗ (−k, n)�2
GCFI (n; K0 ) = ∑K−1 ��.
2∗ (k, n)c∗ (−k, n) + �c∗ (−k, n)�2
n)�
+
c(k,
n)c(−k,
n)
+
c
�c(k,
k=−K
∑ K0

(13)

k=−K0

where it has no terms depending on f . Then, GCF is calculated by substituting Eq. (8) into Eq. (3) or Eq. (5). The termthat changes in the nth direction in Eq. (8) is only A(m, n).
As A(m, n) is the complex amplitude of the envelope, thechange is gradual compared to f . In the case of baseband
demodulation, there is no carrier frequency component in
Eq. (6), and as a result, the power value of c(k, n) is the sameas in Eq. (8).
Next, let us consider the case when GCF is calculated
from the Fourier coefficients obtained by DFT for the real
signal I(m, n) . As the real signal is represented byI(m, n) =]1[{I(m, n) + jQ(m, n)} + {I(m, n) + jQ(m, n)}∗ ,2(9)the Fourier coefficient d(k, n) of I(m, n) is represented byd(k, n) =1{c(k, n) + c∗ (−k, n)}.
2

(10)

Therefore, the power value of d(k, n) is expressed asfollows:

As
K0
∑

k=−K0

2

|c(k, n)| =

K0
∑

k=−K0

|c∗ (−k, n)|2

(14)

and
K0
∑

{c(k, n)c(−k, n) + c∗ (k, n)c∗ (−k, n)}

k=−K0

=2

K0
∑

Re{c(k, n)c(−k, n)},
(15)

k=−K0

Equation (13) can be rewritten as follows:
∑K0

k=−K0

GCFI (n; K0 ) = ∑K−1

�c(k, n)�2 +
2

∑ K0

k=−K0

Re[c(k, n)c(−k, n)]∑K−1n)� + k=−K Re[c(k, n)c(−k, n)]p(−K0 , K0 , n) + x(−K0 , K0 , n)=,p(−K, K − 1, n) + x(−K, K − 1, n)(16)
k=−K �c(k,wherep(k1 , k2 , n) =

k2
∑

k=k1

|c(k, n)|2 ,

(17)

13

184

x(k1 , k2 , n) =Journal of Medical Ultrasonics (2020) 47:179–192
k2∑

Re[c(k, n)c(−k, n)].

(18)

k=k1

In Eq. (16), the cross term represented by Eq. (18) is
added to the denominator and numerator in Eq. (3). As
the term in Eq. (17) is equivalent to the denominator andnumerator in Eq. (3), it is expressed as a GCF term in
the present paper. The term c(k, n)c(−k, n) in Eq. (18) is
calculated byc(k, n)c(−k, n) =

∑
[
] M−1
1
A(m, n)exp
j4𝜋fnT
M2
m=0
[
] M−1
]
[
2mkπ
2mkπ ∑A(m, n)exp j,
exp −j
M m=0
M
(19)

where there is a component of frequency 2f in the time ( nth) direction. Due to this component, GCFI (n) has oscillationof 2f , which does not appear in the GCF value calculated
from the analytic signals. To calculate ­GCFreal equivalent to
GCF by removing this unnecessary cross term, we used thefinite impulse response low pass filter (FIR-LPF) with thecoefficient fLPF (l), (l = −L, … , L) for the nth direction, as in
the following equation:
�∑
��
∑L �
K0
2
d(k,
n
−
l)�
f
(l)
⋅
�
l=−L LPF
k=−K0
GCFreal (n; K0 ) = ∑
�
��
�∑
L
K
f (l) ⋅d(k, n − l)�2
l=−L LPF
k=−K �
��
�∑
∑L �
K0
2
f
d(k,
n
−
l)�
(l)
⋅
�
l=−L LPF
k=−K0
=∑
�
� ∑
�� .
M−1
L
2
f
(l)
⋅
M
s(m,
n
−
l)�
�
LPF
m=0
l=−L
(20)

Figure 1c shows a block diagram for estimating ­GCFrealin the CF estimator unit presented in Fig. 1a. Compared toFig. 1b, the Hilbert transformer at each channel is eliminated, and LPF is added after the denominator and numerator calculations. As LPF is the processing after which thechannel signals are delayed and summed, the number of calculations can be reduced to approximately 2/(channel number) when the computational load of the Hilbert transformfor one channel is comparable to LPF.

Results and discussion
Experimental setup
To verify the proposed method, the RF data were acquired
from 96 channels using the ultrasound diagnostic system
ProSound α10 (Hitachi, Tokyo, Japan). The measurement
object was a multipurpose phantom 403 GS-LE (Gammex,

13

Fig. 2  a Original B-mode image (delay and sum). The blue line A
shows the ROI containing the wire target. The red line B shows the
ROI shifted about 2 mm (1.1°) in the azimuth direction with respect
to the wire target. The yellow line C shows the ROI of the diffuse
scattering medium. b Schematic diagram showing the transmit andreceive aperture used in each ROI

WI, USA). In this experiment, we set M = 96 and K = 48.
The GCF and G­ CFreal values were calculated offline using
the acquired RF data and compared to each other, andB-mode images weighted according to the obtained values
were also constructed. As GCF has better CNR than otherfactors in the case of uniform scattering media [12], a 60R
convex probe (UST-9130, 3.5 MHz, element pitch 0.38°)
mainly used for observation of the abdomen, where there
are many regions of diffuse scattering media, was employedin the experiment. The transmitting center frequency was
5 MHz, and the depth of the transmitting focus position was
105 mm from the probe surface. The sampling frequency of
the RF data was 20 MHz. Figure 2a shows a B-mode image
obtained by the conventional DAS. The image consists of
312 received scanning lines in the range of a viewing angleof 60°. We analyzed the signals of the scanning lines in thethree regions of interest (ROI) shown in the B-mode image.
Region A contained a wire target. Region B was shifted bysix scanning lines (about 2 mm (1.1°) from the wire target)in the azimuth direction with respect to region A. Region C
contained a diffuse scattering medium. All three ROIs wereset at the same depth (42–52 mm). As shown in Fig. 2b, theJournal of Medical Ultrasonics (2020) 47:179–192

transmitting and receiving aperture used for measuring each
ROI were set such that the lateral position of each ROI was
located in the center of the aperture.

Channel data and power spectrum in each ROI
In this section, we discuss the channel data s(m, n) and thepower spectrum |c(k, n)|2 used for the GCF estimation ineach ROI. Figure 3a–c show the channel data s(m, n) in ROIs
A, B, and C after delayed processing. Figure 3d–f show theabsolute values of the analytic signal I(m, n) + jQ(m, n)obtained by the Hilbert transform of s(m, n) for each m in
the n th direction, and their phases are shown in Fig. 3g–i.
Figure 3j–l show the power spectra obtained by DFT ofthe analytic signals in the m th direction for each n, and the
DC component ( k = 0) is shown by the white dashed lines.
The display range ([MIN, MAX]) differs depending on the
ROI, except for Fig. 3g–i. The GCF values were calculated
from |c(k, n)|2 in Fig. 3j–l using Eq. (3). In Fig. 3a, d, g for
ROI A, after the signals from the wire target were received
and delayed, they were aligned in the m th direction as the
position of the receiving focus matched the position of thewire target. Therefore, the power value from the wire target was concentrated on the DC component, and the otherfrequency components were minor, as shown in Fig. 3j. In
Fig. 3b, e, h for ROI B, the position of the wire target andthe receiving focus position do not coincide. Therefore, after
delayed processing, the reception time of the reflected waves
from the wire target was shifted for m , and the phases in them th direction changed. Received signals in ROI C shownin Fig. 3c, f, i tended to be in phase in the m th direction.
However, as the scattered waves from many scatterers interfere, the phases fluctuated, as can be seen in comparison
to Fig. 3g. Therefore, the power spectrum in Fig. 3l was
distributed around DC and was observed to have a wider
bandwidth than that presented in Fig. 3j.

Characteristic of cross termOnce the Fourier coefficient c(k, n) is obtained, the crossterm x(k1 , k2 , n) in Eq. (18), which has to be removed in
Eq. (16), can be calculated. To confirm the necessity of thecross term reduction by LPF and its effect, we comparedmagnitudes of the cross term and the GCF term. As mentioned above, when reducing unnecessary signals by GCF,
K0 should be set to an appropriate value. Here, to investigate
the cross term, two cases of K0 = 1 and K0 = 5 were examined. Figure 4Ia–c show the GCF term p(−1, 1, n) and thecross term x(−1, 1, n) values calculated for K0 = 1 in ROIs
A, B, and C. If GCF is calculated from real signals, the sumof GCF and cross terms are calculated as shown in Eq. (16).
Also, the power spectra by DFT in the nth direction for the185

signals are shown in Fig. 4IIa–c. In ROI A, the peak of thecross term was about − 10 dB with respect to the peak valueof the GCF term (i.e., power of the DC component), while it
was as small as about − 20 dB in ROIs B and C. To illustrate
this difference quantitatively, the ratio of the squared sumof the cross term to the GCF term was calculated from thefollowing equation in each ROI:

∑N2

�x(k , k , n)�2
1 2
�,r(k1 , k2 ) = ∑N2
�p(k1 , k2 , n)�2
�
n=N1 �
n=N1 �

(21)

where N1 and N2 represent the sample numbers at the top
and bottom in each ROI, respectively. The cross term of thenumerator oscillates at 2f , while there is less fluctuation in
the GCF term of the denominator, which is the power valueof the envelope. Therefore, even if the cross term of thenumerator and the GCF term of the denominator have thesame amplitude, the power value of the cross term is halfof that of the GCF term. Thereby, the maximum value ofr(k1 , k2 ) is 0.5. The value of r(−1, 1) in each ROI is shownin Table 1. As it was close to the maximum value of 0.5 in
ROI A, the amplitude value of the cross term was almost thesame as that of the GCF term. r(−1, 1) in ROIs B and C wereclose to each other. Figure 5 shows the GCF term p(−5, 5, n)
and the cross term x(−5, 5, n), and those of power spectra for
K0 = 5 in ROIs A, B, and C. The ratio r(−5, 5) is also presented in Table 1. For ROI A, both amplitudes of the GCF
term and the cross term were equivalent in the case of K0 = 1.
However, in ROI B, the peak value of the cross term in the
power spectrum was about − 30 dB with respect to the GCF
term, as shown in Fig. 5IIb, and it can be seen that r(−5, 5)
is considerably smaller than those corresponding to the other
ROIs. In ROI C, r(−5, 5) is slightly smaller than r(−1, 1).
Considering the above results, we can conclude on the
difference between K0 = 1 and K0 = 5 in each ROI. As canbe seen in Fig. 3j, the power value of the signal in ROI
A, where the wire target was included, is concentrated at
DC ( k = 0) in the power spectrum. Therefore, changing thevalue of K0 from 1 to 5 does not lead to significant changes
in either the GCF term or the cross term. Moreover, assuming that the signal component at k = 0 is dominant, the GCF
term in Eq. (8) and the envelope of the cross term in Eq. (19)have the same amplitude despite the fact that the cross term
oscillates at 2f . In ROI B (Fig. 3k), as the signal component from the wire target hardly exists in the narrow rangeof K0 = 1, the signal from the diffuse scattering medium on
the focal point is dominant. Therefore, r(−1, 1) in ROI B
was equivalent to that in ROI C. In the case of K0 = 5, thefrequency component, which was biased to the positive sideby the reflected wave from the wire target, was included in
the range of [−K0 , K0 ], as shown in Fig. 3k. Therefore, thecross term calculated from multiplication of the positive and

13

186

13

Journal of Medical Ultrasonics (2020) 47:179–192

Journal of Medical Ultrasonics (2020) 47:179–192◂Fig. 3  Delayed RF signals and their analytic signals in each ROI. a–cshow RF signals in ROIs A, B, and C. d–f show amplitude values ofanalytic signals generated from (a)–(c). g–i show phases of analytic
signals generated from (a)–(c). j–l show power spectra obtained by
DFT in the m-th direction on analytic signals generated from (a)–(c)negative frequency components was smaller compared tothat of the GCF term. In ROI C, r(−5, 5) was lower thanr(−1, 1). This is because the observation range of the power
spectrum at K0 = 1 is narrower than the one at K0 = 5, and theinfluence of the DC component (k = 0), where the GCF term
and the envelope of the cross term have the same amplitude,
became dominant.
As mentioned above, the GCF value becomes smaller inthe case of setting K0 in such way that the component of theunnecessary signal is not included in the range of [−K0 , K0 ].
Therefore, in the case of ROI B, K0 should be set to avoid
inclusion of the signal from the wire target. For this reason,

187

K0 = 1 is preferable in this experimental data. Furthermore,
there are few situations where only a strong reflector likea wire target exists in a human body. In general, due to thepresence of similar scatterers around it, the spectrum has abandwidth similar to that shown in ROI C (Fig. 3l). Therefore,
within a narrow range such as K0 = 1, the frequency component is less biased to positive or negative sides. Consequently,
there are few regions where the cross term becomes smaller,
as shown in Fig. 5Ib and IIb. In most of the data obtained fromthe biological object, the cross terms cannot be ignored, andtheir reduction by LPF is essential. If the cross term remains,the ­GCFreal value differs from that of the conventional GCF,
and there is a possibility that speckle noise may be emphasized by the presence of high-frequency components.

Fig. 4  GCF term p(−1, 1, n) andcross term x(−1, 1, n) values
in ROIs Ia A, Ib B, and Ic C.
Power spectra obtained by DFT
in the n-th direction with respect
to ROIs IIa A, IIb B, and IIc C

13

188

Journal of Medical Ultrasonics (2020) 47:179–192

Table 1  Ratio of the squared sum of the cross term to the GCF termr(k1 , k2 )r (− 1,1)r (− 5, 5)

ROI A

ROI B

ROI C

0.482
0.463

0.289
0.002

0.285
0.127

Cross term reduction by FIR‑LPF
Here, the effect of reducing the cross term by applying
FIR-LPF in Eq. (20) is discussed, considering the case ofK0 = 1, which is effective for reducing unnecessary signals. As shown in Eq. (8), the GCF term has only lowfrequency components in the n th direction, because it has
no carrier frequency component. As the center frequencyof the acquired data was approximately 3 MHz, that of thecross term is considered to be distributed around 6 MHz

Fig. 5  GCF term p(−5, 5, n) andcross term x(−5, 5, n) values
in ROIs Ia A, Ib B, and Ic C.
Power spectra obtained by DFT
in the n-th direction with respect
to ROIs IIa A, IIb B, and IIc C

13

according to Eq. (19). For all of the power spectra shownin Figs. 4 and 5, it was confirmed that the GCF term wasdistributed around the DC component, and the cross term
was distributed around 6 MHz. Therefore, the cross term was
reduced using an 8th-order FIR-LPF with a cutoff frequencyof 3 MHz. Figure 6Ia–c show the results of applying LPF to
the numerator values p(−1, 1, n) + x(−1, 1, n) in Eq. (16) in
ROIs A, B, and C. The value of p(−1,[ 1, n) obtained by the]conventional GCF and that of LPF p(−1, 1, n) + x(−1, 1, n)obtained by the proposed method were in agreement for all
ROIs. Although the reduction effect of the cross term wasdifferent depending on the frequency characteristics of LPF,
the differences were minor,[ and the correlation coefficients]of p(−1, 1, n) and LPF p(−1, 1, n) + x(−1, 1, n) were sufficiently high (0.99 or more in each ROI) in the 8th-order
FIR-LPF. Figure 6IIa–c show the power spectra obtained by
DFT in the nth direction in Fig. 6Ia–c. In addition, the power
spectra of the signal p(−1, 1, n) + x(−1, 1, n) before LPF andJournal of Medical Ultrasonics (2020) 47:179–192

189

Fig. 6  Ia–Ic Signals and IIa–
IIc power spectra in whichcross terms have been reduced
by low-pass filtering for signalsof p(−1, 1, n) + x(−1, 1, n). Ia
and IIa show results for ROI A,
Ib and IIb show results for ROI
B, and Ic and IIc show resultsfor ROI C

the frequency characteristic of the 8th-order FIR-LPF are
shown in Fig. 6IIa–c. It can be confirmed that the frequency
component of the cross terms was reduced, and the GCF
term became dominant after applying LPF.

Comparison of ­GCFreal and GCF
Figure 7a, b, c show the values of the proposed G­ CFreal
obtained from Eq. (20) and the conventional GCF obtained
from Eq. (5) at K0 = 1 adopted as outlined in Fig. 6. ­GCFreal
and GCF were close to each other, and the correlation coefficient of them showed 0.99 or more in each ROI. The vicinity of nT = 60𝜇s . corresponds to the time, where the signals from the wire target were received. In ROI B, the value
became almost 0 around nT = 60𝜇s, and it can be concluded

that the output signal including the signal received fromoutside of the focus point was suppressed by weighting of­GCFreal. However, as the G­ CFreal values were too small, theuse of these values as the weight values resulted in the dark
region artifact [22]. Therefore, adjustment of the GCF values
by the power p in Eq. (2) is required to apply weighting.
Figure 7d shows a B-mode image obtained by DAS, andFig. 7e, f show B-mode images weighted by G
­ CFreal and
GCF at p = 0.2. The dynamic range of each image is 80 dB.
Sidelobes were generated in the azimuth direction for wire
targets in Fig. 7d; however, they were reduced in Fig. 7e,f. As a representative example, the intensity profiles in theazimuth direction at the depth of the wire target surroundedby the green square in Fig. 7d are shown in Fig. 7g. Regarding ­GCFreal and GCF, it can be confirmed that the sidelobe
generated around the wire target was reduced, while the

13

190

Journal of Medical Ultrasonics (2020) 47:179–192

Fig. 7  GCFreal and GCF values
calculated in ROIs a A, b
B, and c C. B-mode images
generated by d DAS and thoseweighted by e ­GCFreal and f
GCF. The dynamic range ofeach image is 80 dB. g Lateralamplitude profiles in the depthof the wire target are shown by
the green square in (d)brightness of the wire target was at the same level as in the
DAS. The CNR values were calculated from the following
equation [23] using the brightness values after log compression in the yellow ROI as region 1, and the red ROI as region
2, as shown in Fig. 7d:|𝜇1 − 𝜇2 |||,CNR = √
2
𝜎1 + 𝜎22(22)

where 𝜇1 and 𝜎1 are the average brightness and the standard deviation of region 1, respectively, and 𝜇2 and 𝜎2
are those of region 2. The CNR values in the respective
images are 2.29, 2.65, and 2.64 for DAS, G­ CFreal, and GCF,

13

respectively. Therefore, the CNR improvement effect of the
proposed method was equivalent to that of the conventional
GCF. It was confirmed that the computational complexity
for applying FIR-LPF in the proposed ­GCFreal calculation
was approximately 1/50 compared with that for generating
the analytic signals in all channels in the conventional GCF
calculation so that the computational complexity can be
drastically reduced by the proposed method. The quality of
the B-mode image weighted by the proposed ­GCFreal was
equivalent to that weighted by the conventional GCF. Here,the computational complexity indicates the number of calculations for generating the analytic signals in all channels,
which is only used in the conventional method, or applying
FIR-LPF, which is only used in the proposed method.

Journal of Medical Ultrasonics (2020) 47:179–192

The conventional GCF has the advantage of reducing thesampling frequency when using baseband demodulators.
Reducing the sampling frequency is equivalent to reducing
the number of subsequent GCF value calculations. The conventional GCF is composed of analytic signal generation and
GCF value calculation, whereas the proposed method only
needs ­GCFreal value calculation. The analytic signal generation needs to filter the signals in the time (nth) direction for all
channels. Therefore, the number of multiplications required
per pixel is Lf M , where Lf is the number of filter coefficientsand M is the number of received channels. On the other hand,the processing for the channel signals in the G­ CFreal value
calculation is DFT and power value calculation in the channel (mth) direction. (Therefore,) the number(of multiplications)
required per pixel is 2K0 + 1 M + M = 2 K0 + 1 M , where
K0 = 1 was used
( in the) present paper. In general, since Lf islarger than 2 K0 + 1 , the number of calculations is smallerin the ­GCFreal value calculation than in the analytic signal
generation. Therefore, the number of calculations in the proposed method is small even when compared with the conventional GCF, in which the sampling frequency is reduced bybaseband demodulation.

Conclusions
Adaptive beamforming based on GCF can reduce unnecessary signals without reducing CNR. However, in thismethod, it is necessary to generate analytic signals forthe received signals in individual channels, and, therefore,
the computational complexity increases compared to theconventional DAS beamforming. In the present paper,
we proposed a method to calculate the values equivalent
to that of the conventional GCF method using only real
signals. The proposed method can omit the generation ofanalytic signals without deteriorating the accuracy of the
GCF value. Additional processing is required in terms ofapplying only FIR-LPF after DAS. Regarding this part,the computational complexity can be reduced to approximately 2/(receive channel number) as compared to the conventional GCF method. The proposed method can improve
the feasibility of small and low-end models of ultrasonicdiagnostic apparatuses. There are two factors that increase
the computational complexity of conventional GCF compared to DAS. One is the generation of analytic signal
in each channel, and the other is the part that calculates
the GCF value by DFT in the channel direction for eachpixel. In the present paper, we focused on the former andsucceeded to omit the generation of analytic signals. In
future research, we will consider the latter by developinga method to further reduce the computational complexityof the GCF calculation units.

191

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflicts ofinterest.
Ethical statements This article does not contain any studies withhuman or animal subjects performed by any of the authors.

References
1. Capon J. High-resolution frequency-wavenumber spectrum
analysis. Proc IEEE. 1969;57:1408–18.
2. Synnevag JF, Austeng A, Holm S. Adaptive beamforming
applied to medical ultrasound imaging. IEEE Trans Ultrason
Ferroelectr Freq Control. 2007;54:1606–13.
3. Synnevag JF, Austeng A, Holm S. Benefits of minimum-variance beamforming in medical ultrasound imaging. IEEE TransUltrason Ferroelectr Freq Control. 2009;56:1868–79.
4. Hollman KW, Rigby KW, O’Donnell M. Coherence factorof speckle from a multi-row probe. In: Proceedings of IEEE
ultrasonics symposium. IEEE: Los Alamitos, CA; 1999. pp.
1257–1260.
5. Kanai H, Sato M, Koiwa Y, et al. Transcutaneous measurement and spectrum analysis of hear wall vibrations. IEEE TransUltrason Ferroelectr Freq Control. 1996;43:791–810.
6. Li PC, Li ML. Adaptive imaging using the generalized coherence factor. IEEE Trans Ultrason Ferroelectr Freq Control.
2003;50:128–41.
7. Wang SL, Chang CH, Yang HC, et al. Performance evaluation of coherence-based adaptive imaging using clinicalbreast data. IEEE Trans Ultrason Ferroelectr Freq Control.
2007;54:1669–788.
8. Camacho J, Parrilla M, Fritsch C. Phase coherence imaging. IEEE
Trans Ultrason Ferroelectr Freq Control. 2009;56:958–74.
9. Wang Y, Zheng C, Peng H, et al. An adaptive beamforming
method for ultrasound imaging based on the mean to standard
deviation factor. Ultrasonics. 2018;90:32–41.
10. Burckhardt CB. Speckle in ultrasound B-mode scans. IEEE Trans
Son Ultrason. 1978;SU-25:1–6.
11. Wagner RF, Smith SW, Sandrik JM, et al. Statistics of speckle inultrasound B-scans. IEEE Trans Son Ultrason. 1983;30:156–63.
12. Hverven SM, Rindal OMH, Rodriguez-Molares A, et al. The
influence of speckle statistics on contrast metrics in ultrasoundimaging. In: Proceedings of IEEE Ultrasonics Symposium. IEEE:
Washington, DC; 2017.
13. Tuthill TA, Sperry RH, Parker KJ. Deviations from Rayleigh statistics in ultrasonic speckle. Ultrason Imaging. 1988;10:81–9.
14. Jespersen SK, Wilhjelm JE, Sillesen H. Multi-angle compound
imaging. Ultrason Imaging. 1998;20:81–102.
15. Magnin PA, Von Ramm OT, Thurstone FL. Frequency compounding for speckle contrast reduction in phased array images. Ultrason
Imaging. 1982;4:267–81.
16. Bamber JC, Daft C. Adaptive speckle reduction filter for
log-compressed B-scan images. IEEE Trans Med Imaging.
1996;15:802–13.
17. Karaman M, Kutay MA, Bozdagi G. An adaptive speckle suppression filter for medical ultrasonic imaging. IEEE Trans Med
Imaging. 1995;14:283–92.
18. Hasegawa H, Kanai H. Effect of sub aperture beamforming onphase coherence factor imaging. IEEE Trans Ultrason Ferroelectr
Freq Control. 2014;61:1779–900.
19. Tanter M, Fink M. Ultrafast imaging in biomedical ultrasound.
IEEE Trans Ultrason Ferroelectr Freq Control. 2014;61:102–19.

13

192
20. Yiu BYS, Tsang IKH, Yu ACH. GPU-based beamformer: fast
realization of plane wave compounding and synthetic aperture imaging. IEEE Trans Ultrason Ferroelectr Freq Control.
2011;58:1698–705.
21. Shen CC, Xing YQ, Jeng G. Autocorrelation based generalized
coherence factor for low-complexity adaptive beamforming. Ultrasonics. 2016;72:177–83.
22. Rindal OMH, Rodriguez-Molares A, Austeng A. The dark region
artifact in adaptive ultrasound beamforming. In: Proceedings of13

Journal of Medical Ultrasonics (2020) 47:179–192
IEEE international Ultrasonics Symposium. IEEE: Washington,
DC; 2017.
23. Patterson MS, Foster FS. The improvement and quantitative
assessment of B-mode images produced by an annular array/conehybrid. Ultrason Imaging. 1983;5:195–213.
Publisher’s Note Springer Nature remains neutral with regard tojurisdictional claims in published maps and institutional affiliations.