Vector-radix FFT algorithm: Difference between revisions

Revision as of 18:11, 15 November 2016

Latest revision as of 00:07, 5 July 2025

(25 intermediate revisions by 15 users not shown)

Line 1:

{{Short description|Multidimensional fast Fourier transform algorithm}}

The '''Vectorvector-Radixradix FFT algorithm''', is a multidimensional [[Fastfast Fourier transform]] (FFT) algorithm, which is a generalization of the ordinary [[Cooley–Tukey FFT algorithm]] that divides the transform dimensions by arbitrary radices. It breaks a multidimensional (MD) [[discrete Fourier transform]] (DFT) down into successively smaller MD [[DFT]]sDFTs until, ultimately, only trivial MD [[DFT]]sDFTs need to be evaluated.<ref name="Dudgeon83">{{cite book|last1=Dudgeon|first1=Dan|last2=Russell|first2=Mersereau|title=Multidimensional Digital Signal Processing|date=September 1983|publisher=Prentice Hall|isbn=0136049591|pages=76}}</ref>.<br />

The most common multidimensional [[FFT]] algorithm is row-column algorithm, which means transforming the array first in one index and then in the other, see more in [[FFT]]. Then a radix-2 direct 2-D FFT has been developed<ref name="Rivard77">{{cite journal|last1=Rivard|first1=G.|title=Direct fast Fourier transform of bivariate functions|journal=IEEE Transactions on Acoustics, Speech, and Signal Processing|volume=25|page=250-252|doi=10.1109/TASSP.1977.1162951|url=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=1162951&isnumber=26125}}</ref>, and it can eliminate 25% of the multiplies by the conventional row-column approach. And this algorithm has been extended to rectangular arrays and arbitrary radices<ref name="Harris77">{{cite journal|last1=Harris|first1=D.|last2=McClellan|first2=J.|last3=Chan|first3=D.|last4=Schuessler|first4=H.|title=Vector radix fast Fourier transform|journal=Acoustics, Speech, and Signal Processing, IEEE International Conference on ICASSP '77|volume=2|page=548-551|doi=10.1109/ICASSP.1977.1170349|url=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=1170349&isnumber=26347}}</ref>, which is the general Vector-Radix algorithm. Perhaps it is the simplest non-row-column [[FFT]] algorithm. <br />

'''Vector-radix FFT algorithm''' can reduce the number of complex multiplications significantly, compared to row-vector algorithm. For example, for a <math>N^M</math> time matrix (M dimensions, and size N on each dimension), the number of complex multiples of vector-radix FFT algorithm is <math>\frac{2^M -1}{2^M} N^M log_2 N</math>, meanwhile, for row-column algorithm, it is <math>\frac{M N^M}{2} log_2 N</math>. And generally, even larger savings in multiplies are obtained when this algorithm is operated on larger radices and on higher dimensional arrays<ref name=Harris77/>.

The most common multidimensional [[Fast Fourier transform|FFT]] algorithm is the row-column algorithm, which means transforming the array first in one index and then in the other, see more in [[Fast Fourier transform|FFT]]. Then a radix-2 direct 2-D FFT has been developed,<ref name="Rivard77">{{cite journal|last1=Rivard|first1=G.|title=Direct fast Fourier transform of bivariate functions|journal=IEEE Transactions on Acoustics, Speech, and Signal Processing|volume=25|pageissue=250-2523|pages=250–252|doi=10.1109/TASSP.1977.1162951|urlyear=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=1162951&isnumber=261251977}}</ref>, and it can eliminate 25% of the multiplies byas compared to the conventional row-column approach. And this algorithm has been extended to rectangular arrays and arbitrary radices,<ref name="Harris77">{{cite journalbook|last1=Harris|first1=D.|last2=McClellan|first2=J.|last3=Chan|first3=D.|last4=Schuessler|first4=H.|title=VectorICASSP radix'77. fastIEEE FourierInternational Conference on transform|journal=Acoustics, Speech, and Signal Processing, IEEE|chapter=Vector Internationalradix Conferencefast onFourier ICASSPtransform '77|volume=2|pagepages=548-551548–551|doi=10.1109/ICASSP.1977.1170349|urlyear=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=1170349&isnumber=263471977}}</ref>, which is the general Vectorvector-Radixradix algorithm. Perhaps it is the simplest non-row-column [[FFT]] algorithm. <br />

'''Vector-radix FFT algorithm''' can reduce the number of complex multiplications significantly, compared to row-vector algorithm. For example, for a <math>N^M</math> timeelement matrix (M dimensions, and size N on each dimension), the number of complex multiples of vector-radix FFT algorithm for radix-2 is <math>\frac{2^M -1}{2^M} N^M \log_2 N</math>, meanwhile, for row-column algorithm, it is <math>\frac{M N^M}{ 2} \log_2 N</math>. And generally, even larger savings in multiplies are obtained when this algorithm is operated on larger radices and on higher dimensional arrays.<ref name=Harris77/>.

Overall, the vector-radix algorithm significantly reduces the structural complexity of the traditional DFT having a better indexing scheme, at the expense of a slight increase in arithmetic operations. So this algorithm is widely used for many applications in engineering, science, and mathematics, for example, implementations in image processing,<ref name="Buijs74">{{cite journal|last1=Buijs|first1=H.|last2=Pomerleau|first2=A.|last3=Fournier|first3=M.|last4=Tam|first4=W.|title=Implementation of a fast Fourier transform (FFT) for image processing applications|journal=IEEE Transactions on Acoustics, Speech, and Signal Processing|date=Dec 1974|volume=22|issue=6|pages=420–424|doi=10.1109/TASSP.1974.1162620}}</ref> and high speed FFT processor designing.<ref name="Badar15">{{cite book|last1=Badar|first1=S.|last2=Dandekar|first2=D.|title=2015 International Conference on Industrial Instrumentation and Control (ICIC) |chapter=High speed FFT processor design using radix −⁴ pipelined architecture |pages=1050–1055|doi=10.1109/IIC.2015.7150901|year=2015|isbn=978-1-4799-7165-7|s2cid=11093545 }}</ref>

== 2-D DIT case ==

As with the [[Cooley–Tukey FFT algorithm]], the two dimensional vector-radix FFT is derived by decomposing the regular 2-D [[DFT]] into sums of smaller [[DFT]]'s multiplied by "twiddle" factorfactors.<br />

A decimation-in-time ('''DIT''') algorithm means the decomposition is based on time ___domain <math>x</math>, see more in [[Cooley–Tukey FFT algorithm]].

We suppose the 2-D [[DFT]] is defined

:<math>X(k_1,k_2) = \sum_{n_1=0}^{N_1-1} \sum_{n_2=0}^{N_2-1} x[n_1, n_2] \cdot W_{N_1}^{k_1 n_1} W_{N_2}^{k_2 n_2}, </math>

where <math>k_1 = 0,\dots,N_1-1</math>, and <math>k_2 = 0,\dots,N_2-1</math>, and <math>x[n_1, n_2]</math> is aan <math>N_1 \times N_2</math> matrix, and <math>W_{N}^{k n}W_N = \exp(-j 2\pi /N)</math>.<br />

For simplicity, let us assume that <math>N_1=N_2=N</math>, and the radix-<math>(r\times r)</math>( is such that <math>N/r</math> areis integers).<bran />integer.

Using the change of variables:

* <math>n_i=rp_i+q_i</math>, where <math>p_i=0,\cdotsldots,(N/r)-1; q_i = 0,\cdotsldots,r-1;</math>

* <math>k_i=u_i+v_i N/r</math>, where <math>u_i=0,\cdotsldots,(N/r)-1; v_i = 0,\cdotsldots,r-1;</math>

where <math>i = 1</math> or <math>2</math>, then the two dimensional DFT can be written as:<ref name="Chan92">{{cite journal|last1=Chan|first1=S. C.|last2=Ho|first2=K. L.|title=Split vector-radix fast Fourier transform|journal=IEEE Transactions on Signal Processing|volume=40|issue=8|pages=2029–2039|doi=10.1109/78.150004|bibcode=1992ITSP...40.2029C|year=1992}}</ref>

where <math>i = 1</math> or <math>2</math>.<br />

:<math> X(u_1+v_1 N/r,u_2+v_2 N/r)=\sum_{q_1=0}^{r-1} \sum_{q_2=0}^{r-1} \left[ \sum_{p_1=0}^{N/r-1} \sum_{p_2=0}^{N/r-1} x[rp_1+q_1, rp_1rp_2+q_1q_2] W_{N/r}^{p_1 u_1} W_{N/r}^{p_2 u_2} \right] \cdot W_{N}W_N^{q_1 u_1+q_2 u_2} W_{r}W_r^{q_1 v_1} W_{r}W_r^{q_2 v_2},</math>

[[File:2-D2x2 radix DIT-FFT-butterfly diagram.pngsvg|thumb|400px|One stage "butterfly" for DIT vector-radix 2x2 FFT]]

Then, the two dimensional DFT can be written as:

:<math> X(u_1+v_1 N/r,u_2+v_2 N/r)=\sum_{q_1=0}^{r-1} \sum_{q_2=0}^{r-1} \left[ \sum_{p_1=0}^{N/r-1} \sum_{p_2=0}^{N/r-1} x[rp_1+q_1, rp_1+q_1] W_{N/r}^{p_1 u_1} W_{N/r}^{p_2 u_2} \right] \cdot W_{N}^{q_1 u_1+q_2 u_2} W_{r}^{q_1 v_1} W_{r}^{q_2 v_2},</math>

The equation above defines the basic structure of the 2-D DIT radix-<math>(r\times r)</math> "butterfly". (See 1-D "butterfly" in [[Cooley–Tukey FFT algorithm]])<br />

[[File:2-D DIT-FFT-butterfly.png|thumb|400px|One stage "butterfly" for DIT vector-radix 2x2 FFT]]

When <math>r=2</math>, the equation can be broken into four summations, and this leads to:<ref name=Dudgeon83/>

The equation above defines the basic structure of the 2-D DIT radix-<math>(r\times r)</math> "butterfly". (See 1-D "butterfly" in [[Cooley–Tukey FFT algorithm]])<br />

:<math> X(k_1,k_2) = S_{00}(k_1,k_2) + S_{01}(k_1,k_2) W_{N}W_N^{k_2} +S_{10}(k_1,k_2) W_{N}W_N^{k_1} + S_{11}(k_1,k_2) W_{N}W_N^{k_1+k_2}</math> for <math>0\leq k_1, k_2 <br \frac{N}{2}</math>,

When <math>r=2</math>, the equation can be broken into four summations: one over those samples of x for which both <math>n_1</math> and <math>n_2</math> are even, one for which <math>n_1</math> is even and <math>n_2</math> is odd, one of which <math>n_1</math> is odd and <math>n_2</math> is even, and one for which both <math>n_1</math> and <math>n_2</math> are odd<ref name=Dudgeon83/>, and this leads to:

:<math> X(k_1,k_2) = S_{00}(k_1,k_2) + S_{01}(k_1,k_2) W_{N}^{k_2} +S_{10}(k_1,k_2) W_{N}^{k_1} + S_{11}(k_1,k_2) W_{N}^{k_1+k_2}</math>,<br />

where <math>S_{ij}(k_1,k_2)=\sum_{n_1=0}^{N/2-1} \sum_{n_2=0}^{N/2-1} x[2 n_1 + i, 2 n_2 + j] \cdot W_{N/2}^{n_1 k_1} W_{N/2}^{n_2 k_2}</math>.

The <math>S_{ij}</math> can be viewed as the <math>N/2</math>-dimensional DFT, each over a subset of the original sample:

== 2-D DIF case ==

* <math>S_{00}</math> is the DFT over those samples of <math>x</math> for which both <math>n_1</math> and <math>n_2</math> are even;

Similarly, a decimation-in-frequency ('''DIF''', also called the Sande-Tukey algorithm) algorithm means the decomposition is based on frequency ___domain <math>X</math>, see more in [[Cooley–Tukey FFT algorithm]].<br />

* <math>S_{01}</math> is the DFT over the samples for which <math>n_1</math> is even and <math>n_2</math> is odd;

Using the change of variables:

* <math>n_i=p_i+q_i N/rS_{10}</math>, whereis the DFT over the samples for which <math>p_i=0,\cdots,(Nn_1</r)-1;math> q_iis =odd 0,\cdots,r-1;and <math>n_2</math> is even;

* <math>k_i=r u_i+v_iS_{11}</math>, whereis the DFT over the samples for which both <math>u_i=0,\cdots,(Nn_1</r)-1;math> v_iand = 0,\cdots,r-1;<math>n_2</math> are odd.

where <math>i = 1</math> or <math>2</math>.<br />

And the DFT equation can be written as:

:<math> X(r u_1+v_1,r u_2+v_2)=\sum_{p_1=0}^{N/r-1} \sum_{p_2=0}^{N/r-1} \left[ \sum_{q_1=0}^{r-1} \sum_{q_2=0}^{r-1} x[p_1+q_1 N/r, p_1+q_1 N/r] W_{r}^{q_1 v_1} W_{r}^{q_2 v_2} \right] \cdot W_{N}^{p_1 v_1+p_2 v_2} W_{N/r}^{p_1 u_1} W_{N/r}^{p_2 u_2},</math>

Thanks to the [[List of trigonometric identities#Shifts and periodicity|periodicity of the complex exponential]], we can obtain the following additional identities, valid for <math>0\leq k_1, k_2 < \frac{N}{2}</math>:

:* <math> A_X\biggl(k_1+\frac{11N}{2}(k_1,k_2\biggr) W_= S_{N00}^{(k_1+,k_2}) + A_S_{1301}(k_1,k_2) W_{N}W_N^{k_1+3 k_2} +A_-S_{3110}(k_1,k_2) W_{N}W_N^{3 k_1+k_2} +- A_S_{3311}(k_1,k_2) W_{N}W_N^{3(k_1+k_2)}</math>,<br />;

* <math>X\biggl(k_1,k_2+\frac{N}{2}\biggr) = S_{00}(k_1,k_2) - S_{01}(k_1,k_2) W_N^{k_2} +S_{10}(k_1,k_2) W_N^{k_1} - S_{11}(k_1,k_2) W_N^{k_1+k_2}</math>;

* <math>X\biggl(k_1+\frac{N}{2},k_2+\frac{N}{2}\biggr) = S_{00}(k_1,k_2) - S_{01}(k_1,k_2) W_N^{k_2} -S_{10}(k_1,k_2) W_N^{k_1} + S_{11}(k_1,k_2) W_N^{k_1+k_2}</math>.

== 2-D DIF case ==

Similarly, a decimation-in-frequency ('''DIF''', also called the Sande-TukeySande–Tukey algorithm) algorithm means the decomposition is based on frequency ___domain <math>X</math>, see more in [[Cooley–Tukey FFT algorithm]].<br />

Using the change of variables:

* <math>n_i=p_i+q_i N/r</math>, where <math>p_i=0,\ldots,(N/r)-1; q_i = 0,\ldots,r-1;</math>

* <math>k_i=r u_i+v_i</math>, where <math>u_i=0,\ldots,(N/r)-1; v_i = 0,\ldots,r-1;</math>

where <math>i = 1</math> or <math>2</math>., and the DFT equation can be written as:<brref name=Chan92/>

:<math> X(r u_1+v_1,r u_2+v_2)=\sum_{p_1=0}^{N/r-1} \sum_{p_2=0}^{N/r-1} \left[ \sum_{q_1=0}^{r-1} \sum_{q_2=0}^{r-1} x[p_1+q_1 N/r, p_1p_2+q_1q_2 N/r] W_{r}^{q_1 v_1} W_{r}^{q_2 v_2} \right] \cdot W_{N}^{p_1 v_1+p_2 v_2} W_{N/r}^{p_1 u_1} W_{N/r}^{p_2 u_2},</math>

== Other approaches ==

The [[Splitsplit-radix FFT algorithm]] has been proved to be a useful method for 1-D DFT. And this method has been applied to the vector-radix FFT to obtain a split vector-radix FFT.<ref name=Chan92/><ref name="Pei87">{{cite journalbook|last1=Pei|first1=Soo-Chang|last2=Wu|first2=Ja-Lin|title=SplitICASSP vector'87. radixIEEE 2DInternational fastConference Fourieron transform|journal=Acoustics, Speech, and Signal Processing, IEEE|chapter=Split Internationalvector Conferenceradix on2D ICASSPfast '87.Fourier transform |volume=12|date=April 1987|pagepages=1987-19901987–1990|doi=10.1109/ICASSP.1987.1169345|urls2cid=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=1169345&isnumber=26345}}</ref><ref118173900 name="Chan92">{{cite journal|last1=Chan|first1=S. C.|last2=Ho|first2=K. L.|title=Split vector-radix fast Fourier transform|journal=IEEE Transactions on Signal Processing|volume=40|page=2029-2039|doi=10.1109/78.150004|url=http://ieeexplore.ieee.org/stamp/stamp.jsp?tp=&arnumber=150004&isnumber=3966}}</ref>.

In conventional 2-D vector-radix algorithm, we decompose the incidentindices <math>k_1,k_2</math> into 4 groups:<br />

: <math>

In conventional 2-D vector-radix algorithm, we decompose the incident <math>k_1,k_2</math> into 4 groups:<br />

<math>

\begin{array}{lcl}

X(2 k_1, 2 k_2) & : & \text{even-even} \\

X(2 k_1, 2 k_2 +1) & : & \text{even-odd} \\

X(2 k_1 +1, 2 k_2) & : & \text{odd-even} \\

X(2 k_1+1, 2 k_2+1) & : & \text{odd-odd }

\end{array}

</math>

By the split vector-radix algorithm, the first three groups remain unchanged, the fourth odd-odd group is further decomposed into another four sub-groups, and seven groups in total:

: <br /math>

<math>

\begin{array}{lcl}

X(2 k_1, 2 k_2) & : & \text{even-even} \\

X(2 k_1, 2 k_2 +1) & : & \text{even-odd} \\

X(2 k_1 +1, 2 k_2) & : & \text{odd-even} \\

X(4 k_1+1, 4 k_2+1) & : & \text{odd-odd} \\

X(4 k_1+1, 4 k_2+3) & : & \text{odd-odd} \\

X(4 k_1+3, 4 k_2+1) & : & \text{odd-odd} \\

X(4 k_1+3, 4 k_2+3) & : & \text{odd-odd}

\end{array}

</math>

That means the fourth term in 2-D DIT radix-<math>(2\times 2)</math> equation, <math>S_{11}(k_1,k_2) W_{N}^{k_1+k_2}</math> becomes:<brref name="Wu89">{{cite journal|last1=Wu|first1=H.|last2=Paoloni|first2=F.|title=On the two-dimensional vector split-radix FFT algorithm|journal=IEEE Transactions on Acoustics, Speech, and Signal Processing|date=Aug 1989|volume=37|issue=8|pages=1302–1304|doi=10.1109/29.31283}}</ref>

:<math> A_{11}(k_1,k_2) W_{N}^{k_1+k_2} + A_{13}(k_1,k_2) W_{N}^{k_1+3 k_2} +A_{31}(k_1,k_2) W_{N}^{3 k_1+k_2} + A_{33}(k_1,k_2) W_{N}^{3(k_1+k_2)}</math>,<br />

where: <math> A_{ij11}(k_1,k_2)=\sum_{n_1=0}^{N/4-1} \sum_{n_2=0}W_N^{N/4-1k_1+k_2} x[4 n_1 + iA_{13}(k_1,k_2) 4 n_2 W_N^{k_1+3 j]k_2} \cdot W_+A_{N/431}(k_1,k_2) W_N^{n_13 k_1+k_2} W_+ A_{N/433}(k_1,k_2) W_N^{n_2 3(k_1+k_2)},</math><br />

where <math>A_{ij}(k_1,k_2)=\sum_{n_1=0}^{N/4-1} \sum_{n_2=0}^{N/4-1} x[4 n_1 + i, 4 n_2 + j] \cdot W_{N/4}^{n_1 k_1} W_{N/4}^{n_2 k_2}</math>

The 2-D N by N DFT is then obtained by successive use of the above decomposition, up to the last stage.<br />

It has been shown that the split vector radix algorithm has saved about 30% of the complex multiplications and about the same number of the complex additions for typical <math>1024\times 1024</math> array, compared with the vector-radix algorithm<ref name=Pei87/>.

The 2-D N by N DFT is then obtained by successive use of the above decomposition, up to the last stage.<br />

It has been shown that the split vector radix algorithm has saved about 30% of the complex multiplications and about the same number of the complex additions for typical <math>1024\times 1024</math> array, compared with the vector-radix algorithm.<ref name=Pei87/>.

==References==

{{reflist|30em}}

[[Category:FFT algorithms]]

[[Category:Digital signal processing]]

[[Category:Discrete transforms]]