Discrete Fourier Transform and Inverse Discrete Fourier Transform
To test, it creates an input signal using a Sine wave that has known frequency, amplitude, phase. Later it calculates DFT of the input signal and finds its frequency, amplitude, phase to compare.
Python, 57 lines
Download
Copy to clipboard1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 | # Discrete Fourier Transform (DFT)
# FB - 20141227
import random
import math
import cmath
pi2 = cmath.pi * 2.0
def DFT(fnList):
N = len(fnList)
FmList = []
for m in range(N):
Fm = 0.0
for n in range(N):
Fm += fnList[n] * cmath.exp(- 1j * pi2 * m * n / N)
FmList.append(Fm / N)
return FmList
def InverseDFT(FmList):
N = len(FmList)
fnList = []
for n in range(N):
fn = 0.0
for m in range(N):
fn += FmList[m] * cmath.exp(1j * pi2 * m * n / N)
fnList.append(fn)
return fnList
# TEST
print "Input Sine Wave Signal:"
N = 360 # degrees (Number of samples)
a = float(random.randint(1, 100))
f = float(random.randint(1, 100))
p = float(random.randint(0, 360))
print "frequency = " + str(f)
print "amplitude = " + str(a)
print "phase ang = " + str(p)
print
fnList = []
for n in range(N):
t = float(n) / N * pi2
fn = a * math.sin(f * t + p / 360 * pi2)
fnList.append(fn)
print "DFT Calculation Results:"
FmList = DFT(fnList)
threshold = 0.001
for (i, Fm) in enumerate(FmList):
if abs(Fm) > threshold:
print "frequency = " + str(i)
print "amplitude = " + str(abs(Fm) * 2.0)
p = int(((cmath.phase(Fm) + pi2 + pi2 / 4.0) % pi2) / pi2 * 360 + 0.5)
print "phase ang = " + str(p)
print
### Recreate input signal from DFT results and compare to input signal
##fnList2 = InverseDFT(FmList)
##for n in range(N):
## print fnList[n], fnList2[n].real
|
