99 lines
1.6 KiB
Text
99 lines
1.6 KiB
Text
PROGRAM FFT
|
|
|
|
CONST CNT=8
|
|
|
|
!$DYNAMIC
|
|
DIM REL[0],IMG[0],CMP[0],V[0]
|
|
|
|
BEGIN
|
|
SIG=INT(LOG(CNT)/LOG(2)+0.9999)
|
|
REAL1=2^SIG
|
|
|
|
REAL=REAL1-1
|
|
REAL2=INT(REAL1/2)
|
|
REAL4=INT(REAL1/4)
|
|
REAL3=REAL4+REAL2
|
|
|
|
!$DIM REL[REAL1],IMG[REAL1],CMP[REAL3]
|
|
|
|
FOR I=0 TO CNT-1 DO
|
|
READ(REL[I],IMG[I])
|
|
END FOR
|
|
|
|
DATA(1,0,1,0,1,0,1,0,0,0,0,0,0,0,0,0)
|
|
|
|
SIG2=INT(SIG/2)
|
|
SIG1=SIG-SIG2
|
|
CNT1=2^SIG1
|
|
CNT2=2^SIG2
|
|
|
|
!$DIM V[CNT1-1]
|
|
V[0]=0
|
|
DV=1
|
|
PTR=CNT1
|
|
|
|
FOR J=1 TO SIG1 DO
|
|
HLFPTR=INT(PTR/2)
|
|
PT=CNT1-HLFPTR
|
|
FOR I=HLFPTR TO PT STEP PTR DO
|
|
V[I]=V[I-HLFPTR]+DV
|
|
END FOR
|
|
DV=2*DV
|
|
PTR=HLFPTR
|
|
END FOR
|
|
|
|
K=2*π/REAL1
|
|
|
|
FOR X=0 TO REAL4 DO
|
|
CMP[X]=COS(K*X)
|
|
CMP[REAL2-X]=-CMP[X]
|
|
CMP[REAL2+X]=-CMP[X]
|
|
END FOR
|
|
|
|
PRINT("FFT: BIT REVERSAL")
|
|
|
|
FOR I=0 TO CNT1-1 DO
|
|
IP=I*CNT2
|
|
FOR J=0 TO CNT2-1 DO
|
|
H=IP+J
|
|
G=V[J]*CNT2+V[I]
|
|
IF G>H THEN
|
|
SWAP(REL[G],REL[H])
|
|
SWAP(IMG[G],IMG[H])
|
|
END IF
|
|
END FOR
|
|
END FOR
|
|
|
|
T=1
|
|
FOR STAGE=1 TO SIG DO
|
|
PRINT("STAGE:";STAGE)
|
|
D=INT(REAL2/T)
|
|
FOR II=0 TO T-1 DO
|
|
L=D*II
|
|
LS=L+REAL4
|
|
FOR I=0 TO D-1 DO
|
|
A=2*I*T+II
|
|
B=A+T
|
|
F1=REL[A]
|
|
F2=IMG[A]
|
|
CNT1=CMP[L]*REL[B]
|
|
CNT2=CMP[LS]*IMG[B]
|
|
CNT3=CMP[LS]*REL[B]
|
|
CNT4=CMP[L]*IMG[B]
|
|
REL[A]=F1+CNT1-CNT2
|
|
IMG[A]=F2+CNT3+CNT4
|
|
REL[B]=F1-CNT1+CNT2
|
|
IMG[B]=F2-CNT3-CNT4
|
|
END FOR
|
|
END FOR
|
|
T=2*T
|
|
END FOR
|
|
|
|
PRINT("NUM REAL IMAG")
|
|
FOR I=0 TO REAL DO
|
|
IF ABS(REL[I])<1E-5 THEN REL[I]=0 END IF
|
|
IF ABS(IMG[I])<1E-5 THEN IMG[I]=0 END IF
|
|
PRINT(I;"";)
|
|
WRITE("##.###### ##.######";REL[I];IMG[I])
|
|
END FOR
|
|
END PROGRAM
|