Discrete Fractional Fourier Transform Based on New Nearly Tridiagonal
Commuting Matrices
Soo-Chang Pei, Wen-Liang Hsue, and Jian-Jiun Ding
Department of Electrical Engineering, National Taiwan University, Taipei, Taiwan, R.O.C.
Email address: [email protected]
ABSTRACT
Based on discrete Hermite-Gaussian like functions, a discrete fractional Fourier transform (DFRFT) which provides sample approximations of the continuous fractional Fourier transform was defined and investigated recently. In this paper, we propose a new nearly tridiagonal matrix which commutes with the dis-crete Fourier transform (DFT) matrix. The eigenvectors of the new nearly tridiagonal matrix are shown to be better discrete Hermite-Gaussian like functions than those developed before. Furthermore, by appropriately combining two linearly independ-ent matrices which both commute with the DFT matrix, we de-velop a method to obtain even better discrete Hermite-Gaussian like functions. Then, new versions of DFRFT produce their transform outputs more close to the samples of the continuous fractional Fourier transform, and their application is illustrated.
1. INTRODUCTION
The ath-order continuous fractional Fourier transform (FRT) of
x(t) is defined as [4] , ) , ( ) ( ) (u xt K t u dt Xa
³
f a f (1)where the transform kernel Ka( ut, ) is given by , cot 1 ) , (
D
jS(t2cotD 2tucsc(D) u2cotD) a t u j e K (2)in which өːaӸ˂˅. It is known that the transform kernel )
, ( ut
Ka can also be written as [4]
), ( ) ( ) 2 / exp( ) , ( 0 u t jna u t K n n n a
¦
< < fS
(3) where 2 ) 2 ( ! 2 2 ) ( 4 1 t n n n H t e n t <S
S (4)is the nth-order Hermite-Gaussian function with Hn being the nth -order Hermite polynomial.
The NuN DFT matrix F is defined by 1 , 0 , ) / 1 ( 2 d d e k n N N jNkn kn S F . (5)
In [5], Dickinson and Steiglitz introduced an NuN nearly tridi-agonal matrix S whose nonzero entries are:
. 1 2) -( 0 , 1 1) -( 0 ), 2 cos( 2 1 , 0 0 , 1 , 1 1 , , d d d d N N n n n n n n N n N n n N S S S S S
S
(6)With S defined above, S commutes with F, i.e., SF=FS. There-fore, the DFT matrix F and the above matrix S share a common eigenvector set and we can find the eigenvectors of F from those of the matrix S [3].
Analogous to the spectral expansion of the continuous FRT kernel Ka( ut, ) in (3), and from the fact that the eigenvectors of
S can be used as the discrete Hermite-Gaussian like functions, in
[2], Pei et al. defined the ath-order DFRFT matrix FSa by
° ° ° ¯ °° ° ®
¦
¦
even for , odd for , 2 0 2 2 1 0 2 N e e N e N k T N N Na j T k k ka j N k T k k ka j T a a v v v v v v V VD FS S S S (7)where T denotes the matrix transpose, the matrix
] | | | | [v0 v1 vN2 vN1
V " for odd N and
" | | [v0 v1
V |vN2 |vN] for even N, D is a diagonal ma-trix with its diagonal entries corresponding to the eigenvalues for each column eigenvectors in V, and vkis the kth-order discrete Hermite-Gaussian like function with k zero-crossings and is obtained from the corresponding normalized eigenvector of S. The S-based DFRFT of x can be easily obtained by ya FSax.
2. A NEW NEARLY TRIDIAGONAL
COMMUTING MATRIX T
In [1], GrΏnbaum introduced an exactly tridiagonal matrix com-muting with the centered discrete Fourier transform matrix of even size. Inspired by the work of GrΏnbaum, we propose in this section a novel nearly tridiagonal matrix which commutes with the ordinary DFT matrix of any size, even or odd. Moreover, we will demonstrate that its eigenvectors approximate samples of the continuous Hermite-Gaussian functions better than those of the S matrix. Therefore, we can intuitively expect better per-formance of defining its DFRFT using the new nearly tridiagonal matrix.
IV - 385
Let us define an NuN nearly tridiagonal matrix T whose nonzero entries are (note that the matrix indices are from 0 to
N-1):
. 5 . 0 2) -( 0 , / cos 2 ) 1 ( cos cos 1) -( 0 ), ( cos 1 , 0 0 , 1 , 1 1 , 2 , d d d d N N n n n n n n N n N N n N n N n N n T T T T T
S
S
S
S
(8)Note that except for the two 0.5 entries at the upper-right and lower-left corners, T is tridiagonal, which is similar to the S matrix of (6). Thus we call them nearly tridiagonal. Since T is real and symmetric, T has real and orthogonal eigenvectors. Besides, T has the following important property for this paper.
Property 1: The NͪN matrix T commutes with the NͪN DFT matrix F defined in (5), i.e., TF=FT.
From Property 1, it can be seen that if x is the eigenvector of T corresponding to an eigenvalue of multiplicity 1, then x is also an eigenvector of F. It can be shown that the entries of the eigenvectors of T are solutions of a discrete version of the defin-ing second-order differential equation of the continuous Her-mite-Gaussian functions [3]. Therefore, the eigenvectors of T are discrete Hermite-Gaussian like functions. To motivate our fur-ther discussions, we perform some computer experiments to show that the eigenvectors of T approximate samples of con-tinuous Hermite-Gaussian functions better than those of S.
-2 -1 0 1 2 -0.4 -0.2 0 0.2 0.4 (a) 4th-order -2 -1 0 1 2 -0.4 -0.2 0 0.2 0.4 (b) 6th-order -2 -1 0 1 2 -0.4 -0.2 0 0.2 0.4 (c) 8th-order -2 -1 0 1 2 -0.4 -0.2 0 0.2 0.4 (d) 10th-order
Fig. 1. The continuous Hermite-Gaussian functions (solid line), the discrete Hermite-Gaussian like functions based on S (‘*’) and based on T (‘o’), with N=25. (a) 4th-order: The error norms
of S and T are 0.0719 and 0.0312, respectively. (b) 6th-order:
The error norms of S and T are 0.1427 and 0.0579, respectively. (c) 8th-order: The error norms of S and T are 0.2637 and 0.0959, respectively. (d) 10th-order: The error norms of S and T are 0.4965 and 0.1472, respectively. 0 5 10 15 20 25 0 0.2 0.4 0.6 0.8 1 1.2 1.4 Order E rr o r N o rm
Fig. 2. The error norms of the discrete Hermite-Gaussian like functions based on S (‘*’), and the discrete Hermite-Gaussian like functions based on T (‘o’) of various orders, with N=25.
Computer experiment 1: Fig. 1 (a)-(d) show the 4th, 6th, 8th, and 10th-orders continuous Hermite-Gaussian functions, the crete Hermite-Gaussian like functions based on S, and the dis-crete Hermite-Gaussian like functions based on T, with N=25. The error norms, which are the Euclidean norms of the error vectors between the discrete Hermite-Gaussian like functions based on S (or T) and samples of the continuous Hermite-Gaussian functions, are plotted in Fig. 2 (with N=25). Fig.1 and Fig. 2 both demonstrate that the discrete Hermite-Gaussian like functions based on T are better than those based on S. The error norms of the discrete Hermite-Gaussian like functions based on both S and T tend to increase for higher order ones because of the aliasing effects.
From the definition of T in (8), we can express T in block matrix form as:
1) If N is odd, , 5 . 0 5 . 0 5 . 0 5 . 0 1 1 1 1 1 1 1 » » » ¼ º « « « ¬ ª J JT JAJ Je A T e J e e T T T (9) with e1 being T
]
0
,
,
0
,
1
[
"
of size (N-1)/2, and T1 and A beingthe N21
u
N21 submatrices of T. J is the exchange matrix with ones on the antidiagonal.2) If N is even, , 5 . 0 0 0 5 . 0 5 . 0 0 5 . 0 1 2 1 2 1 1 1 » » » » » ¼ º « « « « « ¬ ª J JT 0 0 Je 0 0 0 0 T e J e e T T T (10)
with e1 being the vector
T
]
0
,
,
0
,
1
[
"
of size (N2 1), and T 2being the (N2 1)u(N2 1) submatrix of T. Then, we have the following property.
Property 2: For the NuN T matrix defined in (8), the trans-formed matrix
» ¼ º « ¬ ª 2 1 M 0 0 M UTU T (11)
is a block diagonal matrix, where U is the NuN unitary symmet-ric matrix defined by
° ° ° ° ° ° ¯ °° ° ° ° ° ® » » » » » » ¼ º « « « « « « ¬ ª » » » » ¼ º « « « « ¬ ª even, is if , 2 0 0 2 2 1 odd, is if , 2 2 1 1 1 1 1 2 2 2 2 2 1 2 1 2 1 2 1 N N N N N N N N N N I 0 J 0 0 0 J 0 I 0 0 0 I J 0 J I 0 0 0 U (12)
with Jq beingʳ the q×q exchange matrix, and M1 and M2 are two
square matrices of sizes
¬
N2 1¼
and¬ ¼
2 1N , respectively.
¬ ¼
xdenotes the largest integer less than or equal to x. Moreover, from (9) and (10), M1 and M2 in (11) are respectively
° ° ° ° ¯ ° ° ° ° ® » » » » » ¼ º « « « « « ¬ ª » » ¼ º « « ¬ ª , if iseven, 0 0 0 1 odd is if , 1 1 1 2 1 2 2 1 2 2 1 1 2 1 1 2 1 1 2 2 N N T T T N N 0 0 T e e AJ T e e M (13) with 1 2 N
0 being the (N2 1)u1 zero vector, and
°¯ ° ® even. is if , odd is if , 1 2 1 1 2 2 2 2 1 2 1 2 1 N N N N N N N J T J A J J T J M (14)
From [6], we know that any symmetric and exactly tridi-agonal matrix with nonzero subditridi-agonal entries has distinct ei-genvalues. From Property 2, it can then be shown that M1 and M2 has distinct eigenvalues with the exception that the zero
ei-genvalue of M1 is of multiplicity two when N is even. We can
also show that most of the even extensions and odd extensions of eigenvectors of M1 and M2, respectively, are eigenvectors of F.
But if N is even, the even extensions of eigenvectors of M1
cor-responding to the zero eigenvalue are not necessarily eigenvec-tors of F. Thus, for N even, we need to develop a method to compute the eigenvectors of F in the even subspace spanned by eigenvectors of T of the zero eigenvalue.
Property 3: If N is even, the two orthogonal eigenvectors
of F in the subspace spanned by even eigenvectors of T of
ei-genvalue zero are
>
@
1 2 1 , 1 , , 1 , 1 , 1 , 1 r N N T e " , where 1 2 N
e is the Nu1 column vector with zero entries except a 1 at the (N2 1)th entry.
3. LINEAR COMBINATIONS OF
MATRICES S AND T
Property 4: If k1 and k2 are any two constants, then k1S+ k2T
commutes with the DFT matrix F, where S and T are defined in (6) and (8), respectively.
From Property 4, we can compute the eigenvectors of DFT matrix F using k1S+k2T. Since k1S+k2T and S(k2 /k1)T
have the same eigenvectors if k1 is nonzero, we discuss in the
following linear combinations of S and T of the form S+kT. If
k>0, we find from computer experiments that the eigenvalues of
S+kT are distinct. Therefore, from Property 4, eigenvectors of S+kT are all eigenvectors of F if k>0. Because the eigenvectors
of both S and T are discrete Hermite-Gaussian like functions, we can expect that eigenvectors of S+kT are also discrete Hermite-Gaussian like functions. We next show through computer ex-periments that eigenvectors of S+kT are new versions of discrete Hermite-Gaussian like functions and, with appropriate choice of
k, these new discrete Hermite-Gaussian like functions
approxi-mate samples of the continuous Hermite-Gaussian functions better than those obtained from both S and T.
Computer experiment 2: To determine the optimal choice of k, we first compute the eigenvectors of S+kT, which are new
versions of discrete Hermite-Gaussian like functions. All of the resulting N eigenvectors are compared with samples of the con-tinuous Hermite-Gaussian functions of the corresponding orders and the total error norms are calculated. For N=25 and N=145, the total error norms are plotted versus various values of k (from
k=0 to k=50 with spacing 1) in Fig. 3(a) and Fig. 3(b),
respec-tively. From these results and other experiments for different values of N (up to 145), we find that the optimal k is approxi-mately 15. 0 10 20 30 40 50 60 70 80 90 100 110 120 130 140 (b) k 0 10 20 30 40 50 4 5 6 7 8 9 10 11 12 (a) k
Fig. 3. Total error norms of discrete Hermite-Gaussian like func-tions based on S+kT. (a) N=25. (b) N=145.
4. DISCRETE FRACTIONAL FOURIER
TRANSFORM BASED ON T OR S+
kT
AND ITS APPLICATION
The DFRFT based on T (or S+kT) is:
° ° ° ¯ °° ° ®
¦
¦
even, for , odd for , 2 0 2 2 1 0 2 N e e N e N r T N N Na j T r r ra j N r T r r ra j T a a u u u u u u U UD FT S S S (15)where U [u0|u1|"|uN2|uN1] for odd N,
] | | | | [u0 u1 uN 2 uN
U " for even N, and ur is the rth-order
discrete Hermite-Gaussian like function with r zero-crossings and is computed from the corresponding normalized eigenvector of T (or S+kT). The performances of the DFRFTs based on S and T (or S+kT) are compared in the following experiment.
Computer experiment 3: We compute the continuous FRT,
and the DFRFTs based on S, T, and S+15T of the following rectangular function
x(t) = 1 when |t|d 17/16, x(t) = 0 elsewhere. (16) The continuous FRT is computed by numerical integration of the definition of FRT in (1). The DFRFTs based on S, T, and S+15T for the samples of x(t) in (16) are computed with sample number
N=64 and sampling interval 1/8. The transform results of x(t) are
plotted in Fig. 4 with transform order a = 0.25. We find that the transform results of the DFRFTs based on T and S+15T are more similar to those of the continuous FRT. Their root mean square errors (RMSE) are obviously less than that of the DFRFT based on S.
(a) Continuous FRT (a=0.25) (b) DFRFT based on S (RMSE=0.0913) -4 -2 0 2 4 -1 0 1 2 -4 -2 0 2 4 -1 0 1 2 (c) DFRFT based on T (d) DFRFT based on S+15T (RMSE=0.0647) (RMSE=0.0526) -4 -2 0 2 4 -1 0 1 2 -4 -2 0 2 4 -1 0 1 2
Fig. 4. Comparing the real parts (solid lines) and the imaginary parts (dashes) of the transform results of the continuous FRT and the DFRFTs based on S, T, and S+15T for a rectangular function (transform order a = 0.25). -200 -100 0 100 200 -2 0 2 4 -200 -100 0 100 200 -2 0 2 4 (a) x[n] (b) y[n]=x[n-56] -2000 -100 0 100 200 1 2 3 4 -2000 -100 0 100 200 1 2 3 4 (c) (d) |X [n]| a |Y [n]| a DFRFT a = 0.3
Fig. 5. If y[n] = x[n-W], after doing the DFRFT based on S+15T, the amplitudes are the same and the distance is reduced to Wcos(aS/2).
We then give an example that uses the DFRFT based on T or S+kT for space-variant pattern recognition. It is known that the continuous FRT has the following property [4]:
W
D
W
D SW D SW cos sin 2 2 2 sin 2 o F u e e u G t f t g a u j j a (17) whereD = aS/2, Fa(u) and Ga(u) are the continuous FRTs of f(t) and g(t), respectively. In other words, if g(t) is the same as f(t) except for the locations, then after doing the FRT, their ampli-tudes are also the same, and the difference of the locations is reduced by multiplying cosD:|Ga(u)| = |Fa(uWcosD)|. (18) Since the DFRFT based on T or S+kT are very similar to the continuous FRT, the properties in (17) and (18) also apply for it with some modification:
> @
n f>
n k@
g> @
D>
D
@
S D S cos sin 2 2 2 sin 2 k R m F e e m G N jNkm a k j a | o , (19) |Ga[m]| | |Fa[mR(kcosD)]|D = aS/2, R( ): rounding operation, (20) where Fa[m] and Ga[m] are the DFRFTs based on T or S+kT for
f[n] and g[n]. It can be shown from the experiments in Fig. 5. In
Figs. 5(a)-(b), x[n] is a signal generated by random variables and
y[n] is a shifting version of x[n]. We do the DFRFT based on
S+15T of order a = 0.3 for x[n] and y[n] and show the
ampli-tudes of the results in Figs. 5(c)-(d). Then we find that |Ya[m]| is very similar to the shifting of |Xa[m]|. From (20), the distance between |Xa[m]| and |Ya[m]| should be
R(56cos(0.3S/2)) = R(49.8964) = 50. (21) From Figs. 5(c)-(d), it can be found that |Ya[m]| is indeed very near to |Xa[m50]|. In fact, their correlation is near to 100%:
% 72 . 99 | ] [ | | ] 50 [ || ] [ |
/
2¦
¦
m a m a a m Y m X m X . (22)Thus we can use the DFRFT based on T or S+kT to do space-variant pattern recognition. In the transform domain, we can use whether there exists an h such that |Ga[m]| = |Fa[mh]| to con-clude whether the two patterns g[n] and f[n] are equivalent. If so, we can use k| h/cosD to estimate the distance between the two patterns in the space domain.
5. REFERENCES
[1] F. A. GrΏnbaum, “The eigenvectors of the Discrete Fourier Transform: a version of the Hermite functions,” J. Math.
Anal. Appl., vol. 88, no. 2, pp.355-363, 1982.
[2] S. C. Pei and M. H. Yeh, “Improved discrete fractional Fou-rier transform,” Opt. Lett., vol. 22, pp. 1047-1049, 1997. [3] C. Candan, M. A. Kutay, and H. M. Ozaktas, “The discrete
fractional Fourier transform,” IEEE Trans. Signal
Process-ing, vol. 48, pp. 1329-1337, May 2000.
[4] H. M. Ozaktas, Z. Zalevsky, and M. A. Kutay, The
Frac-tional Fourier Transform with Applications in Optics and Signal Processing. New York, John Wiley & Sons, 2000.
[5] B. W. Dickinson and K. Steiglitz, “Eigenvectors and func-tions of the discrete Fourier transform,” IEEE Trans. Acoust., Speech, Signal Processing, vol. ASSP-30, pp 25-31, Jan. 1982.
[6] J. H. Wilkinson, The Algebraic Eigenvalue Problem. Oxford, U.K.: Oxford, 1988.
