# Fast DCT implementation

I'm having trouble figuring out how to follow the fast 8x8 DCT algorithm diagrams found in the following two papers:

(1) A Fast Computational Algorithm for the Discrete Cosine Transform by Chen et al.

and

(2) Practical Fast 1-D DCT Algorithms With 11 Multiplications by Loeffler et al.

In particular, the second diagram that shows the algorithm in (2) looks like the following:

The description of the operations in this algorithm are:

There are a few questions I have about this formulation, and I'm not sure where to find the answers:

1. (2) suggests that this algorithm generates a DCT that is scaled by some value $C = \sqrt{2}$. It mentions that this $C$ was chosen arbitrarily to avoid any multiplications in computing the DC coefficient. Really the only requirement is that $C_{DCT} * C_{IDCT} = \frac{4}{N^2}$. So my question is this: What is the scaling factor of the output coefficients using this algorithm? It seems like they are different than the original definition of the DCT, but I don't know by how much (mostly because I don't actually see any relationship between this diagram and the original formulation of the DCT): $$F(k) = \frac{2c(k)}{N}\sum_{n = 0}^{N - 1}f(n)\cos\left(\frac{\left(2n + 1\right)\pi k}{2N}\right)$$ where $c(k) = \frac{1}{\sqrt{2}}$ for $k = 0$ and $c(k) = 1$ for $k \neq 0$.

2. The paper states that performing the IDCT can be done using the exact same algorithm but transforming outputs to inputs and vice versa. First, should the DCT coefficients be ordered in bit-reverse order prior to running them through the IDCT? Second, for the rotation blocks (the squares in the diagram), shouldn't the inverse operation be: \begin{align} O_0 = I_0 \cdot k \cdot \cos\frac{n\pi}{2N} - I_1 \cdot k \cdot \sin\frac{n\pi}{2N} \\ O_1 = I_1 \cdot k \cdot \sin\frac{n\pi}{2N} + I_1 \cdot k \cdot \cos\frac{n\pi}{2N} \\ \end{align} My reasoning is this: The inverse of a rotation by $\theta$ is a rotation by $-\theta$. Hence, we just replace the angle by its inverse and use the identities $\cos(-\theta) = \cos(\theta)$ and $\sin(-\theta) = -\sin(\theta)$. Third, what is the scaling factor of the transformed values after the IDCT? (2) says $\frac{2}{N^2}$, but empirically, this hasn't produced correct results.

3. Suppose after I run the algorithm, I have the result of each lane stored in the values d0 ... d7. Which of the following is correct:

output[0] = d0               or               output[0] = d0
output[4] = d1                                output[1] = d4
output[2] = d2                                output[2] = d2
output[6] = d3                                output[3] = d6
output[7] = d4                                output[4] = d7
output[3] = d5                                output[5] = d3
output[5] = d6                                output[6] = d5
output[1] = d7                                output[7] = d1 

If there are any ways to improve this question, or if I should ask elsewhere, please let me know.

• To practically answer these kind of questions what you actually need is a set of precomputed DCT values and adjusting your implementation until its results conforms to those precomputed ones.... – Fat32 Jan 11 '16 at 13:50
• I have all those questions and more... did you ever figure it out? I found some C implementation that I try to extract stuff form. I'll write something if I find answers. – Pepijn Jun 30 at 13:06

## 1 Answer

Alright, after some days of staring at this problem I hope I can provide a bit of guidance to the next poor soul.

1. Yes the scaling is different. Compared to scipy.fftpack.dct the DC term is $\frac{1}{2}$ and the other terms $\frac{\sqrt{2}}{2}$. But apparently it all nicely cancels out in the inverse transform.
2. The inverse input order is exactly as they come out: bit reversed. Literally as if you flip around the figure and connect the lines. And yes, you are correct that the sin gets a minus. I'm seeing a scaling factor of 8, FWIW.
3. The output order of the inverse is the same as the input order of the forward transform. So $x[n]=\frac{\mathrm{IDCT}(\mathrm{DCT}(x[n]))}{8}$

Also note that there is an error in the graph and it's $\sqrt{2}c6$ for the even-side rotation block.

• Thanks for taking the time to answer these questions! My original motivation for this was in the context of a 2D compression algorithm, so I'm still a bit unclear about the relative ordering of the outputs (I'd like them to be from 0-7 so that I can have them from smallest to largest). They aren't quite bit reversed either: 3 -> 5, 5 -> 6, 7 -> 1 aren't exactly bit-reversals (or did I misunderstand here). – Mokosha Jul 1 at 15:16
• The order matches en.wikipedia.org/wiki/Bit-reversal_permutation so there is that. You can of course reorder them as you like. I would suggest directly taking this into account in the zig-zag step to avoid extra costs. – Pepijn Jul 3 at 19:26