(Comp.sources.hp48)
Item: 138 by kcooke@milton.u.washington.edu [Ken Cooke]
Subj: Machine Code FFT
Date: 12 Jul 1992

Here is some FFT code that I wrote several months ago but never got around
to posting.  I tried to make it fast enough to be useful.  The code is based
on the standard DIT radix-2 FFT, so the input length must be a power of 2.
All intermediate math is done in 15-digit precision, which gives very accurate
results. (a 1024-point FFT matches MATLAB to 12 digits)

I used the following definition of the FFT:
                N-1
         X[k] = Sum( x[n] * exp(-j*2pi*k*n/N) )
                n=0
                      N-1
         x[n] = 1/N Sum( X[k] * exp(j*2pi*k*n/N) )
                      k=0

Quick description of the programs in this UUENCODED directory:

FFT: machine code FFT.  Input is a vector (either real or complex).
Size must be a power of two.  Output is a complex vector.  Includes
error checking.  I believe this is as fast as the HP48 can be at a complex
FFT.  Random 128-pt complex FFT: 7.5 seconds.  About 8 times faster than
system RPL, and 1.5 times faster than the FFT in the HP71 mathpack after
compensating for clock speed differences.  (I did not look at the mathpack
code)

IFFT: inverse FFT.  Could have been combined with the FFT program, but I
hate using an additional argument to tell it which one to do.

ABSV: machine code absolute value.  Input can be a vector or matrix.
Output is a vector or matrix containing the absolute values of the complex
elements.  About 5 times faster than RPL.

ANGLE: same as ABSV, but returns the angle (i.e. phase) according to the
current angle mode.  This is only 2 or 3 times as fast as RPL, because
the mcode arctan is such a bottleneck.

SPEC(trum): User code quickie that does FFT, ABSV, and then plots using the
built-in barplot function.  Note that plot is from 0..2pi as returned by FFT.

DISCLAIMER
While I believe these programs to be accurate and bug-free, I
will assume no responsibility for any Memory Lost, hardware damage, or
poor test scores that result from use of this code.

Please let me know if you find a bug or have suggestions.

Ken Cooke
N7VFE
kcooke@u.washington.edu