Revision history [back]
 | 1 | initial version |
Yes, using the functions in the scipy.signal module. You're probably better off searching the scipy docs and mailing lists for help if you're going to be doing signal processing.
# modified for Sage from an example by Ivo Maljevic
# http://mail.scipy.org/pipermail/scipy-user/attachments/20090917/d2b7a007/attachment-0001.py
# see thread http://mail.scipy.org/pipermail/scipy-user/2009-September/022565.html
import numpy as np
import scipy.signal
t = np.linspace(0,1,1000) # f_s = 1000 Hz
f_1 = 100 # 100 Hz component
f_2 = 300 # 300 Hz component
f = np.linspace(0, 1000,1000)
v = np.sin(2*np.pi*f_1*t) + np.cos(2*np.pi*f_2*t) # f_1 and f_2 freq components added together
V = np.abs(np.fft.fft(v))**2
p = list_plot(zip(f, V),color="blue",plotjoined=True)
p.show()
(b,a) = scipy.signal.butter(12, 200.0/500.0, btype='low') # filter out anything above 200 Hz
v = scipy.signal.lfilter(b,a,v)
V = np.abs(np.fft.fft(v))**2
p2 = list_plot(zip(f,V),color='red',plotjoined=True)
p2.show()
I'd definitely recommend reading up on the different filters before applying them to any important data. N.B. I didn't verify that there were no minor errors in the above (I always have off-by-one errors in this stuff anyway, so my confirmation wouldn't be worth much..)
