JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
35
Computation of the Normalized Detection Threshold for the FFT Filter Bank-Based Summation CFAR Detector Sichun Wang Defence Research and Development Canada-Ottawa Email:
[email protected] Franc ¸ ois Patenaude Communications Research Centre Canada Email:
[email protected] Robert Inkol Defence Research and Development Canada-Ottawa Email:
[email protected] AbstractThe FFT lter bank-based summation CFAR
procedures such as the Newton-Ralphson algorithm. As
detector is widely used for the detection of narrowband
practical implementations of these numerical procedures
signals embedded in wideband noise. The simulation and implementation of this detector involves some problems concerning the reliable computation of the normalized detection
require good lower and upper bounds of
T
for reliable
initialization, several theoretical lower and upper bounds
T
for
presents a comprehensive theoretical treatment of major
However, a self-contained comprehensive treatment of
aspects of the numerical computation of the normalized
the FFT lter bank-based summation CFAR detector has
detection threshold for an AWGN channel model. Equations are derived for the probability of false alarm,
Pf a ,
for
both non-overlapped and overlapped input data and then
have also been obtained [8].
not appeared in the literature and details of relevant mathematical derivations have not been published. This paper presents a comprehensive theoretical treat-
used to compute theoretical upper and lower bounds for
T.
for a given
Pf a
threshold for a given probability of false alarm. This paper
A very useful transformation is
ment of major aspects of the numerical computation of
introduced that guarantees the global quadratic convergence
the normalized detection threshold for an AWGN channel
the detection threshold
of the Newton-Ralphson algorithm in the computation of
T
for overlapped data with an overlap ratio not exceeding
50%.
It is shown that if the product of the number of FFT
model. Formulas for computing
Pf a
as a function of
T
are derived from rst principles. Several theoretical lower
T
bins assigned to a channel for signal power estimation and
and upper bounds for
the number of input data blocks is relatively small, e.g.,
transformation is introduced that guarantees the global
less than
60,
the theoretical normalized detection threshold
can be accurately computed without numerical problems. To handle other cases, good approximations are derived. Index TermsDigital FFT lter bank, Detection and estimation, Constant false alarm rate (CFAR) detection.
are derived and a very useful
quadratic convergence of the Newton-Ralphson algorithm in the computation of the normalized detection threshold for overlapped data. Finally, accurate approximations for
T
are presented. The results presented in this paper are
very useful for the simulation and implementation of the FFT lter bank-based summation CFAR detector.
I. I NTRODUCTION The FFT lter bank-based summation CFAR detector
This paper is organized as follows. Section II introduces the FFT lter bank-based summation CFAR detector.
is an efcient technique for the detection of narrowband
Section III derives the formulas for computing
signals embedded in wideband noise and has important
a given
T
Pf a
for
for overlapped and non-overlapped signal data.
T
applications in civilian spectrum monitoring, electronic
Section IV proves several lower and upper bounds of
warfare, radio astronomy and instrumentation [1]-[7].
for overlapped and non-overlapped input data. Section V
This detector operates by adding estimates of spectral
introduces a very useful transformation that guarantees
power of FFT bins in channels, each of which corresponds
the global quadratic convergence of the Newton-Ralphson
to a group of one or more contiguous FFT bins, and
algorithm while Section VI derives two approximations
comparing the sum of the channel power estimates over
for
multiple input data blocks against a detection threshold,
Throughout this paper, typical numerical examples are
T.
given to illustrate the results.
In its general form, this detector has been the subject
T.
Finally Section VII presents concluding remarks.
of many recent studies [7]-[14]. In particular, closed-form algebraic formulas for computing the probability of false alarm,
Pf a , as a function of T
have been derived [5], [7],
[13]. These results can be used to compute the normalized detection threshold (to be dened later) using numerical
© 2007 ACADEMY PUBLISHER
II. T HE FFT F ILTER BANK -B ASED S UMMATION CFAR D ETECTOR Consider a uniformly sampled band-limited signal embedded in additive white Gaussian noise. Assume that
M
36
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
channels are uniformly distributed across the frequency
detector or simply the summation detector if there is no
range contained within the Nyquist bandwidth and that
risk of confusion. A block diagram is given in Figure 1
K N
for the special case where
FFT bins are assigned to each channel with the FFT bins (N
≤ K)
N
and
K
are equal.
centered within each channel
used to estimate the power contained within the channel.
MK
Consequently, an FFT of length the power levels for the
M
Input data blocks
Detection threshold
is used to compute
channels. For notational
convenience and without loss of generality, assume
2
L stage delay
| |
K−N
w = [w0 , · · · , wM K−1 ]t be a symmetric window of length M K , where the superscript t denotes matrix transposition. Let L ≥ 1 be a positive integer and consider L consecutive overlapping sample vectors, Rl , constructed as follows:
is an even integer. Let
2
| |
S +
S
. . .
Channel 1 detection decision
2
Rl = [rl(1−γ)M K+M K−1 , · · · , rl(1−γ)M K ]t ,
| |
(1)
IFFT
sample of the input data stream and
2
| |
L stage delay
| |2
S
. . .
+
0 ≤ l ≤ L − 1, 0 ≤ γ ≤ 1/2, rn is the n-th γM K is an integer. In practice, the overlap ratio, γ , is often selected to be either 0 or 1/2. These two cases correspond to zero and 50% overlap, respectively. For each l, the two input vectors, Rl and Rl+1 , have γM K samples in common. The vectors, Rl , are windowed by w, resulting in the
where
Channel 2 detection decision
windowed sample vectors:
Xl = [w0 rl(1−γ)M K+M K−1 , · · · , wM K−1 rl(1−γ)M K ]t . The vectors,
Xl ,
Block diagram for the
summation CFAR detector.
F of dimensions M K ×
to yield the FFT lter bank output sample vectors:
III. T HE P ROBABILITY OF FALSE A LARM FOR THE S UMMATION CFAR D ETECTOR
Yl = FXl = [sl,0 , sl,1 , · · · , sl,M K−1 ]t ,
rn
Assume the input data stream
where
complex-valued
F= 1 ··· 1 ··· 1
(2)
1 ··· 2πjl e MK ··· e
of length
M
··· 1 ··· ··· 2πj(M K−1)l MK ··· e ··· ··· 2πj(M K−1)(M K−1) MK ··· e
2πj(M K−1) MK
From each vector,
N
L-block
are then transformed by the inverse
discrete Fourier transform matrix
MK
Fig. 1.
Yl ,
a vector
.
is formed by summing the power from the
zl,k =
N −1 X
|sl,Ik +m | ,
N
and the implementation details of the detector. For brevity, the FFT lter bank-based
L-block
sum-
mation CFAR detector shall be referred to as the
L-
block summation CFAR detector or the summation CFAR
with
(L−1 X
) zl,k ≥ T
.
(4)
l=0 The following three theorems provide the mathematical
T > 0, Pf a
(3)
Ik , Ik + 1, · · ·, Ik + (N − 1) are summed to form the power estimate for the k -th channel for the data block Rl where Ik = kK + K−N 2 . The detection criterion is dened as follows: if for a given detection threshold, T , PL−1 l=0 zl,k ≥ T , a signal is declared to exist in the k th channel. The choice of T is dependent on the desired probability of false alarm, Pf a , the noise spectral density,
© 2007 ACADEMY PUBLISHER
Pf a = Pr
given
FFT bins with indices
sequence
detector is then given by:
Theorem 1. Assume
0 ≤ l ≤ L − 1, 0 ≤ k ≤ M − 1,
is a zero-mean
noise
foundation for computing the threshold
2
m=0
i.e., the power estimates from the
Gaussian
E(rp rq∗ ) = σ 2 δpq , where σ 2 > 0 is the noise variance (noise oor), δpq = 1 if p = q and δpq = 0 if p 6= q . For a given threshold, T , the probability of false alarm, Pf a , for the k -th channel of the L-block summation CFAR
zl = [zl,0 , · · · , zl,M −1 ]t
FFT bins centered within each channel:
white
Pf a =
LN X
L ≥ 2
and
T
for a given
0 < γ ≤ 1/2.
Pf a .
For a
is given by:
Y
−1 λLN l
l=1
(λl − λm )
−
e
T σ 2 λl
,
(5)
1≤m≤LN,m6=l
λl , 1 ≤ l ≤ LN , are the LN distinct positive L × L block matrix H: A B 0 ··· 0 0 BH A B 0 ··· 0 H 0 B A B ··· 0 H= (6) ··· ··· ··· ··· ··· ··· . 0 ··· 0 BH A B H 0 ··· 0 ··· B A
where
eigenvalues of the
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
In (6),
H
LN × LN , 0 is the N × N N × N matrices dened
is of dimensions
zero matrix and
A
and
B
are
respectively by:
τ11 ··· A= τp1 ··· τN 1 τpq
τ12 ··· τp2 ··· τN 2
· · · τ1q ··· ··· · · · τpq ··· ··· · · · τN q
MX K−1
=
wl2
l=0
γ11 ··· γp1 ··· γN 1
B=
γ12 ··· γp2 ··· γN 2
37
··· ··· ··· ··· ···
· · · τ1N ··· ··· · · · τpN ··· ··· · · · τN N
In practical systems, the noise oor,
,
For a given
(7)
l=0
L ≥ 1, N > 1 ³
Amp
m=1 p=1
p−1 X
(8)
T /σ 2 ,
which
is computed using one of the equations (5), (11) or (13).
threshold in this paper. This can be loosely interpreted as single data block.
,
IV. L OWER AND U PPER B OUNDS FOR THE
(9)
D ETECTION T HRESHOLD In practice, the computation of
T
T /σ 2
involves the
numerical solution of one of equations (5), (11) or (13) algorithm. To avoid potential numerical problems and
(10)
and
γ = 0.
For a
ensure fast convergence, a good lower or upper bound for
T
is required to initialize the numerical procedures. This
−
T σ 2 µm
,
(11) A. Overlapped Input Data
Amp , 1 ≤ m ≤ N , 1 ≤ p ≤ L,
Y
The cases for overlapped and non-overlapped signal
the appendix.
e
s!
T.
data are treated separately, with full proofs provided in
are
L
AmL =
−1 µN m
(µm − µl )
,
1≤l≤N,l6=m
L ≥ 2, 0 < γ ≤ 1/2, and let λ1 > λ2 > · · · > λLN be the LN distinct positive eigenvalues of the LN × LN Hermitian matrix H dened by (6). For a given threshold T > 0, let the probability of false 2 alarm be denoted by Pf a . Let T1 (Pf a ) = −λ1 σ ln Pf a and denote the solution of the following equation in x by Tn (Pf a ): Theorem 4. Assume
Amp = AmL × X
Pf a =
(12)
1≤l≤N,l6=m
µl − µm
(kl )!(L − 1)!
γ = 0.
Pf a =
¡ T ¢s λ
s!
s=0
T
e− λ ,
MX K−1 l=0
© 2007 ACADEMY PUBLISHER
,
(15)
For a given
n ≥ 1, let the unique solution of z be denoted by Bn (Pf a ): ¶ µ 2 z n−1 z −z = Pf a . (17) + ··· + e 1+z+ (n − 1)! 2!
Then for
0 < Pf a ≤ 2/e,
(13)
and for
(14)
(16)
the following equation in
Bn (Pf a ) ≤
wl2 .
x σ 2 λm
For any positive integer
where
λ = σ2
−
T1 (Pf a ) ≤ T2 (Pf a ) ≤ · · · ≤ Tl (Pf a ) ≤ Tl+1 (Pf a ) ≤ · · · ≤ TLN (Pf a ) = T.
.
is given by:
L−1 X
(λm − λl )
e
Then
kq , 1 ≤ q ≤ N , are non-negative integers and µq , 1 ≤ q ≤ N , are the N distinct positive eigenvalues of the Hermitian matrix A dened by (7). and
λn−1 m
2 ≤ n ≤ LN.
Note that
L ≥ 1, N = 1
Y
1≤l≤n,l6=m
Γm (k1 , · · · , km−1 , km+1 , · · · , kN ), Γm (k1 , · · · , km−1 , km+1 , · · · , kN ) = ¸kl · Y µl (L + kl − 1)!
n X m=1
k1 +···+km−1 +km+1 +···+kN =L−p
T > 0, Pf a
is obtained by
by
using numerical procedures such as the Newton-Ralphson
Theorem 3. Let
T
the detection threshold for one FFT bin computed from a
´s
T
σ 2 µm
s=0
where the coefcients
σ2
T /(LN σ 2 ), which will be called the normalized detection
· · · γ1N ··· ··· · · · γpN ··· ··· · · · γN N
is given by:
N X L X
dened by:
the detection threshold
section presents theoretical lower and upper bounds for
T > 0, Pf a
Pf a =
Pf a ,
multiplying the estimated noise oor
found that it is more appropriate to focus on the quantity,
2πjl(p − q) . wl wl+(1−γ)M K exp MK
Theorem 2. Assume given
needs to be
In the analysis of the summation CFAR detector, it is
γpq = e−2πj(1−γ)(Ik +q−1) × γM K−1 X
σ2 ,
estimated from the output samples of the FFT lter bank.
2πjl(p − q) , exp MK
γ1q ··· γpq ··· γN q
The proofs of these theorems are relatively lengthy and will be given in the appendix.
n+
√
¸ · 2 2πn , ln ePf a 2
(18)
0 < Pf a < 1, T
≤
λ1 σ 2 BLN (Pf a ).
(19)
38
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
L ≥ 1, N ≥ 2, γ = 0 and let µ1 > N distinct positive eigenvalues of A dened by (7). We have
Theorem 5. Assume
µ 2 > · · · > µN
be the
the Hermitian matrix
2
2
µ1 σ BL (Pf a ) ≤ T ≤ µ1 σ BLN (Pf a ). T1
Moreover, if
T2
and
g 00 (p) < 0
Blackman Window (N=L=5, MK=1024) 1
x: 0.9 0.8
1
0.7
(21)
0.6 Pfa
−1 µN − x m Q e σ2 µm 1≤l≤N,l6=m (µm − µl )
m=1
for
(20)
= (Pf a ) L , N X
g(p) is strictly increasing and (0, 1). Specically, g 0 (p) > 0 and 0 < p < 1.
Theorem 6. The function concave on the interval
−1 µN − x m Q e σ2 µm 1≤l≤N,l6=m (µm − µl )
m=1
is rather
non-trivial and will be given in the appendix:
are respectively the solutions of
the equations (21) and (22) in
N X
g(p)
computed. The proof of the concavity of
B. Non-Overlapped Input Data
0.5
1
= 1 − (1 − Pf a ) L ,
(22)
0.4
then
0.3
with
LT1 ≤ T ≤ LT2 ,
(23)
´ ³ 1 T2 ≤ µ1 σ 2 BN 1 − (1 − Pf a ) L .
(24)
0.2 0.1 0
V.
0
0.002
0.004
0.006
0.008
0.01
p
N UMERICAL C OMPUTATION OF THE N ORMALIZED
D ETECTION T HRESHOLD FOR THE S UMMATION CFAR D ETECTOR
Fig. 2.
Using (5), (11) and (13) and in conjuction with the
e−T /(σ
2
Pf a , as a function of p = N = L = 5, M K = 1024, γ = 0.5, Blackman window.
The probability of false alarm,
λ1 ) ,
lower and upper bounds derived in the previous section, the normalized detection threshold
T /(LN σ 2 )
can be
computed numerically. It has been observed that in gen-
Blackman Window (N=L=5, MK=1024) 1
eral, equation (5) is numerically more unstable and hence more difcult to apply in practice than equation (11). To
0.9
avoid numerical problems with (5), the transformation
function
has proved to be very useful. Dene the
g(p)
0.8
by:
g(p) =
LN X
Y
−1 λLN l
l=1
(λl − λm )
0.7
p
λ1 λl
,
0.6
1≤m≤LN,m6=l
0 < p < 1.
(25)
0.5 0.4
Then the formula (5) can be rewritten as:
µ −T ¶ g e σ2 λ1 = Pf a .
Pfa
p = e
−T σ 2 λ1
0.3 (26)
0.2
T can be obtained by rst solving the equation g(p) = Pf a for p and then using the 2 relationship T = λ1 σ ln(1/p). It can be shown that g(p) is concave on the interval (0, 1), as is demonstrated by Figure 2. For comparison, Pf a is plotted as a function 2 of T /σ using (5) in Figure 3. The concavity of g(p) The detection threshold
is a very desirable property to have in the numerical computation of
T /(LN σ 2 );
it guarantees the fast global
0.1 0
0
20
40
60
80
100
T/σ2 Fig. 3.
The probability of false alarm,
N = L = 5, M K = 1024, γ = 0.5,
Pf a ,
as a function of
T /σ 2 ,
Blackman window.
quadratic convergence of the Newton-Ralphson algorithm provided that the functions
g(p) and g 0 (p) can be reliably
© 2007 ACADEMY PUBLISHER
For
LN
not too large, we have successfully computed
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
T /(LN σ 2 ),
the normalized detection threshold,
L-block
39
for the
3
summation CFAR detector. Typical results are
γ = 0.
plotted in Figures 4-5 for
Note that the windows
used to produce these results have been normalized in the
PM K−1
2.8 Hanning Blackman Blackman−Harris Hamming Nuttall Flat−Top Rectangular
Normalized Threshold T/(LNσ2)
2.6 2.4 2.2
2.6 Normalized Threshold T/(LNσ2)
wl2 = 1. The results for γ = 1/2 are l=0 very similar, but are not shown due to space limitations. sense that
Hanning Blackman Blackman−Harris Hamming Nuttall Flat−Top Rectangular
2.8
2.4 2.2 2 1.8 1.6 1.4
2 1.2
1.8 1 −10 10
1.6
−8
−6
10
−4
10
−2
10
10
0
10
PFA
1.4
Fig. 5.
The normalized detection threshold,
of the probability of false alarm
L = 10, γ = 0.
1.2 1 −10 10
−8
−6
10
−4
10
−2
10
0
10
10
with
T /(LN σ 2 ), as a function N = 8, M K = 1024,
The window functions are normalized.
This implies that the normalized detection threshold,
T /(LN σ 2 ),
PFA
Pf a ,
can be approximately computed via a chi-
square distribution. Using this approach, the following Fig. 4.
T /(LN σ 2 ), as a function Pf a , N = 5, M K = 1024, L = 18,
The normalized detection threshold,
of the probability of false alarm
γ = 0.
The window functions are normalized.
two theorems have been derived for normalized windows
PM K−1
wl2 = 1). l=0 Theorem 7. Assume
(i.e.,
Pf a , T VI. ACCURATE A PPROXIMATIONS TO THE D ETECTION T HRESHOLD For very large
LN ,
there are serious numerical dif-
culties in computing the normalized detection threshold
T /(LN σ 2 )
using (5) or (11). In such cases accurate
approximations to
T /(LN σ 2 )
are necessary. As will be
demonstrated in the appendix in the proofs of Theorems 1-3, the detection decision statistic,
PL−1 l=0
zl,k ,
of the
k-
th channel is distributed as a weighted sum of chi-square random variables. More specically,
L−1 X
zl,k ∼
l=0
PLN
PN
m=1
0 < γ ≤ 1/2,
σ 2 µm 2 2 χ2L (m),
γ = 0, (27)
where
λl , 1 ≤ l ≤ LN ,
are the
eigenvalues of the Hermitian matrix
µm , 1 ≤ m ≤ N ,
are the
of the Hermitian matrix
N A
LN distinct positive H dened by (6) and
distinct positive eigenvalues dened by (7). As pointed
out in [15], one can approximate a weighted sum of chi-square random variables by a random variable of the form
c + dχ2h
where
c
and
d
are constants and
chi-square random variable with
© 2007 ACADEMY PUBLISHER
h
and
γ = 0.
For a given
T ≈
³P
N m=1
σ 2 L N − PN
µ2m
´2
3 m=1 µm ³ ´3 PN 2 m=1 µm h = L ³P ´2 , N 3 m=1 µm
+
PN
µ3m Bh (Pf a ) Pm=1 , N 2 m=1 µm
(28)
µm , 1 ≤ m ≤ N , are the N distinct positive A dened by (7) and bxc denotes the greatest integer less than or equal to x. Theorem 8. Assume L ≥ 2 and 0 < γ ≤ 1/2. For a given Pf a , T is approximately computed by: ´2 ³P LN 2 PLN 3 λ l=1 l λl , + Bh (Pf a ) Pl=1 T ≈ σ 2 LN − PLN LN 2 3 l=1 λl l=1 λl ³ PLN 2 ´3 l=1 λl (29) h = ³P ´2 , LN 3 λ l=1 l
where
eigenvalues of
σ 2 λl 2 2 χ2 (l),
l=1
L ≥ 2
is approximately computed by:
χh
is a
degrees of freedom.
where
λl , 1 ≤ l ≤ LN , are the LN H dened by (6).
eigenvalues of
distinct positive
40
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
The excellent accuracy of approximations (28) and
A PPENDIX
(29) is apparent in Figure 6. This gure shows the true normalized detection thresholds for the cases of
γ = 1/2
γ = 0 and
and the corresponding approximate normalized
detection thresholds plotted as a function of
Pf a
for a
K−N 2 . For
typical example. In this gure's legend, NO and O denote non-overlapped (γ
= 0) and overlapped (γ = 1/2),
zl,k
N −1 X
=
|sl,Ik +m |2 =
respectively. It follows that, the approximations (28 ) and
m=0
(29) can be used to avoid numerical problems in the direct
= ZH l,k Zl,k ,
computation of the normalized detection threshold.
Blackman Window
2
Normalized Threshold T/(LNσ )
12
10
8
N −1 X
and let
sl,Ik +m (sl,Ik +m )∗
m=0 (30)
= [sl,Ik , sl,Ik +1 , · · · , sl,Ik +N −1 ]t = Fk WRl ,
Zl,k
14 True Threshold (NO) True Threshold (O) Upper Bound: Eq.(19) Upper Bound: Eq. (23) Lower Bound: Eq. (20) Lower Bound: Eq. (23) Approximation: Eq. (28) Approximation: Eq. (29)
0 ≤ γ ≤ 1/2 0 ≤ l ≤ L − 1, we have
Proof of Theorems 1-3. Assume
Ik = kK +
(31)
Rl is the l-th sample vector dened by (1), M K × M K diagonal matrix with its diagonal dened by w:
where in (31),
W
is the
elements
w0 0 0 w1 W= ··· ··· 0 ···
6
··· ··· 0 ··· ··· 0 , ··· ··· ··· · · · · · · wM K−1
(32)
Fk is the N × M K matrix consisting of N rows of M K × M K inverse discrete Fourier transform matrix F with row indices Ik , Ik +1, · · ·, Ik +m, · · ·, Ik +N −1:
and
4
the
2
0 −10 10
−8
−6
10
−4
10
−2
10
Fk = 1 ··· 1 ··· 1
0
10
10
PFA Fig. 6.
Upper and lower bounds and true and approximate normalized
detection thresholds as a function of the probability of false alarm
N = L = 5, M K = 1024, γ = 0
or
γ = 1/2,
Pf a ,
··· αIk +N −1
L−1 X
,
The problem of computing the normalized detection detector has been comprehensively treated. Formulas for computing the probability of false alarm,
Pf a ,
as a
Hence
zl,k =
l=0
threshold for the FFT lter bank-based summation CFAR where
Z
L−1 X
H ZH l,k Zl,k = Z Z,
(34)
l=0
is the length-LN zero mean Gaussian random
vector dened by:
were derived from
£ ¤t Z = Zt0,k , Zt1,k , · · · , ZtL−1,k .
rst principles. In addition, accurate approximations for
T
Ik +m
2πj α = exp M K.
VII. C ONCLUSIONS
T,
α
··· α(M K−1)Ik ··· ··· (M K−1)(Ik +m) ··· α ··· ··· · · · α(M K−1)(Ik +N −1)
Blackman window.
where
function of the detection threshold,
(33)
α Ik ···
(35)
that avoid potential numerical problems are presented.
These results are very useful for the simulation and implementation of the FFT lter bank-based summation CFAR detector.
From the identity (4a) of [16], the characteristic function of the quadratic form
φ(t),
are
grateful
to
the
reviewer
for
φ(t) = very
l=0
zl,k = ZH Z,
1 , det |I − jtE (ZZH )|
(36)
helpful
comments and suggestions which greatly improved the
where
readability of this paper. We would also like to thank
LN .
I
denotes the identity matrix of dimensions
LN ×
We shall derive the equations (5), (11) and (13)
φ(t)¡. First,¢it is necessary E ZZH explicitly. In
our colleague Dr. Sreeraman Rajan at DRDC Ottawa for
using the characteristic function
carefully reviewing an earlier version of this paper.
to compute the Hermitian matrix
© 2007 ACADEMY PUBLISHER
denoted by
can be shown to be given by:
VIII. ACKNOWLEDGEMENTS We
PL−1
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
fact,
¡ ¢ E ZZH
¡ H¢ E ZZ = E
¡ ¢ E ZZH
is given by:
41
³ ´ E Z0,k ZH 0,k · · · ³ ´
Z0,k Z1,k ··· Zp,k ··· ZL−1,k ··· ···
E Zp,k ZH ··· 0,k · · · ³ ´ ··· H E ZL−1,k Z0,k ···
Z0,k Z1,k ··· Zp,k ··· ZL−1,k
H
=
³ ´ E Z0,k ZH L−1,k · · · ³ ´
E Zp,k ZH L−1,k · · · ³
E ZL−1,k ZH L−1,k
. ´ (37)
For any
l, m, 0 ≤ l ≤ m ≤ L − 1, ¡ ¢ E Zl,k ZH m,k ³ ´ H = E (Fk WRl ) (Fk WRm ) ¡ ¢ H H = E Fk WRl RH m W Fk ¡ ¢ H = Fk WE Rl RH m WFk ,
simplies to the following tridiagonal
(38)
¡ ¢ E ZZH = A0 A1 0 ··· AH A A 0 0 1 1 ··· ··· ··· ··· 0 ··· 0 AH 1 0 ··· 0 ··· ³ ´ H where A0 = E Z0,k Z0,k and A1 of dimensions N × N .
where
δJ−M K+1 ··· δJ−M K+n , ··· δJ
J = (m − l)(1 − γ)M K . δq = 0, if q¡ 6= 0 and ¢ δq = 1 if q = 0. H It is clear that the matrix E Rl Rm is shift invariant in the sense that it depends only on the value of m − l. Let ¡ ¢ Ap = E Zl,k ZH (40) m,k ¡ ¢ H H = Fk WE Rl Rm WFk , p = m − l, m ≥ l. ¡ H¢ The covariance matrix E ZZ can then be put in the with
Here in (39),
following block Toeplitz form:
¢ ¡ E ZZH = A0 A1 AH A0 1 H A2 AH 1 ··· ··· ··· ··· H A ··· L−3 H AL−2 · · · AH AH L−2 L−1
(41)
A2 A1 A0 ··· ··· AH 1 AH 2 ···
··· A2 A1 ··· ··· A0 AH 1 AH 2
AL−2 ··· ··· ··· ··· A1 A0 AH 1
AL−1 AL−2 AL−3 ··· ··· A2 A1 A0
.
m−l ≥ 2, then |(m−l)(1−γ)M K| ≥ 2(1−γ)M ¡ K ≥¢ 2(1 − 1/2)M K = M K and hence entries in E Rl RH m are all equal to zero. Thus Ap = 0 if p ≥ 2 and H =
0 0 ··· 0 ··· ··· , A0 A1 AH A0 1 ³ ´ = E Z0,k ZH 1,k
= σ 2 Fk W2 FH k . Ik = kK + 2πj M K , we obtain:
With
exp
are
K−N 2 ,
(44)
J = MK
and
α =
2 A0 = σ 2 Fk W2 FH k =σ × 2 (M K−1)Ik · · · wM w12 αIk w02 K−1 α ··· ··· ··· ··· 2 Ik +l w2 ··· ··· w1 α 0 ··· ··· ··· ··· 2 2 Ik +N −1 ··· ··· w0 w1 α 1 ··· 1 ··· ··· ··· −lIk −l(Ik +N −1) α · · · α × ··· ··· ··· α−(M K−1)Ik · · · α−(M K−1)(Ik +N −1) τ11 τ12 · · · τ1q · · · τ1N ··· ··· ··· ··· ··· ··· = σ2 τp1 τp2 · · · τpq · · · τpN ··· ··· ··· ··· ··· ··· τN 1 τN 2 · · · τN q · · · τN N
= σ 2 A, where, for
τpq
(45)
1 ≤ p ≤ q ≤ N, = =
If
© 2007 ACADEMY PUBLISHER
(42)
¡ ¢ A0 = E Z0,k ZH 0,k ¡ ¢ H = Fk WE R0 RH 0 WFk
(39)
··· ··· ··· ··· ···
We now consider the two separate cases γ = 0 and 0 < γ ≤ 1/2. ¡ H¢ 1) γ = 0. In this case, A1 = 0 and E ZZ becomes the following diagonal L × L block matrix ¢ ¡ E ZZH = (43) A0 0 0 ··· 0 0 0 A0 0 0 · · · 0 ··· ··· ··· ··· ··· ··· , 0 ··· 0 0 A0 0 0 · · · 0 · · · 0 A0
where it can be veried that
¡ ¢ 2 E Rl RH m =σ × δJ ··· δJ−m ··· ··· ··· δJ+n−1 ··· δJ+n−1−m ··· ··· ··· δJ+M K−1 · · · δJ+M K−1−m
L×L
block matrix:
MX K−1 l=0 MX K−1
wl2 αl(p−q) wl2 exp
l=0 and
A
is dened by (7).
2πjl(p − q) , MK
(46)
42
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
2)
0 < γ ≤ 1/2.
In this case
and it remains to compute
A0 is computed A1 . In fact,
¡ ¢ E R0 RH = 1 · 0γM K×(1−γ)M K σ2 0(1−γ)M K×(1−γ)M K
by (45)
and
F1k = JI α k ··· αJ(Ik +p−1) ··· αJ(Ik +N −1)
(47)
¸
IγM K×γM K 0(1−γ)M K×γM K
,
0γM K×(1−γ)M K , 0(1−γ)M K×(1−γ)M K and 0(1−γ)M K×γM K are zero matrices of dimensions γM K ×(1−γ)M K, (1−γ)M K ×(1−γ)M K and (1 − γ)M K × γM K respectively and IγM K×γM K is the identity matrix of dimensions γM K ×γM K .
where
Let
w0 0 W1 = ··· 0 w0 0 W2 = ··· 0
0 w1 ··· 0
··· 0 0 ··· , ··· ··· · · · wγM K−1 ··· 0 0 ··· , ··· ··· · · · w(1−γ)M K−1
0 w1 ··· 0
wγM K 0 ··· 0
wγM K+1 ··· 0
W4 = w(1−γ)M K 0 ··· 0
0 w(1−γ)M K+1 ··· 0
¡
··· ··· ··· ···
0 0 ···
where
(49)
·
0 W3
γpq =
,
γM K−1 X
αl(Ik +p−1) ×
(57)
l=0
wl wl+(1−γ)M K α−(l+(1−γ)M K)(Ik +q−1)
.
(52)
¸
=·σ 2 Fk × 0γM K×(1−γ)M K 0(1−γ)M K×(1−γ)M K
= e−2πj(1−γ)(Ik +q−1) × γM K−1 X 2πjl(p − q) . wl wl+(1−γ)M K exp MK
W1 W4 0(1−γ)M K×γM K
Summarizing the preceding calculations, we see that
¡ ¢ E ZZH = σ 2 H where H is dened by (6) with B = 0 if γ = 0. Theorems 1-3 can now be proved using the following identity:
1 2π =
¸
Z
+∞
T p−1 X s=0
where
×FH k .
·Z
¡ T ¢s λ
s!
T > 0, λ > 0
+∞
−∞
¸ e−2πjtx dt dx (1 − jtλ)p
T
e− λ ,
and
p
(58)
is a positive integer. This
identity can be proved using an elaborate but standard
2πj α = exp( M K ), J = (1 − γ)M K
F0k = 1 α Ik ··· ··· 1 αIk +p−1 ··· ··· 1 αIk +N −1
γ1N ··· γpN ··· γN N
(56)
with
= σ Fk × ¸ · 0γM K×(1−γ)M K IγM K×γM K 0(1−γ)M K×(1−γ)M K 0(1−γ)M K×γM K · ¸ W2 0 × FH k 0 W4
Let
··· ··· ··· ··· ···
l=0
¢
W1 0
¡ 1 ¢H Fk · · · γ1q ··· ··· · · · γpq ··· ··· · · · γN q
C = F0k W1 W4 γ11 γ12 ··· ··· γ γp2 = p1 ··· ··· γN 1 γN 2
wM K−1
A1 = E Z0,k ZH 1,k ¡ ¢ H = Fk WE R0 RH 1 WFk 2
(55)
= σ 2 C,
(51)
α ··· (M K−1)(Ik +p−1) . α ··· (M K−1)(Ik +N −1) α
¡ ¢ = E Z0,k ZH 1,k · ¸ 0 W1 W4 = σ 2 Fk FH k 0 0 ¡ ¢ H = σ 2 F0k W1 W4 F1k
(48)
··· 0 0 0 ··· 0 , ··· ··· ··· · · · 0 wM K−1
0
··· ··· ··· ··· ···
(54)
can be further simplied to yield:
A1
(50)
A1
W3 =
We obtain
Then
(M K−1)Ik
contour integral argument from the theory of functions of
and set
one complex variable. Details are omitted due to space
(53)
(γM K−1)Ik
··· α ··· ··· (γM K−1)(Ik +p−1) , ··· α ··· ··· (γM K−1)(Ik +N −1) ··· α
© 2007 ACADEMY PUBLISHER
constraints.
0 < γ ≤ 1/2. H dened by (6) denoted by λ1 >
We rst prove Theorem 1. Assume Let the distinct positive eigenvalues of be arranged in decreasing order and
λ2 > · · · > λLN > 0. From (36) it follows that the PL−1 H characteristic function φ(t) of l=0 zl,k = Z Z is given
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
43
by:
φ(t)
1 (59) = det |I − jtE (ZZH )| 1 1 = QLN = 2 det |I − jtσ 2 H| l=1 (1 − jtσ λl ) LN X
=
Y
λlLN −1
1 . 2λ 1 − jtσ l (λl − λm )
l=1 1≤m≤LN,m6=l This implies that
PL−1 l=0
zl,k = ZH Z
=
LN X σ 2 λl
zl,k = ZH Z ∼
l=0 where
2
l=1
χ22 (l), 1 ≤ l ≤ LN, ∼
m=1
Z
+∞
= =
T LN X
1 2π
Z
L−1 X
χ22 (l),
(60)
Y
T
=
2L
of
degrees of freedom. For a given
given by:
Pf a =
) zl,k > T
Z
is
+∞
=
p(x)dx
(65)
T
l=0
¸ Z +∞ 1 −2πjtx φ(t)e dt dx = 2π −∞ T ¸ ·Z +∞ Z N X L X e−2πjtx Amp +∞ dt dx = 2 p 2π T −∞ (1 − jtσ µm ) m=1 p=1 ´s ³ T p−1 N X L X X 2 σ µm − T e σ2 µm . = Amp s! m=1 p=1 s=0
Z
(61)
+∞
·
This completes the proof of Theorem 2.
¸
−2πjtx
(L−1 X
T > 0, Pf a
It remains to prove Theorem the covariance matrix
dt dx
¡ ¢ E ZZH
3. Since
N = 1, L×L
reduces to a
diagonal matrix with the diagonal elements all equal to
−1 λLN l
σ2
(λl − λm )
−∞
Y
are independent and
−∞
1≤m≤LN,m6=l Z +∞ ·Z +∞ T
(64)
identically distributed chi-square random variables each
φ(t)e−2πjtx dt, x > 0,
l=1
1 2π LN X
is distributed
N X σ 2 µm 2 χ2L (m), 2 m=1
χ22L (m), 1 ≤ m ≤ N,
where
means having identical prob-
φ(t)e
zl,k = ZH Z ∼
l=0
are independent and identi-
+∞
zl,k
variables:
T > 0, Pf a is then computed by: (L−1 ) Z +∞ X = zl,k > T = p(x)dx l=0
l=0
as the following weighted sum of chi-square random
−∞
Z +∞ ·
PL−1
the detection decision statistic
and for a given
Pf a
1 (1−jtσ 2 µm )L
are omitted due to space constraints. It then follows that
σ 2 λl 2 2 χ2 (l) is given by 1 . Let p(x) denote the probability density function 2 1−jtσ λl PL−1 H of the detection decision statistic l=0 zl,k = Z Z. We
1 2π
are dened by
QN
decomposition of the rational function
that the characteristic function of
p(x) =
(63)
(12). The routine but tedious derivations of the fractional
ability density functions. Note that here we used the fact
have
1 , (1 − jtσ 2 µm )p
Amp , 1 ≤ m ≤ N , 1 ≤ p ≤ L,
where
is distributed as the
cally distributed chi-square random variables, each of two degrees of freedom, and
Amp
m=1 p=1
following weighted sum of chi-square random variables:
L−1 X
N X L X
e 1 − jtσ 2 λl
−1 λLN l
l=1
2 from (36) that the characteristic l=0 wl . It follows PL−1 H function φ(t) of l=0 zl,k = Z Z is given by:
×
−2πjtx
(λl − λm )
e
−
φ(t) =
¸ dt dx
T σ 2 λl
,
= (62)
p = 1
was used. This
L−1 X
completes the proof of Theorem 1. To prove Theorem 2, let the distinct positive eigenvalues of the positive denite Hermitian matrix
A
dened
by (7) be arranged in decreasing order and denoted by
µ1 > µ2 > · · · > µN . It follows PL−1 from (36) that the φ(t) of l=0 zl,k = ZH Z is given
characteristic function by:
1 det |I − jtE (ZZH )| 1 = det |I − jtσ 2 H| 1 = QN 2 L m=1 (1 − jtσ µm )
φ(t) =
© 2007 ACADEMY PUBLISHER
1 det |I − jtE (ZZH )| 1 ³ ´L . P MK 1 − jtσ 2 l=0 wl2
(66)
PL−1
H l=0 zl,k = Z Z is distributed as a weighted chi-square random variable:
This implies that
1≤m≤LN,m6=l where the formula (58) with
PM K
zl,k = ZH Z ∼
σ2
PM K−1 l=0
2
l=0
χ22L
where
χ22L ,
(67)
is a chi-square random variable with
degrees of freedom. For a given computed by:
Pf a =
wl2
(L−1 X l=0
) zl,k > T
Z
T > 0, Pf a
2L
is then
+∞
=
p(x)dx T
¸ Z +∞ 1 −2πjtx φ(t)e dt dx = 2π −∞ T # " Z +∞ Z +∞ e−2πjtx 1 = PM K−1 2 L dt dx 2π T wl ) −∞ (1 − jtσ 2 l=0
Z
+∞
·
44
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
¶s
µ =
L−1 X
σ
PMTK−1 2 l=0
−
wl2
e
s!
s=0
σ2
PMTK−1 l=0
The function w2 l
.
(68)
ψn0 (x) = −
This completes the proof of Theorem 3. Proof of Theorem 4. For
fn (x)
by:
fn (x) = Pr
( n X σ 2 λl 2
l=1
n
where
x > 0,
dene the function and
) χ22 (l) > x ,
is an integer in the range
1 ≤ n ≤ LN .
fn+1 (x) = Pr ≥ Pr
( n X σ 2 λl 2
l=1 Also we have
f1 (x) = e
fn (x) =
n X
χ22 (l)
x σ 2 λ1
−
2
l=1
Y
) χ22 (l)
>x
)
>x
= fn (x).
λn−1 m (λm − λl )
e
−
x σ 2 λm
(71)
Let
µ0
be the solution in
gn (x) as follows: ¶ µ xn x2 . + ··· + gn (x) = e−x 1 + x + n! 2!
(72)
xn e−x , gn0 (x) = − n!
(73)
and
xn−1 (x − n)e−x . n!
(74)
gn is a decreasing function on the interval (0, ∞), gn0 (x) < 0 for x > 0. gn is concave on (0, n) 00 and convex on (n, ∞) since gn (x) < 0 on (0, n) and gn00 (x) > 0 on (n, ∞). It can be proved that gn (n) > 1/2 for all n ≥ 1 and gn (n) monotonically decreases to 1/2 as n approaches innity. Let the solution in x to the equation Clearly since
Dene
ψn (x)
µ.
0 < p < 1,
(75)
as the natural logarithm of
gn (x),
n X xl l=0
© 2007 ACADEMY PUBLISHER
ψn (x), it follows
x
(80)
to the equation
This linear equation in
x
is easily solved and
(81)
µ0
is
obtained as:
µ0
= =
We claim that
µ0
ln p − ψn (x0 ) ψn0 (x0 ) ln p − ψn (n) . n+ ψn0 (n)
x0 +
is an upper bound for the solution
(82)
µ
to
ln p ≥ ψn (µ0 ) = ln(gn (µ0 )). Hence p ≥ gn (µ0 ). Since gn (µ) = p, we have gn (µ) ≥ gn (µ0 ). Hence µ ≤ µ0 , since gn (x) is a decreasing function of x. Thus µ0 is an upper bound for µ. We can simplify the expression (82). Assume 0 < p ≤ 2/e. Since gn (n) is a monotone decreasing sequence, gn (n) ≤ g1 (1) = 2/e and therefore ψn (n) ≤ ln(2/e). It follows that
µ0
(83)
ln p − ln(2/e) ln p − ψn (n) ≤n+ =n+ ψn0 (n) ψn0 (n) # µ ¶ " n n3 n2 1 + n + 2! + 3! + · · · + nn! 2 . ln ≤n+ nn pe n! This upper bound can be further simplied. Using the
We now derive an upper bound for
ψn (x) = ln(gn (x)) = −x + ln
(79)
the equation (75). In fact, from (80) and (81), we see that
It can be veried that
be denoted by
ψn00 (x) < 0
ln p = (x − x0 )ψn0 (x0 ) + ψn (x0 ).
dene
= p,
(77)
(x − x0 )ψn0 (x0 ) + ψn (x0 ) ≥ ψn (x), x ∈ (0, ∞).
,
We next prove the inequality (18). For each positive
l!
,
y = (x − x0 )ψn0 (x0 ) + ψn (x0 ).
This proves the chain of inequalities in (16).
l=0
xn n!
that
fn (Tn (Pf a )) = Pf a , 1 ≤ n ≤ LN . Since fn (x) ≤ fn+1 (x), 1 ≤ n ≤ LN − 1, and since fn (x) are decreasing functions of x, we see immediately that Tn (Pf a ) ≤ Tn+1 (Pf a ), 1 ≤ n ≤ LN − 1. Also T1 (Pf a ) = −σ 2 λ1 ln(Pf a ) and TLN (Pf a ) = T .
n X xl
+ ··· +
for all x ∈ (0, ∞). Hence ψn (x) is concave on the positive real axis (0, ∞) and is always below its tangent lines. Let x0 = n. The equation of the tangent line of the curve of ψn (x) at the point (x0 , ψn (x0 )) is given by:
From the denition (15), we see that
gn (x) = e−x
x2 2!
As this tangent line is above the curve of
2 ≤ n ≤ LN.
gn00 (x) =
1+x+
The above expression shows that
1≤l≤n,l6=m
n,
xn n!
(78)
(70)
and
m=1
integer
fact,
Since
probability density functions are supported on the positive
(n+1 X σ 2 λl
shown to be concave on the
i h Pn−1 n 1 )xk − (k−1)! xn−1 n + k=1 ( k! . ψn00 (x) = − ¡ n ¢2 2 n! 1 + x + x2! + · · · + xn!
(69)
chi-square random variables are non-negative (i.e., their real axis), we have:
ψn (x) can be (0, ∞). In
positive real line
l!
.
µ.
that is, (76)
following result for large
n,
gn (n) = (84) ¶ µ n 3 2 n n n ∼ + ··· + + e−n 1 + n + = 1/2, n! 3! 2!
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
45
Proof of Theorem 5. From (64) we see that
the expression
1+n+
n2 2!
n3 3!
+
+ ··· +
nn n!
nn n!
Pf a = Pr
(85)
√ can be shown to be bounded from above by
√
2πn
³ n ´n e
(86)
we see that
Ã
1+n+ ³
≈
exp(n) 2 nn n!
√ ≈
2πn
n2 2!
+
n3 3!
+ ··· +
nn
nn n!
!
= ¡ n ¢n e
n! exp(n) 2nn
T ≥ BL (Pf a ) σ 2 µ1
exp(n)
2nn
=
2πn . 2
Pf a = Pr
(87)
µ≤n+
2πn ln 2
µ
2 pe
¶ , 0 < p ≤ 2/e.
(88)
The inequality (18) then follows directly from (88). Finally we prove (19). Since
Pf a = Pr
(L−1 X
) zl,k > T
where
l=0
= Pr
( LN X σ 2 λl 2
l=1 ( LN X
)
where
χ22LN
³ Pf a ≤
s=0
T σ 2 λ1
s!
χ22LN
e
.
T ≤ BLN (Pf a ), σ 2 λ1
T ≤ σ 2 µ1 BLN (Pf a ).
(95)
in increasing order as follows: (96)
Then we have
(89)
Ly1 ≤
2LN
(90)
L−1 X
zl,k ≤ LyL ,
(97)
l=0 and it follows that
(L−1 X
) zl,k > T
(98)
l=0
≥ Pr {Ly1 > T } = Pr {y1 > T /L} = Pr {zl,k > T /L, 1 ≤ l ≤ L − 1} L−1 Y = Pr {zl,k > T /L} l=0
(91)
" =
which completes the proof of (19).
or
y1 ≤ y2 ≤ y3 ≤ · · · ≤ yL .
Pf a = Pr T σ 2 λ1
2LN
It remains to prove (23). Rearrange the power levels
This implies
© 2007 ACADEMY PUBLISHER
is a chi-square random variable with
zl,k , 0 ≤ l ≤ L − 1,
´s −
(94)
Combining (93) and (95) yields (20).
degrees of freedom, we have
LN −1 X
zl,k > T
l=0
T ≤ BLN (Pf a ) σ 2 µ1
)
is a chi-square random variable of
)
degrees of freedom. This implies that
χ22 (l) > T
σ 2 λ1 2 χ2 (l) > T 2 l=1 ¾ ½ 2 σ λ1 2 χ2LN > T = Pr 2 ´s ³ T LN −1 X σ 2 λ1 − T e σ2 λ1 , = s! s=0
≤ Pr
(93)
) N X σ 2 µm 2 χ2L (m) > T = Pr 2 m=1 ) ( N X σ 2 µ1 χ22L (m) > T ≤ Pr 2 ¾ ½ m=1 σ 2 µ1 2 χ2LN > T = Pr 2 ´s ³ T LN −1 X σ 2 µ1 − T e σ2 µ1 , = s! s=0
µ: √
(92)
σ 2 µ1 BL (Pf a ) ≤ T.
(L−1 X
(
Substituting (87) into (83), we obtain the following upper bound for
or
Similarly from (64) we also see that
√
)
This implies that
n!
´
zl,k > T
N X σ 2 µm 2 χ2L (m) > T = Pr 2 m=1 ¾ ½ 2 σ µ1 2 χ2L (1) > T ≥ Pr 2 ´s ³ L−1 X σ2Tµ1 − T e σ 2 µ1 . = s! s=0
2πn 2 . In
,
)
l=0
(
fact, using Stirling's formula
n! ∼ =
(L−1 X
N X m=1
T
−1 µN − L m Q e σ2 µm 1≤l≤N,l6=m (µm − µl )
#L .
46
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
or
This implies
N X m=1
−1 µN m Q e 1≤l≤N,l6=m (µm − µl )
T − 2L σ µm
≤
1/L Pf a
k6=i (αk − α) αim Q k6=i (αk − αi ) i=1
(99)
T ≥ T1 , L
LT1 ≤ T.
or
n X
(100)
(−1)n−1
i=1
Similarly, we have
(L−1 X
)
If
zl,k > T
m = n − 1,
(101)
n X
l=0
≤ Pr {LyL > T } = Pr {yL > T /L} = 1 − Pr {yL ≤ T /L} = 1 − Pr {zl,k ≤ T /L, 1 ≤ l ≤ L − 1} L−1 Y =1− Pr {zl,k ≤ T /L} =1−
l=0 L−1 Y
=1− 1−
Letting
α
n X
αin−1−m Q
n−1
(−1)
(−1)
βi
T
#L
α)
= 1.
k6=i (αk − αn−1
α→+∞
Q lim
−1 µN − L m Q e σ2 µm 1≤l≤N,l6=m (µm − µl )
k6=i (αk − αn−1
βi lim
i=1
α→+∞
(109)
k6=i (αk − αn−1
α)
α)
= 1.
= (−1)n−1 .
(110)
(111)
Substituting this result into (110) yields:
.
n X
βi = 1.
(112)
i=1
N X
T
Q
m=1
−1 µN − L m e σ2 µm 1≤l≤N,l6=m (µm − µl ) 1/L
≥ 1 − (1 − Pf a )
,
This proves the rst identity in (105). To prove the
T ≤ T2 , L
or
α=0
(102)
T ≤ LT2 .
(103)
in (108). We then obtain the identity
n X
Y
βi αk n−1−m α i=1 i k6=i
Dividing both sides of (113) by
Combining (100) and (103) yields the inequalities of
n X
(23). The inequality (24) is obtained by using (102) and applying the argument in the proof of (19) in Theorem 4. The proof of Theorem 5 is thus completed. To prove Theorem 6, we need to prove two lemmas. Lemma 1. Let
n ≥ 2.
Assume
0 < α1 < α2 < · · ·
0.
(116)
Proof of Lemma 2. We prove this lemma by mathemat-
1 ≤ m ≤ n−1 be an integer and m dene the polynomial Pm by Pm (α) = α . According to Lagrange's Interpolation Theorem, the following identity holds for all
n X
© 2007 ACADEMY PUBLISHER
i=1
where
Proof of Lemma 1. Let
α ∈ (−∞, +∞): Q k6=i (αk − α) = Pm (α), Pm (αi ) Q k6=i (αk − αi ) i=1
1 ≤ m ≤ n−1
remaining identities in (105), assume and let
which in turn yields
Then
(108)
k6=i
Q
This implies
αn
(αk − α) = αm .
go to innity yields the identity
n−1
[1 − Pr {zl,k > T /L}]
m=1
(107)
(108) can be rewritten as
But
N X
Y
βi
i=1
l=0
"
= αm .
This identity can be rewritten as
which in turn yields
Pf a = Pr
Q
n X
(106)
ical induction. First consider the case
n = 2.
We have
C(p, α , α ) = (117) · ¸ ·1 2 ¸ α2 α2 − 1 pα2 /α1 β1 α1 α1 ¸ ¸ · · ¸· α1 α2 α2 < 0, − 1 pα2 /α1 = α1 − α2 α1 α1
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
since
α2 > α1
and
p ∈ (0, 1).
47
∂E ∂p (p, α1 , α2 , · · · , αn ), viewed as a function of p, is strictly on the interval [0, 1]. Using the Pn decreasing βi identity, i=1 α = 0, proved in Lemma 1, we see that
Hence,
On the other hand,
D(p, α1 , α2 ) · ¸h · ¸h i i α2 α2 α2 /α1 pα2 /α2 β2 p β1 + = α2 α ¸ ¸ · · 1¸· α2 α1 α2 p pα2 /α1 − = α1 − α2 α1 − α2 α1 ¸h · i α2 pα2 /α1 − p > 0, (118) = α1 − α2 α2 > α1 and p ∈ (0, 1). This proves Lemma 2 for the case n = 2. Next assume Lemma 2 holds for positive integer n − 1. We show that Lemma 2 must also hold for the integer n. To do this, we reformulate the expression C(p, α1 , · · · , αn ) as follows:
i
∂E (p, α1 , α2 , · · · , αn )|p=1 ∂p ¡ −1 ¢ = p D(p, α1 , α2 , · · · , αn ) p=1 = D(1, α1 , α2 , · · · , αn ) n X βi = 0. = αn α i=1 i
again using the fact that
C(p, α1 , α2 , · · · , αn ) (119) ¸ ¸ · · n−1 X αn αn − 1 pαn /αi βi = α α i i i=1 n−1 X · αn ¸ · αn − αi ¸ αin−1 Q pαn /αi = αi αi k6=i (αi − αk ) i=1 = −αn = −
n−1 X
αin−2 1 Q pαn /αi αi 1≤k≤n−1,k6=i (αi − αk )
i=1 n−1 αn X
αn−1
i=1
³
αn−1 αn /αn−1 p αi
´αn−1 /αi
βi∗ ,
It follows that for all
(127)
and consequently
D(p, α1 , α2 , · · · , αn ) ¶ µ ∂E (p, α1 , α2 , · · · , αn ) > 0. =p ∂p
(128)
By the principle of mathematical induction, Lemma 2
n ≥ 2.
must hold for all integers
Proof of Theorem 6. Applying Lemma 2 to the
λLN < λLN −1 < · · · < λ2 < λ1 , we g 0 (p) > 0 and g 00 (p) for 0 < p < 1.
see immediately that
αin−2 . 1≤k≤n−1,k6=i (αi − αk )
βi∗ = Q
Proof of Theorem 7. Since (c.f. (27))
(120)
L−1 X
Hence
zl,k ∼
l=0
C(p, α1 , α2 , · · · , αn ) = αn αn D(p αn−1 , α1 , α2 , · · · , αn−1 ). − αn−1 D(p, α1 , α2 , · · · , αn−1 ) > 0,
(121)
To prove
Pf a = Pr
we have
D(p, α1 , α2 , · · · , αn ) > 0,
E(p, α1 , α2 , · · · , αn ) =
n X
( = Pr
where
(123)
i=1 It can be veried that
∂E ∂p (p, α1 , α2 , · · · , αn ) = p−1 D(p, α1 , α2 , · · · , αn ), ∂2E ∂ 2 p (p, α1 , α2 , · · · , αn ) = p−2 C(p, α1 , α2 , · · · , αn ).
C(p, α1 , α2 , · · · , αn ) < 0,
) zl,k > T
l=0 N X
)
µm χ22L (m) > 2T /σ 2 (129)
¡ ¢ 1/2 h0 = c32 /c23 , y = 2T /σ 2 − c1 (h0 /c2 ) +h0 and PN c1 = 2L Pm=1 µm = 2LN, N (130) c2 = 2L m=1 µ2m , P c = 2L N µ3 . 3 m=1 m
Here we used the fact that
PN m=1
µm = N , which follows
from the assumption that the window is normalized. We have
³
it follows that
∂2E (p, α1 , α2 , · · · , αn ) ∂2p = p−2 C(p, α1 , α2 , · · · , αn ) < 0. © 2007 ACADEMY PUBLISHER
(124)
(L−1 X
© m=1 ª ≈ Pr χ2h0 > y ,
let
pαn /αi βi .
N X σ 2 µm 2 χ2L (m), 2 m=1
we see from Eq. (4.1) of [15] that
C(p, α1 , α2 , · · · , αn ) = (122) αn αn αn−1 , α1 , α2 , · · · , αn−1 ) < 0. D(p − αn−1
Since
0 < p < 1,
∂E (p, α1 , α2 , · · · , αn ) ∂p ¶ µ ∂E (p, α1 , α2 , · · · , αn ) |p=1 = 0, > ∂p
increasing sequence
where
Since
(126)
h0
=
2L ³ 2L
(125)
³P
PN
´3
PN
´2 = 2L ³P
2 m=1 µm 3 m=1 µm
N m=1
µ2m
N m=1
µ3m
´3 ´2 , (131)
48
JOURNAL OF COMPUTERS, VOL. 2, NO. 6, AUGUST 2007
and
[8] S. Wang, F. Patenaude and R. Inkol, Upper and Lower Bounds for the Threshold of the FFT Filter Bank-Based Summation CFAR
1/2
0
(h /c2 ) ¢3 1/2 ¡PN µ2 m=1 m ¡ ¢2 P 2L PN N µ3 µ2m m=1 m . = Pm=1 = PN N 3 2L m=1 µ2m m=1 µm
Detector, Proceedings of ICASSP 2006, Toulouse, France, vol.3, pp. 289-301, May 2006. [9] S. Wang, R. Inkol, S. Rajan and F. Patenaude, Theory of the FFT Filter Bank-Based Majority and Median CFAR Detectors, DRDC
(132)
Ottawa TR 2007-088, May 2007. [10] S. Wang and R. Inkol, Operating Characteristics of the Wideband FFT Filter Bank
J -out-of-L
CFAR Detector, DRDC Ottawa TR
2003-235, December 2003. [11] R. Inkol and S. Wang, A Comparative Study of FFT-Summation
It follows that
and Polyphase-FFT CFAR Detectors, Canadian Conference on
¡ ¢ 1/2 y = 2T /σ 2 − c1 (h0 /c2 ) + h0 ¡
= 2T /σ 2 − 2LN
¢
³P
PN
2 m=1 µm PN 3 m=1 µm
(133)
N m=1
µ2m
N m=1
µ3m
+ 2L ³P
´3
´2 .
© ª © ª Pf a ≈ Pr χ2h0 > y ≈ Pr χ22h > y ª © ≈ Pr χ22h /2 > y/2 ¶ µ z h−1 z2 , + ··· + = e−z 1 + z + (n − 1)! 2! $ ¡P ¢3 % N
z = y/2
Summation CFAR Detectors: A Comparative Study, Canadian Conference on Electrical and Computer Engineering, pp. 10391044, May 2004. [13] S. Wang and R. Inkol, Theoretical Performance of the FFT Filter Bank Based Summation Detector, DRDC Ottawa TR 2005-153, November 2005.
Consequently
where
Electrical and Computer Engineering, pp. 1175-1178, May 2004. [12] S. Wang and R. Inkol, FFT Filter Bank-Based Majority and
and
h = L ¡P
m=1 N m=1
µ2m
µ3m
¢2
[14] S. Wang, R. Inkol and S. Rajan, A Formula for the Probability of False Alarm for the FFT Filter Bank-Based August 2005.
(134)
¡ ¢ T /σ 2 − LN
2 m=1 µm PN 3 m=1 µm
in Normal Variables, Biometrika, vol.48, pp. 419-426, December 1961. [16] G.L. Turin, The Characteristic Function of Hermitian Quadratic Forms in Complex Normal Variables, Biometrika, vol.47, pp. 199201, June 1960.
. Hence
[17] H. Solomon and M. A. Stephens, Distribution of a Sum of Weighted Chi-Square Variables, Journal of the American Statistical Association, vol.72, no.360, pp. 881-885, December 1977.
³P
PN
CFAR
[15] J. P. Imhof, Computing the Distribution of Quadratic Forms
z = y/2 ≈ Bh (Pf a ), or
J -out-of-L
Detector, Proceedings of MWSCAS 2005, Cincinnati, Ohio, USA,
N m=1
µ2m
N m=1
µ3m
+ L ³P
≈ Bh (Pf a ).
´3 Sichun Wang (DRDC Ottawa and CRC) received the B.Sc. and M.Sc.
´2
degrees, both in Mathematics, from Nankai University in 1983 and 1989, respectively. He received the PhD degree in Mathematics from McMaster University in 1996. He was a postdoctoral fellow at the
(135)
Communications Research Laboratory of McMaster University from February 1996 to August 1996 and a research scientist at Telexis
The approximation (28) then follows immediately from
and Intrinsix Corporations from September 1997 to May 2002. From
(135).
June 2002 to March 2007 he worked as a research consultant for
The proof of Theorem 8, as it is almost identical to that of Theorem 7, is omitted here.
Defence Research and Development Canada. Currently he is a research scientist at the Communications Research Centre (CRC) on secondment to Defence Research and Development Canada. His recent research has focused on the FFT lter bank-based CFAR detection algorithms and
R EFERENCES
turbo codes. Dr. Wang is a senior member of the IEEE.
[1] E. O. Brigham, The Fast Fourier Transform and Its Applications, Prentice-Hall, Inc., NJ, USA, 1988. [2] G. A. Zimmerman and S. Gulkis, Polyphase-Discrete Fourier Transform Spectrum Analysis for the Search for the Extraterrestrial Intelligence Sky Survey, TDA progress report 42-107, Jet Propulsion Laboratory, Pasadena, CA., pp. 141-154, July-September, 1991. [3] G. Lanyi and R. Khan, Tone Detection via Incoherent Averaging of Fourier Transforms to Support the Automated SpacecraftMonitoring Concept, TDA Progress Report 42-129, Jet Propulsion Laboratory, Pasadena, CA., pp. 1-22, May 1997. [4] R. S. Walker, The Detection Performance of FFT Processors for Narrowband Signals, D.R.E.A. Techincal Memorandum 82/A, Defence Research Establishment Atlantic, Research and Development
Franc ¸ ois Patenaude (CRC) received the B.A.Sc. degree from the University of Sherbrooke in 1986 and the M.A.Sc. and Ph.D. degrees from the University of Ottawa in 1990 and 1996, all in Electrical Engineering. His graduate work was conducted in collaboration with the Mobile Satellite Group of the Communications Research Centre (CRC), Ottawa, Canada. In 1995, he joined CRC to work on signal processing applications for communications and for spectrum monitoring. His main research interests include automatic signal recognition, detection and estimation in the spectrum monitoring context, and real-time signal processing.
Branch, Department of National Defence, Canada, February 1982. [5] B. H. Maranda, On the False Alarm Probability for an Overlapped FFT Processor, IEEE Transactions on Aerospace and Electronic Systems, vol.32, no.4, pp. 1452-1456, October 1996. [6] H. C. So, Y. T. Chan, Q. Ma, P. C. Ching, Comparisons of Various
Robert Inkol (DRDC Ottawa) received the B.Sc. and M.A.Sc. degrees
Periodograms for Sinusoidal Detection and Frequency Estimation,
in Applied Physics and Electrical Engineering from the University of
IEEE Transactions on Aerospace and Electronic Systems, vol.35,
Waterloo in 1976 and 1978, respectively. Since 1978 he has been with
no.3, pp. 945-952, July 1999.
Defence Research and Development Canada where he is currently a
[7] F. Patenaude, D. Boudreau and R. Inkol, CFAR Detection Based on
senior scientist. He has made numerous contributions to the application
Windowed and Polyphase FFT Filter Banks for Channel Occupancy
of very large scale integrated circuit technology and digital signal
Measurements, Proceedings of the 19th Biennial Symposium on
processing techniques to electronic warfare systems. He is an adjunct
Communications, Kingston, Ontario, Canada, pp. 339-343, May
professor with the Royal Military College and is a senior member of
1998.
the IEEE.
© 2007 ACADEMY PUBLISHER