home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
HAM Radio 1
/
HamRadio.cdr
/
tech
/
design2
/
fft.pas
< prev
next >
Wrap
Pascal/Delphi Source File
|
1986-07-22
|
11KB
|
379 lines
program fft;
{ Forward & Inverse FFT Program -- Jeff Crawford
Reference: " Add the FFT to Your Box of Design Tools", Robert Cushman,
EDN, September 16, 1981.
Program Calculates Both Forward and Inverse FFT's. Data May be
Entered Manually, Point by Point, Read off Disk, or Calculated by
a Program FFT.In and Inserted in Input Arrays.
Plotting of Both Input & Output Waveforms is Provided Through the
Use of Turbo Graphix Toolbox. At Present, the Largest FFT Capable of
Using the Plotting Package is a 256 Point FFT.
}
{$ITypedef.sys}
{$Igraphix.sys}
{$Ikernel.sys}
{$iwindows.sys}
{$Iaxis.hgh}
{$IPolygon.hgh}
{$IModpoly.hgh}
type vector = array[0..256] of real;
XYdata = record
XX,YY: real;
end;
var Ans,Anz,An : string[1];
k : real;
N : integer;
D1,D2,X1,Y1 : vector;
N3,N4,N6,N7,N8 : integer;
N1,N2 : integer;
a,c,s,t1,t2 : real; { Used in FFt }
i : integer;
file_name : string[10];
write_file : file of XYdata;
read_file : file of XYdata;
file_rec : XYdata;
OK : Boolean;
FI : integer; { Flagged for Forward or Inverse FFT }
PlotTr,PlotTi : PlotArray; { Arrays for Plotting Input }
Graph_title : string[25];
xloc,yloc : integer;
Title : string[15];
PlotFr,PlotFi : PlotArray; { Arrays for Plotting Output }
tinc : real; { Used in FFT_Input Subroutine }
jm : integer; { Provides Plotting for > 100 Points }
label 10,20,30,99;
{$IFFT.In} { Program Automatically Fills Input Arrays }
Procedure PlotG( A,B:PlotArray); { PlotArray Defined as Type in Graphics.SYs}
var ymin,ymax : real; { For Real Components }
ymin1,ymax1 : real; { For Imag Components }
i : integer;
Colornum : integer;
begin
Colornum:= 16; { Sets Color of Graphics Content }
N:= N div jm;
initgraphic;
ymin:= 1.0E20; ymax:= -1.0E20;
for i:= 1 to N do
begin
if A[i,2] > ymax then ymax:= A[i,2];
if A[i,2] < ymin then ymin:= A[i,2];
end;
ymin1:= 1.0E20; ymax1:= -1.0E20;
for i:= 1 to N do
begin
if B[i,2] > ymax1 then ymax1:= B[i,2];
if B[i,2] < ymin1 then ymin1:= B[i,2];
end;
if ymax < ymax1 then ymax:= ymax1;
if ymin > ymin1 then ymin:= ymin1;
Definewindow(1,2,0,Xmaxglb,Ymaxglb);
SetforegroundColor(Colornum);
Defineworld(1,A[1,1],Ymax,A[N,1],Ymin);
Defineheader(1,Graph_title);
setheaderon;
setbackground(0);
Selectworld(1);
Selectwindow(1);
Setlinestyle(0);
Drawborder;
DrawAxis(7,7,8,2,2,10,0,0,false);
DrawPolygon(A,1,N,4,1,5);
SetLinestyle(1);
Drawaxis(7,7,8,2,2,10,0,-1,false); { Must Recall as per Graphix Book p128}
DrawPolygon(B,1,-N,3,1,5); { Draws Imaginary Portion }
SetLinestyle(0);
SetClippingOff;
Drawtext(300,190,1,' C E L L # ' );
xloc:= 50;
yloc:= 25;
Title[1]:= 'M';
Title[2]:= 'A';
Title[3]:= 'G';
Title[4]:= 'N';
Title[5]:= 'I';
Title[6]:= 'T';
Title[7]:= 'U';
Title[8]:= 'D';
Title[9]:= 'E';
for i:= 1 to 9 do
begin
Drawtext(xloc,yloc+i*15,1,Title[i]);
end;
SetClippingOn;
repeat until keypressed;
leavegraphic;
N:= N * jm;
end;
procedure Read_Data;
label 2;
begin
2:clrscr;
gotoXY(10,10);
write('Enter File Name ');
{$I-} readln(file_name) {$I+};
OK:= (IOresult = 0);
if NOT OK then goto 2;
assign(read_file,file_name);
{$I-} Reset(read_file) {$I+}; { Flags if Disk is Non-Existent }
OK:= (IOResult = 0);
if not OK then
begin
clrscr;
gotoXY(10,10);
writeln('File "',file_name,'" Is Not Available ');
writeln;
Delay(1000);
clrscr;
goto 2;
end;
i:= 0;
clrscr;
writeln(' Filename: ',file_name:10);
writeln;
writeln('Index Real Imaginary');
writeln;
while not EOF(read_file) do { If File is Available, Reads From Disk }
begin
read(read_file,file_rec);
with file_rec do
begin
D1[i]:= XX;
D2[i]:= YY;
writeln(' ',i+1,' ',D1[i]:12:6,' ',D2[i]:12:6);
i:= i + 1;
end;
end;
N:= i;
close(read_file);
end;
procedure Write_Data( X,Y: vector);
var Ans1,Ans2: string[1];
label 4;
begin
4: clrscr;
gotoXY(10,10);
write('Write Data to What File ?');
readln(file_name);
assign(write_file,file_name);
{$I-} Reset(write_file) {$I+};
OK:= (Ioresult = 0); { Determines if Write-To File Already Exists }
if OK then
begin
clrscr;
gotoXY(10,10);
write('File Already Exists. . . Overwrite ? Y / N ? ');
readln(Ans2);
Ans2:= upcase(Ans2);
if Ans2 <> 'Y' then goto 4;
end;
assign(write_file,file_name);
rewrite(write_file);
clrscr;
writeln(' Filename: ',file_name:12);
writeln;
writeln('Index Real Imaginary');
writeln;
with file_rec do
begin
for i:= 0 to N-1 do { Writes to Disk File "write_file"}
begin
writeln(' ',i+1,' ',X[i]:12:6,' ',Y[i]:12:6);
XX:= X[i];
YY:= Y[i];
write(write_file,file_rec);
end;
end;
close(write_file);
end;
function XtoY(x,y: real): real;
begin
XtoY:= exp(y*Ln(x));
end;
begin
99: clrscr;
gotoxy(25,10);
write(' Fft Ift "I" or "F" ');
for i:= 0 to 256 do
begin
D1[i]:= 0.0;
D2[i]:= 0.0;
end;
readln(An);
An:= upcase(An);
if An = 'F' then FI:= 1 else FI:= -1;
gotoxy(5,10);
clrscr;
write(' Number of Steps FFT/IFT ');
readln(N);
k:= Trunc(Ln(N)/Ln(2));
clrscr;
writeln;
write(' Read Diskfile = F Manual Input = M ');
readln(Anz);
Anz:= upcase(Anz);
if Anz = 'M' then
begin
clrscr;
writeln;
write(' M - Manual Input F - File Program Entry ');
readln(Ans);
Ans:= upcase(Ans);
if Ans = 'M' then
begin
writeln('Input Data ');
writeln;
writeln(' Data # Real Imag ');
tinc:= N;
tinc:= 1/tinc;
for i:= 0 to N-1 do
begin
write(' ',i+1,' ');
readln(D1[i], D2[i]);
end;
end;
if Ans = 'F' then
begin
Input_Data;
end;
end
else Read_Data;
writeln;
write(' Write Data to Disk Y/N ? ');
readln(Anz);
Anz:= upcase(Anz);
if Anz = 'Y' then Write_Data(D1,D2);
writeln;
write(' Plot of Input Data [real] Y/N ? ');
readln(Anz);
Anz:= upcase (Anz);
if Anz = 'Y' then
begin
clrscr;
gotoxy(10,10);
write('Input Graph Title ');
readln(Graph_title);
if N > 100 then jm:= 4 else jm:= 1;
i:= 1;
while i <= N div jm do
begin
PlotTr[i,1]:= i*jm;
PlotTr[i,2]:= D1[(i-1)*jm];
PlotTi[i,1]:= i*jm;
PlotTi[i,2]:= D2[(i-1)*jm];
i:= i + 1;
end;
PlotG(PlotTr,PlotTi);
end;
N1:= 0;
N2:= N-1;
for N3:= 1 to N2 do
begin
N4:= N;
10: N4:= N4 div 2;
if N1 + N4 > N2 then goto 10;
N1:= N1 - ( N1 div N4 ) * N4 + N4;
if N1 <= N3 then goto 20;
t1:= D1[N3];
D1[N3]:= D1[N1];
D1[N1]:= t1;
t2:= D2[N3];
D2[N3]:= D2[N1];
D2[N1]:= t2;
20: end;
N4:= 1;
30: N6:= 2*N4;
for N3:= 0 to N4 - 1 do
begin
a:= Fi*N3*Pi/N4;
c:= Cos(a);
s:= Sin(a);
N7:= N3;
while N7 <= N-1 do
begin
N8:= N7 + N4;
t1:= c*D1[N8] - s*D2[N8];
t2:= c*D2[N8] + s*D1[N8];
D1[N8]:= D1[N7] - t1;
D2[N8]:= D2[N7] - t2;
D1[N7]:= D1[N7] + t1;
D2[N7]:= D2[N7] + t2;
N7:= N7 + N6;
end;
end;
N4:= N6;
if N4 < N then goto 30;
clrscr;
write(' Output Results ');
if Fi= 1 then writeln('-- Forward Transform');
if Fi= -1 then writeln('-- Inverse Transform');
writeln;
writeln('Index Real Imaginary');
writeln;
if An = 'F' then
begin
for i:= 0 to N-1 do
begin
X1[i]:= D1[i]/N; Y1[i]:= D2[i]/N;
writeln(' ',i+1,' ',X1[i]:12:6,' ',Y1[i]:12:6);
end;
end;
if An = 'I' then
begin
for i:= 0 to N-1 do
begin
X1[i]:= D1[i]; Y1[i]:= D2[i];
writeln(' ',i+1,' ',X1[i]:12:6,' ',Y1[i]:12:6);
end;
end;
gotoxy(60,24);
write('Any Key to Cont');
repeat until keypressed;
clrscr;
gotoxy(10,10);
write(' Write Data to Disk Y/N ? ');
readln(Anz);
Anz:= upcase(Anz);
if Anz = 'Y' then Write_Data(X1,Y1);
writeln;
write(' Plot of Input Data [real] Y/N ? ');
readln(Anz);
Anz:= upcase (Anz);
if Anz = 'Y' then
begin
clrscr;
gotoxy(10,10);
write('Input Graph Title ');
readln(Graph_title);
if N > 100 then jm:= 4 else jm:= 1;
i:= 1;
while i <= N div jm do
begin
PlotFr[i,1]:= i*jm;
PlotFr[i,2]:= X1[(i-1)*jm];
PlotFi[i,1]:= i*jm;
PlotFi[i,2]:= Y1[(i-1)*jm];
i:= i + 1;
end;
PlotG(PlotFr,PlotFi);
end;
goto 99;
end.$