The Danielson-Lanczos FFT method and bit-reversal

The Danielson-Lanczos FFT method and bit-reversal | Khilav Majmudar This note summarizes my understanding of the Danielson-Lanczos method of calculating the Fast Fourier Transform of an array of numbers. I read it in Numerical Recipes1. The method explained there becomes clear, as in most such cases, by expanding out the expressions for small values. The Danielson-Lanczos method is the prototype example of the general class of FFT algorithms. We have a discrete Fourier transform of length \(N\). \begin{equation}\label{main} F_k = \sum_{j=0}^{N-1} e^{2 \pi i j k/N} f_{j} \end{equation} For now, we will only consider cases in which \(N\) is a power of 2 as this is where the approach is most clearly seen. Let us split this according to its even and odd components. \[\begin{align} F_k &= \sum_{j=0}^{N-1} e^{2 \pi i j k/N} f_{j} \nonumber \\ &= \sum_{j=0}^{N/2 - 1} e^{2 \pi i (2j) k/N} f_{2j}+\sum_{j=0}^{N/2 - 1} e^{2 \pi i (2j+1) k/N} f_{2j+1} \nonumber \\ &= \sum_{j=0}^{N/2 - 1} e^{2 \pi ij k/(N/2)} f_{2j}+ W^{k} \sum_{j=0}^{N/2 - 1} e^{2 \pi ij k/(N/2)} f_{2j+1} \nonumber \\ &= F_{k}^{e} + W^{k}F_{k}^{o} \label{eq1} \end{align}\] Here, \(F_k\) is the \(k\)-th component of the Fourier-transformed array. The sum that leads to it has been split into even and odd transforms of length \(N/2\), given by \(F_{k}^{e}\) and \(F_{k}^{o}\) respectively. \(W\) is a constant factor and is given by \(e^{2\pi i/N}\). Now suppose that \(N=2\). Then the sums in the second-last line collapse to one element each and are given by \begin{equation}\label{eq2} F_k = f_0 + W^{k}f_1 \end{equation} This is our DFT. For completeness’ sake, the elements of \(F\) are \[\begin{align*} & F_0=f_0+W^0f_1=f_0+f_1 \\ & F_1=f_0+W^1f_1=f_0+e^{i\pi} f_1=f_0-f_1 \end{align*}\] Let us take \(N=4\). Then we can go further and create four Fourier transforms of length \(N/4\). We start from the second-last line of equation \eqref{eq1} and split each term into even and odd parts. \[\begin{align} F_k &= \sum_{j=0}^{N/2 - 1} e^{2 \pi ij k/(N/2)} f_{2j}+ W^{k} \sum_{j=0}^{N/2 - 1} e^{2 \pi ij k/(N/2)} f_{2j+1} \nonumber \\ &= \sum_{j=0}^{N/4 - 1} e^{2 \pi i (2j) k/(N/2)} f_{2(2j)} + \sum_{j=0}^{N/4 - 1} e^{2 \pi i (2j+1) k/(N/2)} f_{2(2j+1)} \nonumber \\ &\quad + W^{k} \sum_{j=0}^{N/4 - 1} e^{2 \pi i(2j)k/(N/2)} f_{2(2j)+1} \nonumber \\ &\quad + W^{k} \sum_{j=0}^{N/4 - 1} e^{2 \pi i(2j+1)k/(N/2)} f_{2(2j+1)+1} \nonumber \\ &= \sum_{j=0}^{N/4 - 1} e^{2 \pi ij k/(N/4)} f_{4j} + W^{2k} \sum_{j=0}^{N/4 - 1} e^{2 \pi ij k/(N/4)} f_{4j+2} \nonumber \\ &\quad + W^{k} \sum_{j=0}^{N/4 - 1} e^{2 \pi ij k/(N/4)} f_{4j+1} + W^{k} W^{2k} \sum_{j=0}^{N/4 - 1} e^{2 \pi ij k/(N/4)} f_{4j+3} \nonumber \\ &= F_{k}^{ee}+W^{2k}F_{k}^{eo}+W^{k}F_{k}^{oe}+W^{3k}F_{k}^{oo} \nonumber \\ &= f_{0} +W^{2k}f_{2} +W^{k}f_{1} +W^{3k}f_{3} \label{eq3} \end{align}\] Notice how the power of \(W\) is the same as the index of the element of the original array \(f\) multiplying it. Also, \(W^{N}=W^{0}=1\). For \(N=4\), \(W=e^{i\pi/2}\). The complete \(F\) is as follows: \[\begin{align*} & F_0=f_0+f_2+f_1+f_3 \\ & F_1=f_0-f_2+if_1-if_3 \\ & F_2=f_0+f_2-f_1-f_3 \\ & F_3=f_0-f_2-if_1+if_3 \end{align*}\] Let us do this one more time for \(N=8\). \[\begin{align} F_k &= \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j} + W^{4k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+4} \nonumber \\ &\quad + W^{2k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+2} + W^{2k} W^{4k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+6} \nonumber \\ &\quad + W^{k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+1} + W^{k} W^{4k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+5} \nonumber \\ &\quad + W^{k}W^{2k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+3} + W^{k}W^{2k}W^{4k} \sum_{j=0}^{N/8 - 1} e^{2 \pi ij k/(N/8)} f_{8j+7} \nonumber \\ &= F_{k}^{eee}+W^{4k}F_{k}^{eeo}+W^{2k}F_{k}^{eoe}+W^{6k}F_{k}^{eoo} \nonumber \\ &\quad + W^{k}F_{k}^{oee}+W^{5k}F_{k}^{oeo}+W^{3k}F_{k}^{ooe}+W^{7k}F_{k}^{ooo} \nonumber \\ &= f_{0} +W^{4k}f_{4} +W^{2k}f_{2} +W^{6k}f_{6} + W^{k}f_{1} +W^{5k}f_{5} +W^{3k}f_{3} +W^{7k}f_{7} \label{eq4} \end{align}\] Through three decompositions in halves, we have reduced our DFT calculation to 8 one-point transforms. The bookkeeping of which element in \(f\) corresponds to which odd-even subdivision is possible for small \(N\)’s, but it quickly becomes cumbersome. In equation \eqref{eq4} above, it is easy to go from the first line to the third, but to get to line 3 from line 2 is the point of bit-reversal. Take the \(eeoeoeo\ldots oo\) pattern, reverse it and assign 1 to \(o\) and 0 to \(e\). This is the value in binary of the element in \(f\) whose one-point Fourier transform is represented as \(F_{k}^{eeoeoeo\ldots oo}\). For example, \(eoo \mapsto 110 = 6\). Thus, the bit-reversed data array \(f\), in index terms, for \(N=8\), would be \(\{0,4,2,6,1,5,3,7\}\). These bit-reversed elements are then combined in pairs (going backward in our equations) until finally the combined sum of...

The Danielson-Lanczos FFT method and bit-reversal

Related Articles

(no title)

Scientists reverse brain aging, with a nasal spray

AI has torched the market for junior programmers

Al Vigier: Canada's AI strategy shouldn't include secret Palantir bills

The Anatomy of an AI-Native Org