I am trying to write a simple low-pass filter using scipy, but I need help defining the parameters.
I have 3.5 million records in time series data that need to be filtered, and the data is sampled at 1000 Hz.
I am using signal.firwin and signal.lfilter from the scipy library.
The options that I select in the code below do not filter my data at all. Instead, the code below simply creates what graphically looks like the same exact data, with the exception of time phase distortions that shift the graph to the right by a little less than 1000 data points (1 second).
In another program, running the low-pass fir filter using graphical user interface commands produces output that has similar means for each segment of 10 seconds (10,000 data points) but has significantly lower standard deviations, so that we essentially lose noise in this specific data file and replace it with something that maintains an average value, showing long-term trends that are not polluted by higher noise levels. The dialog box for other software parameters contains a check box that allows you to select the number of coefficients so that it is "optimized based on the sample size and sampling frequency." (The mine is 3.5 million samples collected at 1000 Hz, but I need a function that uses these inputs as variables.)
* Can someone show me how to configure the code below so that it removes all frequencies above 0.05 Hz? * I would like to see smooth waves on the graph, and not just temporary distortions the same graph that I get from the code below.
class FilterTheZ0():
def __init__(self,ZSmoothedPylab):
self.n = 1000
self.ZSmoothedPylab=ZSmoothedPylab
self.l = len(ZSmoothedPylab)
self.x = arange(0,self.l)
self.cutoffFreq = 0.05
self.RunLowPassFIR_Filter(self.ZSmoothedPylab, self.n, self.l
, self.x, self.cutoffFreq)
def RunLowPassFIR_Filter(self,data, order, l, x, cutoffFreq):
a = 1
b = signal.firwin(self.n, cutoff = self.cutoffFreq, window = "hamming")
response = signal.lfilter(b,a,data)
responsePylab=p.array(response)
plot(x[10000:20000],data[10000:20000])
plot(x[10000:20000],responsePylab[10000:20000])
show()
return