On matrix factorizations for recursive pruned discrete
cosine transforms
Kuo-Liang Chung
!,*, Wen-Ming Yan1,"
! Department of Information Management and Graduate Program of Information Engineering, National Taiwan University of Science and Technology, No. 43, Section 4, Keelung Road, Taipei, Taiwan 10672, Taiwan, ROC
" Department of Computer Science and Information Engineering, National Taiwan University, Taipei, Taiwan 10764, Taiwan, ROC Received 27 August 1997; received in revised form 6 February 1998
Abstract
This paper presents two matrix factorizations for recursive pruned discrete cosine transforms. Both factorizations have lower computational complexity when compared to the method of El-Sharkawy and Eshmawy (1995). ( 1998 Elsevier Science B.V. All rights reserved.
Zusammenfassung
Dieser Artikel behandelt zwei Methoden zur Matrizenfaktorisierung fu¨r rekursive reduzierte diskrete Cosinustransfor-mationen. Beide Methoden haben einen geringeren Rechenaufwand als die Methode von El-Sharkawy und Eshmawy (1995). ( 1998 Elsevier Science B.V. All rights reserved.
Re´sume´
Cet article pre´sente deux me´thodes de factorisation de matrices pour les transformations discre`tes en cosinus e´lague´es re´cursives. Ces me´thodes pre´sentent toutes deux une complexite´ de calcul infe´rieure a` celle de la me´thode de El-Sharkawy et Eshmawy (1995). ( 1998 Elsevier Science B.V. All rights reserved.
Keywords: Computational complexity; Fast algorithm; Matrix factorization; Recursive pruned discrete cosine transform
1. Introduction
The recursive discrete cosine transforms (DCTs) [1,5] have many important applications in the area of image processing, communication, and digital signal processing. These applications include image
* Corresponding author. E-mail: [email protected]. Prof. Chung was supported in part by the National Science Council of R.O.C. under contracts NSC87-2213-E011-001 and NSC87-2213-E011-003.
1 Prof. Yan was supported in part by the National Science Council of R.O.C. under contract NSC87-2119-M002-006. 0165-1684/98/$19.00 ( 1998 Elsevier Science B.V. All rights reserved.
compression, filtering, feature extraction, teleconferencing, video phones, progressive image transmission, and so on. Due to the energy-packing capabilities, it approaches the statistically optimum Karhunen—Loeve transform for first-order Markov stationary random data [2]. In some applications, e.g. image compression, only the low-frequency DCT components should be kept, thus we only utilize a subset of input points or a subset of output points in order to decrease the computational complexity. This method is referred to as
pruning [6,3].
Recently, El-Sharkawy and Eshmawy [3] presented an efficient matrix factorization for one-dimensional (1-D) recursive pruned DCT (RPDCT) which improves the result of Wang [6] significantly. Later, they [4] extended the proposed method from the 1-D domain to the 2-D domain successfully. Based on the result of [3], the motivation of this research is to present new matrix factorizations for RPDCT in order to gain some computational advantages.
Employing row operations and some addition rules for cosine functions, this paper presents two matrix factorizations for RPDCT. The proposed two factorizations are quite competitive with the result [3], but eliminate all the shift operations in one factorization and all but one in the other factorization.
The remainder of this paper is organized as follows. Section 2 introduces the work of El-Sharkawy and Eshmawy [3]. The proposed two matrix factorizations for RPDCT are presented in Section 3. Some concluding remarks are addressed in Section 4.
2. The work of El-Sharkawy and Eshmawy
For describing the proposed two matrix factorizations more clearly, in this section, the work of El-Sharkawy and Eshmawy [3] is reviewed in detail.
For 0)i)n!1, the Hadamard recurrence is given by h2n(2i)"hn(i) and h2n(2i#1)"2n!1!hn(i) with the initial condition h1(0)"0, where the constraint n"2k must be satisfied. From [6], the pruned DCT matrix is given by CIIn"
C
J(1/n)t0,0 J(1/n)t0,1 2 J(1/n)t0,n~1 J(2/n)t1,0 J(2/n)t1,1 2 J(2/n)t1,n~1 F F 2 F J(2/n)tn~1,0 J(2/n)tn~1,1 2 J(2/n)tn~1,n~1D
,where ti,j"cos(2hn( j)#1)ip/2n"cos2i/ j with /j"(2hn( j)#1)p/4n. For the pruned DCT, we want to compute CIInx, where x is the related input sequence with length n. Let
¹(n)"
C
t0,0 t0,1 2 t0,n~1 t1,0 t1,1 2 t1,n~1 F F 2 F tn~1,0 tn~1,1 2 tn~1,n~1D
,then the following two steps can be used to compute CIInx:
Step 1. Compute y"¹(n)x.
Step 2. Compute diag[J1/n, J2/n, J2/n,2,kJ2/n] y.
The above two steps need O(n2) multiplications and additions since they mainly concern a matrix—vector multiplication.
In [3], El-Sharkawy and Eshmawy presented a faster method to reduce the above time bound to O(n log n) multiplications and additions.
Let P(n) be an n]n permutation matrix formed by reordering the rows of the identity matrix so that the following equality holds: P(n)(0,1,2,2, n!2,n!1)5"(0,2,2,n!2,1,3,2,n!1)5. The previous matrix ¹(n) is reordered by exchanging rows and columns of the matrix, then the resultant matrix is decomposed into four quadrants [3] as shown below:
P(n)¹(n)P(n)T"
C
¹ee(n2) ¹eo(n2) ¹oe(n2) ¹oo(n2)D
, (1) where ¹eeA
n2B
"C
t0,0 t0,2 2 t0,n~2 t2,0 t2,2 2 t2,n~2 F F 2 F tn~2,0 tn~2,2 2 tn~2,n~2D
, ¹eoA
n 2B
"C
t0,1 t0,3 2 t0,n~1 t2,1 t2,3 2 t2,n~1 F F 2 F tn~2,1 tn~2,3 2 tn~2,n~1D
, ¹oeA
n2B
"C
t1,0 t1,2 2 t1,n~2 t3,0 t3,2 2 t3,n~2 F F 2 F tn~1,0 tn~1,2 2 tn~1,n~2D
and ¹ooA
n 2B
"C
t1,1 t1,3 2 t1,n~1 t3,1 t3,3 2 t3,n~1 F F 2 F tn~1,1 tn~1,3 2 tn~1,n~1D
.From the definitions of ti,j and the Hadamard recurrence, we have
ti,2j"cos(2hn(2j)#1)i2np"cos(2hn@2(j)#1)i p
2n"cos i/@j, where/@j"(2hn@2(j)#1)p/2n. Similarly, we have
ti,2j`1"cos(2hn(2j#1)#1)i2np"cos(2(n!1!hn@2(j))#1)i p 2n
"cos ip!(2hn@2(j))#1)ip
2n"(!1)*cos i/@j. Then it follows that
t2i,2j"t2i,2j`1"cos 2i/@j
and
t2i`1,2j"!t2i`1,2j`1"cos(2i#1)/@j.
It is not hard to verify that ¹ee(n/2)"¹eo(n/2)"¹(n/2) and ¹oo(n/2)"!¹oe(n/2). For analyzing the computational complexity, let m"n/2. We then have
¹ee(m)"
C
1 1 2 1
cos 2/@0 cos 2/@1 2 cos 2/@m~1
F F 2 F
and
¹oe(m)"
C
cos/@0 cos/@1 2 cos/@m~1 cos 3/@0 cos 3/@1 2 cos 3/@m~1
F F 2 F
cos(2m!1)/@0 cos(2m!1)/@1 2 cos(2m!1)/@m~1
D
.The matrix ¹oe(m) has the following matrix factorization: ¹oe(m)"A(m)¹ee(m)B(m), where A(m)"
C
1 0 0 0 2 0 !1 2 0 0 2 0 1 !2 2 0 2 0 !1 2 !2 2 2 0 F F F F 2 F !1 2 !2 2 2 2D
mCm and B(m)"diag[cos/@0,cos/@1,cos/@2,2,cos/@m~1].From Eq. (1) and ¹oe(m)"A(m)¹ee(m)B(m), it is rather straightforward that ¹(2m)"P(2m)T
C
¹ee (m) ¹eo(m) ¹oe(m) ¹oo(m)D
P(2m)"P(2m)TC
¹(m) ¹(m) A(m)¹(m)B(m) !A(m)¹(m)B(m)D
P(2m).We now want to compute y"¹(2m)x, where x is the input sequence with length 2m. First, we compute
"
C
1 2D
"P(2m)x. Then, we computeC
z1 z2D
"C
¹(m) ¹(m) A(m)¹(m)B(m) !A(m)¹(m)B(m)DC
1 2D
"C
¹(m)(1#2) A(m)¹(m)B(m)(1!2)D
. Finally, we compute y"P(2m)TC
z1 z2D
.Equivalently, the procedure [3] consisting of the following six steps for computing y"¹(2m)x is shown below:
Step 1: Compute "P(2m)x.
Step 2: Compute w1"1#2 and w2"1!2. Step 3: Compute w@2"B(m)w2.
Step 4: Compute z1"¹(m)w1 and z@2"¹(m)w@2. Step 5: Compute z2"A(m)z@2.
Step 6: Compute y"P(2m)Tz@, where z@"
C
z1In the next section, two new matrix factorizations for ¹oe(m)"A(m)¹ee(m)B(m) are presented. For the first proposed factorization, the computational complexity required in Step 5 of the above procedure can be reduced from ((m!1) shift operations and (m!1) subtractions) to (one shift and (m!1) subtractions), while preserving the same bound in the remaining steps. For the second proposed factorization, the computational complexity required in Step 5 in the above procedure can be reduced to (m!1) additions, while still preserving the same bound in the remaining steps.
3. Two new matrix factorizations
Let ¹oe(m)"
C
r0 r1 F rm~1D
and ¹ee(m)"C
r@0 r@1 F r@m~1D
, whererk"(cos(2k#1)/@0, cos(2k#1)/@1,2, cos(2k#1)/@m~1), r@k"(cos2k/@0,cos2k/@1,2,cos2k/@m~1),
for k"0,1,2, m!1. From /@j"(2hm(j)#1)p/4m, it follows that r@m"(cos2m/@0,cos2m/@1,2, cos 2m/@m~1)"(0,0,2,0,0).
We add the (k!1)th row vector in ¹oe(m) to the kth row vector for 1)k)m!1, then we have
rk~1#rk"(cos(2k!1)/@0#cos(2k#1)/@0, cos(2k!1)/@1#cos(2k#1)/@1,2,
cos(2k!1)/@m~1#cos(2k#1)/@m~1). From the trigonometric addition identity:
cos(2k!1)/@j#cos(2k#1)/@j"2cos/@jcos2k/@j"ajcos2k/@j, it follows that
rk~1#rk"(a0 cos 2k/@0, a1 cos 2k/@1,2, am~1 cos 2k/@m~1),
whereaj"2cos/@j.
By the same argument, we add the kth row vector in ¹ee(m) to the (k#1)th row vector for 0)k)m!2. From the trigonometric addition identity:
cos 2k/@j#cos(2k#2)/@j"2cos/@jcos(2k#1)/@j"ajcos(2k#1)/@j, we have r@k#r@k`1"(a0cos(2k#1)/@0,a1cos(2k#1)/@1,2, am~1cos(2k#1)/@m~1). We then have
C
2r0 r0#r1 F rm~2#rm~1D
"C
a0 a1 2 am~1a0cos2/@0 a1cos2/@1 2 am~1cos2/@m~1
F F 2 F
"
C
1 1 2 1cos 2/@0 cos 2/@1 2 cos 2/@m~1
F F 2 F
cos(2m!2)/@0 cos(2m!2)/@1 2 cos(2m!2)/@m~1
DC
a0a1 2
am~1
D
"¹ee(m)D(m)and from r@m"(0,0,2,0,0), it yields
C
r@0#r@1 F r@m~2#r@m~1 r@m~1#r@mD
"C
r@0#r@1 F r@m~2#r@m~1 r@m~1D
"C
a0cos/@0 a1cos/@1 2 am~1cos/@m~1 a0cos3/@0 a1cos3/@1 2 am~1cos3/@m~1F F 2 F
a0cos(2m!1)/@0 a1cos(2m!1)/@1 2 am~1cos(2m!1)/@m~1
D
"C
cos/@0 cos/@1 2 cos/@m~1 cos 3/@0 cos 3/@1 2 cos 3/@m~1F F 2 F
cos(2m!1)/@0 cos(2m!1)/@1 2 cos(2m!1)/@m~1
DC
a0 a1 2 am~1D
"¹oe(m)D(m), where D(m)"diag[a0,a1,2,am~1]"diag[2cos/@0,2cos/@1,2,2cos/@m~1].In summary, we have the following two important equations:
C
2r0 r0#r1 F rm~2#rm~1D
"¹ee(m)D(m) (2) andC
r@0#r@1 F r@m~2#r@m~1 r@m~1D
"¹oe(m)D(m). (3)Let ¸(m)"
C
2 0 0 0 2 0 1 1 0 0 2 0 0 1 1 0 2 0 F F F F F F 0 2 0 1 1 0 0 2 0 0 1 1D
mCm and º(m)"C
1 1 0 0 2 0 0 1 1 0 2 0 F F F F F F 0 2 0 1 1 0 0 2 0 0 1 1 0 2 0 0 0 1D
mCm .Eqs. (2) and (3) can be rewritten as
¸ (m)¹oe(m)"¸(m)
C
r0 r1 F rm~1D
"C
2r0 r0#r1 F rm~2#rm~1D
"¹ee(m)D(m) and º (m)¹ee(m)"º(m)C
r@0 F r@m~2 r@m~1D
"C
r@0#r@1 F r@m~2#r@m~1 r@m~1D
"¹oe(m)D(m).As a result, two new factorizations for ¹oe(m) are derived as shown below:
¹oe(m)"¸~1(m)¹ee(m)D(m), (4)
¹oe(m)"º(m)¹ee(m)D~1(m). (5)
From Eq. (4), computing Step 5 is equal to solving the bidiagonal linear system ¸(m)z2"z@2, and it takes one shift operation and m!1 subtractions. From Eq. (5), it takes m!1 additions to compute Step 5.
Recall that the computational complexity required in Step 5 of the above procedure mentioned in Section 2 is (m!1) shift operations and (m!1) subtractions. The computational complexity required in Step 5 of any one of the proposed two methods is lower than the that of in the above procedure, while preserving the same bound in remaining steps. Table 1 illustrates the computational complexity required in Step 5 among the three methods.
In Table 1, if we assume that the computational load of one subtraction is equal to that of one addition, the second method is the fastest among the three methods.
Table 1
Computational complexity comparison for Step 5 Shift operations Subtrac-tions Additions [3] m!1 m!1 0 First method 1 m!1 0 Second method 0 0 m!1
4. Concluding remarks
We have presented two matrix factorizations to improve some steps in the procedure of El-Sharkawy and Eshmawy [3] for recursive pruned discrete cosine transforms. How to design efficient matrix factorizations to handle the case nO2k is our future research topic. Instead of using trigonometric identities to derive the proposed matrix factorizations, how to just change some multiplicators to derive the same results is another interesting research issue.
Acknowledgements
The authors are indebted to the two reviewers and Prof. Murat Kunt for making some valuable suggestions and corrections that lead to the improved version of the paper.
References
[1] N. Ahmed, K.R. Rao, Orthogonal Transforms in Digital Signal Processing, Springer, New York, 1975.
[2] R.J. Clack, Relation between the Karhunen—Loeve and cosine transform, Proc. IEE. Part F 128 (November 1981) 359—360. [3] M. El-Sharkawy, W. Eshmawy, A fast recursive pruned discrete cosine transform, Digital Signal Process. 5 (1995) 140—148. [4] M. El-Sharkawy, W. Eshmawy, A fast 8]8 pruned DCT algorithm, Digital Signal Process. 6 (1996) 145—154.
[5] K.R. Rao, P. Yip, Discrete Cosine Transform: Algorithms, Advantages, and Applications, Academic Press, New York, 1990. [6] Z. Wang, Pruning the fast discrete cosine transform, IEEE Trans. Commun. 39 (1991) 640—643.