projects
/
hooke.git
/ commitdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
| commitdiff |
tree
raw
|
patch
|
inline
| side by side (parent:
224004c
)
Got Parseval's checks working.
author
W. Trevor King
<wking@drexel.edu>
Sat, 4 Oct 2008 13:43:54 +0000
(09:43 -0400)
committer
W. Trevor King
<wking@drexel.edu>
Sat, 4 Oct 2008 13:43:54 +0000
(09:43 -0400)
FFT_tools.py
patch
|
blob
|
history
diff --git
a/FFT_tools.py
b/FFT_tools.py
index 4abb0a7523ad47b3479c0910fd41a4f35d0fe2d5..b6b7b353ecf12dd25c9ed1da8be136ce50d223d8 100644
(file)
--- a/
FFT_tools.py
+++ b/
FFT_tools.py
@@
-54,17
+54,18
@@
def _test_rfft(xs, Xs) :
j = complex(0,1)
n = len(xs)
Xa = []
j = complex(0,1)
n = len(xs)
Xa = []
- for k in range(
len(Xs)
) :
+ for k in range(
n
) :
Xa.append(sum([x*exp(-j*2*pi*k*m/n) for x,m in zip(xs,range(n))]))
Xa.append(sum([x*exp(-j*2*pi*k*m/n) for x,m in zip(xs,range(n))]))
- assert (Xs[k]-Xa[k])/abs(Xa[k]) < 1e-6, \
- "rfft mismatch on element %d: %g != %g, relative error %g" \
- % (k, Xs[k], Xa[k], (Xs[k]-Xa[k])/abs(Xa[k]))
+ if k < len(Xs):
+ assert (Xs[k]-Xa[k])/abs(Xa[k]) < 1e-6, \
+ "rfft mismatch on element %d: %g != %g, relative error %g" \
+ % (k, Xs[k], Xa[k], (Xs[k]-Xa[k])/abs(Xa[k]))
# Which should satisfy the discrete form of Parseval's theorem
# n-1 n-1
# SUM |x_m|^2 = 1/n SUM |X_k|^2.
# m=0 k=0
timeSum = sum([abs(x)**2 for x in xs])
# Which should satisfy the discrete form of Parseval's theorem
# n-1 n-1
# SUM |x_m|^2 = 1/n SUM |X_k|^2.
# m=0 k=0
timeSum = sum([abs(x)**2 for x in xs])
- freqSum = sum([abs(X)**2 for X in Xa])
+ freqSum = sum([abs(X)**2 for X in Xa])
assert abs(freqSum/float(n) - timeSum)/timeSum < 1e-6, \
"Mismatch on Parseval's, %g != 1/%d * %g" % (timeSum, n, freqSum)
assert abs(freqSum/float(n) - timeSum)/timeSum < 1e-6, \
"Mismatch on Parseval's, %g != 1/%d * %g" % (timeSum, n, freqSum)
@@
-118,23
+119,31
@@
def unitary_rfft(data, freq=1.0) :
freq_axis = linspace(0, freq/2, nsamps/2+1)
return (freq_axis, trans)
freq_axis = linspace(0, freq/2, nsamps/2+1)
return (freq_axis, trans)
-def _test_unitary_rfft_parsevals():
+def _test_unitary_rfft_parsevals(
xs, freq, freqs, Xs
):
print "Test unitary rfft on Parseval's theorem"
print "Test unitary rfft on Parseval's theorem"
- xs = [1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1]
- dt = pi
- freqs,Xs = unitary_rfft(xs, 1.0/dt)
# Which should satisfy the discretized integral form of Parseval's theorem
# n-1 n-1
# SUM |x_m|^2 dt = SUM |X_k|^2 df
# m=0 k=0
# Which should satisfy the discretized integral form of Parseval's theorem
# n-1 n-1
# SUM |x_m|^2 dt = SUM |X_k|^2 df
# m=0 k=0
+ dt = 1.0/freq
df = freqs[1]-freqs[0]
assert (df - 1/(len(xs)*dt))/df < 1e-6, \
"Mismatch in spacing, %g != 1/(%d*%g)" % (df, len(xs), dt)
df = freqs[1]-freqs[0]
assert (df - 1/(len(xs)*dt))/df < 1e-6, \
"Mismatch in spacing, %g != 1/(%d*%g)" % (df, len(xs), dt)
+ Xa = list(Xs)
+ for k in range(len(Xs)-1,1,-1) :
+ Xa.append(Xa[k])
+ assert len(xs) == len(Xa), "Length mismatch %d != %d" % (len(xs), len(Xa))
lhs = sum([abs(x)**2 for x in xs]) * dt
lhs = sum([abs(x)**2 for x in xs]) * dt
- rhs = sum([abs(X)**2 for X in X
s
]) * df
- assert abs(lhs - rhs)/lhs < 1e-
6
, "Mismatch on Parseval's, %g != %g" \
+ rhs = sum([abs(X)**2 for X in X
a
]) * df
+ assert abs(lhs - rhs)/lhs < 1e-
4
, "Mismatch on Parseval's, %g != %g" \
% (lhs, rhs)
% (lhs, rhs)
+def _test_unitary_rfft_parsevals_suite():
+ xs = [1,2,3,1,2,3,1,2,3,1,2,3,1,2,3,1]
+ dt = pi
+ freqs,Xs = unitary_rfft(xs, 1.0/dt)
+ _test_unitary_rfft_parsevals(xs, 1.0/dt, freqs, Xs)
+
def _rect(t) :
if abs(t) < 0.5 :
return 1
def _rect(t) :
if abs(t) < 0.5 :
return 1
@@
-556,7
+565,7
@@
def _test_unitary_avg_power_spectrum_sin_suite() :
def test() :
_test_rfft_suite()
def test() :
_test_rfft_suite()
- _test_unitary_rfft_parsevals()
+ _test_unitary_rfft_parsevals
_suite
()
_test_unitary_rfft_rect_suite()
_test_unitary_rfft_gaussian_suite()
_test_unitary_power_spectrum_sin_suite()
_test_unitary_rfft_rect_suite()
_test_unitary_rfft_gaussian_suite()
_test_unitary_power_spectrum_sin_suite()