101 lines
1.9 KiB
Plaintext
101 lines
1.9 KiB
Plaintext
cnt = 8
|
|
sig = int(log(cnt) /log(2) +0.9999)
|
|
|
|
pi = 3.14159265
|
|
real1 = 2^sig
|
|
|
|
real = real1 -1
|
|
real2 = int(real1 / 2)
|
|
real4 = int(real1 / 4)
|
|
real3 = real4 +real2
|
|
|
|
dim rel(real1)
|
|
dim img(real1)
|
|
dim cmp(real3)
|
|
|
|
for i = 0 to cnt -1
|
|
read rel(i)
|
|
read img(i)
|
|
next i
|
|
|
|
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
|
|
hlfPtr = int(ptr / 2)
|
|
pt = cnt1 - hlfPtr
|
|
for i = hlfPtr to pt step ptr
|
|
v(i) = v(i -hlfPtr) + dv
|
|
next i
|
|
dv = dv + dv
|
|
ptr = hlfPtr
|
|
next j
|
|
|
|
k = 2 *pi /real1
|
|
|
|
for x = 0 to real4
|
|
cmp(x) = cos(k *x)
|
|
cmp(real2 - x) = 0 - cmp(x)
|
|
cmp(real2 + x) = 0 - cmp(x)
|
|
next x
|
|
|
|
print "fft: bit reversal"
|
|
|
|
for i = 0 to cnt1 -1
|
|
ip = i *cnt2
|
|
for j = 0 to cnt2 -1
|
|
h = ip +j
|
|
g = v(j) *cnt2 +v(i)
|
|
if g >h then
|
|
temp = rel(g)
|
|
rel(g) = rel(h)
|
|
rel(h) = temp
|
|
temp = img(g)
|
|
img(g) = img(h)
|
|
img(h) = temp
|
|
end if
|
|
next j
|
|
next i
|
|
|
|
t = 1
|
|
for stage = 1 to sig
|
|
print " stage:- "; stage
|
|
d = int(real2 / t)
|
|
for ii = 0 to t -1
|
|
l = d *ii
|
|
ls = l +real4
|
|
for i = 0 to d -1
|
|
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
|
|
next i
|
|
next ii
|
|
t = t +t
|
|
next stage
|
|
|
|
print " Num real imag"
|
|
for i = 0 to real
|
|
if abs(rel(i)) <10^-5 then rel(i) = 0
|
|
if abs(img(i)) <10^-5 then img(i) = 0
|
|
print " "; i;" ";using("##.#",rel(i));" ";img(i)
|
|
next i
|
|
end
|