想当年，俺毕业论文也是做的DFFT，用的是PC机，用BASIC编程，绘一张DFT的信息图要一天！呵呵
下面这个是从网上找到的，自己看
{$N+,R-}
unit fftunit;
interface
type real = extended;
type realArray = array [0..1024] of real;
var fr, fi : realarray;
procedure fft(var fr, fi : realarray; lnN : integer; sign : integer);
procedure Init(lnN : integer);
const lnN = 10;
var n : integer;
divN : real;
implementation
function power2(i : integer) : integer;
var n : integer;
begin
n := 1;
for i := 1 to i do n := n*2;
power2:= n;
end;
{From "Image Processing" by Jan Teuber}
procedure fft(var fr, fi : realarray; lnN : integer; sign : integer);
var t, nd2, nm1 : integer;
i,j,k,l, le, le1, ip : integer;
label 1,2,3;
var s, ur, ur1, ui, wr, wi, tr, ti, dv : real;
begin
n := power2(lnN);
nd2 := N div 2;
nm1 := n-1;
j := 1;
for i := 1 to N-1 do
begin
if i >= j then goto 1;
s := fr; fr := fr[j]; fr[j] := s;
s := fi; fi := fi[j]; fi[j] := s;
1: k := nd2;
2: if k>=j then goto 3;
j := j - k;
k := k div 2;
goto 2;
3: j := j + k;
end;
For l := 1 to lnN do
begin
le := Power2(l);
le1 := le div 2;
ur := 1;
ui := 0;
wr := cos(pi/le1);
wi := sign*sin(pi/le1);
for j:=1 to le1 do
begin
i := j;
while (i<=N) do
begin
ip := i + le1;
tr := fr[ip]*ur - fi[ip]*ui;
ti := fr[ip]*ui + fi[ip]*ur;
fr[ip] := fr - tr;
fi[ip] := fi - ti;
fr := fr + tr;
fi := fi + ti;
i := i + le;
end;
ur1 := ur*wr - ui*wi;
ui := ur*wi + ui*wr;
ur := ur1;
end;
end; {l - l沰ke}
for i := 1 to N do
begin
fr := fr*divN;
fi := fi*divN;
end;
end;
procedure Init(lnN : integer);
var i : integer;
begin
n := power2(lnN);
for i := 0 to n do
begin
fi := 0;
fr := sin(i/n*4*pi)+0.2*sin(i/n*40*pi) + 0.1*random;
{random*2 +1;}
end;
divN := 1/sqrt(n*1.0);
end;
end.
init(lnN);
grafik;
draw(fr,1);
wtg;
setcolor(0);
{erase: }
draw(fr,1);
fft(fr, fi, LnN, -1);
{remove high frequencies}
{Try different values of 'keep'}
{keep := 5;}
keep := 5;
for t := 1+ keep to N-keep do
begin
fr[t] := 0;
fi[t] := 0;
end;
fft(fr, fi, lnN, 1);
setcolor(green);
draw(fr, 1);
wtg;
end.
