Note that there are some explanatory texts on larger screens.

plurals
  1. POApplying time-variant filter in Python
    primarykey
    data
    text
    <p>I'm attempting to apply a bandpass filter with time-varying cutoff frequencies to a signal, using Python. The routine I am currently using partitions my signal into equal-length time segments, then for each segment I apply a filter with time-specific parameters, before merging the signal back together. The parameters are based on pre-existing estimates.</p> <p>The problem I seem to be having is that there are "ripples" at the edge of each time-segment that appear after the filter has been applied. This causes discontinuities in my signal, which interfere with my post-filtering data analysis.</p> <p>I was hoping someone could inform me whether there are any existing routines for implementing filters with time-varying parameters in Python? Alternatively, advice on how I might get around this problem would be much appreciated.</p> <p><strong>EDIT</strong> -- example of what I want to do is added below.</p> <p>Let's say I have a signal x(t). I want to filter the first half with a bandpass filter with parameters (100,200) Hz. The second half I want to filter with parameters (140, 240) Hz. I iterate over x(t), applying my filter to each half, then recombine the results. Some example code might look like:<br/><br/></p> <pre><code>outputArray = np.empty(len(x)) segmentSize = len(x) / 2 filtParams = [(100, 200), (140, 240)] for i in range(2): tempData = x[i*segmentSize:(i+1)*segmentSize] tempFiltered = bandPassFilter(tempData, parameters=filtParams[i]) outputArray[i*segmentSize:(i+1)*segmentSize] = tempFiltered </code></pre> <p>(To save space let's assume I have a function which performs bandpass filtering).</p> <p>As you can see, the data segments do not overlap and are simply "pasted" back together in the new array.</p> <p><strong>EDIT 2</strong> -- some sample code of my problem @H.D.</p> <p>First of all, thanks for your significant input thus far. The audiolazy package looks like a great tool.</p> <p>I thought it would be a bit more useful if I describe my goals in further detail. As I have posted <a href="https://stackoverflow.com/questions/17302240/filter-design-and-frequency-extraction-in-python">elsewhere</a>, I am attempting to extract the <a href="http://www.scholarpedia.org/article/Hilbert-Huang_transform#Instantaneous_frequency_and_the_Hilbert_transform" rel="nofollow noreferrer">instantaneous frequency</a> (IF) of a signal, using the Hilbert transform. My data contains significant noise but I have a good estimate of the bandwidth where my IF signal lies. A problem I have come up against, however, is that the IF is often nonstationary. Using a "static" filter approach I am often therefore required to use a broad bandpass region, to ensure all frequencies are captured.</p> <p>The following code demonstrates the effect of increasing the filter bandwidth on an IF signal. It includes a signal generating function, an implementation of a bandpass filter using the scipy.signal package, and a method to extract the IF of the resultant filtered signal.</p> <pre><code>from audiolazy import * import scipy.signal as sig import numpy as np from pylab import * def sineGenerator( ts, f, rate, noiseLevel=None ): """generate a sine tone with time, frequency, sample rate and noise specified""" fs = np.ones(len(ts)) * f y = np.sin(2*np.pi*fs*ts) if noiseLevel: y = y + np.random.randn(len(y))/float(noiseLevel) return y def bandPassFilter( y, passFreqs, rate, order ): """STATIC bandpass filter using scipy.signal Butterworth filter""" nyquist = rate / 2.0 Wn = np.array([passFreqs[0]/nyquist, passFreqs[1]/nyquist]) z, p, k = sig.butter(order, Wn, btype='bandpass', output='zpk') b, a = sig.zpk2tf(z, p, k) return sig.lfilter(b, a, y) if __name__ == '__main__': rate = 1e4 ts = np.arange(0, 10, 1/rate) # CHANGING THE FILTER AFFECTS THE LEVEL OF NOISE ys = sineGenerator(ts, 600.0, 1e4, noiseLevel=1.0) # a 600Hz signal with noise filts = [[500, 700], [550, 650], [580, 620]] for f in filts: tempFilt = bandPassFilter( ys, f, rate, order=2 ) tempFreq = instantaneousFrequency( tempFilt, rate ) plot( ts[1:], tempFreq, alpha=.7, label=str(f).strip('[]') ) ylim( 500, 750 ) xlabel( 'time' ) ylabel( 'instantaneous frequency (Hz)' ) legend(frameon=False) title('changing filter passband and instantaneous frequency') savefig('changingPassBand.png') </code></pre> <p><img src="https://i.stack.imgur.com/aSTU3.png" alt="changing passband"></p> <p>There is a single frequency component in the signal (at 600Hz). The legend shows the filter parameters used in each case. Using a narrower "static" filter gives a "cleaner" output. But how narrow my filter can be is limited by what the frequencies are. For instance, consider the following signal with two frequency components (one at 600Hz, another at 650Hz).</p> <p><img src="https://i.stack.imgur.com/82Iel.png" alt="varying frequency"></p> <p>In this example I have been forced to use a broader bandpass filter, which has resulted in extra noise creeping in to the IF data.</p> <p><strong>My idea is that by using a time varying filter, I can "optimise" the filter for my signal at certain time increments.</strong> So for the above signal I might want to filter around 580-620Hz for the first 5 seconds, then 630-670Hz for the next 5 seconds. Essentially I want to minimise noise in the final IF signal.</p> <p>Based on the example code you sent I have written a function that uses audiolazy to implement a static Butterworth filter on a signal.</p> <pre><code>def audioLazyFilter( y, rate, Wp, Ws ): """implement a Butterworth filter using audiolazy""" s, Hz = sHz(rate) Wp = Wp * Hz # Bandpass range in rad/sample Ws = Ws * Hz # Bandstop range in rad/sample order, new_wp_divpi = sig.buttord(Wp/np.pi, Ws/np.pi, gpass=dB10(.6), gstop=dB10(.1)) ssfilt = sig.butter(order, new_wp_divpi, btype='bandpass') filt_butter = ZFilter(ssfilt[0].tolist(), ssfilt[1].tolist()) return list(filt_butter(y)) </code></pre> <p>The IF data obtained using this filter in conjunction with the Hilbert transform routine compares well to those obtained using scipy.signal:</p> <pre><code>AL_filtered = audioLazyFilter( ys, rate, np.array([580, 620]), np.array([570, 630]) ) SP_filtered = bandPassFilter( ys, [580, 620], rate, order=2 ) plot(ts[1:], instantaneousFrequency( SP_filtered, 1e4 ), alpha=.75, label='scipy.signal Butterworth filter') plot(ts[1:], instantaneousFrequency( AL_filtered, 1e4 ), 'r', alpha=.75, label='audiolazy Butterworth filter') </code></pre> <p><img src="https://i.stack.imgur.com/JMRzm.png" alt="compare audiolazy with scipy.signal"></p> <p>My question is now can I combine the audiolazy Butterworth routine with the time-varying properties you described in your original posts?</p>
    singulars
    1. This table or related slice is empty.
    plurals
    1. This table or related slice is empty.
    1. This table or related slice is empty.
 

Querying!

 
Guidance

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload