C. H. Wei and N. C. Chen, “Sigma-delta modulation adaptive digital filter,” in Proc. ISCAS (Espoo, Finland), June 1988, pp. 523-526.
P. W. Wong and R. M. Gray, “FIR filters with sigma-delta modu- lation encoding,” IEEE Trans. Acoust., Speech, Signal Processing,
vol. 38, pp. 979-990, June 1990.
T. Saramaki, T. Karema, T. Ritoniemi, and H. Tenhunen, “Multi- plier-free decimator algorithms for superresolution oversampled con- verters,” in Proc. ISCAS (New Orleans, LA), May 1990, pp. 3275- 3278.
R. Roncella and R. Saletti, “A VLSI systolic adder for digital filter- ing of delta-modulated signals,” IEEE Trans. Acoust., Speech, Sig- nal Processing, vol. 37, pp. 749-754, May 1989.
G. Picchi, G. Prati, and R. Raheli, “Incremental adaptive equalizers using ternary
+
1, 0, - 1 tap coefficients and/or delta modulation in- put encoding,” Ann. Telecommun., vol. 41, pp. 260-268, May-June1986.
R. Raheli and L. E. Franks, “Seimsystolic incremental realization of FIR digital filters using ternary arithmetic,” IEEE Trans. Acoust., Speech, Signal Processing, vol. 37, pp. 1231-1240, Aug. 1989.
H. h o s e and Y . Yasuda, “A unity bit coding method by negative feedback,” Proc. IEEE, vol. 51, pp. 1524-1535, Nov. 1963.
J. C. Candy, Y. C. Ching, and D. S . Alexander, “Using triangularly weighted interpolation to get 13-bit PCM from a sigma-delta modu- lator,” IEEE Trans. Commun., vol. COM-24, pp. 1268-1275, Nov. 1976.
J. J. Paulos, G. T. Brauns, M . B. Steer, and S . H. Ardalan, “Im- proved signal-to-noise ratio using trilevel delta-sigma modulation,” in Proc. ISCAS (Philadelphia, PA), May 1987, pp. 463-466.
R. M . Gray, “Spectral analysis of quantization noise in a single-loop sigma-delta modulator with dc input,” IEEE Trans. Commun., vol. 37, pp. 588-599, June 1989.
P. W. Wong and R. M. Gray, “Two-stage sigma-delta modulation,”
IEEE Trans. Acoust., Speech, Signal Processing, vol. 38, pp. 1937-
1952, Nov. 1990.
W. Chou, P. W. Wong, and R. M. Gray, “Multistage sigma-delta modulation,” IEEE Trans. Inform. Theory, vol. 35, pp. 784-796,
July 1989.
N. He, F. Kuhlmann, and A. Buzo, “Double-loop sigma-delta mod- ulation with dc input,” IEEE Trans. Commun., vol. 38, pp. 487-495, Apr. 1990.
R. M. Gray, “Quantization noise spectra,” IEEE Trans. Inform. Theory, vol. 36, pp. 1220-1244, Nov. 1990.
R. M. Gray, W. Chou, and P. W. Wong, “Quantization noise in single-loop sigma-delta modulation with sinusoidal inputs,” IEEE Trans. Commun., vol. 37, pp. 956-968, Sept. 1989.
P. W. Wong and R. M. Gray, “Sigma-delta modulation with i.i.d. Gaussian inputs,” IEEE Trans. Inform. Theory, vol. 36, pp. 784-
798, July 1990.
W. R. Bennett, “Spectra of quantized signals,” Bell Syst. Tech. J . , vol. 27, pp. 446-472, 1948.
IEEE TRANSACTIONS ON SIGNAL PROCESSING, VOL. 40, NO. 6 , J U N E 1992 I. INTRODUCTION
Since the introduction of the discrete cosine transform (DCT) [ l ] , it has found a number of applications in image and speech processing. It has been shown that DCT performs very close to the statistically optimal Karhunen-Lobve transform (KLT) for a large number of signal classes [ 11, [2]. Since it is used mostly in various signal processing applications, the DCT processor designed with real-time computation capabilities is required urgently. Some DCT architectures were proposed to meet this requirement, such as dis- tributed arithmetic DCT [3], the discrete Fourier-cosine transform chip [4], and the concurrent CORDIC DCT [ 5 ] .
Recently, Ersoy [6] proposed a very efficient two-stage FFT al- gorithm which is based on the linear vector transformation of sines and cosines into new basis functions using the Mobius inversion of number theory. Such an approach has been proved valid also for the computation of discrete Hartley transform (DHT) [7].
In this correspondence, based on the idea of [6], a novel two- stage algorithm for DCT and IDCT is derived. In this approach, only additions/subtractions are required for the preprocessing stage, while the postprocessing stage consists of circular convolution or correlation operations which can be efficiently implemented by using parallel structures such as semisystolic arrays [8] or the num- ber theoretic transforms (NTT’s) technique [ 9 ] .
11. TWO-STAGE ALGORITHM FOR IDCT
The DCT and IDCT of an N-point real sequence n(k) are, re- spectively, defined as (where N is assumed to be powers of two):
f o r 0 5 k , n 5 N - 1, where
e(n) =
[:/a,
i f n = Ootherwise. (3)
Since (1) can be realized simply by transposing the flow graph for (2), and since e(n) is a normalization constant, it is sufficient for our derivation to consider the IDCT-like equation
A Novel Two-Stage Algorithm for DCT and IDCT
N - Ix(k) =
X(n)
cosn = O
Ja-Ling Wu, Shyh-Huei Hsu, and Wei-Jou Duh
Abstract-In this correspondence, on the basis of Mobius function, a two-stage algorithm for the discrete cosine transform (DCT) and the IDCT is proposed. In this approach, the DCT matrix is factorized into the preprocessing and postprocessing matrices. The preprocessing ma- a circular convolutionlcorrelation matrix.
NOW let
xl(k) = x(2k)
trix has elements of values 1 and - 1, and the postprocessing matrix is x , ( k
+
N/2) = x(N - (2k+
l))k = 0, 1,
. .
’ , N / 2Manuscript received July 22, 1991; revised January 2, 1992. The authors are with the Department of Computer Science and Infor-
then (4) can be rewritten as mation Engineering, National Taiwan University, Taipei, Taiwan, Repub-
lic of China.
IEEE Log Number 9107634.
N - I x , ( k ) = n = O X(n) cos . , N - l . (4) 1 1053-587X/92$03.00 0 1992 IEEE
I
IEEE TRANSACTIONS ON SIGNAL PROCESSING, VOL. 40, NO. 6, JUNE 1992
N - I ( r ( 2 N
+2F
+
1 ) n x l ( k+
N / 2 ) =X(n)
cos n = O (7) k = 0, 1,. .
* , N / 2 - 1. With k' = k+
N / 2 , (7) reduces to N - I X l ( k ' ) =c
X(n) cos " = OFrom (6) and (8), it follows that N - I
x l ( k ) = X(n) cos
t(4k2
k = 0, 1,. . .
, N - 1 n = O(9)
and x ( k ) can be obtained from x I ( k ) by performing an appropriate
permutation. Now let us consider the following new sequence: N - I
~ 2 ( k ) =
c
X(n)
COS ( & n / 2 N ) , k =0,
1,. .
. , 4N - 1.n = O
(10)
It is clear that the relationship between xI ( k ) and x2 ( k ) can be
expressed as
~ l ( k ) = ~ 2 ( 4 k
+
I), k = 0 , 1 ,. . .
, N - 1. (11) In other words, x, (k) or x ( k ) can be obtained from x2 ( k ) . So, let us pay attention to the computation of (10).From [ 6 ] , on the basis of Mobius inversion formula, the cosine function can be represented as
M
cos ( 2 ? m k / N ) = b ( m 1 , N ) u ( m , n k / N
+
1 / 4 ) ( 1 2 )where m l equals 1, 5 , 9,
. . .
, M
and M = N - 3 . Substituting ( 1 2 ) into (IO), we obtainm l = I M x * ( k ) = b ( m l , 4 N ) h ( m l k modulo 4 N ) ( 1 3 ) m l = I where N - I h ( f ) = X ( n ) u ( l n / 4 N
+
l / 4 ) , k = 0, 1,. . .
, 4 N - 1 . n = O ( 1 4 )Expression ( 1 4 ) is the preprocessing and ( 1 3 ) is the postprocess-
ing equations of IDCT, respectively. Indices k and m l can be per- muted a s follows:
k + a" modulo 4N
ml + a'"'' modulo 4 N (15)
where a = - 3 or ( 4 N - 3 ) (since N is a power of two) and s =
k
1. Then, ( 1 3 ) becomesN - I
.. .
x2(ak' modulo 4 N ) = b(a, 4 N ) h ( ~ ~ ' + ~ ' " ' modulo 4 N ) ,
m i = O
k' = 0 , 1, 2 ,
. .
. , N - 1 . ( 1 6 )From (16), it follows that the postprocessing stage of IDCT can be obtained from the N-point circular convolution (for s = - 1) or correlation ( f o r s = 1) of the output sequence of the preprocessing stage and the sequence b ( m I , 4 N ) .
1611 1 1 1 1 1 1 1 1 1 1 1 -1 -1 -1 -1 -1 1 -1 -1 1 1 -1 -1 1 1 1 - 1 - 1 - 1 1 1 1 1 -1 1 -1 1 -1 1 -1 1 - 1 1 1 . 1 1 - 1 1 1 1 -1 -1 1 1 -1 -1 1 -1 -1 1 -1 -1 1 -1 0.8297 0.2105 -0.0400 -0.1901 -0.0763 -0.0020 0.1401 0.1280 0.1280 0.2105 0.8297 -0.0400 -0.1901 -0.0763 -0.0020 0.1401 0.1401 0.1280 0.2105 0.8297 -0.0400 -0.1901 -0.0763 -0.0020 -0.0020 0.1401 0.1280 0.2105 0.8297 -0.0400 -0.1901 -0.0763 -0.0763 -0.0020 0.1401 0.1280 0.2105 0.8297 -0.0400 -0.1901 -0.1901 -0.0763 -0,0020 0.1401 0.1280 0.2105 0.8297 -0.0400 -0.0400 -0.1901 -0.0763 -0.0020 0.1401 0.1280 0.2105 0.8297 0.8297 -0.0400 -0.1901 -0.0763 -0.0020 0.1401 0.1280 0.2105
Fig. 1. The preprocessing and the postprocessing of an 8-point IDCT. The mapping used in ( 5 ) is not new [ 101 and the derivation from
( 1 2 ) to (16) is similar to that in [ 6 ] ; however, the Mobius function
cannot be applied to the computation of DCT without integrating them together. Therefore, the proposed two-stage Mobius function based DCTlIDCT algorithm is believed to be new.
111. DISCUSSIONS
For simplicity, let [DCT ( N ) ] and [IDCT ( N ) ] denote, without
including e(n) and 2 / N , the transform matrices of N-point DCT and IDCT, respectively. From previous derivations, we have
[IDCT ( N ) ] = [POST ( N ) ] [PRE ( N ) ] (17) where [PRE ( N ) ] and [POST ( N ) ] denote the preprocessing and
postprocessing matrices of IDCT, respectively. Fig. 1 depicts the values of [PRE ( 8 ) ] and [POST (8)]. It follows that
[DCT ( N ) ] = [IDCT ( N ) ] ' = [PRE ( N ) ] ' [ P O S T ( N ) ] ' (18) where the superscript t denotes the operation of transposition. Note that [POST ( N ) ] ' also be a circular convolutionlcorrelation matrix.
On the basis of (17) and ( 1 8 ) a convolver-based DCTlIDCT archi-
tecture can be constructed as shown in Fig. 2 , where the convolver
can be implemented by distributed arithmetics technique, systolic array, or NTT processors.
As pointed out by one of the referees, some convolution-based DCT algorithms have been recently derived [ 1 I ] , [ 121. The length variation of the algorithm derived in [ 1 I ] introduces different mod- ules and irregular structures. After appropriate precomputation, the algorithm in [I21 is more efficient than the proposed one. However,
the preprocessing stage of the proposed algorithm (although some extra additions are required) is itself a meaningful transformation, which has been successfully applied to extract the features of sim- ilar patterns [ 131.
1612 Time Reverse
rn
tramform mic
D c r transform Domain t t I I+
p w t+
+
: Signal flowI
b(m ,4N)1
+
: Der signal flowFig. 2. The convolver-based architecture for DCT/IDCT POST
convolver)
+:aconstantcoef.
+-
PE+
.
REFERENCES
[l] N. Ahmed, T. Natarajan, and K. R. Rao, “Discrete cosine trans- form,” IEEE Trans. Cornput., vol. C-25, pp. 90-93, Jan. 1974. [2] M. Hamidi and J. Pearl, “Comparison of the cosine and Fourier
transforms of Markov-1 signals,” IEEE Trans. Acoust., Speech, Sig- nal Processing, vol. ASSP-24, pp. 428-429, Oct. 1976.
I
IEEE TRANSACTIONS ON SIGNAL PROCESSING, VOL. 40, NO. 6, JUNE 1992 M. T. Sun, L. Wu, and M. L. Liou, “A concurrent architecture for VLSI implementation of discrete cosine transform,” IEEE Trans. Circuirs S y s t . , vol. CAS-34, pp. 992-994, Aug. 1987.
M. Vetterli and A. Ligtenberg, “A discrete Fourier-cosine transform chip,” IEEE J. Select. Area Commun., vol. SAC-4, pp. 49-61, Jan.
1986.
J.-L. Wu and W.-J. Duh, “A novel concurrent architecture to imple- ment discrete cosine transform based on index partitations,” Inr. J.
Electron., vol. 68, no. 2, pp. 165-174, 1990.
0. K. Ersoy, “Two stage representation of DFT and its applica- tions,” IEEE Trans. Acoust., Speech, Signal Processing, vol. ASSP-
35, pp. 825-831, June 1987.
J.-L. Wu and C.-Y. Hsu, “Comment on ‘A two stage representation of DFT and its applications’,’’ IEEE Trans. Acoust., Speech, Signal
Processing, vol. ASSP-36, pp. 1687-1688, Oct. 1988.
0. K. Ersoy, “Semisystolic array implementation of circular, skew- circular, and linear convolutions,” IEEE Trans. Cornput., vol. C-34, pp. 190-194, Feb. 1985.
J. H. McClellan and C. M. Rader, Number Theory in Digital Signal
Processing.
M. J. Narasimha and A. M. Peterson, “On the computation of the DCT,” IEEE Trans. Comrnun., vol. COM-26, no. 6, pp. 934-936,
June 1978.
W. Li, “A new algorithm to compute the DCT and its inverse,” IEEE Trans. Signalfrocessing, vol. 39, no. 6, pp. 1305-1313, June 1991.
P. Duhamel and H. H’Mida, “New 2” DCT algorithms suitable for VLSI implementation,” in Proc. ICASSP, 1987, pp. 1805-1809.
J.-L. Wu and W.-J. Duh, “Feature extraction capability of some dis- crete transforms,” in Proc. IEEE Inr. Symp. Circuits Syst. (ISCAS)
(Singapore), June 1991, pp. 2649-2652.
Englewood Cliffs, NJ: Prentice-Hall, 1979.