4,494
edits
Fast Fourier transform: Difference between revisions
Jump to navigation
Jump to search
m
Fast Fourier transform (view source)
Revision as of 21:06, 10 September 2022
, 21:06, 10 September 2022Text replacement - "<source" to "<syntaxhighlight"
(→Simple recursive implementation: otherwise produces text garbage on incorrect input) |
m (Text replacement - "<source" to "<syntaxhighlight") |
||
| (One intermediate revision by the same user not shown) | |||
| Line 5: | Line 5: | ||
== Implementations == | == Implementations == | ||
=== APLX === | === APLX === | ||
This FFT code is implemented with the [[wikipedia:Cooley–Tukey FFT algorithm|Cooley–Tukey FFT algorithm]] by dividing the transform into two pieces of size < | This FFT code is implemented with the [[wikipedia:Cooley–Tukey FFT algorithm|Cooley–Tukey FFT algorithm]] by dividing the transform into two pieces of size <syntaxhighlight lang=apl inline>N÷2</syntaxhighlight> at each step. It will run under [[APLX]]. | ||
This is as given in Robert J. Korsan's article in APL Congress 1973, p 259-268, with just line labels and a few comments added. | This is as given in Robert J. Korsan's article in APL Congress 1973, p 259-268, with just line labels and a few comments added. | ||
* X and Z are two-row [[Matrix|matrices]] representing the input and output real and imaginary data. The data length must be < | * X and Z are two-row [[Matrix|matrices]] representing the input and output real and imaginary data. The data length must be <syntaxhighlight lang=apl inline>2*N</syntaxhighlight> (N integer), and the algorithm will cope with varying N, unlike many DSP versions which are for fixed N. | ||
* Zero frequency is at < | * Zero frequency is at <syntaxhighlight lang=apl inline>Z[1;]</syntaxhighlight>, maximum frequency in the middle; from there to <syntaxhighlight lang=apl inline>¯1↑[1] Z</syntaxhighlight> are negative frequencies. i.e. for an input Gaussian they transform a 'bath-tub' to a 'bath-tub'. | ||
* This is an elegant algorithm, and works by transforming the input data into an array of 2×2 [[wikipedia:Butterfly diagram|FFT Butterflies]]. | * This is an elegant algorithm, and works by transforming the input data into an array of 2×2 [[wikipedia:Butterfly diagram|FFT Butterflies]]. | ||
< | <syntaxhighlight lang=apl> | ||
Z←fft X;N;R;M;L;P;Q;S;T;O | Z←fft X;N;R;M;L;P;Q;S;T;O | ||
⍝ | ⍝ | ||
| Line 30: | Line 30: | ||
→loop | →loop | ||
done:Z←⍉(N,2)⍴(+⌿X),[O-0.5]-⌿X | done:Z←⍉(N,2)⍴(+⌿X),[O-0.5]-⌿X | ||
</ | </syntaxhighlight> | ||
=== Simple recursive implementation === | === Simple recursive implementation === | ||
< | <syntaxhighlight lang=apl> | ||
fft←{ | fft←{ | ||
1>K←2÷⍨M←≢⍵:⍵ | 1>K←2÷⍨M←≢⍵:⍵ | ||
| Line 42: | Line 42: | ||
(odd+T),odd-T←even×exp | (odd+T),odd-T←even×exp | ||
} | } | ||
</ | </syntaxhighlight> | ||
{{Works in|[[Dyalog APL]]}} | {{Works in|[[Dyalog APL]]}} | ||
Sample usage: | Sample usage: | ||
| Line 50: | Line 50: | ||
Inverse FFT can be defined for testing: | Inverse FFT can be defined for testing: | ||
< | <syntaxhighlight lang=apl> | ||
ifft←{(≢⍵)÷⍨+fft+⍵} | ifft←{(≢⍵)÷⍨+fft+⍵} | ||
test←{⌈/(10○⊢)(⍵-ifft fft ⍵)} | test←{⌈/(10○⊢)(⍵-ifft fft ⍵)} | ||
test 1 1 1 1 0 0 0 0 | test 1 1 1 1 0 0 0 0 | ||
7.850462E¯17 | 7.850462E¯17 | ||
</ | </syntaxhighlight> | ||
2-dimensional FFT and inverse 2D FFT: | 2-dimensional FFT and inverse 2D FFT: | ||
< | <syntaxhighlight lang=apl> | ||
fft2D←{ | fft2D←{ | ||
∨/0≠1|2⍟⍴⍵:'Matrix dimensions must be powers of 2' | ∨/0≠1|2⍟⍴⍵:'Matrix dimensions must be powers of 2' | ||
| Line 64: | Line 64: | ||
} | } | ||
ifft2D←{(≢∊⍵)÷⍨+fft2D+⍵} | ifft2D←{(≢∊⍵)÷⍨+fft2D+⍵} | ||
</ | </syntaxhighlight> | ||
Sample usage: | Sample usage: | ||
| Line 76: | Line 76: | ||
=== Dyalog APL === | === Dyalog APL === | ||
FFT appears in [[dfns.dws]], a [[workspace]] supplied with [[Dyalog APL]], in the context of fast multi-digit multiplication<ref>dfns.dws: [http://dfns.dyalog.com/n_xtimes.htm xtimes] — Fast multi-digit product using FFT</ref>. Extracted from there, it is there defined as: | FFT appears in [[dfns.dws]], a [[workspace]] supplied with [[Dyalog APL]], in the context of fast multi-digit multiplication<ref>dfns.dws: [http://dfns.dyalog.com/n_xtimes.htm xtimes] — Fast multi-digit product using FFT</ref>. Extracted from there, it is there defined as: | ||
< | <syntaxhighlight lang=apl> | ||
roots←{×\1,1↓(⍵÷2)⍴¯1*2÷⍵} | roots←{×\1,1↓(⍵÷2)⍴¯1*2÷⍵} | ||
cube←{⍵⍴⍨2⍴⍨2⍟⍴⍵} | cube←{⍵⍴⍨2⍴⍨2⍟⍴⍵} | ||
floop←{(⊣/⍺)∇⍣(×m)⊢(+⌿⍵),[m-0.5]⍺×[⍳m←≢⍴⍺]-⌿⍵} | floop←{(⊣/⍺)∇⍣(×m)⊢(+⌿⍵),[m-0.5]⍺×[⍳m←≢⍴⍺]-⌿⍵} | ||
FFT←{,(cube roots⍴⍵)floop cube ⍵} | FFT←{,(cube roots⍴⍵)floop cube ⍵} | ||
</ | </syntaxhighlight> | ||
== References == | == References == | ||