U.S. patent application number 10/643164 was filed with the patent office on 2005-05-12 for methods and apparatus for fast fourier transforms.
Invention is credited to Greene, Jonathan E..
Application Number | 20050102342 10/643164 |
Document ID | / |
Family ID | 27739057 |
Filed Date | 2005-05-12 |
United States Patent
Application |
20050102342 |
Kind Code |
A1 |
Greene, Jonathan E. |
May 12, 2005 |
Methods and apparatus for fast fourier transforms
Abstract
A system for calculating fast Fourier transforms includes a
non-final stage calculating element for repetitively performing
in-place butterfly calculations for n-1 stages as well as a final
stage calculating element for performing a final stage of butterfly
calculations. The final stage calculation includes a first loop and
a second loop. The first loop performs a portion of the final stage
butterfly calculations and includes control logic to perform groups
of butterfly calculations and to store the butterfly calculation
outputs in a shuffled order in place of the inputs to result in a
correct ordering of transform outputs. The second loop performs a
remaining portion of the final stage butterfly calculations and
includes control logic to perform butterfly calculations and to
store the butterfly calculation outputs in a shuffled order in
place of the inputs to result in a correct ordering of transform
outputs.
Inventors: |
Greene, Jonathan E.; (Great
Barrington, MA) |
Correspondence
Address: |
NUTTER MCCLENNEN & FISH LLP
WORLD TRADE CENTER WEST
155 SEAPORT BOULEVARD
BOSTON
MA
02210-2604
US
|
Family ID: |
27739057 |
Appl. No.: |
10/643164 |
Filed: |
August 18, 2003 |
Related U.S. Patent Documents
|
|
|
|
|
|
Application
Number |
Filing Date |
Patent Number |
|
|
10643164 |
Aug 18, 2003 |
|
|
|
09728469 |
Nov 30, 2000 |
|
|
|
6609140 |
|
|
|
|
60168027 |
Nov 30, 1999 |
|
|
|
60192639 |
Mar 27, 2000 |
|
|
|
Current U.S.
Class: |
708/409 |
Current CPC
Class: |
G06F 17/142
20130101 |
Class at
Publication: |
708/409 |
International
Class: |
G06F 015/00 |
Claims
What is claimed is:
1-46. (canceled)
47. A system for performing a fast Fourier transform on N ordered
inputs in n stages comprising: a non-final stage calculating means
for repetitively performing in-place butterfly calculations for n-1
stages; a final stage calculating means for performing a final
stage of butterfly calculations including: a first loop means for
performing a portion of the final stage butterfly calculations, the
first loop means performing the set of butterfly calculations, and
storing butterfly calculation outputs in shuffled order in place of
the selected inputs to result in a correct ordering of transform
outputs; and a second loop means for performing a remaining portion
of the final stage butterfly calculations, the second loop means
performing two sets of butterfly calculations, and storing
butterfly calculation outputs from a first one of the two sets of
butterfly calculations in shuffled order in place of the inputs
selected for a second one of the two sets of butterfly calculations
and storing butterfly calculation outputs from the second one of
the two sets of butterfly calculations in shuffled order in place
of the inputs selected for the first one of the two sets of
butterfly calculations to result in a correct ordering of transform
outputs.
48. The system of claim 47, wherein the final stage calculating
means performs all butterfly calculations as radix-4 butterflies
having four inputs and four outputs.
49. The system of claim 48, wherein N is a power of two.
50. The system of claim 49, wherein the non-final stage calculating
means performs a first stage of radix-8 butterfly calculations
followed by n-2 stages of radix-4 butterfly calculations.
51. The system of claim 48, wherein the non-final and final stage
calculating means include a four-fold SIMD processor for performing
four radix-4 butterfly calculations at a time.
52. A method for performing a fast Fourier transform on N ordered
inputs in n stages comprising: performing non-final stage
calculations by repetitively performing in-place butterfly
calculations for n-1 stages; performing final stage calculations by
performing a final stage of butterfly calculations in a first loop
for performing a portion of the final stage butterfly calculations
and in a second loop for performing a remaining portion of the
final stage butterfly calculations; wherein each of the butterfly
calculations in the first loop and the second loop includes storing
butterfly calculation outputs in shuffled order in place of
selected inputs to result in a correct ordering of transform
outputs.
53. The method of claim 52, wherein the final stage butterfly
calculations are all performed as radix-4 butterflies having four
inputs and four outputs.
54. The method of claim 53, further comprising storing twiddle
factors for application in the butterfly calculations in groups of
four, each group having an index and the groups being stored in bit
reversed order based on the index.
55. A method for performing a fast Fourier transform on N ordered
inputs in n stages comprising: performing non-final stage
calculations by repetitively performing in-place butterfly
calculations for n-1 stages; performing final stage calculations by
performing a final stage of butterfly calculations wherein
butterfly calculation outputs are stored in shuffled order in place
of selected inputs to result in a correct ordering of transform
outputs.
56. The method of claim 55, wherein the final stage butterfly
calculations are all performed as radix-4 butterflies having four
inputs and four outputs.
57. The method of claim 56, further comprising storing twiddle
factors for application in the butterfly calculations in groups of
four, each group having an index and the groups being stored in bit
reversed order based on the index.
Description
RELATED APPLICATIONS
[0001] This application claims priority from U.S. patent
application Ser. Nos. 60/168,027, filed Nov. 30, 1999, and Ser. No.
60/192,639, filed Mar. 27, 2000, all of which are hereby
incorporated by reference.
BACKGROUND
[0002] Fast Fourier transforms (FFTs) are a set of algorithms for
computing discrete Fourier transforms (DFTs). The discrete Fourier
transform is the central computation in many spectral analysis
problems that are encountered in signal processing. For example, in
speech recognition, spectral analysis is usually the first, or at
least an early, step in processing an acoustic signal to discern
words or word components. In sonar systems, sophisticated spectrum
analysis is required for the location of surface vessels and
submarines. In radar systems, obtaining target velocity information
generally involves the measurement of spectra. Fast Fourier
transforms can in some cases reduce the computational time required
to compute discrete Fourier transforms by a factor of 100 or
more.
[0003] In general, FFTs operate by breaking an original N-point
sequence into shorter sequences, the DFTs of which can be combined
to give the DFT of the original N-point sequence. Because a direct
evaluation of an N-point DFT requires (N-1) complex multiplications
and N(N-1) complex additions, breaking down the calculation into
smaller units results in a significant time savings, even in light
of the further steps required to combine the results of the smaller
DFT calculations.
[0004] The basic calculating operation in most FFT algorithms is
the so-called "butterfly." A butterfly is a complex
adder-subtractor operation that has the same number of inputs as
outputs. The number of inputs, which corresponds to the level
ofdecomposition of the DFT, is sometimes referred to as the
"radix." Thus, a radix-4 calculation will have four inputs to, and
four outputs from, each butterfly. Taking the example of a
256-point DFT calculation, the DFT may be broken into 64 4-point
DFTs which may be combined using 64 (N/4) butterflies in each of 4
(log.sub.4(N)) stages, making a total of 4 stages in the FFT
including the 4-point DFT calculations. Combining multiplications
are also applied at each butterfly. The multipliers employed are
referred to as "twiddle factors" or sometimes as phase or rotation
factors.
[0005] One important aspect of the butterfly structure is that,
because the operation has the same number of inputs as outputs,
intermediate stage FFT calculations require no additional memory
allocation (except of course for the twiddle factors) during
calculation because the outputs of the butterfly can be written
back to memory at the same location that the inputs were taken
from. An FFT algorithm that uses the same locations to store both
the input and output sequences is called an in-place algorithm.
[0006] The order of the inputs to, and outputs from, an FFT
calculation is important because typical FFT algorithms result in
the outputs being in "bit-reversed" order from the inputs. While
algorithms for reordering the outputs to be in natural order (or
reordering the inputs to be in bit reversed order) are known in the
art, additional processing steps must be employed to perform the
reordering. Because FFTs are often calculated many at a time,
processing time could be meaningfully reduced if an improved
ordering scheme could be applied.
[0007] In addition, the processing of FFT algorithms can be
improved by employing vector processors. Traditionally, many
high-performance computational devices for processing spectra have
included a combination of a single microprocessor controlling
system functions with off-chip devices, such as DSP (digital signal
processor) farms or custom ASICs (application specific integrated
circuits) to perform specialized computations. More recently, array
or vector processors have been employed to address high-bandwidth
data processing and algorithmic-intensive computations in a single
chip. In one common type of array processor, the SIMD ("single
instruction, multiple data) processor, a single instruction
operates on a vector having multiple values rather than on a single
value in a register. One exemplary SIMD processor is the
PowerPC.TM. G4.TM. processor having AltiVec.TM. technology. FFT
algorithms, however, must be improved upon to take full advantage
of the promise offered by vector processing systems.
SUMMARY OF THE INVENTION
[0008] The methods and apparatus of the invention improve on
existing FFT calculations by providing a multistage FFT calculation
in which the final stage is characterized by two processing loops
that store the outputs of butterfly calculations in a shuffled
order that results in the FFT outputs being correctly ordered with
no need to perform an additional bit-reversal ordering pass. In one
embodiment, the invention provides a system for performing a fast
Fourier transform on N ordered inputs in n stages. The system
includes a non-final stage calculating section that repetitively
performs in-place butterfly calculations for n-1 stages, as well as
a final stage calculating section that performs a final stage of
butterfly calculations.
[0009] The final stage calculating section executes a first loop
and a second loop. The first loop performs a portion of the final
stage butterfly calculations by iterating on a table of first loop
index values consisting of values that bit-reverse into themselves.
The first loop also executes control logic to select inputs for
groups of butterfly calculations based on the first loop index
values, to perform the groups of butterfly calculations, and to
store the butterfly calculation outputs in shuffled order in place
of the selected inputs to result in a correct ordering of transform
outputs.
[0010] The second loop performs the remaining portion of the final
stage butterfly calculations by iterating on a table of second loop
index value pairs, where each pair consists of two values that
bit-reverse into each other. The second loop executes control logic
to select inputs for two groups of butterfly calculations based on
the two second loop index pair values respectively, to perform the
two groups of butterfly calculations, and to store the butterfly
calculation outputs from one of the two groups butterfly
calculations in shuffled order in place of the inputs selected for
the other of the two groups of butterfly calculations and vice
versa.
[0011] In more specific embodiments, N is a power of two and in the
final calculating stage, calculations are carried out as radix-4
butterflies. The system of the invention can also be a computer
system having a four-fold SIMD processor wherein the butterfly
calculations are carried out as four simultaneous radix-4
butterflies.
[0012] Where radix-4 butterflies are performed in groups of four, a
specific shuffling order can be employed by representing the inputs
and outputs for the four radix-4 butterflies as 4.times.4 matrices.
In this scenario, the first loop iterates through a list of first
loop index values between 0 and N/16-1 that bit reverse into
themselves. Control logic executed by the first loop selects four
groups of four consecutive inputs by transforming the first loop
index value into four input indices by multiplying the first loop
index value by four, and successively adding N/4 to result in four
input indices. Each group of four consecutive inputs is then
selected beginning with one input index. Four radix-4 butterfly
calculations, one calculation for each group of four inputs, are
performed and the outputs are stored in place of the inputs in
shuffled order. The shuffled order results from a 4.times.4 matrix
transposition and subsequent swapping of two inner columns.
[0013] Following this same scenario, the second loop iterates
through a list of second loop index pair values which includes
pairs of values between 0 and N/16-1 that bit-reverse into each
other. In the second loop, two sets of four groups of four
consecutive inputs are selected by transforming each value in the
second loop index pair into four input indices and selecting four
consecutive inputs at each input index as in the first loop. Two
sets of four radix-4 butterfly calculations are then performed and
the outputs are stored in place of each other (one set's outputs
over the other set's inputs and vice versa) using the same
shuffling of order used in the first loop.
BRIEF DESCRIPTION OF THE DRAWINGS
[0014] The invention will be more fully understood from the
following detailed description taken in conjunction with the
accompanying drawings, in which like reference numerals designate
like parts throughout the figures, and wherein:
[0015] FIG. 1 illustrates a fast Fourier transform method of the
invention;
[0016] FIG. 2 is a flow chart for a first processing loop in a
final calculation stage of the fast Fourier transform method of
FIG. 1;
[0017] FIG. 3 illustrates a shuffling of elements useful in the
processing loop of FIG. 2; and
[0018] FIG. 4 is a flow chart for a second processing loop in a
final calculation stage of the fast Fourier transform method of
FIG. 1.
DETAILED DESCRIPTION OF THE INVENTION
[0019] FIG. 1 illustrates an exemplary Fast Fourier Transform (FFT)
method of the invention for calculating the Fourier Transform of a
complex sequence of length N. While the example of FIG. 1 is
described for a 256-point FFT (N=256) with a radix-4 (R=4)
decimation-in-time decomposition, a person of ordinary skill in the
art will recognize that the invention may be applied to other than
this specific calculation. For example, the method of the invention
can readily be employed where N is a power of two and N is at least
16. In addition the method may be employed for forward or inverse
FFT calculations.
[0020] According to FIG. 1, the inputs 12, consisting of N complex
"points" to the FFT calculation are provided from computer memory.
For SIMD processing, it can be beneficial to store complex inputs
12 in two separate floating point arrays: Re[N], containing the
real portion of the sequence, and Im[N], containing the imaginary
portion of the sequence. In one example, SIMD registers can be
viewed as arrays or vectors containing four floating point values,
and SIMD operations can be viewed as four independent operations.
In addition, many SIMD architectures provide instructions for
loading and storing SIMD registers from consecutive memory
locations. Thus, storing inputs 12 in separate floating point
arrays having the real and imaginary portions of the sequence
allows for straightforward reading of the real and imaginary values
into registers where complex additions, subtractions and
multiplications can be performed in groups of four without first
having to arrange the elements within the SIMD registers. In one
embodiment, the method of the invention is carried out on a
processor having a four-fold SIMD architecture capable of
performing four identical floating point operations at one time and
having eight vector registers, each having four floating point
values, for a total of 32 floating point registers.
[0021] Applied using radix-4 decomposition and a four-fold SIMD
architecture, the FFT methodology of the invention uses in-place
radix-4 decimation-in-time (DIT) butterflies 14 for each stage 16
except for the last stage 18. It is worth noting that if N is an
odd power of two, for example, if N is 128 (2.sup.7) or 512
(2.sup.9), the first stage 16 can employ radix-8 butterflies and
the following stages can employ radix-4 butterflies. Where N does
not naturally occur as a power of two, the method of the invention
can still be employed by adding additional zero elements to the
sequence to cause N (which, as noted above, is the length of the
sequence) to be a power of two. Where N=256, each of four stages
will have 64 radix-4 butterfly 14 calculations. When four
butterflies 14 are performed simultaneously on a four-fold SIMD
processor, each stage 16 can loop 20 through the four butterfly 14
calculations 16 times to complete the stage.
[0022] For each radix-4 butterfly 14, 4 complex elements are
accessed from the input sequence (I.sub.0, I.sub.1, I.sub.2,
I.sub.3) and mathematically combined with 4 complex exponentials
(the twiddle factors, denoted as W) to produce 4 complex output
elements (O.sub.0, O.sub.1, O.sub.2, O.sub.3). In a given stage 16
of the FFT, every element in the sequence participates in exactly
one butterfly 14, and the entire sequence is read once, and written
once in-place, per stage. In writing the output elements from each
butterfly 14 in-place, the two central output elements (See O.sub.1
and O.sub.2) in each butterfly are interchanged so as to generate a
bit-reversed output for each successive stage. In operation on a
four-fold SIMD processor, 32 floating point registers are required
for inputs to four simultaneous butterfly 14 calculations as each
of the four radix-4 butterfly operations has 4 complex inputs, each
consisting of real and imaginary floating point elements.
[0023] In the butterfly 14 notation illustrated in FIG. 1 (where
one radix-4 butterfly calculation is shown), each open circle
represents an adder-subtractor with the sum of the two inputs
appearing as the upper output from the open circle and the
difference between the inputs appearing as the lower output to the
open circle. The arrows represent the complex multiplications using
the twiddle factors (denoted as W). Further details of the
operation of butterfly evaluations and their role in FFT algorithms
in general are readily available in the art, and Rabiner et al.,
THEORY AND APPLICATION OF DIGITAL SIGNAL PROCESSING, Prentice Hall,
Inc., Chapter 6 "Spectrum Analysis and the Fast Fourier Transform,"
pp. 356 to 437 (1975), is hereby incorporated by reference for this
purpose.
[0024] As is known in the art, in each stage 16 prior to the last,
the indices (zero based) of elements that participate in a given
radix-4 butterfly 14 are divisible by 4, and the same butterfly
calculation (i.e., application of the same twiddle factors) is
performed on consecutive elements in groups that are a multiple of
4 elements in size. In each stage 16 prior to the last, elements in
consecutive memory locations participate in independent
butterflies. In this same way, known algorithms for loading input
and twiddle values from memory for SISD (Single Instruction Single
Data) processor systems can also be applied on SIMD processor
systems.
[0025] In the final stage 18, radix-4 butterflies are performed on
adjacent elements (i.e., consecutive indices). If the same
calculation methodology is applied to the final stage 18 that is
employed in earlier stages 16, the output sequence will be in an
entirely bit-reversed order, requiring an additional processing
loop to reorder the data. In the final stage 18 of the invention,
the computation is performed in a first loop 22 and a second loop
24 with a particular ordering of the output elements that results
in output sequences 26 being written in natural order without the
need for a separate reordering procedure.
[0026] Beginning the first loop 22, illustrated in more detail in
FIG. 2, the lowest index value for the input elements to the first
butterfly is read from a lookup table 28. For calculating an FFT
for a complex sequence having N elements, the lookup table contains
index values that are numbers from 0 to N/16-1 that bit-reverse
into themselves, multiplied by four. For example, where N=256, the
numbers from 0 to N/16-1 that bit-reverse into themselves are 0
(0000), 6 (0110), 9 (1001), and 15 (1111). The values inthe lookup
table would then be these numbers multiplied by four, or 0, 24, 36,
and 60. These values may be stored as 0, 24, 36 and 60, or they may
be stored as the bit-reversal numbers themselves as it is
relatively inexpensive (in terms of computer resources) to multiply
by four, especially where the processor provides instructions for
shifting values left by some number of bits (and thereby
multiplying by powers of two). In addition, because the first value
in the table will generally be 0, the 0 in the lookup table can be
replaced with the number of values stored in the table. This
information can be valuable in processing the first loop 22 where
the FFT method allows for N to vary.
[0027] In the first pass through the first loop 22, the lowest
value in the lookup table is used to select the inputs for the
first four radix-4 butterflies 30. For example, if the index read
from the table is the value j, the indices for the four input
elements to the first radix-4 butterfly will elements j, j+1, j+2,
and j+3. For the other three radix-4 butterflies that are
simultaneously calculated on the four-fold SIMD processor, the
lowest indices for the input elements to the butterfly calculation
are M=N/4 indices apart. Accordingly, the inputs to the second
butterfly will have indices j+M, j+1+M, j+2+M, and j+3+M; the third
butterfly j+2M, j+1+2M, j+2+2M, and j+3+2M; and the fourth
butterfly j+3M, j+1+3M, j+2+3M, and j+3+3M. Where N=256 and the
index value is 0, the indices for the 16 complex inputs are 0, 1,
2, 3, 64, 65, 66, 67, 128, 129, 130, 131, 192, 193, 194, and 195.
Four radix-4 butterflies are then calculated 32.
[0028] The sixteen complex outputs from the first loop 22 of four
radix-4 butterflies are written back to memory in-place 36, but in
a shuffled order 34. As illustrated in FIG. 3, the 16 complex input
values corresponding to j, j+1, . . . j+3+3M can be represented as
a.sub.0, a.sub.1, a.sub.2, a.sub.3, b.sub.0, b.sub.1, b.sub.2,
b.sub.3, c.sub.0, c.sub.1, c.sub.2, c.sub.3, d.sub.0, d.sub.1,
d.sub.2, and d.sub.3. Using this convention, the inputs to the
first butterfly are a.sub.0, a.sub.1, a.sub.2, and a.sub.3; the
inputs to the second butterfly are b.sub.0, b.sub.1, b.sub.2, and
b.sub.3; and so on. The corresponding 16 complex outputs can be
denominated as a.sub.0', a.sub.1', a.sub.2', a.sub.3', b.sub.0',
b.sub.1', b.sub.2', b.sub.3', c.sub.0', c.sub.1', c.sub.2',
c.sub.3', d.sub.0', d.sub.1', d.sub.2', and d.sub.3' and the output
40 can be represented in matrix form.
[0029] The shuffling can be illustrated, as shown in FIG. 3, as two
matrix operations where the outputs, illustrated as matrix 40, are
first transposed 42, or flipped about a diagonal through elements
a.sub.0', b.sub.1', c.sub.2', and d.sub.3', as illustrated in
intermediate matrix 44. In the second matrix operation step, the
two internal columns are swapped 46 to result in final matrix 48.
The 16 complex outputs are written to memory in the order of this
final matrix 48. That is, output element a.sub.0' is written to the
same location in memory that input element a.sub.0 was retrieved
from, output element c.sub.0' is written to the same location in
memory that input element a.sub.1 was retrieved from, output
element b.sub.0' is written to the same location in memory that
input element a.sub.2 was retrieved from, and so on.
[0030] Using a four-fold SIMD processor, the shuffling can be
implemented by performing modified 4.times.4 element transposition
on the 16 real and 16 imaginary input values corresponding to a
four butterfly set prior to performing the SIMD butterfly. Using
the notation above, the real portion of the 16 input values can be
loaded into vector registers (each having four floating point
values) v.sub.0 to v.sub.3, and the imaginary portion in vector
registers v.sub.4 to v.sub.7 as follows:
1 v.sub.0: (a.sub.0r, a.sub.1r, a.sub.2r, a.sub.3r) v.sub.1:
(b.sub.0r, b.sub.1r, b.sub.2r, b.sub.3r) v.sub.2: (c.sub.0r,
c.sub.1r, c.sub.2r, c.sub.3r) v.sub.3: (d.sub.0r, d.sub.1r,
d.sub.2r, d.sub.3r) v.sub.4: (a.sub.0i, a.sub.1i, a.sub.2i,
a.sub.3i) v.sub.5: (b.sub.0i, b.sub.1i, b.sub.2i, b.sub.3i)
v.sub.6: (c.sub.0i, c.sub.1i, c.sub.2i, c.sub.3i) v.sub.7:
(d.sub.0i, d.sub.1i, d.sub.2i, d.sub.3i)
[0031] Performing the 4.times.4 transposition and swapping the
middle columns produces:
2 v.sub.0: (a.sub.0r, c.sub.0r, b.sub.0r, d.sub.0r) v.sub.1:
(a.sub.1r, c.sub.1r, b.sub.1r, d.sub.1r) v.sub.2: (a.sub.2r,
c.sub.2r, b.sub.2r, d.sub.2r) v.sub.3: (a.sub.3r, c.sub.3r,
b.sub.3r, d.sub.3r) v.sub.4: (a.sub.0i, c.sub.0i, b.sub.0i,
d.sub.0i) v.sub.5: (a.sub.1i, c.sub.1i, b.sub.1i, d.sub.1i)
v.sub.6: (a.sub.2i, c.sub.2i, b.sub.2i, d.sub.2i) v.sub.7:
(a.sub.3i, c.sub.3i, b.sub.3i, d.sub.3i)
[0032] Four butterflies can be performed on the SIMD processor
using these vectors as inputs, and the corresponding outputs, which
will now be in the correct order, can be stored to consecutive
memory locations.
[0033] Referring again to FIG. 2, if after performing the four
butterflies and shuffling the output and writing it to memory, more
bit-reversal numbers remain in the lookup table, then the method
loops 38 back and selects inputs 30, calculates butterflies 32,
shuffles the output elements 34 and writes the outputs to memory 36
for another group of four butterflies based on the next value in
the lookup table. For N=256 and using radix-4 butterflies, four
loops of four butterflies will be calculated in the first loop 22.
When all of the elements in the lookup table have been used,
processing continues to the second loop 24 of the final stage
18.
[0034] The second loop 24, illustrated in more detail in FIG. 4,
has a number of steps similar to those in the first loop 22,
however, the second loop performs eight radix-4 butterflies
corresponding a pair of index values read from a bit-reversal pair
lookup table during each pass of the loop. The bit-reversal pair
table includes all of the pairs of values between 0 and N/16-1 that
bit-reverse into each other--as opposed to the bit-reversal table
used in the first loop 22 that contained values that bit-reversed
into themselves. For example, for N=256, the table would contain
the pairs: 1, 8; 2, 4; 3, 12; 5, 10; 7, 14; and 11, 13. Thus, for
this example, the bit-reversal pair lookup table used for the
second loop 24 includes six pairs of numbers between 0 and 15 that
bit-reverse into each other. As one example, 3 (0011) bit-reverses
to 12 (1100), so 3 and 12 form one bit-reversal pair in the
bit-reversal pair lookup table. Because the bit-reversal pair table
can always be arranged so that the first number in the first pair
is unity, the first value in the table can be replaced with the
number of index values contained in the table. For the example
table described above, an N=256 bit-reversal pair lookup table
could be stored as 12, 8; 2, 4; 3, 12; 5, 10; 7, 14; and 11, 13. In
this way, the processing logic performing the second loop 24 could
use 1 as the first value for selecting inputs to the first group of
four butterflies, but would also know exactly how many values were
contained in the bit-reversal pair lookup table and thus how many
times the second loop would have to be executed to complete the
second loop calculations. In addition, as with the first loop 22
bit-reversal table, it is the bit-reversal index value multiplied
by four that will ultimately be used to select complex inputs for
butterfly calculations, and the index values may either be stored
as their actual bit-reversal pair values, or with these values
multiplied by four.
[0035] The second loop 24 begins with the selection of the first
index pair from the bit-reversal pair lookup table 50. Using the
same scheme as for selecting 16 complex inputs 30 for four radix-4
butterflies in the first loop 22, the first index value from the
index pair is used to select 16 complex inputs 52 for four radix-4
butterflies in the second loop 24. That is, the first index value
is multiplied by four (if it is not stored multiplied by four),
four consecutive elements are read from memory starting with this
index value, and four consecutive elements are read from memory
starting with index values that are M=N/4 elements apart. Using the
example of the first index value from the first index pair from the
bit-reversal pair lookup table being equal to one and N=256, the
indices for the 16 complex inputs to the first group of four
radix-4 butterflies are 4, 5, 6, 7, 68, 69, 70, 71, 132, 133, 134,
135, 196, 197, 198, and 199. The complex inputs are read from
memory locations corresponding to these index values, and
computation of the first group of four radix-4 butterflies 54 can
begin.
[0036] The second index value of the index pair selected from the
bit-reversal pair lookup table is then used to select 16 complex
inputs 56 for a second group of four radix-4 butterflies. Using the
same example from above, where the first index pair is 1, 8 (and
thus the second index value is 8), the 16 indices used to read
complex inputs from memory are 32, 33, 34, 35, 96, 97, 98, 99, 160,
161, 162, 163, 224, 225, 226, 227 and 228. The complex inputs are
read from memory locations corresponding to these index values, and
computation of the second group of four radix-4 butterflies 60 can
begin.
[0037] The sixteen complex outputs from the first group of four
radix-4 butterflies are stored, in an order shuffled from the
inputs in the manner illustrated in FIG. 3, in place of the inputs
for the second group of four radix-4. butterflies in memory 58.
Similarly, the sixteen complex outputs from the second group of
four radix-four butterflies are stored in shuffled order in memory
in place of the sixteen complex inputs for the first group of four
radix-4 butterflies 62. Temporary storage can be minimized in
performing these reading, calculating and storing operations if the
second sixteen element input set is loaded from memory 56 after the
four butterflies have been calculated for the first group 54 of
sixteen elements. The sixteen output elements from the first group
of four butterflies can then be stored immediately to the memory
locations corresponding to the sixteen input elements for the
second group of four butterflies 58.
[0038] The second loop 24 continues until all of the bit-reversal
pairs have been exhausted 64. Continuing with the N=256 example,
the second loop 24 will iterate through 6 bit-reversal pairs, with
each iteration performing two groups of four radix-4 butterflies,
for a total of 48 butterflies performed in the second loop. As
noted above, 16 butterfly calculations were performed in this
example for the first loop 22, so, for the first and second loops
combined, the required 64 radix-4 butterfly calculations for the
final stage 18 have been performed, and the outputs have been
written to memory in the required order so that a separate
bit-reversal pass need not be performed and the FFT calculation may
end 66.
[0039] The method of the invention can also benefit from storing
the twiddle factors in a way that compliments the calculation of
butterflies in the final stage of the FFT transform method of the
invention. Typically, N/4 twiddle sets (each having four twiddle
factors) are stored in a table in bit-reversed order (the
bit-reversal being done on the indices of the twiddle sets from 0
to N/4-1). Each stage of the FFT uses more twiddle sets than are
used in the prior stage. In the example used above, with N=256 and
using radix-4 butterflies, each stage uses four times as many
twiddle factors as the stage before. In the final stage of the FFT
method of the invention, access to twiddle factors stored in the
typical manner is effectively done in a random access way.
[0040] The twiddle factors can be ordered in memory, however, in a
manner that allows them to be accessed much more efficiently during
the final stage. The twiddle factors can be organized into groups
of four sets of four factors, much like a 4.times.4 floating point
matrix. These groups are then ordered in a bit-reversed order over
indices of 0 to N/16-1. When the twiddle factors are ordered in
this way, they are ordered so that the same indices that are used
to select inputs for the four-fold SIMD radix-4 butterflies can
also be used to select the twiddle factors from the table, simply
by multiplying the index by 16 instead of by four. The group of
twiddle factors can then be loaded consecutively into four SIMD
registers (i.e., four floating point values into each of four
vector registers).
[0041] In stages prior to the final stage, the same twiddle table
can be accessed sequentially from its start in each stage. The next
sequential twiddle set in this context is simply the next set of
twiddle factors from the first column of twiddle factors in each
4.times.4 group. Alternatively, a second twiddle table may be
provided for the stages prior to the final stage. The second
twiddle table can be smaller, one fourth of the size of the table
used for the final stage, and can be organized like a traditional
bit-reversed twiddle table for a N/4 FFT.
[0042] The FFT method of the invention can be implemented in any
manner in which its logic can practically be carried out. The
method can be embodied in software that can run on a variety of
computers known in the art, including without limitation, personal
computers, workstation computers, and supercomputers known in the
art. In one preferred embodiment, the FFT method of the invention
is implemented in software and executed on a computer having one or
more vector processors such as the PowerPC.TM. G4.TM. processor
having AltiVec.TM. technology. An example of computer code useful
for programming such a computer to carry out the method of the
invention is included in U.S. patent application Ser. No.
60/168,027, from which this application claims priority, and which
is incorporated by reference into this disclosure. In addition, the
logic of the FFT method of the invention can be implemented in a
special purpose circuit, such as a circuit including ASICs having
the specific logic employed herein deployed within the circuit, or
in a system that combines such special purpose circuits with a
general purpose computer with software.
[0043] One skilled in the art will appreciate further features and
advantages of the invention based on the above-described
embodiments. Accordingly, the invention is not to be limited by
what has been particularly shown and described, except as indicated
by the appended claims. All publications and references cited
herein are expressly incorporated herein by reference in their
entirety.
* * * * *