pwtools.signal.FIRFilter

class pwtools.signal.FIRFilter(cutoff, nyq, ntaps=None, ripple=None, width=None, window='hamming', mode='lowpass')[source]

Bases: object

Build and apply a digital FIR filter (low-, high-, band-pass, band-stop). Uses firwin() and in some cases kaiserord().

Doc strings stolen from scipy.signal.

Notes

To plot the frequency response (the frequency bands), use:

>>> f = Filter(...)
>>> plot(f.w, abs(f.h))

Examples

# Example for using a digital filter.
#
# Aliasing
# --------
# Suppose you have a signal with 50 + 80 or 120 Hz and a Nyquist freq of 100.
# The 80 Hz part can be filtered out by using a lowpass with e.g. cutoff=70 or
# a bandpass with cutoff=[10,70] or so.
#
# Because the Nyquist frequency is 100 Hz, the 120 Hz signal is aliased back
# (folded back at 100 Hz) to 80 Hz by sampling the signal and shows then up as
# a peak in FFT. This can be also taken care of by using a filter in time
# domain on the signal before FFT, but make sure that the filter cutoff
# frequencies are such that the aliased peak is excluded (i.e. smaller then 80
# Hz). As such, the aliased 50+120 signal behaves exactly like a 50+80 signal.
#
# Note that in general, you don't know to which frequency aliases have been put
# and just using a bandpass around you desired frequency band won't help. The
# only solution in this case is to avoid aliasing in the first place :)


import numpy as np
from pwtools import mpl
from scipy.signal import hann
from scipy.fftpack import fft
from pwtools.signal import fftsample, FIRFilter, pad_zeros
pi = np.pi
plt = mpl.plt

plots = mpl.prepare_plots(['freq', 'filt_pad', 'filt_nopad'])
nyq = 100 # Hz
df = 1.0  # Hz
dt, nstep = fftsample(nyq, df, mode='f')
t = np.linspace(0, 1, int(nstep))
filt1 = FIRFilter(cutoff=[10,50], nyq=nyq, mode='bandpass', ripple=60,
                  width=10)
filt2 = FIRFilter(cutoff=[10,50], nyq=nyq, mode='bandpass', ntaps=100,
                  window='hamming')
plots['freq'].ax.plot(filt1.w, abs(filt1.h), label='filt1')
plots['freq'].ax.plot(filt2.w, abs(filt2.h), label='filt2')
plots['freq'].ax.legend()

for pad in [True,False]:
    x = np.sin(2*pi*20*t) + np.sin(2*pi*80*t)
    if pad:
        x = pad_zeros(x, nadd=len(x))
        pl = plots['filt_pad']
    else:
        pl = plots['filt_nopad']
    f = np.fft.fftfreq(len(x), dt)
    sl = slice(0, len(x)//2, None)
    win = hann(len(x))
    pl.ax.plot(f[sl], np.abs(fft(x)[sl]), label='fft(x)')
    pl.ax.plot(f[sl], np.abs(fft(filt1(x))[sl]),     label='fft(filt1(x))')
    pl.ax.plot(f[sl], np.abs(fft(filt1(win*x))[sl]), label='fft(filt1(hann*x))')
    pl.ax.plot(f[sl], np.abs(fft(filt2(win*x))[sl]), label='fft(filt2(hann*x))')
    pl.ax.set_title('zero pad = %s' %pad)
    pl.ax.legend()

plt.show()

References

__init__(cutoff, nyq, ntaps=None, ripple=None, width=None, window='hamming', mode='lowpass')[source]
Parameters:

  • cutoff (float or 1D array_like) – Cutoff frequency of filter (expressed in the same units as nyq) OR an array of cutoff frequencies (that is, band edges). In the latter case, the frequencies in cutoff should be positive and monotonically increasing between 0 and nyq. The values 0 and nyq must not be included in cutoff.

  • nyq (float) – Nyquist frequency [Hz]. Each frequency in cutoff must be between 0 and nyq.

  • ntaps (int) – Length of the filter (number of coefficients, i.e. the filter order + 1). ntaps must be even if a passband includes the Nyquist frequency. Use either ntaps or ripple + width for a Kaiser window.

  • ripple (float) – Positive number specifying maximum ripple in passband (dB) and minimum ripple in stopband. Large values (like 1000) remove the “rippling” in the pass band almost completely and make frequency response almost “square” (if width is small) but also lead to a large number of filter coeffs (ntaps).

  • width (float) – Width of transition region (same unit as nyq, e.g. Hz).

  • window (string or tuple of string and parameter values) – Desired window to use. See scipy.signal.get_window for a list of windows and required parameters. Default is “hamming”. Ignored if width and ripple given b/c then kaiserord is used to build a Kaiser window.

  • mode (str) – ‘lowpass’, ‘highpass’, ‘bandpass’, ‘bandstop’

__call__(x, axis=-1)[source]

Apply filter to signal.

Parameters:

  • x (1d array)

  • axis (int)

Methods