This page explains how the inverse discrete cosine transform(IDCT) is implemented. It follows the proposals of Arai, Agiu and Nakajima* as well as of Tseng and Miller** completed by the extensions made by Feig***.
*Y.Arai,T.Agui,M.Nakajima, A Fast DCT-SQ Scheme for Images,
Trans. of the IEICE.E 71(11):1095(Nov.1988)
**B.D.Tseng,W.C.Miller, On Computing the Discrete Cosine Transform, IEEE Trans. on Computers, C-27(10):966-8(Oct. 1978)
***E.Feig,E.Linzer, Discrete Cosine Transform Algorithms for Image Data Compression, Proceedings Electronic Imaging '90 East, pp.84-7, Boston, MA(Oct.29 - Nov.1,1990)
also published in:W.B.Pennebaker,J.L.Mitchell,JPEG Still Image Compression Standard, Van Nostrand Reinhold, 1993
The starting point is the formula for the 1 dimensional DCT(1D-DCT):
8 means that the DCT was built by 8 coefficients.
First the 1D-DCT is applied to all rows of the DCT matrix:
Then the 1D-DCT is applied to the columns of the the result of this operation:
By substitution we get the formula of the 2D-DCT:
which proofs: The 2D-DCT can be regarded as a double application of the 1D-DCT. Therefore we first regard the 1D-DCT.
By definition of:
we can rewrite the formula:
For some reason we compute :
Therefore the following is valid:
If we constitute a sequence of elements
f(k), k = 0...15 with:
we can rewrite the formula:
Because the 16 point discrete Fourier transform is defined by
the following relation is valid:
That is the basic idea: Instead of performing an inverse DCT an inverse DFT is performed. To obtain the DFT values we first multiply the DCT values with . It seems as this were very inefficient. But bear in mind that the last operation before the IDCT is the quantization. That means every value is to be multiplied with a certain constant depending on it's position in DCT-matrix. To obtain the DFT values we simply multiply the quantization constants with . As a result the quantization and DCT-DFT transformation is performed in one step. The IDFT actually computes 16 values. But we can simply ignore the higher 8 values (They are equal to the first 8 values but in reverse order).
To obtain the IDFT equations we first bring the DFT equations in matrix form. If we define a line of the original values by:
a line of the scaled transform values
and a function:
we can establish a Matrix
which transforms the values according to equation
cos computes only
0, 1, -1 an the positive and
negative values of the following 3 constants:
Therefore the matrix can be rewritten:
To get the matrix equation for IDFT the matrix
The inversion gives:
L can be factored:
So the IDFT can be performed by 5 steps:
M a simplification can be applied:
If you count the operations you'll get 29 additions and 5 multiplications.
The ordinary way is to exploit the horizontal-vertical properties of the DCT (equations
First the matrix
L is applied to all rows. Then the result is rotated by 90° to
the right and the
L is applied once again. Because the 1D-IDFT is applied 16 times
this results in
29 * 16 = 464 additions and
16 * 5 = 80
multiplications. The known
libjpeg.a uses the horizontal-vertical properties of the DCT wereby some
special cases are regarded.
But Feig*** describes an 2D-algorithm with less operations.
The original matrix is re-packed into a one-line-matrix
concatenating all lines:
Then a matrix
KK must be found which transforms
FF in a
XX containing all retransformed lines in one line:
KK is the tensor or Kronecker product of
To explain what it means regard the matrices
The tensor or Kronecker product is:
By means of the following substitution:
KK can be represented as:
Because of some distributive properties of the tensor product this is the same as:
Now we can decide how to compute the 3 tensor products: according to equation
by exploiting the horizontal-vertical-properties of the tensor product (equations
According to Feig*** the current application uses the horizontal-vertical-property for the 1st and last tensor
product. This requires
16 * 26 = 416 additions.
Note that the rotation by 90° isn't a real operation! It is simply a permutation of the indexes.
The middle tensor product is computed according to equation
Note that every
0 represents a 8x8 zero matrix!
The lines 1,2,4 and 8 are multiplications with the original Matrix
We already know this can be performed by 3 additions and 5 multiplications per line:
The lines 3 and 6 are treated the same way extended by a multiplication with
This requires 7 multiplications, 3 additions and 2 shifts to the left.
To transform the values of lines 5 and 7, we build a matrix
containing the lines 5 and 7 without the elements 5 and 7:
M by discarding both, lines 5 and
7 as well as columns 5 and 7:
N contains only the values which are both, in lines and columns
5 and 7 of matrix
The resulting matrix
can be computed by:
This gives equations like:
They can be computed similarly to equations
(6) by 3 additions
and 3 multiplications.
To compute the 5th and the 7the value of lines 5 and 7 the values
are to be multiplied with the tensor product of
If we compute the matrix we get:
But because of:
this results in:
which can be factored:
Using the relation: we get:
This implies 2 multiplications 10 additions and 2 shifts to the left.
If you count the operations for the tensor product of
M you'll get:
|shifts to the left|
|shifts to the left|
|5 and 7||6*3+10||6*3+2||2||28||20||2|
If we add the 416 additions for the first and the last tensor product we get: 462 additions, 54 multiplications and 6 shifts to the left.